source: flexpart.git/src/getfields_mpi.f90 @ 5f9d14a

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

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