source: flexpart.git/src/readpartpositions_mpi.f90 @ 32b49c3

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 32b49c3 was 32b49c3, checked in by Espen Sollum ATMOS <eso@…>, 8 years ago

Parallel version can now save/restart simulations with IPOUT/IPIN

  • Property mode set to 100644
File size: 8.4 KB
Line 
1!**********************************************************************
2! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010         *
3! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
4! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
5!                                                                     *
6! This file is part of FLEXPART.                                      *
7!                                                                     *
8! FLEXPART is free software: you can redistribute it and/or modify    *
9! it under the terms of the GNU General Public License as published by*
10! the Free Software Foundation, either version 3 of the License, or   *
11! (at your option) any later version.                                 *
12!                                                                     *
13! FLEXPART is distributed in the hope that it will be useful,         *
14! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
15! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
16! GNU General Public License for more details.                        *
17!                                                                     *
18! You should have received a copy of the GNU General Public License   *
19! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
20!**********************************************************************
21
22subroutine readpartpositions
23
24  !*****************************************************************************
25  !                                                                            *
26  !   This routine opens the particle dump file and reads all the particle     *
27  !   positions from a previous run to initialize the current run.             *
28  !                                                                            *
29  !                                                                            *
30  !     Author: A. Stohl                                                       *
31  !                                                                            *
32  !     24 March 2000                                                          *
33  !                                                                            *
34  !   CHANGES                                                                  *
35  !     12/2014 eso: MPI version                                               *
36  !                  Root process reads positions and distributes the data     *
37  !                  :TODO: The above..                                        *
38  !*****************************************************************************
39  !                                                                            *
40  ! Variables:                                                                 *
41  !                                                                            *
42  !*****************************************************************************
43
44  use par_mod
45  use com_mod
46  use random_mod, only: ran1
47  use mpi_mod !, only: mp_seed
48
49  implicit none
50
51  integer :: ibdatein,ibtimein,nspecin,itimein,numpointin,i,j,ix,ip
52  integer :: id1,id2,it1,it2
53  integer :: addone,numparticlecount_all,numpart_all,lbnd,ubnd
54  real :: xlonin,ylatin,topo,hmixi,pvi,qvi,rhoi,tri,tti
55  character :: specin*7
56  real(kind=dp) :: julin,julpartin,juldate
57
58  integer :: idummy = -8
59
60  ! These variables are allocated at the root process for all particles in file.
61  real,dimension(maxpart) :: xtra1_all,ytra1_all,ztra1_all
62  real,dimension(maxpart,maxspec) :: xmass1_all
63  integer,dimension(maxpart) :: npoint_all,itramem_all
64
65  ! Different seed for each process
66  idummy=idummy+mp_seed
67
68  numparticlecount=0
69  numparticlecount_all=0
70  numpart_all=0
71
72  ! Open header file of dumped particle data
73  ! Each MPI process sequentially access file (just in case)
74  !*********************************************************
75
76  do ip=0, mp_partgroup_np-1
77    call mpif_mpi_barrier
78    open(unitpartin,file=path(2)(1:length(2))//'header', &
79         form='unformatted',err=998)
80
81    read(unitpartin) ibdatein,ibtimein
82    read(unitpartin)
83    read(unitpartin)
84
85    read(unitpartin)
86    read(unitpartin)
87    read(unitpartin) nspecin
88    nspecin=nspecin/3
89    if ((ldirect.eq.1).and.(nspec.ne.nspecin)) goto 997
90
91    do i=1,nspecin
92      read(unitpartin)
93      read(unitpartin)
94      read(unitpartin) j,specin
95      if ((ldirect.eq.1).and.(species(i)(1:7).ne.specin)) goto 996
96    end do
97   
98    read(unitpartin) numpointin
99    if (numpointin.ne.numpoint) then ! goto 995
100
101! eso 2016: moved this warning here to avoid out-of-block goto
102!995   write(*,*) ' #### FLEXPART MODEL WARNING IN READPARTPOSITIONS#### '
103      write(*,*) ' #### FLEXPART MODEL WARNING IN READPARTPOSITIONS#### '
104      write(*,*) ' #### NUMBER OF RELEASE LOCATIONS DOES NOT     #### '
105      write(*,*) ' #### AGREE WITH CURRENT SETTINGS!             #### '
106!  stop
107      goto 999
108    end if
109
110999 continue
111    do i=1,numpointin
112      read(unitpartin)
113      read(unitpartin)
114      read(unitpartin)
115      read(unitpartin)
116      do j=1,nspec
117        read(unitpartin)
118        read(unitpartin)
119        read(unitpartin)
120      end do
121    end do
122    read(unitpartin)
123    read(unitpartin)
124
125    do ix=0,numxgrid-1
126      read(unitpartin)
127    end do
128    close(unitpartin)
129
130
131  ! Open data file of dumped particle data
132  ! All processes read the whole file
133  !***************************************
134
135  open(unitpartin,file=path(2)(1:length(2))//'partposit_end', &
136         form='unformatted',err=998)
137   
138100 read(unitpartin,end=99) itimein
139    i=0
140200 i=i+1
141    read(unitpartin) npoint_all(i),xlonin,ylatin,ztra1_all(i),itramem_all(i), &
142         topo,pvi,qvi,rhoi,hmixi,tri,tti,(xmass1_all(i,j),j=1,nspec)
143   
144    if (xlonin.eq.-9999.9) goto 100
145    xtra1_all(i)=(xlonin-xlon0)/dx
146    ytra1_all(i)=(ylatin-ylat0)/dy
147    numparticlecount_all=max(numparticlecount_all,npoint(i))
148    goto 200
149
15099  numpart_all=i-1
151   
152    close(unitpartin)
153
154  end do
155
156
157  ! Distribute particles among processes
158  !**************************************
159  lbnd=1
160  ubnd=0
161
162  do ip=0, mp_partgroup_np-1
163   
164! Extra particle distributed in case remainder in the division
165    if (ip.lt.mod(numpart_all,mp_partgroup_np)) then
166      addone=1
167    else
168      addone=0
169    end if
170
171    if (ip==0) then
172      ubnd=ubnd + numpart_all/mp_partgroup_np + addone
173    else
174      ubnd=lbnd + numpart_all/mp_partgroup_np + addone - 1
175    end if
176   
177    if (ip==mp_pid) then
178     
179      numparticlecount=numparticlecount_all/mp_partgroup_np+addone
180      numpart=numpart_all/mp_partgroup_np+addone
181
182      xtra1(1:numpart) = xtra1_all(lbnd:ubnd)
183      ytra1(1:numpart) = ytra1_all(lbnd:ubnd)
184      ztra1(1:numpart) = ztra1_all(lbnd:ubnd)
185      xmass1(1:numpart,:) = xmass1_all(lbnd:ubnd,:)
186      npoint(1:numpart) = npoint_all(lbnd:ubnd)
187      itramem(1:numpart) = itramem_all(lbnd:ubnd)
188
189    end if
190
191    lbnd=ubnd+1
192
193  end do
194
195  julin=juldate(ibdatein,ibtimein)+real(itimein,kind=dp)/86400._dp
196  if (abs(julin-bdate).gt.1.e-5) goto 994
197  do i=1,numpart
198    julpartin=juldate(ibdatein,ibtimein)+ &
199         real(itramem(i),kind=dp)/86400._dp
200    nclass(i)=min(int(ran1(idummy)*real(nclassunc))+1, &
201         nclassunc)
202    idt(i)=mintime
203    itra1(i)=0
204    itramem(i)=nint((julpartin-bdate)*86400.)
205    itrasplit(i)=ldirect*itsplit
206  end do
207
208  return
209
210
211994   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
212  write(*,*) ' #### ENDING TIME OF PREVIOUS MODEL RUN DOES   #### '
213  write(*,*) ' #### NOT AGREE WITH STARTING TIME OF THIS RUN.#### '
214  call caldate(julin,id1,it1)
215  call caldate(bdate,id2,it2)
216  write(*,*) 'julin: ',julin,id1,it1
217  write(*,*) 'bdate: ',bdate,id2,it2
218  stop
219
220996   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
221  write(*,*) ' #### SPECIES NAMES TO BE READ IN DO NOT       #### '
222  write(*,*) ' #### AGREE WITH CURRENT SETTINGS!             #### '
223  stop
224
225997   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
226  write(*,*) ' #### THE NUMBER OF SPECIES TO BE READ IN DOES #### '
227  write(*,*) ' #### NOT AGREE WITH CURRENT SETTINGS!         #### '
228  stop
229
230998   write(*,*) ' #### FLEXPART MODEL ERROR!   THE FILE         #### '
231  write(*,*) ' #### '//path(2)(1:length(2))//'grid'//' #### '
232  write(*,*) ' #### CANNOT BE OPENED. IF A FILE WITH THIS    #### '
233  write(*,*) ' #### NAME ALREADY EXISTS, DELETE IT AND START #### '
234  write(*,*) ' #### THE PROGRAM AGAIN.                       #### '
235  stop
236
237end subroutine readpartpositions
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG