<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=Windows-1252">
</head>
<body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space;">
Hi Paul,
<div><br>
</div>
<div>Thank you for your quick answer, this is really helpful.</div>
<div><br>
</div>
<div>I still have a question concerning the computation of a gradient.</div>
<div><br>
</div>
<div>Why did you construct the diffusion operator with the helmholtz routine? Why not using directly the routines gradm1 ?</div>
<div><br>
</div>
<div><br>
</div>
<div>By the way I am not sure to understand the differences between gradm1 (navier5.f), wgradm1 (navier1.f) , gradm11 (navier2.f) and gradm11ts (navier2.f).</div>
<div>I guess that the w letter refers to the weak formulation, but I don’t understand in which cases we need to choose a formulation instead of others.</div>
<div><br>
</div>
<div>In summary, if I want to modify the temperature/passive scalar equations to add a term that contains a gradient of the solved temperature/PS variable, for example, which “gradm” routine do I need to use?</div>
<div><br>
</div>
<div>Thank you for your help !</div>
<div><br>
</div>
<div><br>
<div>
<div style="color: rgb(0, 0, 0); letter-spacing: normal; orphans: auto; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; widows: auto; word-spacing: 0px; -webkit-text-stroke-width: 0px; word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space;">
<div>Best regards,</div>
<div>Emmanuel</div>
</div>
<br class="Apple-interchange-newline">
<br class="Apple-interchange-newline">
</div>
<br>
<div>
<div>On 29 May 2014, at 20:20, <a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a> wrote:</div>
<br class="Apple-interchange-newline">
<blockquote type="cite"><br>
<br>
Hi Emmanuel,<br>
<br>
axhelm effects the weak Laplacian plus the mass matrix.<br>
<br>
Specifically, if you wanted to solve:<br>
<br>
(1) -nabla . [ h1 * nabla ] u + h2 * u = f<br>
<br>
with h1, h2, and f all functions of X=(x,y,z),<br>
the discrete form would be:<br>
<br>
<br>
(2) A(h1) u + B(h2) u = B f<br>
<br>
or<br>
<br>
(3) H u = B f<br>
<br>
where A is the (h1-dependent stiffness matrix and<br>
B(h2) is the h2-dependent mass matrix.<br>
<br>
For the SEM, B is diagonal, so B(h2) = h2*B, where<br>
we view h2 as a diagonal matrix.<br>
<br>
In Nek5000, a matrix-vector product w = Hu is effected via<br>
<br>
call axhelm (w,u,h1,h2,imesh,isd)<br>
<br>
with w,u,h1,h2 being arrarys. Here, imesh indicates<br>
if you are on the fluid mesh or thermal mesh. The fluid<br>
mesh is a subset of the thermal mesh in the case of conjugate heat transfer between a fluid and a solid.<br>
isd indicates which spatial direction is considered<br>
in axisymmetric applications where u is a component<br>
of a vector field.<br>
<br>
Notice the sign of the Laplacian in Eq. (1) -- that may<br>
help resolve your second question? (There is a sign<br>
change that ariss from the standard integration-by-parts<br>
stepin the weak form.)<br>
<br>
Paul<br>
<br>
<br>
On Thu, 29 May 2014, <a href="mailto:nek5000-users@lists.mcs.anl.gov">nek5000-users@lists.mcs.anl.gov</a> wrote:<br>
<br>
<blockquote type="cite">Dear nek5000 team,<br>
<br>
Could you help me with the 2 questions below? Thank you very much !<br>
<br>
1) Concerning the routine “axhelm” in hmholtz.f. If I am right, one<br>
possibility of this routine can be the discretisation of the diffusion<br>
operator. However it is not clear to me if axhelm performs:<br>
<br>
a) vdiff x Laplacian(phi)<br>
<br>
or<br>
<br>
b) grad ( vdiff x grad (phi) )<br>
<br>
phi being the variable solved for (T or PS for example).<br>
<br>
So, is it actually a) or b) ?<br>
<br>
<br>
2) In the routine wlaplacian in navier1.f, I am not sure to understand the<br>
consequence of the line:<br>
<br>
call sub2 (out,wrk,ntot)<br>
<br>
after the lines:<br>
<br>
call bcneusc(out,1) call axhelm(wrk,a,diff,h2,imesh,1)<br>
<br>
If I am right, this will change the sign in front of the diffusion operator,<br>
right? Because in makeq.f, I am not sure that the line call<br>
add2(bq(1,1,1,1,ifield-1),w1,ntot)<br>
<br>
will actually assemble the RHS of the transport equation : - (vtrans) x U<br>
grad (phi) + grad ( vdiff x grad (phi) )<br>
<br>
<br>
Thank you for your help !<br>
<br>
Best regards, Emmanuel<br>
<br>
</blockquote>
_______________________________________________<br>
Nek5000-users mailing list<br>
<a href="mailto:Nek5000-users@lists.mcs.anl.gov">Nek5000-users@lists.mcs.anl.gov</a><br>
https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users<br>
</blockquote>
</div>
<br>
</div>
</body>
</html>