[Nek5000-users] Calculation of Divergence Term
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Tue Nov 15 12:39:08 CST 2011
Hi Paul,
If I'm not mistaken, aren't the variable names ux, uy, and uz already
taken in the common blocks in NEKUSE? If so, wouldn't the declaration
of the myvel common block violate those names? Just wanted to make
sure Alireza does not run into any pesky compilation errors that may
be difficult to spot.
Josh
On Nov 15, 2011, at 9:30 AM, "nek5000-users at lists.mcs.anl.gov"
<nek5000-users at lists.mcs.anl.gov> wrote:
>
>
> Alireza,
>
> nek supports a Robin boundary condition - aka Newton law of cooling.
>
> The standard form for this is:
>
> -k dT/dn = h*(T-T_inf)
>
> where k = vdiff(.,.,.,.,ifld) is the conductivity inside the medium
>
> h = heat transfer coefficient
>
> T_inf = temperature outside the domain
>
> As a user, you can trigger this BC by using the "c " (convective)
> boundary condition specification. Lower case "c" indicates that
> the user-supplied quantities are specified in subroutine userbc().
> These quantities are:
>
> T_inf
> h
>
> and passed through the NEKUSE common block. I believe that the
> variable names are "tinf" and "h" ... but am not at a terminal
> where I can verify this at the moment. The conductivity is
> stored in vdiff() and may be variable or constant.
>
> In your case, you need something like the following.
>
>
> subroutine userbc(ix,iy,iz,iside,eg)
> include 'SIZE'
> include 'TOTAL'
> include 'NEKUSE'
>
> integer e,eg
>
> common /myvel/ ux(lx1,ly1,lz1,lelt),uy(lx1,ly1,lz1,lelt)
> $ , uz(lx1,ly1,lz1,lelt)
>
> e = gllel(eg)
>
> uvx = ux(ix,iy,iz,e) + vx(ix,iy,iz,e)
> uvy = uy(....
> uvz = uz ....
>
> if (iside.eq.1.or.iside.eq.3) then
> j1=ix
> j2=iz
> elseif (iside.eq.2.or.iside.eq.4) then
> j1=iy
> j2=iz
> else
> j1=ix
> j2=iy
> endif
>
> c = uvx*unx(j1,j2,iside,e)+uvy*uny(j1,j2,iside,e)+uvz*unz(j1,j2,iside,e) (this line too long)
>
> tinf = 0
> h = c/vdiff(ix,iy,iz,e) ! be sure to check the sign here
>
> return
> end
>
> ...
>
> hth,
>
> Paul
>
>
>
>
>
>
>
>
> ----- Original Message -----
> From: nek5000-users at lists.mcs.anl.gov
> To: nek5000-users at lists.mcs.anl.gov
> Sent: Monday, November 14, 2011 10:38:57 PM
> Subject: Re: [Nek5000-users] Calculation of Divergence Term
>
> Hi Alireza,
>
> If you can cast your advection-diffusion equation for c in terms of new advection vector w = u + V then you get a standard mixed/Robin BC
>
> dc/dn = w.n * c
>
> that can be prescribed with existing Nek5000 infrastructure by specifying the appropriate variables from NEKUSE in userbc such as flux,hc,tinf.
>
> Let me know if you need more details.
> Aleks
>
>
> ----- Original Message -----
> From: nek5000-users at lists.mcs.anl.gov
> To: nek5000-users at lists.mcs.anl.gov
> Sent: Monday, November 14, 2011 8:33:57 PM
> Subject: Re: [Nek5000-users] Calculation of Divergence Term
>
> Hi Neks,
>
> Considering the advection-diffusion equation discussed in my previous
> messages, I want to impose a mixed type of boundary condition for the
> concentration c. In general, on the boundaries, this BC takes the
> following form as,
>
> dc/dn = (u.n + V.n)*c
>
> where n is the normal vector to the boundary, u is the velocity field,
> and V is the given vector field. I wonder if you have any idea how to
> impose this mixed Dirichlet-Neumann BC in the usr file.
>
> Regards,
> Alireza Karimi
>
>
>
>
>
>
>
> Quoting nek5000-users at lists.mcs.anl.gov:
>
>>
>> No problem... sorry for the delay in getting you what you needed.
>>
>>
>>
>> On Sat, 12 Nov 2011, nek5000-users at lists.mcs.anl.gov wrote:
>>
>>> Hi,
>>>
>>> Thanks a lot for your hints and I appreciate your help. Sorry for
>>> interrupting when you were so busy.
>>>
>>> Regards,
>>> Alireza Karimi
>>>
>>>
>>>
>>>
>>> Quoting nek5000-users at lists.mcs.anl.gov:
>>>
>>>>
>>>> Dear Alireza,
>>>>
>>>> You would add this to qvol, since it is the convection-diffusion
>>>> equation -- so, the routine is userq.
>>>>
>>>> You should precompute the term in userchk, store the result in
>>>> a common block, and then dereference in userq. Expanding on
>>>> Josh's suggestion below you would have something like
>>>>
>>>> c-----------------------------------------------------------------------
>>>> subroutine userchk
>>>> :
>>>> :
>>>>
>>>> common /mydiv/ divcv(lx1,ly1,lz1,lelt)
>>>> common /scrns/ dcv1dx(lx1,ly1,lz1,lelt)
>>>> , dcv2dy(lx1,ly1,lz1,lelt)
>>>> , dcv3dz(lx1,ly1,lz1,lelt)
>>>> , work1 (lx1,ly1,lz1,lelt)
>>>> , work2 (lx1,ly1,lz1,lelt)
>>>>
>>>>
>>>> then his calls... be sure to define n
>>>>
>>>> return
>>>> end
>>>> c-----------------------------------------------------------------------
>>>> subroutine userq (ix,iy,iz,eg)
>>>> include 'SIZE'
>>>> include 'TOTAL'
>>>> include 'NEKUSE'
>>>>
>>>> common /mydiv/ divcv(lx1,ly1,lz1,lelt)
>>>>
>>>> integer e,f,eg
>>>> e = gllel(eg)
>>>>
>>>> qvol = divcv(ix,iy,iz,e) Be sure to check the sign of this term
>>>>
>>>> return
>>>> end
>>>> c-----------------------------------------------------------------------
>>>>
>>>>
>>>>
>>>> It's not a priori clear that you need to dealias this term.
>>>> Dealiasing is needed in certain circumstance for stability -- it
>>>> does not materially affect accuracy. Without further context
>>>> I cannot tell you whether you need to dealias or not.
>>>>
>>>> The approach suggested by Josh is perfectly adequate for the
>>>> standard (non-dealiased) approach.
>>>>
>>>>
>>>> Paul
>>>>
>>>>
>>>> On Fri, 11 Nov 2011, nek5000-users at lists.mcs.anl.gov wrote:
>>>>
>>>>> Hi,
>>>>>
>>>>> I'm sending this message again, because I got no answer for my
>>>>> last email. I have some questions regarding the additional term
>>>>> given in my previous message. I don't know where and in which
>>>>> subroutine I should add this term to the right hand side of the
>>>>> advection-diffusion equation. Also, could you please provide me
>>>>> with the best option I can take to dealias this term (if needed
>>>>> at all)? I wonder if you could help me in this regard.
>>>>>
>>>>> Best Regards,
>>>>> Alireza Karimi
>>>>>
>>>>>
>>>>>
>>>>> Quoting nek5000-users at lists.mcs.anl.gov:
>>>>>
>>>>>> Alireza,
>>>>>>
>>>>>> For this, I'll assume that you have split V into V1, V2, V3, and that
>>>>>> they have the structure
>>>>>>
>>>>>> V1(lx1,ly1,lz1,lelv) (and likewise for the other two)
>>>>>>
>>>>>>
>>>>>> First, to multiply c*V1 (V2, V3), you can use col2. It would look like
>>>>>>
>>>>>> call col2(V1,c,n)
>>>>>> call col2(V2,c,n)
>>>>>> call col2(V3,c,n)
>>>>>>
>>>>>> where n = nx1*ny1*nz1*nelv. A couple of things to note here:
>>>>>>
>>>>>> 1) col2 will overwrite V1 (V2, V3). Effectively, the operation is
>>>>>> V1 = c*V1. If this is not desirable, then you could use col3 with
>>>>>> something like
>>>>>>
>>>>>> call col3(cV1,V1,c,n)
>>>>>>
>>>>>> This will store it in cV1 and leave V1 untouched.
>>>>>>
>>>>>> 2) You could run into aliasing issues by multiplying c and V
>>>>>> directly, and so it may be wise to dealias the product before taking
>>>>>> the derivative. Developers, maybe you could suggest the best way of
>>>>>> going about doing this in Nek???
>>>>>>
>>>>>> Moving on, to take the derivatives, there are a few choices. One way
>>>>>> would be do to something like
>>>>>>
>>>>>> call gradm1(dcV1dx,work1 ,work2 ,cV1)
>>>>>> call gradm1(work1 ,dcV2dy,work2 ,cV2)
>>>>>> call gradm1(work1 ,work2 ,dcV3dz,cV3)
>>>>>>
>>>>>> call add4(divcV, dcV1dx, dcV2dy, dcV3dz, n)
>>>>>>
>>>>>> which would give you div(cV) in the array divcV. Note that dcV1dx,
>>>>>> dcV2dy, dcV3dz, work1, work2, and divcV all should have the same
>>>>>> storage structure as V1, c, etc.
>>>>>>
>>>>>> This is probably the simplest way that I can think to do it, although
>>>>>> it is somewhat inefficient since we are computing, for instance,
>>>>>> d(cV1)/dy unnecessarily. To make it more efficient, you can manually
>>>>>> compute the derivatives individually (take a look at the gradm1
>>>>>> subroutine to get an idea of how to do this, located in navier5.f in
>>>>>> the source)
>>>>>>
>>>>>> Hope this helps, and developers, please correct me if I did something
>>>>>> horribly wrong here :).
>>>>>>
>>>>>> Josh
>>>>>>
>>>>>> On Fri, Nov 4, 2011 at 4:59 PM, <nek5000-users at lists.mcs.anl.gov> wrote:
>>>>>>> Hi,
>>>>>>>
>>>>>>> I want to solve an advection-diffusion equation in conjunction with the
>>>>>>> momentum and mass conservation equations. But my equation has an extra
>>>>>>> divergence term which I don't know how to calculate it in Nekton. So,
>>>>>>> denoting the scalar concentration field with c and velocity
>>>>>>> field with u,
>>>>>>> the equation turns out to be:
>>>>>>>
>>>>>>> dc/dt + u.grad c = Laplacian c - div (c V)
>>>>>>>
>>>>>>> which by d/dt I mean partial derivative with respect to time. V
>>>>>>> is a known
>>>>>>> vector field calculated from the vorticity tensor. I wonder if you could
>>>>>>> show me how to evaluate this additional term (the last term) in the user
>>>>>>> file.
>>>>>>>
>>>>>>>
>>>>>>> Best Regards,
>>>>>>> Alireza Karimi
>>>>>>>
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> Nek5000-users mailing list
>>>>>>> Nek5000-users at lists.mcs.anl.gov
>>>>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> Josh Camp
>>>>>>
>>>>>> "All that is necessary for the triumph of evil is that good men do
>>>>>> nothing" -- Edmund Burke
>>>>>> _______________________________________________
>>>>>> Nek5000-users mailing list
>>>>>> Nek5000-users at lists.mcs.anl.gov
>>>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>>>>
>>>>>
>>>>>
>>>>>
>>>>> _______________________________________________
>>>>> Nek5000-users mailing list
>>>>> Nek5000-users at lists.mcs.anl.gov
>>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>> _______________________________________________
>>>> Nek5000-users mailing list
>>>> Nek5000-users at lists.mcs.anl.gov
>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>>
>>>
>>> _______________________________________________
>>> Nek5000-users mailing list
>>> Nek5000-users at lists.mcs.anl.gov
>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>> _______________________________________________
>> Nek5000-users mailing list
>> Nek5000-users at lists.mcs.anl.gov
>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>
>
>
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
More information about the Nek5000-users
mailing list