source: flexpart.git/src/getfields_mpi.f90 @ 3481cc1

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

move license from headers to a different file

  • Property mode set to 100644
File size: 10.9 KB
Line 
1subroutine getfields(itime,nstop,memstat,metdata_format)
2!                       i     o       o
3!*****************************************************************************
4!                                                                            *
5!  This subroutine manages the 3 data fields to be kept in memory.           *
6!  During the first time step of petterssen it has to be fulfilled that the  *
7!  first data field must have |wftime|<itime, i.e. the absolute value of     *
8!  wftime must be smaller than the absolute value of the current time in [s].*
9!  The other 2 fields are the next in time after the first one.              *
10!  Pointers (memind) are used, because otherwise one would have to resort the*
11!  wind fields, which costs a lot of computing time. Here only the pointers  *
12!  are resorted.                                                             *
13!                                                                            *
14!     Author: A. Stohl                                                       *
15!                                                                            *
16!     29 April 1994                                                          *
17!                                                                            *
18!*****************************************************************************
19!  Changes, Bernd C. Krueger, Feb. 2001:
20!   Variables tth,qvh,tthn,qvhn (on eta coordinates) in common block.
21!   Function of nstop extended.
22!
23!  eso 2014:
24!    MPI version.
25!    If running with number of processes >= mpi_mod::read_grp_min,
26!    only one process (mpi_mod::lmpreader=.true.) does the actual reading, while
27!    the rest call this routine only to update memind, memstat etc.
28!
29!    If mpi_mod::lmp_sync=.true., uses 3 fields instead of 2, to allow reading
30!    the newest in the background.
31!
32!    Return memstat, which is the sum of
33!   
34!    memstat=16:     one new field to be read
35!    memstat=32:     two new fields must be read
36!    memstat=1,2,3:  position(s) in array to read next field
37!    memstat=0:      no new fields to be read
38!         
39!   Unified ECMWF and GFS builds                                             
40!   Marian Harustak, 12.5.2017                                               
41!     - Added passing of metdata_format as it was needed by called routines 
42!         
43!*****************************************************************************
44!                                                                            *
45! Variables:                                                                 *
46! lwindinterval [s]    time difference between the two wind fields read in   *
47! indj                 indicates the number of the wind field to be read in  *
48! indmin               remembers the number of wind fields already treated   *
49! memind(2[3])         pointer, on which place the wind fields are stored    *
50! memtime(2[3]) [s]    times of the wind fields, which are kept in memory    *
51! itime [s]            current time since start date of trajectory calcu-    *
52!                      lation                                                *
53! nstop                > 0, if trajectory has to be terminated               *
54! memstat              keep track of change of status for field pointers     *
55! nx,ny,nuvz,nwz       field dimensions in x,y and z direction               *
56! uu(0:nxmax,0:nymax,nuvzmax,numwfmem)   wind components in x-direction [m/s]*
57! vv(0:nxmax,0:nymax,nuvzmax,numwfmem)   wind components in y-direction [m/s]*
58! ww(0:nxmax,0:nymax,nwzmax,numwfmem)    wind components in z-direction      *
59!                                          [deltaeta/s]                      *
60! tt(0:nxmax,0:nymax,nuvzmax,numwfmem)   temperature [K]                     *
61! ps(0:nxmax,0:nymax,numwfmem)           surface pressure [Pa]               *
62! metdata_format     format of metdata (ecmwf/gfs)                           *
63!                                                                            *
64! Constants:                                                                 *
65! idiffmax             maximum allowable time difference between 2 wind      *
66!                      fields                                                *
67!                                                                            *
68!*****************************************************************************
69
70  use par_mod
71  use com_mod
72  use mpi_mod, only: lmpreader,lmp_use_reader,lmp_sync
73  use class_gribfile
74
75  implicit none
76
77  integer :: metdata_format
78  integer :: indj,itime,nstop,memaux,mindread
79! mindread: index of where to read 3rd field
80  integer, intent(out) :: memstat
81! keep track of 3rd field index. updated when new fields are read
82  integer, save :: mind3=0
83
84  real :: uuh(0:nxmax-1,0:nymax-1,nuvzmax)
85  real :: vvh(0:nxmax-1,0:nymax-1,nuvzmax)
86  real :: pvh(0:nxmax-1,0:nymax-1,nuvzmax)
87  real :: wwh(0:nxmax-1,0:nymax-1,nwzmax)
88  real :: uuhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
89  real :: vvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
90  real :: pvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
91  real :: wwhn(0:nxmaxn-1,0:nymaxn-1,nwzmax,maxnests)
92
93  integer :: indmin = 1
94
95
96! Check, if wind fields are available for the current time step
97!**************************************************************
98
99  nstop=0
100  memstat=0
101
102  if ((ldirect*wftime(1).gt.ldirect*itime).or. &
103       (ldirect*wftime(numbwf).lt.ldirect*itime)) then
104    write(*,*) 'FLEXPART WARNING: NO WIND FIELDS ARE AVAILABLE.'
105    write(*,*) 'A TRAJECTORY HAS TO BE TERMINATED.'
106    nstop=4
107    return
108  endif
109
110
111  if ((ldirect*memtime(1).le.ldirect*itime).and. &
112       (ldirect*memtime(2).gt.ldirect*itime)) then
113
114! The right wind fields are already in memory -> don't do anything
115!*****************************************************************
116    memstat=0
117    continue
118
119  else if ((ldirect*memtime(2).le.ldirect*itime).and. &
120       (memtime(2).ne.999999999)) then
121
122
123! Current time is after 2nd wind field
124! -> Resort wind field pointers, so that current time is between 1st and 2nd
125!***************************************************************************
126
127! 2 fields in memory
128!*******************
129    if (lmp_sync) then
130      memaux=memind(1)
131      memind(1)=memind(2)
132      memind(2)=memaux
133      memtime(1)=memtime(2)
134      memstat=16+memind(2)
135      mindread=memind(2)
136
137    else
138! 3 fields in memory
139! Note that the reader process never needs to keep 3 fields in memory,
140! (2 is enough) so it is possible to save some memory
141!*********************************************************************
142      if (mind3.eq.32.or.mind3.eq.19) then
143        if (lmpreader) then
144          memind(1)=2
145          memind(2)=3
146          memind(3)=3
147        else
148          memind(1)=2
149          memind(2)=3
150          memind(3)=1
151        end if
152        memstat=17
153      else if (mind3.eq.17) then
154        if (lmpreader) then
155          memind(1)=3
156          memind(2)=1
157          memind(3)=1
158        else
159          memind(1)=3
160          memind(2)=1
161          memind(3)=2
162        end if
163        memstat=18
164      else if (mind3.eq.18) then
165        if (lmpreader) then
166          memind(1)=1
167          memind(2)=2
168          memind(3)=2
169        else
170          memind(1)=1
171          memind(2)=2
172          memind(3)=3
173        end if
174        memstat=19
175      else
176        write(*,*) '#### getfields_mpi.f90> ERROR: ', &
177             & 'unknown mind3, exiting.', mind3,' ####'
178        stop
179      end if
180      mind3=memstat
181      memtime(1)=memtime(2)
182      mindread=memind(3)
183    end if
184
185
186! Read a new wind field and store it on place memind(2)
187! or memind(3) for the 3-field read-ahead version
188!******************************************************
189    do indj=indmin,numbwf-1
190      if (ldirect*wftime(indj+1).gt.ldirect*itime) then
191        if (lmpreader.or..not. lmp_use_reader) then
192          if (metdata_format.eq.GRIBFILE_CENTRE_ECMWF) then
193            call readwind_ecmwf(indj+1,mindread,uuh,vvh,wwh)
194          else
195            call readwind_gfs(indj+1,mindread,uuh,vvh,wwh)
196          end if
197          call readwind_nests(indj+1,mindread,uuhn,vvhn,wwhn)
198          call calcpar(mindread,uuh,vvh,pvh,metdata_format)
199          call calcpar_nests(mindread,uuhn,vvhn,pvhn,metdata_format)
200          if (metdata_format.eq.GRIBFILE_CENTRE_ECMWF) then
201            call verttransform_ecmwf(mindread,uuh,vvh,wwh,pvh)
202          else
203            call verttransform_gfs(mindread,uuh,vvh,wwh,pvh)
204          end if
205          call verttransform_nests(mindread,uuhn,vvhn,wwhn,pvhn)
206        end if
207        memtime(2)=wftime(indj+1)
208        nstop = 1
209        goto 40
210      endif
211    end do
21240  indmin=indj
213
214   if (WETBKDEP.and.(lmpreader.or.(.not.lmp_use_reader.and.lroot))) then
215        call writeprecip(itime,memind(1))
216      endif
217
218  else
219
220! No wind fields, which can be used, are currently in memory
221! -> read both wind fields
222!***********************************************************
223    memstat=32
224
225    do indj=indmin,numbwf-1
226      if ((ldirect*wftime(indj).le.ldirect*itime).and. &
227           (ldirect*wftime(indj+1).gt.ldirect*itime)) then
228        memind(1)=1
229        if (lmpreader.or..not.lmp_use_reader) then
230          if (metdata_format.eq.GRIBFILE_CENTRE_ECMWF) then
231            call readwind_ecmwf(indj,memind(1),uuh,vvh,wwh)
232          else
233            call readwind_gfs(indj,memind(1),uuh,vvh,wwh)
234          end if
235          call readwind_nests(indj,memind(1),uuhn,vvhn,wwhn)
236          call calcpar(memind(1),uuh,vvh,pvh,metdata_format)
237          call calcpar_nests(memind(1),uuhn,vvhn,pvhn,metdata_format)
238          if (metdata_format.eq.GRIBFILE_CENTRE_ECMWF) then
239            call verttransform_ecmwf(memind(1),uuh,vvh,wwh,pvh)
240          else
241            call verttransform_gfs(memind(1),uuh,vvh,wwh,pvh)
242          end if
243          call verttransform_nests(memind(1),uuhn,vvhn,wwhn,pvhn)
244        end if
245        memtime(1)=wftime(indj)
246        memind(2)=2
247        if (lmpreader.or..not.lmp_use_reader) then
248          if (metdata_format.eq.GRIBFILE_CENTRE_ECMWF) then
249            call readwind_ecmwf(indj+1,memind(2),uuh,vvh,wwh)
250          else
251            call readwind_gfs(indj+1,memind(2),uuh,vvh,wwh)
252          end if
253          call readwind_nests(indj+1,memind(2),uuhn,vvhn,wwhn)
254          call calcpar(memind(2),uuh,vvh,pvh,metdata_format)
255          call calcpar_nests(memind(2),uuhn,vvhn,pvhn,metdata_format)
256          if (metdata_format.eq.GRIBFILE_CENTRE_ECMWF) then
257            call verttransform_ecmwf(memind(2),uuh,vvh,wwh,pvh)
258          else
259            call verttransform_gfs(memind(2),uuh,vvh,wwh,pvh)
260          end if
261          call verttransform_nests(memind(2),uuhn,vvhn,wwhn,pvhn)
262        end if
263        memtime(2)=wftime(indj+1)
264! DEV: not used
265        if (.not.lmp_sync) memind(3)=3 ! indicate position of next read
266        nstop = 1
267        goto 60
268      endif
269    end do
27060  indmin=indj
271
272!   if (WETBKDEP.and.lroot) then
273    if (WETBKDEP.and.(lmpreader.or.(.not.lmp_use_reader.and.lroot))) then
274      call writeprecip(itime,memind(1))
275    endif
276
277  endif
278
279  lwindinterv=abs(memtime(2)-memtime(1))
280
281  if (lwindinterv.gt.idiffmax) nstop=3
282
283end subroutine getfields
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG