Spanwise Average

266 views
Skip to first unread message

Dewu Yang

unread,
Aug 31, 2023, 7:24:02 AM8/31/23
to Nek5000
Dear all,
   I'm a beginner in nek5000, and I want to use the planar_average_r() function to achieve averaging along the z-directionI'm using an older version. Can you please guide me on how to do this? Currently, I've obtained some values as output, but the total count of these values doesn't seem to match the total number of two-dimensional plane nodes. I have assigned nelx, nely, and nelz as the numbers of elements in the x-direction, y-direction, and z-direction, respectively (N-1, where N is the number of nodes). My goal is to obtain a three-dimensional flow field and calculate the average velocity along the z-direction (spanwise). Could you please advise me on how to modify the approach? I would greatly appreciate it if you could provide me with some suggestions!

Below is a portion of my .usr code

c-----------------------------------------------------------------------
      subroutine userchk
      INCLUDE 'SIZE'
      INCLUDE 'TOTAL'
      INCLUDE 'ZPER'
c     INCLUDE 'RESTART'
c     Output vorticity
      parameter(lt=lx1*ly1*lz1*lelv)
      common /myexact/ud(lt),vd(lt)
      common /myenrgy/  vx2(lx1,ly1,lz1,lelv),
     &                  vy2(lx1,ly1,lz1,lelv),
     &                  vz2(lx1,ly1,lz1,lelv)
      common /myplane/  uavg(lx1*lelx),
     &                  vavg(lx1*lelx),
     &                  wavg(lx1*lelx),
     &                  w1(lx1*lelx),
     &                  w2(lx1*lelx),
     &                  xx(lx1*lelx),
     &                  yy(lx1*lelx)
      character*80 filename(9999)
      nelx = 144      ! Number of elements in x,y, and z directions.
      nely = 19    ! NOTE, this may vary from one mesh to the next.
      nelz = 9        !
      write(6,*) ' Num :=', lelx, lely, lelz
      ntot   = nx1*ny1*nz1*nelv
      if (nid .eq. 0) then
        open(unit=199,file='file.list',form='formatted',status='old')
        read(199,*) nfiles
        read(199,'(A80)') (filename(i),i=1,nfiles)
        close(199)
      endif
      call bcast(nfiles,isize)
      call bcast(filename,nfiles*80)
      do j = 1, nfiles
         call load_fld(filename(j))
      enddo
      do i=1,ntot
         vx2(i,1,1,1) =  vx(i,1,1,1)**2
         vy2(i,1,1,1) =  vy(i,1,1,1)**2
         vz2(i,1,1,1) =  vz(i,1,1,1)**2
      enddo
      call planar_average_r(uavg,vx,w1,w2)
      call planar_average_r(vavg,vy,w1,w2)
      call planar_average_r(wavg,vz,w1,w2)
      if(nid.eq.0) then
        open(unit=56,file='avgeng.dat')
        m = nx1*nelx
        do i=1,m
                write(56,3) uavg(i),vavg(i)
3           format(1p115e17.9)
        enddo
        close(56)
      endif
      call exitt
      return
      end
c-----------------------------------------------------------------------
  subroutine planar_average_r(ua,u,w1,w2)
c
c     Compute s-t planar average of quantity u()
c
      include 'SIZE'
      include 'GEOM'
      include 'PARALLEL'
      include 'WZ'
      include 'ZPER'
c
      real ua(nx1,nelx),u(nx1,ny1,nx1,nelv),w1(nx1,nelx),w2(nx1,nelx)
      integer e,eg,ex,ey,ez
c
      nx = nx1*nelx
      call rzero(ua,nx)
      call rzero(w1,nx)
c
      do e=1,nelt
c
         eg = lglel(e)
         call get_exyz(ex,ey,ez,eg,nelx,nely,nelz)
c
         do i=1,nx1
         do k=1,nz1
         do j=1,ny1
            zz = (1.-zgm1(i,1))/2.  ! = 1 for i=1, = 0 for k=nx1
            aa = zz*area(j,k,4,e) + (1-zz)*area(j,k,2,e)  ! wgtd jacobian
            w1(i,ex) = w1(i,ex) + aa
            ua(i,ex) = ua(i,ex) + aa*u(i,j,k,e)
         enddo
         enddo
         enddo
      enddo
c
      call gop(ua,w2,'+  ',nx)
      call gop(w1,w2,'+  ',nx)
c
      do i=1,nx
         ua(i,1) = ua(i,1) / w1(i,1)   ! Normalize
      enddo

      return
      end
c-----------------------------------------------------------------------
Best,
Dewu
Reply all
Reply to author
Forward
0 new messages