source: flexpart.git/src/readpartpositions_mpi.f90

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bug
Last change on this file was 92fab65, checked in by Ignacio Pisso <ip@…>, 4 years ago

add SPDX-License-Identifier to all .f90 files

  • Property mode set to 100644
File size: 7.1 KB
Line 
1! SPDX-FileCopyrightText: FLEXPART 1998-2019, see flexpart_license.txt
2! SPDX-License-Identifier: GPL-3.0-or-later
3
4subroutine readpartpositions
5
6  !*****************************************************************************
7  !                                                                            *
8  !   This routine opens the particle dump file and reads all the particle     *
9  !   positions from a previous run to initialize the current run.             *
10  !                                                                            *
11  !                                                                            *
12  !     Author: A. Stohl                                                       *
13  !                                                                            *
14  !     24 March 2000                                                          *
15  !                                                                            *
16  !   CHANGES                                                                  *
17  !     12/2014 eso: MPI version                                               *
18  !                  Root process reads positions and distributes the data     *
19  !                  :TODO: The above..                                        *
20  !*****************************************************************************
21  !                                                                            *
22  ! Variables:                                                                 *
23  !                                                                            *
24  !*****************************************************************************
25
26  use par_mod
27  use com_mod
28  use random_mod, only: ran1
29  use mpi_mod !, only: mp_seed
30
31  implicit none
32
33  integer :: ibdatein,ibtimein,nspecin,itimein,numpointin,i,j,ix,ip
34  integer :: id1,id2,it1,it2
35  integer :: addone,numparticlecount_all,numpart_all,lbnd,ubnd
36  real :: xlonin,ylatin,topo,hmixi,pvi,qvi,rhoi,tri,tti
37  character :: specin*7
38  real(kind=dp) :: julin,julpartin,juldate
39
40  integer :: idummy = -8
41
42  ! These variables are allocated at the root process for all particles in file.
43  real,dimension(maxpart) :: xtra1_all,ytra1_all,ztra1_all
44  real,dimension(maxpart,maxspec) :: xmass1_all
45  integer,dimension(maxpart) :: npoint_all,itramem_all
46
47  ! Different seed for each process
48  idummy=idummy+mp_seed
49
50  numparticlecount=0
51  numparticlecount_all=0
52  numpart_all=0
53
54  ! Open header file of dumped particle data
55  ! Each MPI process sequentially access file (just in case)
56  !*********************************************************
57
58  do ip=0, mp_partgroup_np-1
59    call mpif_mpi_barrier
60    open(unitpartin,file=path(2)(1:length(2))//'header', &
61         form='unformatted',err=998)
62
63    read(unitpartin) ibdatein,ibtimein
64    read(unitpartin)
65    read(unitpartin)
66
67    read(unitpartin)
68    read(unitpartin)
69    read(unitpartin) nspecin
70    nspecin=nspecin/3
71    if ((ldirect.eq.1).and.(nspec.ne.nspecin)) goto 997
72
73    do i=1,nspecin
74      read(unitpartin)
75      read(unitpartin)
76      read(unitpartin) j,specin
77      if ((ldirect.eq.1).and.(species(i)(1:7).ne.specin)) goto 996
78    end do
79   
80    read(unitpartin) numpointin
81    if (numpointin.ne.numpoint) then ! goto 995
82
83! eso 2016: moved this warning here to avoid out-of-block goto
84!995   write(*,*) ' #### FLEXPART MODEL WARNING IN READPARTPOSITIONS#### '
85      write(*,*) ' #### FLEXPART MODEL WARNING IN READPARTPOSITIONS#### '
86      write(*,*) ' #### NUMBER OF RELEASE LOCATIONS DOES NOT     #### '
87      write(*,*) ' #### AGREE WITH CURRENT SETTINGS!             #### '
88!  stop
89      goto 999
90    end if
91
92999 continue
93    do i=1,numpointin
94      read(unitpartin)
95      read(unitpartin)
96      read(unitpartin)
97      read(unitpartin)
98      do j=1,nspec
99        read(unitpartin)
100        read(unitpartin)
101        read(unitpartin)
102      end do
103    end do
104    read(unitpartin)
105    read(unitpartin)
106
107    do ix=0,numxgrid-1
108      read(unitpartin)
109    end do
110    close(unitpartin)
111
112
113  ! Open data file of dumped particle data
114  ! All processes read the whole file
115  !***************************************
116
117    open(unitpartin,file=path(2)(1:length(2))//'partposit_end', &
118         form='unformatted',err=998)
119   
120100 read(unitpartin,end=99) itimein
121    i=0
122200 i=i+1
123    read(unitpartin) npoint_all(i),xlonin,ylatin,ztra1_all(i),itramem_all(i), &
124         topo,pvi,qvi,rhoi,hmixi,tri,tti,(xmass1_all(i,j),j=1,nspec)
125   
126    if (xlonin.eq.-9999.9) goto 100
127    xtra1_all(i)=(xlonin-xlon0)/dx
128    ytra1_all(i)=(ylatin-ylat0)/dy
129    numparticlecount_all=max(numparticlecount_all,npoint(i))
130    goto 200
131
13299  numpart_all=i-1
133   
134    close(unitpartin)
135
136  end do
137
138
139  ! Distribute particles among processes
140  !**************************************
141  lbnd=1
142  ubnd=0
143
144  do ip=0, mp_partgroup_np-1
145   
146! Extra particle distributed in case remainder in the division
147    if (ip.lt.mod(numpart_all,mp_partgroup_np)) then
148      addone=1
149    else
150      addone=0
151    end if
152
153    if (ip==0) then
154      ubnd=ubnd + numpart_all/mp_partgroup_np + addone
155    else
156      ubnd=lbnd + numpart_all/mp_partgroup_np + addone - 1
157    end if
158   
159    if (ip==mp_pid) then
160     
161      numparticlecount=numparticlecount_all/mp_partgroup_np+addone
162      numpart=numpart_all/mp_partgroup_np+addone
163
164      xtra1(1:numpart) = xtra1_all(lbnd:ubnd)
165      ytra1(1:numpart) = ytra1_all(lbnd:ubnd)
166      ztra1(1:numpart) = ztra1_all(lbnd:ubnd)
167      xmass1(1:numpart,:) = xmass1_all(lbnd:ubnd,:)
168      npoint(1:numpart) = npoint_all(lbnd:ubnd)
169      itramem(1:numpart) = itramem_all(lbnd:ubnd)
170
171    end if
172
173    lbnd=ubnd+1
174
175  end do
176
177  julin=juldate(ibdatein,ibtimein)+real(itimein,kind=dp)/86400._dp
178  if (abs(julin-bdate).gt.1.e-5) goto 994
179  do i=1,numpart
180    julpartin=juldate(ibdatein,ibtimein)+ &
181         real(itramem(i),kind=dp)/86400._dp
182    nclass(i)=min(int(ran1(idummy)*real(nclassunc))+1, &
183         nclassunc)
184    idt(i)=mintime
185    itra1(i)=0
186    itramem(i)=nint((julpartin-bdate)*86400.)
187    itrasplit(i)=ldirect*itsplit
188  end do
189
190  return
191
192
193994   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
194  write(*,*) ' #### ENDING TIME OF PREVIOUS MODEL RUN DOES   #### '
195  write(*,*) ' #### NOT AGREE WITH STARTING TIME OF THIS RUN.#### '
196  call caldate(julin,id1,it1)
197  call caldate(bdate,id2,it2)
198  write(*,*) 'julin: ',julin,id1,it1
199  write(*,*) 'bdate: ',bdate,id2,it2
200  stop
201
202996   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
203  write(*,*) ' #### SPECIES NAMES TO BE READ IN DO NOT       #### '
204  write(*,*) ' #### AGREE WITH CURRENT SETTINGS!             #### '
205  stop
206
207997   write(*,*) ' #### FLEXPART MODEL ERROR IN READPARTPOSITIONS#### '
208  write(*,*) ' #### THE NUMBER OF SPECIES TO BE READ IN DOES #### '
209  write(*,*) ' #### NOT AGREE WITH CURRENT SETTINGS!         #### '
210  stop
211
212998   write(*,*) ' #### FLEXPART MODEL ERROR!   THE FILE         #### '
213  write(*,*) ' #### '//path(2)(1:length(2))//'grid'//' #### '
214  write(*,*) ' #### CANNOT BE OPENED. IF A FILE WITH THIS    #### '
215  write(*,*) ' #### NAME ALREADY EXISTS, DELETE IT AND START #### '
216  write(*,*) ' #### THE PROGRAM AGAIN.                       #### '
217  stop
218
219end subroutine readpartpositions
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG