[Nek5000-users] About the iobj in the drag calculation
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Thu May 29 19:24:24 CDT 2014
Hi,Tanmoy:
Thank you for your reply. Yes, the x and z are the homogenous directions, so vy_avg and vz_avg should be 0 in this way. I have another question on the output.
I want to output the instaneous velocity field in the x y z vx vy vz format, so that I can do the energy spectrum and two-point correlation in the Matlab(in Nek5000, it's not so convenient to output Ek and two-point correlation directly, in my opinion). So I search the mailist and find a subroutine named 'sem2lex' which can convert us(ix,iy,iz,le) to ul(ix,ie,iy,je,iz,ke), my question is: how can I output this with a different timestep interval?I need more sample point of the instaneous velocity than the means.
Xianbei
At 2014-05-30 01:00:15, nek5000-users at lists.mcs.anl.gov wrote:
Hi Xianbi,
You are right in defining the Reynolds Stresses.
However, if you look carefully, this is an example of fully developed turbulent Channel. Hence the statistics is independant of the x direction (streamwise). Also, the statistics is taken independant of z direction (no wall in z). x and z both have periodic BC. Hence before calculating R_uu, R_vv etc, we just do the spatial planar averaging in the xz plane. Since for homogeneous flow, that should equal an ensemble average.
You can modify the turbChannel to a non-homogeneous flow, in which case you can define Reynolds stresses just like you said, and they would spatially varying.
Thanks,
Tanmoy
On Thu, May 29, 2014 at 5:22 AM, <nek5000-users at lists.mcs.anl.gov> wrote:
Hi,Paul:
I'm confused about part of the subroutine userchk in turbChannel.usr
if(icalld.eq.0) then
call rzero(uavg,n)
call rzero(urms,n)
call rzero(vrms,n)
call rzero(wrms,n)
call rzero(uvms,n)
atime = 0.
timel = time
call planar_average_s(yy ,ym1 ,w1,w2)
icalld = 1
endif
dtime = time - timel
atime = atime + dtime
if (atime.ne.0. .and. dtime.ne.0.) then
beta = dtime/atime
alpha = 1.-beta
ifverbose = .false.
call avg1(uavg,vx ,alpha,beta,n,'uavg',ifverbose)
call avg2(urms,vx ,alpha,beta,n,'urms',ifverbose)
call avg2(vrms,vy ,alpha,beta,n,'vrms',ifverbose)
call avg2(wrms,vz ,alpha,beta,n,'wrms',ifverbose)
call avg3(uvms,vx,vy,alpha,beta,n,'uvmm',ifverbose)
dragx_avg = alpha*dragx_avg + beta*0.5*(dragx(1)+dragx(2))
! averaging over statistical homogeneous directions (r-t)
call planar_average_s(uavg_pl,uavg,w1,w2)
call planar_average_s(urms_pl,urms,w1,w2)
call planar_average_s(vrms_pl,vrms,w1,w2)
call planar_average_s(wrms_pl,wrms,w1,w2)
call planar_average_s(uvms_pl,uvms,w1,w2)
! average over half the channel height
m = ny1*nely
do i=1,ny1*nely/2
uavg_pl(i) = 0.5 * (uavg_pl(i) + uavg_pl(m-i+1))
urms_pl(i) = 0.5 * (urms_pl(i) + urms_pl(m-i+1))
vrms_pl(i) = 0.5 * (vrms_pl(i) + vrms_pl(m-i+1))
wrms_pl(i) = 0.5 * (wrms_pl(i) + wrms_pl(m-i+1))
enddo
endif
tw = dragx_avg/A_w + 1.e-50
u_tau = sqrt(tw/rho)
Re_tau = u_tau*delta/dnu
diff = (tw-ffx_new)/tw
if(nid.eq.0) write(6,1) istep,time,e2,ffx_new,diff
1 format(i6,1p4e17.9,' err2')
! write statistics to file
iostep_avg = param(68)
if(nid.eq.0 .and. istep.gt.0 .and.
& mod(istep,iostep_avg).eq.0) then
write(6,*) 'Dumping statistics ...'
open(unit=56,file='reystresses.dat')
write(56,'(A,1pe14.7)') '#time = ', time
write(56,'(A)')
& '# y y+ R_uu R_vv R_ww R_uv'
open(unit=57,file='means.dat')
write(57,'(A,1pe14.7)') '#time = ', time
write(57,'(A)')
& '# y y+ Umean'
m = ny1*nely/2
do i=1,m
write(56,3) yy(i)+1
& ,(yy(i)+1)*Re_tau
& , (urms_pl(i)-(uavg_pl(i))**2)/u_tau**2
& , vrms_pl(i)/u_tau**2
& , wrms_pl(i)/u_tau**2
& , uvms_pl(i)/u_tau**2
write(57,3) yy(i) + 1.
& , (yy(i)+1.)*Re_tau
& , uavg_pl(i)/u_tau
3 format(1p15e17.9)
enddo
close(56)
close(57)
endif
return
end
I look into the code of navier5.f and find that avg1, avg2 returns the time average of vx and vx**2 in this example, so if we want to calculate the Reynolds stress, the formulation should be :
R_uu=urms-uavg**2
also, if normalized,
R_uu=(urms-uavg**2)/u_tau**2
so are the other components
I think this is the right way instead of that in the example
At 2014-05-29 06:54:10, nek5000-users at lists.mcs.anl.gov
wrote:
>
>Hi Xianbei,
>
>You can define more objects by picking your
>own surface discriminators (which will
>be dependent on your particular application).
>
>If you need more than 4 objects, you can increase maxobj
>in SIZE.
>
>Typically, dragx(0) would be the sum of the drag on all
>the objects --- this is sometimes useful.
>
>If you have (say) a 4-sided channel, you can just define
>one object instead of 4, and set it to be the collection
>of all faces for which cbc(f,e,1)='W '.
>
>Paul
>
>
>
>On Thu, 29 May 2014, nek5000-users at lists.mcs.anl.gov wrote:
>
>> Hi,all:
> I read Paul's reply in this mail https://lists.mcs.anl.gov/mailman/htdig/nek5000-users/2014-May/002756.html
> And got to know something about objects. As you say, the object is a collection of faces, here , I have a question about the iboj, how could I know which face it represent with iboj? As in turbChannel.usr,
> dragx_avg = alpha*dragx_avg + beta*0.5*(dragx(1)+dragx(2))
> I think it means to do the time-average of the average dragx on planey1 and planey2, as
> 0.5*(dragx(1)+dragx(2))
> indicates, while how to calculate the dragx on the planex1, planex2, planez1 and planez2? for the maxobj in this example is only 4 and dragx is defined as
> dragx(0:maxobj)
>
>Thank you
>Xianbei
>_______________________________________________
>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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20140530/e8176541/attachment-0001.html>
More information about the Nek5000-users
mailing list