Re: particle tracking on a DA

26 Sep 2007

      Hi Richard
  I am working on Lagrangian particles as well(using DA objects).
  the way that I am implementing it,is using DA for domain decompositioning.
  but the problem will appears later on.
  problems:

  1-How do you want to store your particles(using link-list or arrays)

  2-how are you going to take care of particle communications?

  right now,by using DA object you can get rid of domain decomposition.
  I did not use the master/slave approach.for my case,I am sending the particles using the point to point communication. if a particle goes out of a sub-domain,it should be send to its neighbor.
  if you want to use this approach,you don't know who are the neighbors using DA object.
  so I wrote simple code to find out who are the neighboring sub-domains(I attached the code).

  I am still working on it, but it would be great if there was partcile object (using link list) in petsc as well.

  Mehdi

Richard Katz <[email protected]> wrote:
  Hi All,   

  Is anyone aware of PETSc particle tracking code?  I'm looking for something that will do Lagrangian particle tracking in parallel for 2D fluid flow on a regular cartesian mesh.

  If this doesn't exist, I may try to create it---quick and dirty style.  My thinking was to a master-slaves parallelism: keep a complete directory of the particle locations on node 0 then send the entire directory to each node after each timestep to do an update.  Each node updates any points that fall into it's subdomain, flags those particle entries, then returns the entire directory to node 0.  Node 0 merges the updates into the master directory.  

  This isn't scalable but it is simpler to program.  Does anyone have any suggestions about strategy, implementation, PETSc tools that would facilitate this, etc?

  Cheers
  Rich

    Richard Foa Katz     NSF IRFP Postdoc
  http://www.damtp.cam.ac.uk/user/rfk22

---------------------------------
Catch up on fall's hot new shows on Yahoo! TV.  Watch previews, get listings, and more!

#include "include/finclude/petscdef.h"

c*********************************************************************

      subroutine findneighbours  
      use sharedvars
      implicit none

      integer :: ierr

      integer :: i

      xsdomains(rank+1)=xs
      call MPI_ALLGATHER(xs,1,MPI_INTEGER, 
     &                   xsdomains        ,1,MPI_INTEGER, 
     &                   MPI_COMM_WORLD,ierr)

      CHKERRQ(ierr) 

      call MPI_ALLGATHER(xe,1,MPI_INTEGER, 
     &                   xedomains        ,1,MPI_INTEGER, 
     &                   MPI_COMM_WORLD,ierr)
      CHKERRQ(ierr) 

      call MPI_ALLGATHER(ys,1,MPI_INTEGER, 
     &                   ysdomains        ,1,MPI_INTEGER, 
     &                   MPI_COMM_WORLD,ierr)
      CHKERRQ(ierr) 

      call MPI_ALLGATHER(ye,1,MPI_INTEGER, 
     &                   yedomains        ,1,MPI_INTEGER, 
     &                   MPI_COMM_WORLD,ierr)
      CHKERRQ(ierr) 

      call MPI_ALLGATHER(zs,1,MPI_INTEGER, 
     &                   zsdomains        ,1,MPI_INTEGER, 
     &                   MPI_COMM_WORLD,ierr)
      CHKERRQ(ierr) 

      call MPI_ALLGATHER(ze,1,MPI_INTEGER, 
     &                   zedomains        ,1,MPI_INTEGER, 
     &                   MPI_COMM_WORLD,ierr)
      CHKERRQ(ierr)

      write(logfile,*)'xs',xsdomains
      write(logfile,*)'xe',xedomains

      write(logfile,*)'ys',ysdomains
      write(logfile,*)'ye',yedomains

      write(logfile,*)'zs',zsdomains
      write(logfile,*)'ze',zedomains

      xleftdomainrank=-10
      yleftdomainrank=-10
      zleftdomainrank=-10

      xrightdomainrank=-10
      yrightdomainrank=-10
      zrightdomainrank=-10

      do i=1,domainsize

       if(xs==1)then

        xleftdomainrank=-10

       endif

       if(ys==1)then

        yleftdomainrank=-10

       endif

       if(zs==1)then

        zleftdomainrank=-10

       endif

       if(xe==nx)then

        xrightdomainrank=-10

       endif

       if(ye==ny)then

        yrightdomainrank=-10

       endif

       if(ze==nz)then

        zrightdomainrank=-10

       endif

      enddo

      do i=1,domainsize

       if((xs-1)==xedomains(i)  .and. 
     &        ys==ysdomains(i)  .and. ye==yedomains(i) .and.
     &        zs==zsdomains(i)  .and. ze==zedomains(i)
     & )then

        xleftdomainrank=i-1

       endif

       if((ys-1)==yedomains(i)  .and. 
     &        xs==xsdomains(i)  .and. xe==xedomains(i) .and.
     &        zs==zsdomains(i)  .and. ze==zedomains(i)
     & )then

        yleftdomainrank=i-1

       endif

       if((zs-1)==zedomains(i)  .and. 
     &        xs==xsdomains(i)  .and. xe==xedomains(i) .and.
     &        ys==ysdomains(i)  .and. ye==yedomains(i)
     & )then

        zleftdomainrank=i-1

       endif

      enddo

      do i=1,domainsize

       if((xe+1)==xsdomains(i) .and.
     &        ys==ysdomains(i) .and. ye==yedomains(i) .and.
     &        zs==zsdomains(i) .and. ze==zedomains(i)
     & )then

        xrightdomainrank=i-1

       endif

       if((ye+1)==ysdomains(i) .and.
     &        xs==xsdomains(i) .and. xe==xedomains(i) .and.
     &        zs==zsdomains(i) .and. ze==zedomains(i)
     & )then

        yrightdomainrank=i-1

       endif

       if((ze+1)==zsdomains(i) .and.
     &        xs==xsdomains(i) .and. xe==xedomains(i) .and.
     &        ys==ysdomains(i) .and. ye==yedomains(i)
     & )then

        zrightdomainrank=i-1

       endif

      enddo

      return
      end 
c*********************************************************************

Re: particle tracking on a DA

Mehdi Bostandoost