source: flexpart.git/src/FLEXPART_MPI.f90 @ df967a99

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

Parallel version exits with error message if trying to run backwards

  • Property mode set to 100644
File size: 14.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
22program flexpart
23
24  !*****************************************************************************
25  !                                                                            *
26  !     This is the Lagrangian Particle Dispersion Model FLEXPART.             *
27  !     The main program manages the reading of model run specifications, etc. *
28  !     All actual computing is done within subroutine timemanager.            *
29  !                                                                            *
30  !     Author: A. Stohl                                                       *
31  !                                                                            *
32  !     18 May 1996                                                            *
33  !                                                                            *
34  !*****************************************************************************
35  !                                                                            *
36  ! Variables:                                                                 *
37  !                                                                            *
38  ! Constants:                                                                 *
39  !                                                                            *
40  !*****************************************************************************
41
42  use point_mod
43  use par_mod
44  use com_mod
45  use conv_mod
46  use mpi_mod
47  use netcdf_output_mod, only: writeheader_netcdf
48  use random_mod, only: gasdev1
49
50  implicit none
51
52  integer :: i,j,ix,jy,inest
53  integer :: idummy = -320
54  character(len=256) :: inline_options  !pathfile, flexversion, arg2
55
56
57! Initialize mpi
58!*********************
59  call mpif_init
60
61  if (mp_measure_time) call mpif_mtime('flexpart',0)
62
63! Initialize arrays in com_mod
64!*****************************
65  call com_mod_allocate(maxpart_mpi)
66
67  ! Generate a large number of random numbers
68  !******************************************
69
70  ! eso: Different seed for each MPI process
71  idummy=idummy+mp_seed
72
73  do i=1,maxrand-1,2
74    call gasdev1(idummy,rannumb(i),rannumb(i+1))
75  end do
76  call gasdev1(idummy,rannumb(maxrand),rannumb(maxrand-1))
77
78  ! FLEXPART version string
79  flexversion='Ver. 10 Beta MPI (2015-05-01)'
80  verbosity=0
81
82  ! Read the pathnames where input/output files are stored
83  !*******************************************************
84
85  inline_options='none'
86  select case (iargc())
87  case (2)
88    call getarg(1,arg1)
89    pathfile=arg1
90    call getarg(2,arg2)
91    inline_options=arg2
92  case (1)
93    call getarg(1,arg1)
94    pathfile=arg1
95    if (arg1(1:1).eq.'-') then
96      write(pathfile,'(a11)') './pathnames'
97      inline_options=arg1
98    endif
99  case (0)
100    write(pathfile,'(a11)') './pathnames'
101  end select
102 
103  if (lroot) then
104  ! Print the GPL License statement
105  !*******************************************************
106    print*,'Welcome to FLEXPART ', trim(flexversion)
107    print*,'FLEXPART is free software released under the GNU General Public License.'
108  endif
109 
110  if (inline_options(1:1).eq.'-') then
111    if (trim(inline_options).eq.'-v'.or.trim(inline_options).eq.'-v1') then
112      print*, 'Verbose mode 1: display detailed information during run'
113      verbosity=1
114    endif
115    if (trim(inline_options).eq.'-v2') then
116      print*, 'Verbose mode 2: display more detailed information during run'
117      verbosity=2
118    endif
119    if (trim(inline_options).eq.'-i') then
120      print*, 'Info mode: provide detailed run specific information and stop'
121      verbosity=1
122      info_flag=1
123    endif
124    if (trim(inline_options).eq.'-i2') then
125      print*, 'Info mode: provide more detailed run specific information and stop'
126      verbosity=2
127      info_flag=1
128    endif
129  endif
130           
131  if (verbosity.gt.0) then
132    write(*,*) 'call readpaths'
133  endif
134  call readpaths(pathfile)
135 
136  if (verbosity.gt.1) then !show clock info
137     !print*,'length(4)',length(4)
138     !count=0,count_rate=1000
139    call system_clock(count_clock0, count_rate, count_max)
140     !WRITE(*,*) 'SYSTEM_CLOCK',count, count_rate, count_max
141     !WRITE(*,*) 'SYSTEM_CLOCK, count_clock0', count_clock0
142     !WRITE(*,*) 'SYSTEM_CLOCK, count_rate', count_rate
143     !WRITE(*,*) 'SYSTEM_CLOCK, count_max', count_max
144  endif
145
146  ! Read the user specifications for the current model run
147  !*******************************************************
148
149  if (verbosity.gt.0) then
150    write(*,*) 'call readcommand'
151  endif
152  call readcommand
153  if (verbosity.gt.0 .and. lroot) then
154    write(*,*) '    ldirect=', ldirect
155    write(*,*) '    ibdate,ibtime=',ibdate,ibtime
156    write(*,*) '    iedate,ietime=', iedate,ietime
157    if (verbosity.gt.1) then   
158      CALL SYSTEM_CLOCK(count_clock, count_rate, count_max)
159      write(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
160    endif
161  endif
162
163  ! Exit if trying to run backwards
164  if (ldirect.le.0) then
165    write(*,FMT='(80("#"))')
166    write(*,*) '#### FLEXPART_MPI> ERROR: ', &
167         & 'MPI version not (yet) working with backward runs. ',&
168         & 'Use the serial version instead.'
169    write(*,FMT='(80("#"))')
170    stop
171  end if
172
173
174
175! Read the age classes to be used
176!********************************
177  if (verbosity.gt.0 .and. lroot) then
178    write(*,*) 'call readageclasses'
179  endif
180  call readageclasses
181
182  if (verbosity.gt.1 .and. lroot) then   
183    call system_clock(count_clock, count_rate, count_max)
184    write(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
185  endif     
186
187  ! Read, which wind fields are available within the modelling period
188  !******************************************************************
189
190  if (verbosity.gt.0 .and. lroot) then
191    write(*,*) 'call readavailable'
192  endif 
193  call readavailable
194!end if
195
196! Read the model grid specifications,
197! both for the mother domain and eventual nests
198!**********************************************
199
200  if (verbosity.gt.0 .and. lroot) then
201    write(*,*) 'call gridcheck'
202  endif
203
204  call gridcheck
205
206  if (verbosity.gt.1 .and. lroot) then   
207    call system_clock(count_clock, count_rate, count_max)
208    write(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
209  endif     
210
211
212  if (verbosity.gt.0 .and. lroot) then
213    write(*,*) 'call gridcheck_nests'
214  endif 
215  call gridcheck_nests
216
217
218! Read the output grid specifications
219!************************************
220
221  if (verbosity.gt.0 .and. lroot) then
222    WRITE(*,*) 'call readoutgrid'
223  endif
224
225  call readoutgrid
226
227  if (nested_output.eq.1) then
228    call readoutgrid_nest
229    if (verbosity.gt.0.and.lroot) then
230      WRITE(*,*) '# readoutgrid_nest'
231    endif
232  endif
233
234  ! Read the receptor points for which extra concentrations are to be calculated
235  !*****************************************************************************
236
237  if (verbosity.eq.1 .and. lroot) then
238    print*,'call readreceptors'
239  endif
240  call readreceptors
241
242  ! Read the physico-chemical species property table
243  !*************************************************
244  !SEC: now only needed SPECIES are read in readreleases.f
245  !call readspecies
246
247
248  ! Read the landuse inventory
249  !***************************
250
251  if (verbosity.gt.0 .and. lroot) then
252    print*,'call readlanduse'
253  endif
254  call readlanduse
255
256
257! Assign fractional cover of landuse classes to each ECMWF grid point
258!********************************************************************
259
260  if (verbosity.gt.0 .and. lroot) then
261    print*,'call assignland'
262  endif
263  call assignland
264
265  ! Read the coordinates of the release locations
266  !**********************************************
267
268  if (verbosity.gt.0 .and. lroot) then
269    print*,'call readreleases'
270  endif
271  call readreleases
272
273
274! Read and compute surface resistances to dry deposition of gases
275!****************************************************************
276
277  if (verbosity.gt.0 .and. lroot) then
278    print*,'call readdepo'
279  endif
280  call readdepo
281
282  ! Convert the release point coordinates from geografical to grid coordinates
283  !***************************************************************************
284
285  call coordtrafo 
286  if (verbosity.gt.0 .and. lroot) then
287    print*,'call coordtrafo'
288  endif
289
290
291  ! Initialize all particles to non-existent
292  !*****************************************
293
294  if (verbosity.gt.0 .and. lroot) then
295    print*,'Initialize all particles to non-existent'
296  endif
297
298  do j=1, size(itra1) ! maxpart_mpi
299    itra1(j)=-999999999
300  end do
301
302  ! For continuation of previous run, read in particle positions
303  !*************************************************************
304
305  if (ipin.eq.1) then
306    if (verbosity.gt.0 .and. lroot) then
307      print*,'call readpartpositions'
308    endif
309    call readpartpositions
310  else
311    if (verbosity.gt.0 .and. lroot) then
312      print*,'numpart=0, numparticlecount=0'
313    endif
314    numpart=0
315    numparticlecount=0
316  endif
317
318
319  ! Calculate volume, surface area, etc., of all output grid cells
320  ! Allocate fluxes and OHfield if necessary
321  !***************************************************************
322
323
324  if (verbosity.gt.0 .and. lroot) then
325    print*,'call outgrid_init'
326  endif
327  call outgrid_init
328  if (nested_output.eq.1) call outgrid_init_nest
329
330  ! Read the OH field
331  !******************
332
333  if (OHREA.eqv..true.) then
334    if (verbosity.gt.0 .and. lroot) then
335      print*,'call readOHfield'
336    endif
337    call readOHfield
338  endif
339
340  !! testing !!
341  ! open(999,file=trim(path(1))//'OH_FIELDS/jscalar_50N.txt',action='write',status='new')
342  ! open(998,file=trim(path(1))//'OH_FIELDS/jscalar_50S.txt',action='write',status='new')
343
344
345  ! Write basic information on the simulation to a file "header"
346  ! and open files that are to be kept open throughout the simulation
347  !******************************************************************
348
349  if (mp_measure_time) call mpif_mtime('iotime',0)
350  if (lroot) then ! MPI: this part root process only
351
352  if (lnetcdfout.eq.1) then
353    call writeheader_netcdf(lnest=.false.)
354  else
355    call writeheader
356  end if
357
358  if (nested_output.eq.1) then
359    if (lnetcdfout.eq.1) then
360      call writeheader_netcdf(lnest=.true.)
361    else
362      call writeheader_nest
363    endif
364  endif
365
366!
367    if (verbosity.gt.0) then
368      print*,'call writeheader'
369    endif
370
371    call writeheader
372! FLEXPART 9.2 ticket ?? write header in ASCII format
373    call writeheader_txt
374!if (nested_output.eq.1) call writeheader_nest
375    if (nested_output.eq.1.and.surf_only.ne.1) call writeheader_nest
376    if (nested_output.eq.1.and.surf_only.eq.1) call writeheader_nest_surf
377    if (nested_output.ne.1.and.surf_only.eq.1) call writeheader_surf
378  end if ! (mpif_pid == 0)
379
380  if (mp_measure_time) call mpif_mtime('iotime',0)
381
382  !open(unitdates,file=path(2)(1:length(2))//'dates')
383
384  if (verbosity.gt.0 .and. lroot) then
385    print*,'call openreceptors'
386  endif
387  call openreceptors
388  if ((iout.eq.4).or.(iout.eq.5)) call openouttraj
389
390  ! Releases can only start and end at discrete times (multiples of lsynctime)
391  !***************************************************************************
392
393  if (verbosity.gt.0 .and. lroot) then
394    print*,'discretize release times'
395  endif
396  do i=1,numpoint
397    ireleasestart(i)=nint(real(ireleasestart(i))/real(lsynctime))*lsynctime
398    ireleaseend(i)=nint(real(ireleaseend(i))/real(lsynctime))*lsynctime
399  end do
400
401  ! Initialize cloud-base mass fluxes for the convection scheme
402  !************************************************************
403
404  if (verbosity.gt.0 .and. lroot) then
405    print*,'Initialize cloud-base mass fluxes for the convection scheme'
406  endif
407
408  do jy=0,nymin1
409    do ix=0,nxmin1
410      cbaseflux(ix,jy)=0.
411    end do
412  end do
413  do inest=1,numbnests
414    do jy=0,nyn(inest)-1
415      do ix=0,nxn(inest)-1
416        cbasefluxn(ix,jy,inest)=0.
417      end do
418    end do
419  end do
420
421
422! Calculate particle trajectories
423!********************************
424
425  if (verbosity.gt.0.and. lroot) then
426    if (verbosity.gt.1) then   
427      CALL SYSTEM_CLOCK(count_clock, count_rate, count_max)
428      WRITE(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
429    endif
430    print*,'call timemanager'
431  endif
432  if (info_flag.eq.1) then
433    print*, 'info only mode (stop)'   
434    call mpif_finalize
435    stop
436  endif
437   
438
439
440  call timemanager
441
442! NIK 16.02.2005
443  if (lroot) then
444    call MPI_Reduce(MPI_IN_PLACE, tot_blc_count, 1, mp_pp, MPI_SUM, id_root, &
445         & mp_comm_used, mp_ierr)
446    call MPI_Reduce(MPI_IN_PLACE, tot_inc_count, 1, mp_pp, MPI_SUM, id_root, &
447         & mp_comm_used, mp_ierr)
448  else
449    if (mp_partgroup_pid.ge.0) then ! Skip for readwind process
450      call MPI_Reduce(tot_blc_count, tot_blc_count, 1, mp_pp, MPI_SUM, id_root, &
451           & mp_comm_used, mp_ierr)
452      call MPI_Reduce(tot_inc_count, tot_inc_count, 1, mp_pp, MPI_SUM, id_root, &
453           & mp_comm_used, mp_ierr)
454    end if
455  end if
456
457  if (lroot) then
458    write(*,*) '**********************************************'
459    write(*,*) 'Total number of occurences of below-cloud scavenging', tot_blc_count
460    write(*,*) 'Total number of occurences of in-cloud    scavenging', tot_inc_count
461    write(*,*) '**********************************************'
462
463    write(*,*) 'CONGRATULATIONS: YOU HAVE SUCCESSFULLY COMPLETED A FLE&
464         &XPART MODEL RUN!'
465  end if
466
467  if (mp_measure_time) call mpif_mtime('flexpart',1)
468
469  call mpif_finalize
470
471end program flexpart
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG