[Nek5000-users] Two-way coupling and mpi calls to userf()

[nek5000-users at lists.mcs.anl.gov]

Hi Neks,


I trust you are all well,


I'm having some trouble with my two-way LPT code, in particular using information from the usrchk() subroutine (where all my LPTing is done) and userf() (where I want to apply a force to all meshpoints based on the number of particles in the surrounding cell).


My issue is that the code doesn't count the number of particles in each cell properly, even though the xmc, xxc etc. arrays store the correct minimum and maximum values for the cell boundaries as in the usrchk() code I use these to calculate a full volume of the mesh and it is correct.


I feel like my understanding on how the mpi works is lacking, so any advice on that would also be helpful - I understand that subroutine userf (ix,iy,iz,ieg) is called once per local element?


So far my 2W code in the userf() looks like this:


c     !=== TWO WAY COUPLING ===!
      if(TWOWAYON .EQ. 1 .AND. istep .GT. 1) then

      real twowayforcex(lx1,ly1,lz1,lelv)
      real twowayforcey(lx1,ly1,lz1,lelv)
      real twowayforcez(lx1,ly1,lz1,lelv)
      real xmc(lx1,ly1,lz1,lelv)
      real xxc(lx1,ly1,lz1,lelv)
      real ymc(lx1,ly1,lz1,lelv)
      real yxc(lx1,ly1,lz1,lelv)
      real zmc(lx1,ly1,lz1,lelv)
      real zxc(lx1,ly1,lz1,lelv)
      real vol(lx1,ly1,lz1,lelv)
      real partpostw(ldim,lpart)
      real partftw(ldim,lpart)
      integer npartcell(lx1,ly1,lz1,lelv)
      integer TWOWAYON

      integer iel

      common /twowayr/ twowayforcex
     $               , twowayforcey
     $               , twowayforcez
     $               , xmc, xxc
     $               , ymc, yxc
     $               , zmc, zxc, vol
     $               , partpostw, partftw

     common /twowayi/ npartcell
     $              , TWOWAYON

iel = gllel(ieg)

if(istep .GT. 2) then
 do jp = 1,300000
            if(partpostw(1,jp) .GE. xmc(ix,iy,iz,iel) .AND.
     $         partpostw(1,jp) .LE. xxc(ix,iy,iz,iel) .AND.
     $         partpostw(2,jp) .GE. ymc(ix,iy,iz,iel) .AND.
     $         partpostw(2,jp) .LE. yxc(ix,iy,iz,iel) .AND.
     $         partpostw(3,jp) .GE. zmc(ix,iy,iz,iel) .AND.
     $         partpostw(3,jp) .LE. zxc(ix,iy,iz,iel)) then

               npartcell(ix,iy,iz,iel) = npartcell(ix,iy,iz,iel) + 1

               ffx = ffx + partftw(1,jp) ! Still need to average these
               ffy = ffy + partftw(2,jp)
               ffz = ffz + partftw(3,jp)

             endif
   enddo

        else
          ffx = 0.0
          ffz = 0.0
          ffy = 0.0
        endif

      endif

Thankyou for reading this, any help would be greatly appreciated!


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/nek5000-users/attachments/20160513/93ae1444/attachment.html>