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

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

Fixed a bug with number of occurances of wet scavenging written to stdout

  • 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    write(*,*) '#### Use the serial version instead.'
169    write(*,FMT='(80("#"))')
170    ! call mpif_finalize
171    ! stop
172  end if
173
174
175
176! Read the age classes to be used
177!********************************
178  if (verbosity.gt.0 .and. lroot) then
179    write(*,*) 'call readageclasses'
180  endif
181  call readageclasses
182
183  if (verbosity.gt.1 .and. lroot) then   
184    call system_clock(count_clock, count_rate, count_max)
185    write(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
186  endif     
187
188  ! Read, which wind fields are available within the modelling period
189  !******************************************************************
190
191  if (verbosity.gt.0 .and. lroot) then
192    write(*,*) 'call readavailable'
193  endif 
194  call readavailable
195!end if
196
197! Read the model grid specifications,
198! both for the mother domain and eventual nests
199!**********************************************
200
201  if (verbosity.gt.0 .and. lroot) then
202    write(*,*) 'call gridcheck'
203  endif
204
205  call gridcheck
206
207  if (verbosity.gt.1 .and. lroot) then   
208    call system_clock(count_clock, count_rate, count_max)
209    write(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
210  endif     
211
212
213  if (verbosity.gt.0 .and. lroot) then
214    write(*,*) 'call gridcheck_nests'
215  endif 
216  call gridcheck_nests
217
218
219! Read the output grid specifications
220!************************************
221
222  if (verbosity.gt.0 .and. lroot) then
223    WRITE(*,*) 'call readoutgrid'
224  endif
225
226  call readoutgrid
227
228  if (nested_output.eq.1) then
229    call readoutgrid_nest
230    if (verbosity.gt.0.and.lroot) then
231      WRITE(*,*) '# readoutgrid_nest'
232    endif
233  endif
234
235  ! Read the receptor points for which extra concentrations are to be calculated
236  !*****************************************************************************
237
238  if (verbosity.eq.1 .and. lroot) then
239    print*,'call readreceptors'
240  endif
241  call readreceptors
242
243  ! Read the physico-chemical species property table
244  !*************************************************
245  !SEC: now only needed SPECIES are read in readreleases.f
246  !call readspecies
247
248
249  ! Read the landuse inventory
250  !***************************
251
252  if (verbosity.gt.0 .and. lroot) then
253    print*,'call readlanduse'
254  endif
255  call readlanduse
256
257
258! Assign fractional cover of landuse classes to each ECMWF grid point
259!********************************************************************
260
261  if (verbosity.gt.0 .and. lroot) then
262    print*,'call assignland'
263  endif
264  call assignland
265
266  ! Read the coordinates of the release locations
267  !**********************************************
268
269  if (verbosity.gt.0 .and. lroot) then
270    print*,'call readreleases'
271  endif
272  call readreleases
273
274
275! Read and compute surface resistances to dry deposition of gases
276!****************************************************************
277
278  if (verbosity.gt.0 .and. lroot) then
279    print*,'call readdepo'
280  endif
281  call readdepo
282
283  ! Convert the release point coordinates from geografical to grid coordinates
284  !***************************************************************************
285
286  call coordtrafo 
287  if (verbosity.gt.0 .and. lroot) then
288    print*,'call coordtrafo'
289  endif
290
291
292  ! Initialize all particles to non-existent
293  !*****************************************
294
295  if (verbosity.gt.0 .and. lroot) then
296    print*,'Initialize all particles to non-existent'
297  endif
298
299  do j=1, size(itra1) ! maxpart_mpi
300    itra1(j)=-999999999
301  end do
302
303  ! For continuation of previous run, read in particle positions
304  !*************************************************************
305
306  if (ipin.eq.1) then
307    if (verbosity.gt.0 .and. lroot) then
308      print*,'call readpartpositions'
309    endif
310    call readpartpositions
311  else
312    if (verbosity.gt.0 .and. lroot) then
313      print*,'numpart=0, numparticlecount=0'
314    endif
315    numpart=0
316    numparticlecount=0
317  endif
318
319
320  ! Calculate volume, surface area, etc., of all output grid cells
321  ! Allocate fluxes and OHfield if necessary
322  !***************************************************************
323
324
325  if (verbosity.gt.0 .and. lroot) then
326    print*,'call outgrid_init'
327  endif
328  call outgrid_init
329  if (nested_output.eq.1) call outgrid_init_nest
330
331  ! Read the OH field
332  !******************
333
334  if (OHREA.eqv..true.) then
335    if (verbosity.gt.0 .and. lroot) then
336      print*,'call readOHfield'
337    endif
338    call readOHfield
339  endif
340
341  !! testing !!
342  ! open(999,file=trim(path(1))//'OH_FIELDS/jscalar_50N.txt',action='write',status='new')
343  ! open(998,file=trim(path(1))//'OH_FIELDS/jscalar_50S.txt',action='write',status='new')
344
345
346  ! Write basic information on the simulation to a file "header"
347  ! and open files that are to be kept open throughout the simulation
348  !******************************************************************
349
350  if (mp_measure_time) call mpif_mtime('iotime',0)
351  if (lroot) then ! MPI: this part root process only
352
353  if (lnetcdfout.eq.1) then
354    call writeheader_netcdf(lnest=.false.)
355  else
356    call writeheader
357  end if
358
359  if (nested_output.eq.1) then
360    if (lnetcdfout.eq.1) then
361      call writeheader_netcdf(lnest=.true.)
362    else
363      call writeheader_nest
364    endif
365  endif
366
367!
368    if (verbosity.gt.0) then
369      print*,'call writeheader'
370    endif
371
372    call writeheader
373! FLEXPART 9.2 ticket ?? write header in ASCII format
374    call writeheader_txt
375!if (nested_output.eq.1) call writeheader_nest
376    if (nested_output.eq.1.and.surf_only.ne.1) call writeheader_nest
377    if (nested_output.eq.1.and.surf_only.eq.1) call writeheader_nest_surf
378    if (nested_output.ne.1.and.surf_only.eq.1) call writeheader_surf
379  end if ! (mpif_pid == 0)
380
381  if (mp_measure_time) call mpif_mtime('iotime',0)
382
383  !open(unitdates,file=path(2)(1:length(2))//'dates')
384
385  if (verbosity.gt.0 .and. lroot) then
386    print*,'call openreceptors'
387  endif
388  call openreceptors
389  if ((iout.eq.4).or.(iout.eq.5)) call openouttraj
390
391  ! Releases can only start and end at discrete times (multiples of lsynctime)
392  !***************************************************************************
393
394  if (verbosity.gt.0 .and. lroot) then
395    print*,'discretize release times'
396  endif
397  do i=1,numpoint
398    ireleasestart(i)=nint(real(ireleasestart(i))/real(lsynctime))*lsynctime
399    ireleaseend(i)=nint(real(ireleaseend(i))/real(lsynctime))*lsynctime
400  end do
401
402  ! Initialize cloud-base mass fluxes for the convection scheme
403  !************************************************************
404
405  if (verbosity.gt.0 .and. lroot) then
406    print*,'Initialize cloud-base mass fluxes for the convection scheme'
407  endif
408
409  do jy=0,nymin1
410    do ix=0,nxmin1
411      cbaseflux(ix,jy)=0.
412    end do
413  end do
414  do inest=1,numbnests
415    do jy=0,nyn(inest)-1
416      do ix=0,nxn(inest)-1
417        cbasefluxn(ix,jy,inest)=0.
418      end do
419    end do
420  end do
421
422
423! Calculate particle trajectories
424!********************************
425
426  if (verbosity.gt.0.and. lroot) then
427    if (verbosity.gt.1) then   
428      CALL SYSTEM_CLOCK(count_clock, count_rate, count_max)
429      WRITE(*,*) 'SYSTEM_CLOCK',(count_clock - count_clock0)/real(count_rate) !, count_rate, count_max
430    endif
431    print*,'call timemanager'
432  endif
433  if (info_flag.eq.1) then
434    print*, 'info only mode (stop)'   
435    call mpif_finalize
436    stop
437  endif
438
439
440  call timemanager
441
442
443! NIK 16.02.2005
444  if (lroot) then
445    call MPI_Reduce(MPI_IN_PLACE, tot_blc_count, 1, MPI_INTEGER, MPI_SUM, id_root, &
446         & mp_comm_used, mp_ierr)
447    call MPI_Reduce(MPI_IN_PLACE, tot_inc_count, 1, MPI_INTEGER, MPI_SUM, id_root, &
448         & mp_comm_used, mp_ierr)
449  else
450    if (mp_partgroup_pid.ge.0) then ! Skip for readwind process
451      call MPI_Reduce(tot_blc_count, tot_blc_count, 1, MPI_INTEGER, MPI_SUM, id_root, &
452           & mp_comm_used, mp_ierr)
453      call MPI_Reduce(tot_inc_count, tot_inc_count, 1, MPI_INTEGER, MPI_SUM, id_root, &
454           & mp_comm_used, mp_ierr)
455    end if
456  end if
457
458  if (lroot) then
459    write(*,*) '**********************************************'
460    write(*,*) 'Total number of occurences of below-cloud scavenging', tot_blc_count
461    write(*,*) 'Total number of occurences of in-cloud    scavenging', tot_inc_count
462    write(*,*) '**********************************************'
463
464    write(*,*) 'CONGRATULATIONS: YOU HAVE SUCCESSFULLY COMPLETED A FLE&
465         &XPART MODEL RUN!'
466  end if
467
468  if (mp_measure_time) call mpif_mtime('flexpart',1)
469
470  call mpif_finalize
471
472end program flexpart
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG