source: flexpart.git/src/FLEXPART.f90 @ 8a65cb0

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

Added code, makefile for dev branch

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