Adjustments to Poisson treatment

STEllAR-GROUP · Feb 22, 2013 · 58b0a27 · 58b0a27
1 parent 75c928e
commit 58b0a27
Show file tree

Hide file tree

Showing 2 changed files with 66 additions and 26 deletions.
diff --git a/examples/gtcx/gtcx_hpx/server/loop.F90 b/examples/gtcx/gtcx_hpx/server/loop.F90
@@ -910,7 +910,8 @@ subroutine poisson(iflag,hpx4_bti,&
               partd_comm,nproc_partd,myrank_partd,&
               toroidal_comm,nproc_toroidal,myrank_toroidal,&
               left_pe,right_pe,&
-              toroidal_domain_location,particle_domain_location&
+              toroidal_domain_location,particle_domain_location,&
+              nindex,indexp,ring,mindex& 
                   )
   !use global_parameters
   !use field_array
@@ -953,7 +954,13 @@ subroutine poisson(iflag,hpx4_bti,&
   integer  :: left_pe,right_pe
   integer  :: toroidal_domain_location,particle_domain_location
 
-  integer iflag,i,it,ij,j,k,n,iteration,mring,mindex,mtest,ierr
+  integer iflag
+
+  integer  :: mindex
+  integer,dimension(mgrid,mzeta) :: nindex
+  integer,dimension(mindex,mgrid,mzeta) :: indexp
+  real(kind=wp),dimension(mindex,mgrid,mzeta) :: ring
+
     end subroutine poisson
     subroutine pushe(icycle,irke,hpx4_bti,&
              t_gids, p_gids,&
@@ -1304,8 +1311,16 @@ end subroutine chargee
   real(kind=wp),dimension(:,:,:),allocatable :: ptracked
   integer,dimension(:),allocatable :: ntrackp
 
+! Poisson thread safety
+  integer,dimension(:,:),allocatable :: nindex
+  integer,dimension(:,:,:),allocatable :: indexp
+  real(kind=wp),dimension(:,:,:),allocatable :: ring
+  integer  :: mindex
+
+
 ! local variables
   integer i
+  integer mring,mtest,nmem
 
   mype = hpx4_mype
   numberpe = hpx4_numberpe
@@ -1470,6 +1485,24 @@ end subroutine chargee
               left_pe,right_pe,&
               toroidal_domain_location,particle_domain_location) ! }}}
 
+  ! Poisson thread safety {{{
+  ! For thread safety in poisson, allocate this now
+  ! number of gyro-ring
+  mring=2
+
+  ! number of summation: maximum is 32*mring+1
+  mindex=32*mring+1
+
+  allocate(indexp(mindex,mgrid,mzeta),ring(mindex,mgrid,mzeta),&
+          nindex(mgrid,mzeta),STAT=mtest)
+  if (mtest /= 0) then
+!   hjw
+    nmem = (2*(mindex*mgrid*mzeta))+(mgrid*mzeta)
+!   write(0,*)mype,'*** Cannot allocate indexp: mtest=',mtest
+    write(0,*)mype,'*** indexp: Allocate Error: ',nmem, ' words mtest= ',mtest
+  endif
+  ! }}}
+
 ! main time loop
   do istep=1,mstep
      if ( mype .eq. 0 ) print*,' step ', istep
@@ -1815,7 +1848,6 @@ end subroutine chargee
                  ) ! }}}
 
 ! solve GK Poisson equation using adiabatic electron
-#if 0
          call poisson(0,hpx4_bti,& ! {{{
              t_gids, p_gids,&
 ! global parameters
@@ -1845,9 +1877,9 @@ end subroutine chargee
               partd_comm,nproc_partd,myrank_partd,&
               toroidal_comm,nproc_toroidal,myrank_toroidal,&
               left_pe,right_pe,&
-              toroidal_domain_location,particle_domain_location&
+              toroidal_domain_location,particle_domain_location,&
+              nindex,indexp,ring,mindex&
                   )   ! }}}
-#endif
 
         do ihybrid=1,nhybrid
 ! smooth potential
@@ -2136,7 +2168,6 @@ end subroutine chargee
                  ) ! }}}
 
 ! solve GK Poisson equation using non-adiabatic electron
-#if 0
          call poisson(1,hpx4_bti,& ! {{{
              t_gids, p_gids,&
 ! global parameters
@@ -2166,9 +2197,9 @@ end subroutine chargee
               partd_comm,nproc_partd,myrank_partd,&
               toroidal_comm,nproc_toroidal,myrank_toroidal,&
               left_pe,right_pe,&
-              toroidal_domain_location,particle_domain_location&
+              toroidal_domain_location,particle_domain_location,&
+              nindex,indexp,ring,mindex&
                   )   ! }}}
-#endif
         enddo
      enddo
 

diff --git a/examples/gtcx/gtcx_hpx/server/poisson.F90 b/examples/gtcx/gtcx_hpx/server/poisson.F90
@@ -27,7 +27,8 @@ subroutine poisson(iflag,hpx4_bti,&
               partd_comm,nproc_partd,myrank_partd,&
               toroidal_comm,nproc_toroidal,myrank_toroidal,&
               left_pe,right_pe,&
-              toroidal_domain_location,particle_domain_location&
+              toroidal_domain_location,particle_domain_location,&
+              nindex,indexp,ring,mindex&
                   )
   !use global_parameters
   !use field_array
@@ -130,10 +131,18 @@ end subroutine poisson_initial
   integer  :: left_pe,right_pe
   integer  :: toroidal_domain_location,particle_domain_location
 
-  integer iflag,i,it,ij,j,k,n,iteration,mring,mindex,mtest,ierr
-  integer,dimension(:,:),allocatable :: nindex
-  integer,dimension(:,:,:),allocatable :: indexp
-  real(kind=wp),dimension(:,:,:),allocatable :: ring
+! Poisson thread safety
+  integer  :: mindex
+  integer,dimension(mgrid,mzeta) :: nindex
+  integer,dimension(mindex,mgrid,mzeta) :: indexp
+  real(kind=wp),dimension(mindex,mgrid,mzeta) :: ring
+
+
+! Local variables
+  integer iflag,i,it,ij,j,k,n,iteration,mring,mtest,ierr
+!  integer,dimension(:,:),allocatable :: nindex
+!  integer,dimension(:,:,:),allocatable :: indexp
+!  real(kind=wp),dimension(:,:,:),allocatable :: ring
   real(kind=wp) gamma,tmp,prms,perr(mgrid)
   real(kind=wp) ptilde(mgrid),phitmp(mgrid),dentmp(mgrid)
 
@@ -143,13 +152,14 @@ end subroutine poisson_initial
   integer :: nmem
 ! hjw
 
-  save nindex,indexp,ring
+  integer :: master_flag
+!  save nindex,indexp,ring
 
 ! number of gyro-ring
   mring=2
 
 ! number of summation: maximum is 32*mring+1
-  mindex=32*mring+1
+  !mindex=32*mring+1
 !  mindex=53
 
 ! gamma=0.75: max. resolution for k=0.577
@@ -158,17 +168,16 @@ end subroutine poisson_initial
 
 ! initialize poisson solver
   if(istep==1 .and. irk==1 .and. iflag==0)then
-     allocate(indexp(mindex,mgrid,mzeta),ring(mindex,mgrid,mzeta),&
-          nindex(mgrid,mzeta),STAT=mtest)
-     if (mtest /= 0) then
-! hjw
-     nmem = (2*(mindex*mgrid*mzeta))+(mgrid*mzeta)
-!     write(0,*)mype,'*** Cannot allocate indexp: mtest=',mtest
-      write(0,*)mype,'*** indexp: Allocate Error: ',nmem, ' words mtest= ',mtest
-! hjw
-        !call MPI_ABORT(MPI_COMM_WORLD,1,ierr)
-     endif
-
+!     allocate(indexp(mindex,mgrid,mzeta),ring(mindex,mgrid,mzeta),&
+!          nindex(mgrid,mzeta),STAT=mtest)
+!     if (mtest /= 0) then
+!! hjw
+!     nmem = (2*(mindex*mgrid*mzeta))+(mgrid*mzeta)
+!!     write(0,*)mype,'*** Cannot allocate indexp: mtest=',mtest
+!      write(0,*)mype,'*** indexp: Allocate Error: ',nmem, ' words mtest= ',mtest
+!! hjw
+!        !call MPI_ABORT(MPI_COMM_WORLD,1,ierr)
+!     endif
 
 ! initialize
      call poisson_initial(mring,mindex,nindex,indexp,ring,&