source: flexpart.git/src/FLEXPART_MPI.f90 @ 92fab65

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