[Nek5000-users] Calculation of Divergence Term
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Fri Nov 18 10:26:04 CST 2011
Hi,
Thanks for your time. In the system which I'm simulating, the base
flow is zero and the flow field is induced by the gradient of the
concentration. So, Reynolds and Peclet numbers cannot be well defined.
But their order of magnitude based on the average value of velocity
is: Re~0.5 and Pe~10. So, the flow is not turbulent. Also, I'm using
3rd order time integration.
I'm not sure about the resolving of the boundary layers. My domain is
a box with the size of 1x1x1 and I have 2 elements in each direction
(8 elements at all) with the degree of polynomial equal to 11. I
wonder if I have to use finer grid with/without nonuniform spacing
close to the boundaries.
Regards,
Alireza Karimi
Quoting nek5000-users at lists.mcs.anl.gov:
>
> Hi Alireza,
>
> What Peclet and Reynolds numbers are you considering? How well
> resolved are your boundary layers? Is the background flow turbulent?
> What order of timestepping are you using (as set by TORDER)?
> The concentration transport is not evaluated in conservation form,
> so drift is possible, particularly if you are not fully resolved.
>
> I see no reason for fixing a boundary value unless you are
> interested in a steady state solution, in which case the solution
> is indeterminant up to a constant. For such situations in nek we
> typically fix the mean to be zero. This should not be an issue,
> however, for a time-marching approach to steady state.
>
> Let me know if you'd like me to dig deeper into the issue you are
> having.
>
> Paul
>
>
>
> On Thu, 17 Nov 2011, nek5000-users at lists.mcs.anl.gov wrote:
>
>> Hi,
>>
>> Thanks a lot for your helpful hints and comments. I have another
>> problem because of boundary conditions. The special type of BC,
>> given in my previous message, is equivalent to no flux for the
>> concentration c and it must applied on all boundaries. So, the
>> value of the integral of the concentration over the whole volume
>> should be conserved and does not change with time. But after some
>> time, c reaches very high values which causes the Helmholtz solver
>> to fail. In some numerical methods, when we deal with the no flux
>> BC, a single node on the boundary should be excluded and be
>> assigned a specific value. I'm not sure if it is the case here. I
>> wonder if you have any idea how to implement this BC and keep the
>> whole concentration constant.
>>
>> Best Regards,
>> Alireza Karimi
>>
>>
>>
>>
>> Quoting nek5000-users at lists.mcs.anl.gov:
>>
>>>
>>> Yes... sorry... you have to change the names.
>>>
>>> Thanks,
>>>
>>> Paul
>>>
>>> PS - also, noted Alek's earlier email - I think "hc" is the ht
>>> trnsfer coef.
>>>
>>>
>>>
>>> On Tue, 15 Nov 2011, nek5000-users at lists.mcs.anl.gov wrote:
>>>
>>>> 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
>>>> _______________________________________________
>>>> 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
>>>
>>
>>
>>
>> --
>> Alireza Karimi
>> PhD Candidate
>> Department of Engineering Science and Mechanics
>> Virginia Polytechnic Institute and State University
>> Blacksburg, VA 24061
>> _______________________________________________
>> 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
>
--
Alireza Karimi
PhD Candidate
Department of Engineering Science and Mechanics
Virginia Polytechnic Institute and State University
Blacksburg, VA 24061
More information about the Nek5000-users
mailing list