[Nek5000-users] Stream function
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Thu Jul 12 13:47:10 CDT 2018
Hi Philipp,
Thanks very much for the userchk routine.
I have been testing this with the Blasius example in Nek. I am hoping to implement this with my geometry (a 2-D backward facing step with the Blasius profile at the inlet upstream from the step).
For the Blasius example in Nek, I used VisIt and checked rhs and psi from the stfc0.f* file, but psi does not look physical.
1. What is the purpose of using v2mask? I wonder whether I should specify a different mask...
2. Including declarations from 'INPUT', 'SOLN', and 'MASS' produced a long list of errors (please see below). To run the code, I commented three lines related to including 'INPUT', 'SOLN', and 'MASS'so that these are not included, but I wonder whether that would create problems?
Any advice would be appreciated.
Thanks very much.
Hannah
INPUT:14:16: Error: Symbol ‘param’ at (1) already has basic type of REAL
INPUT:22:27: Error: Symbol ‘param’ at (1) is already in a COMMON block
INPUT:25:20: Error: Symbol ‘matype’ at (1) already has basic type of INTEGER
INPUT:30:28: Error: Symbol ‘matype’ at (1) is already in a COMMON block
INPUT:34:26: Error: Symbol ‘if3d’ at (1) already has basic type of LOGICAL
INPUT:50:26: Error: Symbol ‘if3d’ at (1) is already in a COMMON block
INPUT:66:27: Error: Symbol ‘ifnav’ at (1) already has basic type of LOGICAL
INPUT:69:27: Error: Symbol ‘hcode’ at (1) already has basic type of CHARACTER
INPUT:70:27: Error: Symbol ‘ocode’ at (1) already has basic type of CHARACTER
INPUT:71:27: Error: Symbol ‘drivc’ at (1) already has basic type of CHARACTER
INPUT:72:26: Error: Symbol ‘rstv’ at (1) already has basic type of CHARACTER
INPUT:73:28: Error: Symbol ‘textsw’ at (1) already has basic type of CHARACTER
INPUT:74:27: Error: Symbol ‘initc’ at (1) already has basic type of CHARACTER
INPUT:75:27: Error: Symbol ‘hcode’ at (1) is already in a COMMON block
INPUT:77:29: Error: Symbol ‘turbmod’ at (1) already has basic type of CHARACTER
INPUT:80:28: Error: Symbol ‘reafle’ at (1) already has basic type of CHARACTER
INPUT:81:28: Error: Symbol ‘reafle’ at (1) is already in a COMMON block
INPUT:83:29: Error: Symbol ‘session’ at (1) already has basic type of CHARACTER
INPUT:84:29: Error: Symbol ‘session’ at (1) is already in a COMMON block
INPUT:86:20: Error: Symbol ‘cr_re2’ at (1) already has basic type of INTEGER
INPUT:87:33: Error: Symbol ‘cr_re2’ at (1) is already in a COMMON block
INPUT:89:24: Error: Symbol ‘re2off_b’ at (1) already has basic type of INTEGER
INPUT:90:31: Error: Symbol ‘re2off_b’ at (1) is already in a COMMON block
INPUT:94:13: Error: Symbol ‘xc’ at (1) already has basic type of REAL
INPUT:98:24: Error: Symbol ‘xc’ at (1) is already in a COMMON block
INPUT:100:20: Error: Symbol ‘igroup’ at (1) already has basic type of INTEGER
INPUT:101:28: Error: Symbol ‘igroup’ at (1) is already in a COMMON block
INPUT:103:28: Error: Symbol ‘ccurve’ at (1) already has basic type of CHARACTER
INPUT:104:25: Error: Symbol ‘cbc’ at (1) already has basic type of CHARACTER
INPUT:105:25: Error: Symbol ‘cbc’ at (1) is already in a COMMON block
INPUT:107:19: Error: Symbol ‘ieact’ at (1) already has basic type of INTEGER
INPUT:108:27: Error: Symbol ‘ieact’ at (1) is already in a COMMON block
INPUT:112:20: Error: Symbol ‘numsts’ at (1) already has basic type of INTEGER
INPUT:113:26: Error: PARAMETER attribute of ‘numsts’ conflicts with PARAMETER attribute at (1)
INPUT:115:20: Error: Symbol ‘numflu’ at (1) already has basic type of INTEGER
INPUT:118:29: Error: Symbol ‘numflu’ at (1) is already in a COMMON block
INPUT:121:17: Error: Symbol ‘bcf’ at (1) already has basic type of INTEGER
INPUT:122:26: Error: Symbol ‘bcf’ at (1) is already in a COMMON block
INPUT:124:24: Error: Symbol ‘bctyps’ at (1) already has basic type of CHARACTER
INPUT:125:29: Error: Symbol ‘bctyps’ at (1) is already in a COMMON block
SOLN:4:18: Error: Symbol ‘lvt1’ at (1) already has basic type of INTEGER
SOLN:5:41: Error: PARAMETER attribute of ‘lvt1’ conflicts with PARAMETER attribute at (1)
SOLN:6:41: Error: PARAMETER attribute of ‘lvt2’ conflicts with PARAMETER attribute at (1)
SOLN:7:45: Error: PARAMETER attribute of ‘lbt1’ conflicts with PARAMETER attribute at (1)
SOLN:8:45: Error: PARAMETER attribute of ‘lbt2’ conflicts with PARAMETER attribute at (1)
SOLN:10:42: Error: PARAMETER attribute of ‘lorder2’ conflicts with PARAMETER attribute at (1)
SOLN:14:13: Error: Symbol ‘bq’ at (1) already has basic type of REAL
SOLN:15:22: Error: Symbol ‘bq’ at (1) is already in a COMMON block
SOLN:18:16: Error: Symbol ‘vxlag’ at (1) already has basic type of REAL
SOLN:33:13: Error: Symbol ‘vx’ at (1) already has basic type of REAL
SOLN:52:27: Error: Symbol ‘vxlag’ at (1) is already in a COMMON block
SOLN:59:13: Error: Symbol ‘bx’ at (1) already has basic type of REAL
SOLN:78:13: Error: Symbol ‘bx’ at (1) is already in a COMMON block
SOLN:82:18: Error: Symbol ‘nu_star’ at (1) already has basic type of REAL
SOLN:83:29: Error: Symbol ‘nu_star’ at (1) is already in a COMMON block
SOLN:85:13: Error: Symbol ‘pr’ at (1) already has basic type of REAL
SOLN:86:22: Error: Symbol ‘pr’ at (1) is already in a COMMON block
SOLN:88:14: Error: Symbol ‘qtl’ at (1) already has basic type of REAL
SOLN:89:25: Error: Symbol ‘qtl’ at (1) is already in a COMMON block
SOLN:91:15: Error: Symbol ‘p0th’ at (1) already has basic type of REAL
SOLN:92:27: Error: Symbol ‘p0th’ at (1) is already in a COMMON block
SOLN:94:18: Error: Symbol ‘v1mask’ at (1) already has basic type of REAL
SOLN:106:28: Error: Symbol ‘v1mask’ at (1) is already in a COMMON block
SOLN:111:15: Error: Symbol ‘vxp’ at (1) already has basic type of REAL
SOLN:133:25: Error: Symbol ‘vxp’ at (1) is already in a COMMON block
SOLN:138:16: Error: Symbol ‘jp’ at (1) already has basic type of INTEGER
SOLN:139:25: Error: Symbol ‘jp’ at (1) is already in a COMMON block
MASS:4:14: Error: Symbol ‘bm1’ at (1) already has basic type of REAL
MASS:11:23: Error: Symbol ‘bm1’ at (1) is already in a COMMON block
/home/Hannah/NekExamples/blasius/c3.f:360:25:
call compute_f (f1,w ,t )
1
Warning: Actual argument contains too few elements for dummy argument ‘w’ (1/4) at (1)
makefile:119: recipe for target 'nek5000' failed
make: *** [nek5000] Error 1
________________________________________
From: Nek5000-users [nek5000-users-bounces at lists.mcs.anl.gov] on behalf of nek5000-users at lists.mcs.anl.gov [nek5000-users at lists.mcs.anl.gov]
Sent: Thursday, July 12, 2018 3:37 AM
To: nek5000-users at lists.mcs.anl.gov
Subject: Re: [Nek5000-users] Stream function
Ok, so can you please just try the attached code in userchck? It should
work (at least does for us) in 2D. Essentially, the steps are:
1. copmute z-vorticity in a new array vorticity
2. construct the rhs by multiplying with -bm1
3. setup the Helmholtz problem and tolerances
4. solve the Helmholtz problem and storing the results in a new array phi
5. write out everything in two files.
In this code you do not use the temperature at all, which I guess is
fine for a 2D case (which I guess is what you look at).
best regards,
Philipp
On 2018-07-10 22:33, Philipp Schlatter wrote:
> Dear Hannah,
> sorry that I only see the post now. Are you still having problems with
> the stream function, or does it work now?
> Best regards,
> Philipp
>
> On 2018-07-02 16:40, nek5000-users at lists.mcs.anl.gov wrote:
>> Hi JC,
>>
>> Yes, the mv_wall example computes the streamfunction, but it applies
>> tmask, tmult that I could not find where they are defined in the
>> file. I was experimenting with different masks.
>>
>> Do you know how postnek compute the streamfunction?
>>
>> I think it is computed in the subroutine strfct(psi) in postnek2.f,
>> but I am not sure what happens in the routine?
>>
>> Thanks,
>>
>> Hannah
>> ________________________________________
>> From: Nek5000-users [nek5000-users-bounces at lists.mcs.anl.gov] on
>> behalf of nek5000-users at lists.mcs.anl.gov
>> [nek5000-users at lists.mcs.anl.gov]
>> Sent: Friday, June 29, 2018 3:47 AM
>> To: nek5000-users at lists.mcs.anl.gov
>> Subject: Re: [Nek5000-users] Stream function
>>
>> Hi Hannah,
>>
>> I tried to implement the streamfunction as well (for a 2D cylinder flow)
>> and it would appear that, doing naïvely, the tmask is set to zero
>> everywhere. Hence, when rhs * tmask is computed in the hmholtz
>> subroutine, rhs is set to zero everywhere and it then makes sense that
>> the returned streamfunction is zero everywhere as well.
>>
>> Not sure how to do fix this properly though. Any help would be
>> appreciated indeed :)
>>
>>
>> ++
>> JC
>>
>>
>> On 06/28/2018 09:41 PM, nek5000-users at lists.mcs.anl.gov wrote:
>>> Hi Philipp,
>>>
>>> Thanks very much for your comments. Yes, the delta* works only for
>>> the parallel plate.
>>>
>>> The streamfunction computation with the Blasius example is still not
>>> working. I suspect I am not using the correct mask (vmult /tmult,
>>> imsh as well).
>>>
>>> Using the Blasius example, v1mask, vmult (imsh = 1), I dumped the
>>> streamfuntion to the temperature field. However, it only contains
>>> negative values and does not look physical (I pasted userchk below).
>>>
>>> Another Nek Example, “mv_wall”, computes the streamfunction from the
>>> Possison equation using vorticity. It uses tmask, tmult, and imsh = 2.
>>>
>>> I tried implementing what was done in the “mv_wall” example in the
>>> Blasius example, but the computed streamfunction is zero everywhere.
>>>
>>> Any ideas what I might be doing wrong?
>>>
>>> Any advice would be appreciated!
>>>
>>> Thanks again,
>>>
>>> Hannah
>>>
>>> c-----------------------------------------------------------------------
>>> subroutine userchk
>>> include 'SIZE'
>>> include 'TOTAL'
>>>
>>> common /myblas/ ub(lx1,ly1,lz1,lelt),vb(lx1,ly1,lz1,lelt)
>>>
>>> common /blasiusr/ u(lx1,ly1,lz1,lelv),v(lx1,ly1,lz1,lelv)
>>> common /blasiusg/ bin(lx1,ly1,lz1,lelv)
>>> common /scrns/ dum(lx1*ly1*lz1*lelv)
>>>
>>> common /xtream/ psi(lx1*ly1*lz1*lelt)
>>> $ , rhs(lx1*ly1*lz1*lelt)
>>> $ , h1 (lx1*ly1*lz1*lelt)
>>> $ , h2 (lx1*ly1*lz1*lelt)
>>>
>>> parameter (lt=lx1*ly1*lz1*lelt)
>>> common /scrns/ vort(lt,3),w1(lt),w2(lt)
>>>
>>> integer bin
>>> real normu(50),normv(50)
>>>
>>> real xmax,xmin,dx
>>> save xmax,xmin,dx
>>>
>>> n = nx1*ny1*nz1*nelv
>>> ifheat=.false.
>>> call comp_vort3(vort , w1, w2, vx, vy, vz)
>>> call copy (t,vort,n) ! Vorticity --> T
>>> ifto = .true. ! Dump vorticity as T
>>> c . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
>>> c Compute streamfunction from vorticity
>>> call col3 (rhs,t,bm1,n) ! B*omega
>>> call rone (h1,n)
>>> call rzero(h2,n)
>>>
>>> ifield = 1
>>> tol = param(22)
>>> imsh = 1
>>> call hmholtz('psi ',psi,rhs,h1,h2,v1mask,vmult,imsh,
>>> $ tol,200,1)
>>> call copy (t,psi,n) ! Psi --> T (for output)
>>> c. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
>>>
>>> nfldt = 1 + npscal
>>>
>>> if (istep.eq.0) call outpost2(vx,vy,vz,pr,t,nfldt,' ')
>>>
>>> do i=1,n
>>> dum(i) = 1.-vx(i,1,1,1)
>>> enddo
>>> tot_deficit = glsc2(dum,bm1,n)
>>> xmin = glmin(xm1,n)
>>> xmax = glmax(xm1,n)
>>> delta_star = tot_deficit/(xmax-xmin)
>>> tmin = glmin(t,n)
>>> tmax = glmax(t,n)
>>> vortmin = glmin(vort,n)
>>> vortmax = glmax(vort,n)
>>> if (nid.eq.0) write(6,1) istep,time,delta_star,vortmin,vortmax
>>> if (nid.eq.0) write(6,1) istep,tmin,tmax
>>> 1 format(i8,1p4e13.5,' blasius: delta*')
>>>
>>> return
>>> end
>>> c-----------------------------------------------------------------------
>>>
>>>
>>> ________________________________________
>>> From: Nek5000-users [nek5000-users-bounces at lists.mcs.anl.gov] on
>>> behalf of nek5000-users at lists.mcs.anl.gov
>>> [nek5000-users at lists.mcs.anl.gov]
>>> Sent: Tuesday, June 26, 2018 1:22 PM
>>> To: nek5000-users at lists.mcs.anl.gov
>>> Subject: Re: [Nek5000-users] Stream function
>>>
>>> Hi,
>>> To correctly compute the stream function, the rhs of the Poisson
>>> equation needs to be the (negative) vorticity. Therefore, the line
>>>
>>> call col2(rhs,bm1,n)
>>>
>>> is wrong, and should be something like
>>>
>>> call col3(rhs, vorticity, bm1, n)
>>>
>>> where vorticity should be the computed vorticity using comp_vort.
>>>
>>> I guess your calculation of delta* only works for parallel boundary
>>> layers.
>>>
>>> Philipp
>>>
>>>
>>>
>>> On 2018-06-26 18:35, nek5000-users at lists.mcs.anl.gov wrote:
>>>> Hi Phillip,
>>>>
>>>> I am using the “Blasius” example test outputting the streamfunction
>>>> as a
>>>> scalar.
>>>>
>>>> I included my code below. I am not sure what I am doing wrong.
>>>>
>>>> In userchk, I added “compute streamfunction from vorticity” from the
>>>> psi_omega example, but it does not compute the psi properly… (max_psi =
>>>> min_psi = 0 in every step)
>>>>
>>>> c-----------------------------------------------------------------------
>>>>
>>>>
>>>> subroutine userchk
>>>>
>>>> include 'SIZE'
>>>>
>>>> include 'TOTAL'
>>>>
>>>> common /myblas/ ub(lx1,ly1,lz1,lelt),vb(lx1,ly1,lz1,lelt)
>>>>
>>>> common /blasiusr/ u(lx1,ly1,lz1,lelv),v(lx1,ly1,lz1,lelv)
>>>>
>>>> common /blasiusg/ bin(lx1,ly1,lz1,lelv)
>>>>
>>>> common /scrns/ dum(lx1*ly1*lz1*lelv)
>>>>
>>>> common /xtream/ psi(lx1*ly1*lz1*lelt)
>>>>
>>>> $ , rhs(lx1*ly1*lz1*lelt)
>>>>
>>>> $ , h1 (lx1*ly1*lz1*lelt)
>>>>
>>>> $ , h2 (lx1*ly1*lz1*lelt)
>>>>
>>>> integer bin
>>>>
>>>> real normu(50),normv(50)
>>>>
>>>> real xmax,xmin,dx
>>>>
>>>> save xmax,xmin,dx
>>>>
>>>> n = nx1*ny1*nz1*nelv
>>>>
>>>> c. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
>>>>
>>>> c Compute streamfunction from vorticity
>>>>
>>>> call col2(rhs,bm1,n)
>>>>
>>>> call rone (h1 ,n)
>>>>
>>>> call rzero(h2 ,n)
>>>>
>>>> tol = param(22)
>>>>
>>>> imsh = 1
>>>>
>>>> call hmholtz('psi ',psi,rhs,h1,h2,v1mask,vmult,imsh,tol,200,1)
>>>>
>>>> c. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
>>>>
>>>> if (istep.eq.0) call outpost(vx,vy,vz,pr,psi,' ')
>>>>
>>>> c Compute delta*
>>>>
>>>> do i=1,n
>>>>
>>>> dum(i) = 1.-vx(i,1,1,1)
>>>>
>>>> enddo
>>>>
>>>> tot_deficit = glsc2(dum,bm1,n)
>>>>
>>>> xmin = glmin(xm1,n)
>>>>
>>>> xmax = glmax(xm1,n)
>>>>
>>>> delta_star = tot_deficit/(xmax-xmin)
>>>>
>>>> psimin = glmin(psi,n)
>>>>
>>>> psimax = glmax(psi,n)
>>>>
>>>> if (nid.eq.0) write(6,1) istep,time,delta_star,psimin,psimax
>>>>
>>>> 1 format(i8,1p4e13.5,' blasius: delta*')
>>>>
>>>> return
>>>>
>>>> end
>>>>
>>>>
>>>>
>>>> _______________________________________________
>>>> 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