Context Navigation

FLEXPART.f90 @ e9e0f06

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bug

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format