July 2018 - Nek5000-users - mailman3.cels.anl.gov

Planar average in z-direction problem
by nek5000-users＠lists.mcs.anl.gov 24 Jul '18

24 Jul '18

Hi, NeksI wonder why planar_average_z in navier5.f shows wrong result.here I attached 2 images.First image shows the z-directional grids points from zm1(I wrote a code to extract z coordinate at y = 0 as follows:======================================================insnum = nidinsuni = insnum + 1000WRITE(insnam,98) 'inst_field_', insnumOPEN(unit=insuni,file=insnam insnum = nid insuni = insnum + 1000 WRITE(insnam,98) 'inst_field_', insnum Open(unit=insuni,file=insnam) do i = 1,n if (ym1(i,1,1,1).eq.0.0)then WRITE(insuni,5) xm1(i,1,1,1) % , zm1(i,1,1,1) % , vx(i,1,1,1) % , vy(i,1,1,1) % , vz(i,1,1,1) % , phig(i,1,1,1) endif enddo CLOSE(insuni)======================================================Then, I can collect lots of field data from each processor and sort data with ascending order using matlab.So, I got the result as 1st image.But when I use plannar_average_zI got the result as shown in 2nd image.1. I really want to modify the problems, so could you give me your idea on this problem?2. Additionally, could you explain 'zgm1' and 'area'?I noticed that, planar_average_r or s ... in turbChannel case has also 'area' but, those case are used to calculate ondifferent direction.But they used same input parameters.What I understand about area is, surface area of element.But I thought that different direction needs different indexing======================================================= subroutine planar_average_z(ua,u,w1,w2)cc Compute r-s planar average of quantity u()c include 'SIZE' include 'GEOM' include 'PARALLEL' include 'WZ' include 'ZPER'c real ua(lz1,nelz),u(lx1*ly1,lz1,nelv),w1(lz1,nelz),w2(lz1,nelz) integer e,eg,ezc melxy = nelx*nelyc nz = lz1*nelz call rzero(ua,nz) call rzero(w1,nz)c do e=1,neltc eg = lglel(e)c ez = 1 + (eg-1)/melxy call get_exyz(ex,ey,ez,eg,nelx,nely,nelz)c do k=1,lz1 do i=1,lx1*ly1 zz = (1.-zgm1(k,3))/2. ! = 1 for k=1, = 0 for k=lz1 aa = zz*area(i,1,5,e) + (1-zz)*area(i,1,6,e) ! wgtd jacobian w1(k,ez) = w1(k,ez) + aa ua(k,ez) = ua(k,ez) + aa*u(i,k,e) enddo enddo enddoc call gop(ua,w2,'+ ',nz) call gop(w1,w2,'+ ',nz)c do i=1,nz ua(i,1) = ua(i,1) / w1(i,1) ! Normalize enddoc return end======================================================Regards,Hyunduk. ------------------------------<wbr style="color: rgb(34, 34, 34); font-family: arial, sans-serif; font-size: 14px; background-color: rgb(255, 255, 255);">------------------------------------------------------------------------------------ HyundukSeo(Mr.) IntegratedPh.D.course,Schoo<wbr style="color: rgb(34, 34, 34); font-family: arial, sans-serif; font-size: 14px; background-color: rgb(255, 255, 255);">PusanNationalUniversity ExperimentalThermo-FluidsMec<wbr style="color: rgb(34, 34, 34); font-family: arial, sans-serif; font-size: 14px; background-color: rgb(255, 255, 255);">ENsys)Laboratory http://fluid.pusan.ac.kr T.:+82-51-510-1536M.:+82-10<wbr style="color: rgb(34, 34, 34); font-family: arial, sans-serif; font-size: 14px; background-color: rgb(255, 255, 255);">-2964-0994

1 0

Re: [Nek5000-users] Using an external file for forcing in 'ffx'
by nek5000-users＠lists.mcs.anl.gov 23 Jul '18

23 Jul '18

Hey Philipp, Thanks this helps a lot! Couple of questions, 1) So you are saying the interpolation can be done in userdat2 at the start? My forcing field does not change with time. 2) So essentially you are saying the routine in red should be just, ffx = new_array(ix,iy,iz,iel) ? Does 'iel' by default designate the local index? Or do I have to declare 'iel' before? One of the problems I was having doing it this way is that I didn't get a flowfield when I tried to use more than one processor, that is, parallelizing it. However it returned me a flowfield when it was with just one processor. That is in series. I am uploading some files if you want to take a further look. https://drive.google.com/drive/folders/1i_TUfUePTHK2He9hGyRoFcC_1emeJhXU?us… Thanks, Saikat Saikat Mukherjee, PhD Student, Paul Research Group - http://www.me.vt.edu/mpaul/ Engineering Science and Mechanics, Virginia Tech. On Mon, Jul 23, 2018 at 5:56 PM, Saikat Mukherjee <msaikat6(a)vt.edu> wrote: > Hey Philipp, > > Thanks this helps a lot! Couple of questions, > > 1) So you are saying the interpolation can be done in userdat2 at the > start? My forcing field does not change with time. > > 2) So essentially you are saying the routine in red should be just, ffx = > new_array(ix,iy,iz,iel) ? Does 'iel' by default designate the local index? > Or do I have to declare 'iel' before? > > One of the problems I was having doing it this way is that I didn't get a > flowfield when I tried to use more than one processor, that is, > parallelizing it. However it returned me a flowfield when it was with just > one processor. That is in series. > > I am uploading some files if you want to take a further look. > > > Thanks, > Saikat > > > > On Mon, Jul 23, 2018 at 3:54 PM, <nek5000-users(a)lists.mcs.anl.gov> wrote: > >> Hi, >> - userdat2 looks ok. Just some comments: use real instead of double >> precision as we compile with -r8. Then you could define your parameter >> right in the beginning and just use that one, instead of hardcoded 651. But >> these are details, it should work, and just read and broadcast your data to >> all processors. >> >> - in userf: again, use real instead. Also, you might read the FD >> coordinates in userdat2 as well. You can do as you do, but that is most >> inefficient (reading every time step from every processor). >> >> - you do not need to reshape the array, you can just define my_array with >> the proper size in userf. (my_array(31,21)). In fact, you do not even need >> that as you can simply pass my_array to pwl_intepr_2d. >> >> - the whole interpolation can be done once in the beginning anyway, or do >> you change the forcing field as a function of time? >> >> - you do not need separate fo and new_array (see above) >> >> - then the red bit: You have now interpolated your forcing on the gll >> points. So ffx simply needs to be the entry of this forcing field at the >> current ix,iy,iz,e. No need for a loop over the (local) elements, and no >> need for local-to-global mapping. just take the local index iel. Everything >> is local, as your interpolation was local as well (i.e. only interpolated >> to those elements that a specific processor had because you used xm1 and >> ym1. >> >> Sorry if I am misunderstanding. It is not easy to debug code without a >> running case. >> >> Hope this helps, >> philipp >> >> On 2018-07-23 16:57, nek5000-users(a)lists.mcs.anl.gov wrote: >> >>> Hello, >>> >>> Following up on this question again. Would appreciate any input. >>> >>> Thanks, >>> Saikat >>> >>> >>> >>> On Mon, Jul 16, 2018 at 10:17 AM, <nek5000-users(a)lists.mcs.anl.gov >>> <mailto:[email protected]>> wrote: >>> >>> Hey everyone, >>> >>> Just following up on this question. Could really use some suggestion. >>> >>> Thanks, >>> Saikat >>> >>> Saikat Mukherjee, >>> PhD Student, >>> Paul Research Group - http://www.me.vt.edu/mpaul/ >>> Engineering Science and Mechanics, >>> Virginia Tech. >>> >>> On Wed, Jul 11, 2018 at 11:34 AM, Saikat Mukherjee <msaikat6(a)vt.edu >>> <mailto:[email protected]>> wrote: >>> >>> Hey All, >>> >>> I have been stuck in this problem for a while now but I think it >>> should be simple and probably you can tell where I am going >>> wrong. >>> >>> - So I have a forcing function on a finite difference X-Y grid, >>> which I am reading in as a text file in the 'userdat2' >>> subroutine. >>> >>> - I am then parallelizing it using a previous suggestion by >>> Paul, here is a snippet of the routine. >>> >>> ------------------------------------------------------------ >>> -------------------------------------------- >>> >>> subroutine usrdat2 >>> >>> >>> include 'SIZE' >>> >>> include 'TOTAL' >>> >>> double precision :: my_array(651) >>> >>> integer m >>> >>> common /myarray/ my_array >>> >>> >>> parameter (my_array_size=651) >>> >>> >>> real work(my_array_size) >>> >>> >>> >>> n = my_array_size >>> >>> >>> call rzero(my_array,n) >>> >>> call rzero(work,n) >>> >>> >>> if (nid.eq.0) then >>> >>> open (unit=49,file='fint41.txt') >>> >>> do k = 1,n >>> >>> read(49,*) my_array(k) >>> >>> enddo >>> >>> close(49) >>> >>> endif >>> >>> >>> call gop(my_array,work,'+ ',n) ! Sum over all >>> processors >>> >>> >>> return >>> >>> end >>> >>> ------------------------------------------------------------ >>> --------------------------------------------- >>> >>> - Next I am interpolating this forcing function to the GLL grid >>> points in the userf routine, by a piecewise linear interpolant, >>> here's how, >>> >>> ------------------------------------------------------------ >>> ---------------------------------------------- >>> >>> subroutine userf (ix,iy,iz,ieg) >>> >>> >>> include 'SIZE' >>> >>> include 'TOTAL' >>> >>> include 'NEKUSE' >>> >>> EQUIVALENCE (prandtl,param(2)) >>> >>> >>> double precision :: my_array(651) >>> >>> double precision :: shaped_array(31,21) >>> >>> real new_array(nx1,ny1,nz1,nelv) >>> >>> double precision :: xin(21),yin(31) >>> >>> double precision :: xn(nx1*ny1*nelv),yn(nx1*ny1*nelv) >>> >>> double precision :: fo(nx1*ny1*nelv) >>> >>> integer i,j,k,l,p,q >>> >>> integer e,eg >>> >>> common /myarray/ my_array >>> >>> >>> open(16,file='xin.txt') >>> >>> do p=1,21 >>> >>> read(16,*) xin(p) >>> >>> end do >>> >>> close(16) >>> >>> >>> open(12,file='yin.txt') >>> >>> do q=1,31 >>> >>> read(12,*) yin(q) >>> >>> end do >>> >>> close(12) >>> >>> >>> shaped_array = reshape(my_array,(/31,21/)) >>> >>> >>> do i=1,nx1*ny1*nelv >>> >>> xn(i) = xm1(i,1,1,1) >>> >>> yn(i) = ym1(i,1,1,1) >>> >>> end do >>> >>> call >>> pwl_interp_2d(31,21,yin,xin,shaped_array,nx1*ny1*nelv,yn,xn,fo) >>> >>> >>> new_array=reshape(fo,(/nx1,ny1,nz1,nelv/)) >>> >>> >>> do e=1,nelv >>> >>> eg = lglel(e) >>> >>> ffx = new_array(ix,iy,iz,eg) >>> >>> end do >>> >>> >>> ffy = 0.0 >>> >>> ffz = 0.0 >>> >>> >>> return >>> >>> end >>> >>> >>> >>> ------------------------------------------------------------ >>> ------------------------------------------------------------ >>> --------------- >>> >>> This is a test domain I made, with 6 elements. It seems to be >>> working when I am using only 1 core, i.e., in series. In >>> parallel, it seems to be giving me a different flowfield. My >>> actual domain has 252 elements which I intend to parallelize >>> with 48 cores. So would like to ask you how to proceed. I think >>> I am close but I am slightly confused about part of the code >>> that I have made red. >>> >>> When I am trying this in my actual domain, I get an error, >>> "Segmentation fault occurred." I think the key here is >>> parallelization and identifying the elements correctly which I >>> am not sure I am doing right. Thanks for reading this lengthy >>> email. Will be eager to hear from you. >>> >>> Thanks, >>> Saikat >>> >>> >>> >>> >>> >>> >>> >>> _______________________________________________ >>> Nek5000-users mailing list >>> Nek5000-users(a)lists.mcs.anl.gov <mailto:[email protected] >>> s.anl.gov> >>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users >>> <https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users> >>> >>> >>> >>> >>> _______________________________________________ >>> Nek5000-users mailing list >>> Nek5000-users(a)lists.mcs.anl.gov >>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users >>> >>> _______________________________________________ >> Nek5000-users mailing list >> Nek5000-users(a)lists.mcs.anl.gov >> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users >> > >

1 0

Re: [Nek5000-users] Using an external file for forcing in 'ffx'
by nek5000-users＠lists.mcs.anl.gov 23 Jul '18

23 Jul '18

Hey everyone, Just following up on this question. Could really use some suggestion. Thanks, Saikat Saikat Mukherjee, PhD Student, Paul Research Group - http://www.me.vt.edu/mpaul/ Engineering Science and Mechanics, Virginia Tech. On Wed, Jul 11, 2018 at 11:34 AM, Saikat Mukherjee <msaikat6(a)vt.edu> wrote: > Hey All, > > I have been stuck in this problem for a while now but I think it should be > simple and probably you can tell where I am going wrong. > > - So I have a forcing function on a finite difference X-Y grid, which I am > reading in as a text file in the 'userdat2' subroutine. > > - I am then parallelizing it using a previous suggestion by Paul, here is > a snippet of the routine. > > ------------------------------------------------------------ > -------------------------------------------- > > subroutine usrdat2 > > > include 'SIZE' > > include 'TOTAL' > > double precision :: my_array(651) > > integer m > > common /myarray/ my_array > > > parameter (my_array_size=651) > > > real work(my_array_size) > > > > n = my_array_size > > > call rzero(my_array,n) > > call rzero(work,n) > > > if (nid.eq.0) then > > open (unit=49,file='fint41.txt') > > do k = 1,n > > read(49,*) my_array(k) > > enddo > > close(49) > > endif > > > call gop(my_array,work,'+ ',n) ! Sum over all processors > > > return > > end > ------------------------------------------------------------ > --------------------------------------------- > > - Next I am interpolating this forcing function to the GLL grid points in > the userf routine, by a piecewise linear interpolant, here's how, > > ------------------------------------------------------------ > ---------------------------------------------- > > subroutine userf (ix,iy,iz,ieg) > > > include 'SIZE' > > include 'TOTAL' > > include 'NEKUSE' > > EQUIVALENCE (prandtl,param(2)) > > > double precision :: my_array(651) > > double precision :: shaped_array(31,21) > > real new_array(nx1,ny1,nz1,nelv) > > double precision :: xin(21),yin(31) > > double precision :: xn(nx1*ny1*nelv),yn(nx1*ny1*nelv) > > double precision :: fo(nx1*ny1*nelv) > > integer i,j,k,l,p,q > > integer e,eg > > common /myarray/ my_array > > > open(16,file='xin.txt') > > do p=1,21 > > read(16,*) xin(p) > > end do > > close(16) > > > open(12,file='yin.txt') > > do q=1,31 > > read(12,*) yin(q) > > end do > > close(12) > > > shaped_array = reshape(my_array,(/31,21/)) > > > do i=1,nx1*ny1*nelv > > xn(i) = xm1(i,1,1,1) > > yn(i) = ym1(i,1,1,1) > > end do > > call pwl_interp_2d(31,21,yin,xin,shaped_array,nx1*ny1*nelv,yn,xn, > fo) > > new_array=reshape(fo,(/nx1,ny1,nz1,nelv/)) > > > do e=1,nelv > > eg = lglel(e) > > ffx = new_array(ix,iy,iz,eg) > > end do > > > ffy = 0.0 > > ffz = 0.0 > > > return > > end > > > ------------------------------------------------------------ > ------------------------------------------------------------ > --------------- > > This is a test domain I made, with 6 elements. It seems to be working when > I am using only 1 core, i.e., in series. In parallel, it seems to be giving > me a different flowfield. My actual domain has 252 elements which I intend > to parallelize with 48 cores. So would like to ask you how to proceed. I > think I am close but I am slightly confused about part of the code that I > have made red. > > When I am trying this in my actual domain, I get an error, "Segmentation > fault occurred." I think the key here is parallelization and identifying > the elements correctly which I am not sure I am doing right. Thanks for > reading this lengthy email. Will be eager to hear from you. > > Thanks, > Saikat > > > > > >

1 2

Dirichlet-Neumann BC
by nek5000-users＠lists.mcs.anl.gov 23 Jul '18

23 Jul '18

Dear Neks, I have been looking into the different boundary conditions available by default in bdry.f. However, I haven't seen any mixed conditions. I would like to impose a no slip condition and specify the shear stress at the boundary. The boundary is flat and aligned with the mesh. Moreover it seems that shear stress boundary conditions are applicable only to the Pn-Pn-2 formulation. Is there a way I can specify the wall normal streamwise velocity gradient u_y as well ? Best regards, Armand. -- *Armand Shams* Master Thesis Intern Département aérodynamique, aéroélasticité, acoustique MSAT Tél: +33 6 99 51 15 84 ONERA - The French Aerospace Lab - Centre de Meudon 8, rue des Vertugadins - 92190 MEUDON Nous suivre sur : www.onera.fr <https://www.onera.fr> | Twitter <http://www.twitter.com/@onera_fr> | LinkedIn <http://www.linkedin.com/company/onera> | Facebook <http://www.facebook.fr/thefrenchaerospacelab> Avertissement/disclaimer https://www.onera.fr/en/emails-terms

1 2

Run time error
by nek5000-users＠lists.mcs.anl.gov 23 Jul '18

23 Jul '18

Hi Neks, I am facing an error while running the simulation which I never had before. after 5001 time step the CFL number disappears and it gives an error stops running. (It is not CFL blowing. I mean if CFL blows it will show something like CFL : ******* right? In my case it is not like this. The CFL dis appears it is not even showing the CFL letters. ) Also it saves my output files unto there. Any one had the same problem before? Thank you Sijo

1 4

mesh converter ICEM CFD to NEK5000
by nek5000-users＠lists.mcs.anl.gov 22 Jul '18

22 Jul '18

Hi Neks, I am trying to convert 2D ICEM CFD mesh to NEK5000. Could you tell me, which is the best way to convert grid? My geometry is a 2D cylinder and 2D ellipse with a rectangular domain for the simulation. -- Alok Mishra Research Scholar Aerospace Department IIT Kanpur

1 5

Mesh deformation 3D mathematical doubt
by nek5000-users＠lists.mcs.anl.gov 17 Jul '18

17 Jul '18

Hi Neks, I had made a question before regarding the mesh deformation using Laplace equation here I continue with the same question but in different way. Since I am using Laplace equation (Elliptic equation) to make the blending function, it is necessary to bound the equation with boundary conditions right? But in my case the right walls and left walls should also move according to the velocity of the object. So I think that, the Dirichlet BCs are applied only on the inlet, up, down and outlet. And in homogenous boundary condition should be applied on the right and left walls. It is just a thought from my mind. If I am right, can I apply a mixed boundary condition on this domain to the Laplace equation from the available source code ? Thanks in advance Sijo GEORGE

1 0

Mesh Deformation in 3D
by nek5000-users＠lists.mcs.anl.gov 17 Jul '18

17 Jul '18

Hi Neks, I have come up with a new problem regarding the mesh deformation. I was trying to make a rigid motion near the walls of a 3D cylinder . I took the advantage of ocyl2.usrfile to make that. But I am getting a jacobian error when I am runnning the simulation. I slightly changed the boundary condition part in the ocyl2.usr nface = 2*ndim do e=1,nelv do f=1,nface c Set Dirichlet for mesh velocity on all non-interior boundaries if ( (cbc(f,e,1).ne.'E ' . or. (cbc(f,e,1).ne.'P ' ) ) $ call facev(msk,e,f,z0,nx1,ny1,nz1) c Set inhomogeneous Dirichlet data on cylinder if (cbc(f,e,1).eq.'mv ') then call fcaver(xavg,xm1,e,f) if (xavg.gt.-1.2) call facev(tmp,e,f,m1,nx1,ny1,nz1) endif enddo enddo Since I use a 3D domain the right and left walls are Periodic bcs. So I also added one more condition (in bold letters above) to make the smoothing function. I feel it is correct but it throwing error once it starts running. I made the 3D geometry by using n2to3 on ocyl.rea file and gave PeriodicBoundary conditions. Could you please tell me how to solve this in 3D Thank you Sijo GEORGE

1 0

Re: [Nek5000-users] LES Application to Natural Circulation Problem
by nek5000-users＠lists.mcs.anl.gov 17 Jul '18

17 Jul '18

Nek5000 implements two regularization methods (filtering = explicit or hpfRT) which you can apply in a DNS _and_ LES context. In the LES context you can think of it like a SGS model stabilizing the under-resolved computation. There is a turbChannel example you can experiment with. Note, we found that classical Smagorinsky cannot challenge the simple use of regularization while filtered/VMS Smagorinsky versions (not avaiable yet) perform similiar. Cheers, Stefan -----Original message----- > From:[email protected] <nek5000-users(a)lists.mcs.anl.gov> > Sent: Tuesday 17th July 2018 0:50 > To: nek5000-users(a)lists.mcs.anl.gov > Subject: [Nek5000-users] LES Application to Natural Circulation Problem > > Good afternoon Nek5000, > Where can I find more information or an example on Neks LES model and how to properly make use of its capabilities. Am I correct in understanding the method used for solving the > periodic hill example uses DNS? > I have modified the periodic hill example to make my natural circulation problem and want to use of the LES dynamic Smagorinsky model to reduce the simulations computational cost. > Specifically my questions are, which parameters are needed and what guidelines are to be used when inputting these parameters? > extrapolation = OIFS? > filtering = none > filterWeight = ? > filterCutoffRatio = 0.9 vs 0.65 > ? > ? > ? > Thanks, > Mario > _______________________________________________ > Nek5000-users mailing list > Nek5000-users(a)lists.mcs.anl.gov > https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users

1 0

Turbulent kinetic energy in TurbInflow example
by nek5000-users＠lists.mcs.anl.gov 16 Jul '18

16 Jul '18

Dear Neks, In TurbInflow example there is a usr file for postprocessing where TKE is calculated. It says: ccc avg all the avgturbInflow0.f0* files call auto_averager(fname1) ccc copy uk=vx**2, vk=vy**2, wk=vz**2 call col3(uk,vx,vx,lt) call col3(vk,vy,vy,lt) call col3(wk,vz,vz,lt) ccc avg all the rmsturbInflow0.f0* files call auto_averager(fname2) ccc u'^2,v'^2,w'^2 = rms - avg**2 call sub2(vx,uk,lt) call sub2(vy,vk,lt) call sub2(vz,wk,lt) ccc Average in thex direction call my_x_avg(uk,vx,gs_avg_hndl,nelx,nely,nelz,ifld,rcmask) call my_x_avg(vk,vy,gs_avg_hndl,nelx,nely,nelz,ifld,rcmask) call my_x_avg(wk,vz,gs_avg_hndl,nelx,nely,nelz,ifld,rcmask) ccc TKE = 0.5*(u'^2+v'^2+w'^2) call add4(tk,uk,vk,wk,lt) call cmult(tk,0.5,lt) ccc sqrt uk,vk,wk to get urms, vrms and wrms call vsqrt(uk,lt) call vsqrt(vk,lt) call vsqrt(wk,lt) ccc output urms,vrms and wrms along with the TKE ccc TKE is in temperature field ifto = .true. call outpost(uk,vk,wk,pr,tk,'avx') As far as I know, the sub2 subroutine works like sub2(A,B,lt) A=A+B That means, it saves the sum in the first argument. The code says: call sub2(vx,uk,lt) call sub2(vy,vk,lt) call sub2(vz,wk,lt) Meaning that u'^2,v'^2,w'^2 is being stored in vx, vy, and vz. So should the code say:.. ? call add4(tk,vx,vy,vz,lt) instead of call add4(tk,uk,vk,wk,lt) Or there is something I'm missing? Thank you all. Juan Pablo.

1 0