Changeset 1228ef7 in flexpart.git
- Timestamp:
- Aug 2, 2021, 12:11:23 PM (3 years ago)
- Branches:
- dev
- Parents:
- 6cb0801
- Location:
- src
- Files:
-
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
src/makefile
r71f2128 r1228ef7 90 90 ## OPTIMIZATION LEVEL 91 91 O_LEV = 2 # [0,1,2,3,g,s,fast] 92 O_LEV_DBG = g# [0,g]92 O_LEV_DBG = 0 # [0,g] 93 93 94 94 ## LIBRARIES … … 100 100 101 101 LDFLAGS = $(FFLAGS) -L$(LIBPATH1) -Wl,-rpath,$(LIBPATH1) $(LIBS) #-L$(LIBPATH2) 102 LDDEBUG = $(DBGFLAGS) -L$(LIBPATH1) $(LIBS) #-L$(LIBPATH2)102 LDDEBUG = $(DBGFLAGS) -L$(LIBPATH1) -Wl,-rpath,$(LIBPATH1) $(LIBS) #-L$(LIBPATH2) 103 103 104 104 MODOBJS = \ -
src/mpi_mod.f90
ra803521 r1228ef7 178 178 implicit none 179 179 180 integer :: i,j,s ,addmaxpart=0180 integer :: i,j,s 181 181 182 182 ! Each process gets an ID (mp_pid) in the range 0,..,mp_np-1 … … 486 486 487 487 integer,intent(in) :: num_part 488 integer :: i ,jj, addone488 integer :: i 489 489 490 490 ! Exit if too many particles … … 589 589 integer, intent(in) :: itime 590 590 real :: pmin,z 591 integer :: i,jj,nn, num_part=1,m,imin, num_trans 592 logical :: first_iter 591 integer :: i,jj,nn,m,imin, num_trans 593 592 integer,allocatable,dimension(:) :: idx_arr 594 593 real,allocatable,dimension(:) :: sorted ! TODO: we don't really need this … … 613 612 end do 614 613 615 ! Do not rebalance particles for ipout=3 616 if ( ipout.eq.3) return614 ! Do not rebalance particles for ipout=3 or mquasilag=1 615 if ((ipout.eq.3).or.(mquasilag.eq.1)) return 617 616 618 617 ! For each successive element in index array, see if a lower value exists -
src/partoutput_short_mpi.f90
r92fab65 r1228ef7 13 13 ! * 14 14 ! 12/2014 eso: Version for MPI * 15 ! P rocesses sequentially access and append data to file*16 ! NB: Do not use yet!*15 ! Particle positions are sent to root process for output * 16 ! * 17 17 !***************************************************************************** 18 18 ! * … … 28 28 29 29 real(kind=dp) :: jul 30 integer :: itime,i,j,jjjjmmdd,ihmmss,numshortout,numshortall 30 integer, dimension(:), allocatable :: numshorts,displs 31 integer :: itime,i,j,jjjjmmdd,ihmmss,numshortout,numshortall,numshortmpi 31 32 integer :: ix,jy,ixp,jyp 32 33 real :: xlon,ylat,zlim,dt1,dt2,dtt,ddx,ddy,rddx,rddy,p1,p2,p3,p4,topo … … 35 36 integer(kind=2) :: idump(3,maxpart) 36 37 integer :: i4dump(maxpart) 37 character(LEN=8) :: file_stat='OLD' 38 integer(kind=2),dimension(:,:),allocatable :: idump_all(:,:) 39 integer,dimension(:), allocatable :: i4dump_all(:) 40 ! character(LEN=8) :: file_stat='OLD' 41 character(LEN=8) :: file_stat='REPLACE' 38 42 39 ! MPI root process creates the file, other processes append to it 40 if (lroot) file_stat='REPLACE' 43 ! This is not needed, in this version only root process writes the file 44 ! if (lroot) then 45 ! file_stat='REPLACE' 46 ! end if 47 48 ! Array to gather numshortout from all processes 49 allocate(numshorts(mp_partgroup_np), displs(mp_partgroup_np)) 41 50 42 51 ! Determine current calendar date, needed for the file name … … 115 124 i4dump(numshortout)=npoint(i) 116 125 endif 117 118 126 endif 119 127 end do 120 128 121 129 130 131 ! Get total number of particles from all processes 132 !************************************************ 133 call MPI_Allgather(numshortout, 1, MPI_INTEGER, numshorts, 1, MPI_INTEGER, & 134 mp_comm_used, mp_ierr) 135 136 numshortmpi = sum(numshorts(:)) 137 138 139 ! Gather all data at root process 140 !******************************** 141 allocate(idump_all(3,numshortmpi), i4dump_all(numshortmpi)) 142 displs(1)=0 143 do i=2,mp_partgroup_np 144 displs(i)=displs(i-1)+numshorts(i-1) 145 end do 146 147 call MPI_gatherv(i4dump, numshortout, MPI_INTEGER, i4dump_all, numshorts(:), & 148 & displs, MPI_INTEGER, id_root, mp_comm_used, mp_ierr) 149 displs = displs*3 150 call MPI_gatherv(idump, 3*numshortout, MPI_INTEGER2, idump_all, 3*numshorts(:), & 151 & displs, MPI_INTEGER2, id_root, mp_comm_used, mp_ierr) 152 122 153 ! Open output file and write the output 123 154 !************************************** 124 155 125 open(unitshortpart,file=path(2)(1:length(2))//'shortposit_'//adate// & 126 atime,form='unformatted',status=file_stat,position='append') 156 if (lroot) then ! MPI root process only 157 open(unitshortpart,file=path(2)(1:length(2))//'shortposit_'//adate// & 158 atime,form='unformatted',status=file_stat,position='append') 159 write(unitshortpart) itime 160 write(unitshortpart) numshortmpi 161 write(unitshortpart) & 162 (i4dump_all(i),(idump_all(j,i),j=1,3),i=1,numshortmpi) 163 close(unitshortpart) 164 end if 127 165 128 ! Write current time to file 129 !*************************** 130 131 if (lroot) write(unitshortpart) itime ! MPI root process only 132 ! :TODO: get total numshortout (MPI reduction), add MPI barrier, open file 133 ! sequentially below 134 write(unitshortpart) numshortout 135 write(unitshortpart) & 136 (i4dump(i),(idump(j,i),j=1,3),i=1,numshortout) 137 138 139 write(*,*) numshortout,numshortall 140 141 close(unitshortpart) 166 deallocate(idump_all, i4dump_all) 142 167 143 168 end subroutine partoutput_short -
src/releaseparticles_mpi.f90
r92fab65 r1228ef7 51 51 real(kind=dp) :: juldate,julmonday,jul,jullocal,juldiff 52 52 real,parameter :: eps=nxmax/3.e5,eps2=1.e-6 53 integer :: mind2 53 integer :: mind2, numpartcount_mpi 54 54 ! mind2 eso: pointer to 2nd windfield in memory 55 56 55 integer :: idummy = -7 57 56 !save idummy,xmasssave … … 60 59 logical :: first_call=.true. 61 60 62 ! Use different seed for each process.63 !******************************************* *********************************61 ! Use different random seed for each process 62 !******************************************* 64 63 if (first_call) then 65 64 idummy=idummy+mp_seed … … 68 67 69 68 mind2=memind(2) 69 70 ! For mquasilag=1, assign unique particle ID across processes 71 numpartcount_mpi=mp_partid-mp_partgroup_np+1 70 72 71 73 ! Determine the actual date and time in Greenwich (i.e., UTC + correction for daylight savings time) … … 210 212 nclassunc) 211 213 numparticlecount=numparticlecount+1 214 ! Use a stride equal to number of processes for the MPI version 215 numpartcount_mpi=numpartcount_mpi+mp_partgroup_np 212 216 if (mquasilag.eq.0) then 213 217 npoint(ipart)=i 214 218 else 215 npoint(ipart)=numpart iclecount219 npoint(ipart)=numpartcount_mpi 216 220 endif 217 221 idt(ipart)=mintime ! first time step -
src/timemanager_mpi.f90
rb1f28c3 r1228ef7 448 448 endif 449 449 450 ! :TODO: MPI output of particle positions; each process sequentially451 ! access the same file452 450 if ((mquasilag.eq.1).and.(itime.eq.(loutstart+loutend)/2)) & 453 call partoutput_short(itime) 451 call partoutput_short(itime) ! dump particle positions in extremely compressed format 454 452 455 453 … … 557 555 if ((iout.eq.4).or.(iout.eq.5)) call plumetraj(itime) 558 556 if ((iflux.eq.1).and.(lnetcdfout.eq.0)) call fluxoutput(itime) 557 #ifdef USE_NCF 559 558 if ((iflux.eq.1).and.(lnetcdfout.eq.1)) call fluxoutput_netcdf(itime) 559 #endif 560 560 if (mp_measure_time) call mpif_mtime('iotime',1) 561 561
Note: See TracChangeset
for help on using the changeset viewer.