Context Navigation

start_em.F

Last change on this file was 11, checked in by aslmd, 14 years ago
spiga@svn-planeto:ajoute le modele meso-echelle martien
File size: 40.1 KB

Line
1	!-------------------------------------------------------------------
2
3	SUBROUTINE start_domain_em ( grid, allowed_to_read &
4	! Actual arguments generated from Registry
5	# include "em_dummy_new_args.inc"
6	!
7	)
8
9	USE module_domain
10	USE module_dm
11	! USE module_io_domain
12	USE module_state_description
13	USE module_model_constants
14	USE module_bc
15	USE module_bc_em
16	! USE module_timing
17	USE module_configure
18	USE module_tiles
19
20	USE module_physics_init
21	#ifdef WRF_CHEM
22	USE module_aerosols_sorgam, only: sum_pm_sorgam
23	USE module_mosaic_driver, only: sum_pm_mosaic
24	#endif
25
26	#ifdef DM_PARALLEL
27	USE module_dm
28	#endif
29
30	!!debug
31	!USE module_compute_geop
32
33	USE module_model_constants
34	IMPLICIT NONE
35	! Input data.
36	TYPE (domain) :: grid
37
38	LOGICAL , INTENT(IN) :: allowed_to_read
39
40	! Definitions of dummy arguments to this routine (generated from Registry).
41	# include "em_dummy_new_decl.inc"
42
43	! Structure that contains run-time configuration (namelist) data for domain
44	TYPE (grid_config_rec_type) :: config_flags
45
46	! Local data
47	INTEGER :: &
48	ids, ide, jds, jde, kds, kde, &
49	ims, ime, jms, jme, kms, kme, &
50	ips, ipe, jps, jpe, kps, kpe, &
51	its, ite, jts, jte, kts, kte, &
52	ij,i,j,k,ii,jj,kk,loop,error,l
53
54	INTEGER :: i_m
55
56	REAL :: p00, t00, a, p_surf, pd_surf, tiso, temp
57	#ifdef WRF_CHEM
58	REAL RGASUNIV ! universal gas constant [ J/mol-K ]
59	PARAMETER ( RGASUNIV = 8.314510 )
60	REAL,DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: &
61	z_at_w,convfac
62	REAL :: tempfac
63	#endif
64
65	REAL :: qvf1, qvf2, qvf
66	REAL :: MPDT
67	REAL :: spongeweight
68	LOGICAL :: first_trip_for_this_domain, start_of_simulation
69
70	REAL :: lat1 , lat2 , lat3 , lat4
71	REAL :: lon1 , lon2 , lon3 , lon4
72	INTEGER :: num_points_lat_lon , iloc , jloc
73	CHARACTER (LEN=132) :: message
74
75	! Needed by some comm layers, grid%e.g. RSL. If needed, nmm_data_calls.inc is
76	! generated from the registry. The definition of REGISTER_I1 allows
77	! I1 data to be communicated in this routine if necessary.
78	#ifdef DM_PARALLEL
79	# include "em_data_calls.inc"
80	#endif
81	CALL get_ijk_from_grid ( grid , &
82	ids, ide, jds, jde, kds, kde, &
83	ims, ime, jms, jme, kms, kme, &
84	ips, ipe, jps, jpe, kps, kpe )
85
86	kts = kps ; kte = kpe ! note that tile is entire patch
87	its = ips ; ite = ipe ! note that tile is entire patch
88	jts = jps ; jte = jpe ! note that tile is entire patch
89
90	CALL model_to_grid_config_rec ( grid%id , model_config_rec , config_flags )
91
92	IF ( ( MOD (ide-ids,config_flags%parent_grid_ratio) .NE. 0 ) .OR. &
93	( MOD (jde-jds,config_flags%parent_grid_ratio) .NE. 0 ) ) THEN
94	WRITE(message, FMT='("Nested dimensions are illegal for domain ",I2,": Both &
95	&MOD(",I4,"-",I1,",",I2,") and MOD(",I4,"-",I1,",",I2,") must = 0" )') &
96	grid%id,ide,ids,config_flags%parent_grid_ratio,jde,jds,config_flags%parent_grid_ratio
97	CALL wrf_error_fatal ( message )
98	END IF
99
100	! here we check to see if the boundary conditions are set properly
101
102	CALL boundary_condition_check( config_flags, bdyzone, error, grid%id )
103
104	!kludge - need to stop CG from resetting precip and phys tendencies to zero
105	! when we are in here due to a nest being spawned, we want to still
106	! recompute the base state, but that is about it
107	! This is temporary and will need to be changed when grid%itimestep is removed.
108
109	IF ( grid%itimestep .EQ. 0 ) THEN
110	first_trip_for_this_domain = .TRUE.
111	ELSE
112	first_trip_for_this_domain = .FALSE.
113	END IF
114
115	IF ( .not. ( config_flags%restart .or. grid%moved ) ) THEN
116	grid%itimestep=0
117	ENDIF
118
119	IF ( config_flags%restart .or. grid%moved ) THEN
120	first_trip_for_this_domain = .TRUE.
121	ENDIF
122
123	IF (config_flags%specified) THEN
124	!
125	! Arrays for specified boundary conditions
126	! wig: Add a combined exponential+linear weight on the mother boundaries
127	! following code changes by Ruby Leung. For the nested grid, there
128	! appears to be some problems when a sponge is used. The points where
129	! processors meet have problematic values.
130
131	DO loop = grid%spec_zone + 1, grid%spec_zone + grid%relax_zone
132	grid%fcx(loop) = 0.1 / grid%dt * (grid%spec_zone + grid%relax_zone - loop) / (grid%relax_zone - 1)
133	grid%gcx(loop) = 1.0 / grid%dt / 50. * (grid%spec_zone + grid%relax_zone - loop) / (grid%relax_zone - 1)
134	! spongeweight=exp(-(loop-2)/3.)
135	! grid%fcx(loop) = grid%fcx(loop)*spongeweight
136	! grid%gcx(loop) = grid%gcx(loop)*spongeweight
137	ENDDO
138
139	ELSE IF (config_flags%nested) THEN
140	!
141	! Arrays for specified boundary conditions
142
143	DO loop = grid%spec_zone + 1, grid%spec_zone + grid%relax_zone
144	grid%fcx(loop) = 0.1 / grid%dt * (grid%spec_zone + grid%relax_zone - loop) / (grid%relax_zone - 1)
145	grid%gcx(loop) = 1.0 / grid%dt / 50. * (grid%spec_zone + grid%relax_zone - loop) / (grid%relax_zone - 1)
146	! spongeweight=exp(-(loop-2)/3.)
147	! grid%fcx(loop) = grid%fcx(loop)*spongeweight
148	! grid%gcx(loop) = grid%gcx(loop)*spongeweight
149	! grid%fcx(loop) = 0.
150	! grid%gcx(loop) = 0.
151	ENDDO
152
153	grid%dtbc = 0.
154
155	ENDIF
156
157	! IF ( grid%id .NE. 1 ) THEN
158
159	! Every time a domain starts or every time a domain moves, this routine is called. We want
160	! the center (middle) lat/lon of the grid for the metacode. The lat/lon values are
161	! defined at mass points. Depending on the even/odd points in the SN and WE directions,
162	! we end up with the middle point as either 1 point or an average of either 2 or 4 points.
163	! Add to this, the need to make sure that we are on the correct patch to retrieve the
164	! value of the lat/lon, AND that the lat/lons (for an average) may not all be on the same
165	! patch. Once we find the correct value for lat lon, we need to keep it around on all patches,
166	! which is where the wrf_dm_min_real calls come in.
167
168	IF ( ( MOD(ide,2) .EQ. 0 ) .AND. ( MOD(jde,2) .EQ. 0 ) ) THEN
169	num_points_lat_lon = 1
170	iloc = ide/2
171	jloc = jde/2
172	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
173	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
174	lat1 = grid%xlat (iloc,jloc)
175	lon1 = grid%xlong(iloc,jloc)
176	ELSE
177	lat1 = 99999.
178	lon1 = 99999.
179	END IF
180	lat1 = wrf_dm_min_real ( lat1 )
181	lon1 = wrf_dm_min_real ( lon1 )
182	CALL nl_set_cen_lat ( grid%id , lat1 )
183	CALL nl_set_cen_lon ( grid%id , lon1 )
184	ELSE IF ( ( MOD(ide,2) .NE. 0 ) .AND. ( MOD(jde,2) .EQ. 0 ) ) THEN
185	num_points_lat_lon = 2
186	iloc = (ide-1)/2
187	jloc = jde /2
188	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
189	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
190	lat1 = grid%xlat (iloc,jloc)
191	lon1 = grid%xlong(iloc,jloc)
192	ELSE
193	lat1 = 99999.
194	lon1 = 99999.
195	END IF
196	lat1 = wrf_dm_min_real ( lat1 )
197	lon1 = wrf_dm_min_real ( lon1 )
198
199	iloc = (ide+1)/2
200	jloc = jde /2
201	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
202	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
203	lat2 = grid%xlat (iloc,jloc)
204	lon2 = grid%xlong(iloc,jloc)
205	ELSE
206	lat2 = 99999.
207	lon2 = 99999.
208	END IF
209	lat2 = wrf_dm_min_real ( lat2 )
210	lon2 = wrf_dm_min_real ( lon2 )
211
212	CALL nl_set_cen_lat ( grid%id , ( lat1 + lat2 ) * 0.50 )
213	CALL nl_set_cen_lon ( grid%id , ( lon1 + lon2 ) * 0.50 )
214	ELSE IF ( ( MOD(ide,2) .EQ. 0 ) .AND. ( MOD(jde,2) .NE. 0 ) ) THEN
215	num_points_lat_lon = 2
216	iloc = ide /2
217	jloc = (jde-1)/2
218	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
219	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
220	lat1 = grid%xlat (iloc,jloc)
221	lon1 = grid%xlong(iloc,jloc)
222	ELSE
223	lat1 = 99999.
224	lon1 = 99999.
225	END IF
226	lat1 = wrf_dm_min_real ( lat1 )
227	lon1 = wrf_dm_min_real ( lon1 )
228
229	iloc = ide /2
230	jloc = (jde+1)/2
231	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
232	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
233	lat2 = grid%xlat (iloc,jloc)
234	lon2 = grid%xlong(iloc,jloc)
235	ELSE
236	lat2 = 99999.
237	lon2 = 99999.
238	END IF
239	lat2 = wrf_dm_min_real ( lat2 )
240	lon2 = wrf_dm_min_real ( lon2 )
241
242	CALL nl_set_cen_lat ( grid%id , ( lat1 + lat2 ) * 0.50 )
243	CALL nl_set_cen_lon ( grid%id , ( lon1 + lon2 ) * 0.50 )
244	ELSE IF ( ( MOD(ide,2) .NE. 0 ) .AND. ( MOD(jde,2) .NE. 0 ) ) THEN
245	num_points_lat_lon = 4
246	iloc = (ide-1)/2
247	jloc = (jde-1)/2
248	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
249	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
250	lat1 = grid%xlat (iloc,jloc)
251	lon1 = grid%xlong(iloc,jloc)
252	ELSE
253	lat1 = 99999.
254	lon1 = 99999.
255	END IF
256	lat1 = wrf_dm_min_real ( lat1 )
257	lon1 = wrf_dm_min_real ( lon1 )
258
259	iloc = (ide+1)/2
260	jloc = (jde-1)/2
261	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
262	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
263	lat2 = grid%xlat (iloc,jloc)
264	lon2 = grid%xlong(iloc,jloc)
265	ELSE
266	lat2 = 99999.
267	lon2 = 99999.
268	END IF
269	lat2 = wrf_dm_min_real ( lat2 )
270	lon2 = wrf_dm_min_real ( lon2 )
271
272	iloc = (ide-1)/2
273	jloc = (jde+1)/2
274	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
275	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
276	lat3 = grid%xlat (iloc,jloc)
277	lon3 = grid%xlong(iloc,jloc)
278	ELSE
279	lat3 = 99999.
280	lon3 = 99999.
281	END IF
282	lat3 = wrf_dm_min_real ( lat3 )
283	lon3 = wrf_dm_min_real ( lon3 )
284
285	iloc = (ide+1)/2
286	jloc = (jde+1)/2
287	IF ( ( ips .LE. iloc ) .AND. ( ipe .GE. iloc ) .AND. &
288	( jps .LE. jloc ) .AND. ( jpe .GE. jloc ) ) THEN
289	lat4 = grid%xlat (iloc,jloc)
290	lon4 = grid%xlong(iloc,jloc)
291	ELSE
292	lat4 = 99999.
293	lon4 = 99999.
294	END IF
295	lat4 = wrf_dm_min_real ( lat4 )
296	lon4 = wrf_dm_min_real ( lon4 )
297
298	CALL nl_set_cen_lat ( grid%id , ( lat1 + lat2 + lat3 + lat4 ) * 0.25 )
299	CALL nl_set_cen_lon ( grid%id , ( lon1 + lon2 + lon3 + lon4 ) * 0.25 )
300	END IF
301	! END IF
302
303	IF ( .NOT. config_flags%restart .AND. &
304	(( config_flags%input_from_hires ) .OR. ( config_flags%input_from_file ))) THEN
305
306	IF ( config_flags%map_proj .EQ. 0 ) THEN
307	CALL wrf_error_fatal ( 'start_domain: Idealized case cannot have a separate nested input file' )
308	END IF
309
310	CALL nl_get_base_pres ( 1 , p00 )
311	CALL nl_get_base_temp ( 1 , t00 )
312	CALL nl_get_base_lapse ( 1 , a )
313	CALL nl_get_tiso ( 1 , tiso )
314
315	! Base state potential temperature and inverse density (alpha = 1/rho) from
316	! the half eta levels and the base-profile surface pressure. Compute 1/rho
317	! from equation of state. The potential temperature is a perturbation from t0.
318
319	DO j = jts, MIN(jte,jde-1)
320	DO i = its, MIN(ite,ide-1)
321
322	! Base state pressure is a function of eta level and terrain, only, plus
323	! the hand full of constants: p00 (sea level pressure, Pa), t00 (sea level
324	! temperature, K), and A (temperature difference, from 1000 mb to 300 mb, K).
325
326	p_surf = p00 * EXP ( -t00/a + ( (t00/a)*2 - 2.ggrid%ht(i,j)/a/r_d ) *0.5 )
327
328	DO k = 1, kte-1
329	grid%em_pb(i,k,j) = grid%em_znu(k)*(p_surf - grid%p_top) + grid%p_top
330	! grid%em_t_init(i,k,j) = (t00 + ALOG(grid%em_pb(i,k,j)/p00))(p00/grid%em_pb(i,k,j))**(r_d/cp) - t0
331	temp = MAX ( tiso, t00 + A*LOG(grid%em_pb(i,k,j)/p00))
332	grid%em_t_init(i,k,j) = temp(p00/grid%em_pb(i,k,j))*(r_d/cp) - t0
333	grid%em_alb(i,k,j) = (r_d/p1000mb)(grid%em_t_init(i,k,j)+t0)(grid%em_pb(i,k,j)/p1000mb)**cvpm
334	END DO
335
336	! Base state mu is defined as base state surface pressure minus grid%p_top
337
338	grid%em_mub(i,j) = p_surf - grid%p_top
339
340	! Integrate base geopotential, starting at terrain elevation. This assures that
341	! the base state is in exact hydrostatic balance with respect to the model equations.
342	! This field is on full levels.
343
344	grid%em_phb(i,1,j) = grid%ht(i,j) * g
345	DO k = 2,kte
346	grid%em_phb(i,k,j) = grid%em_phb(i,k-1,j) - grid%em_dnw(k-1)grid%em_mub(i,j)grid%em_alb(i,k-1,j)
347	END DO
348	END DO
349	END DO
350
351	ENDIF
352
353	IF(.not.config_flags%restart)THEN
354
355	! if this is for a nested domain, the defined/interpolated fields are the _2
356
357	IF ( first_trip_for_this_domain ) THEN
358
359	! data that is expected to be zero must be explicitly initialized as such
360	grid%h_diabatic = 0.
361
362	DO j = jts,min(jte,jde-1)
363	DO k = kts,kte-1
364	DO i = its, min(ite,ide-1)
365	IF ( grid%imask_nostag(i,j) .EQ. 1 ) THEN
366	grid%em_t_1(i,k,j)=grid%em_t_2(i,k,j)
367	ENDIF
368	ENDDO
369	ENDDO
370	ENDDO
371
372	DO j = jts,min(jte,jde-1)
373	DO i = its, min(ite,ide-1)
374	IF ( grid%imask_nostag(i,j) .EQ. 1 ) THEN
375	grid%em_mu_1(i,j)=grid%em_mu_2(i,j)
376	ENDIF
377	ENDDO
378	ENDDO
379	END IF
380
381	! reconstitute base-state fields
382
383	IF(config_flags%max_dom .EQ. 1)THEN
384	! with single domain, grid%em_t_init from wrfinput is OK to use
385	DO j = jts,min(jte,jde-1)
386	DO k = kts,kte-1
387	DO i = its, min(ite,ide-1)
388	IF ( grid%imask_nostag(i,j) .EQ. 1 ) THEN
389	grid%em_pb(i,k,j) = grid%em_znu(k)*grid%em_mub(i,j)+grid%p_top
390	grid%em_alb(i,k,j) = (r_d/p1000mb)(grid%em_t_init(i,k,j)+t0)(grid%em_pb(i,k,j)/p1000mb)**cvpm
391	ENDIF
392	ENDDO
393	ENDDO
394	ENDDO
395	ELSE
396	! with nests, grid%em_t_init generally needs recomputations (since it is not interpolated)
397	DO j = jts,min(jte,jde-1)
398	DO k = kts,kte-1
399	DO i = its, min(ite,ide-1)
400	IF ( grid%imask_nostag(i,j) .EQ. 1 ) THEN
401	grid%em_pb(i,k,j) = grid%em_znu(k)*grid%em_mub(i,j)+grid%p_top
402	grid%em_alb(i,k,j) = -grid%em_rdnw(k)*(grid%em_phb(i,k+1,j)-grid%em_phb(i,k,j))/grid%em_mub(i,j)
403	grid%em_t_init(i,k,j) = grid%em_alb(i,k,j)(p1000mb/r_d)/((grid%em_pb(i,k,j)/p1000mb)*cvpm) - t0
404	ENDIF
405	ENDDO
406	ENDDO
407	ENDDO
408	ENDIF
409
410	DO j = jts,min(jte,jde-1)
411
412	k = kte-1
413	DO i = its, min(ite,ide-1)
414	IF ( grid%imask_nostag(i,j) .EQ. 1 ) THEN
415	qvf1 = 0.5*(moist(i,k,j,P_QV)+moist(i,k,j,P_QV))
416	qvf2 = 1./(1.+qvf1)
417	qvf1 = qvf1*qvf2
418	grid%em_p(i,k,j) = - 0.5(grid%em_mu_1(i,j)+qvf1grid%em_mub(i,j))/grid%em_rdnw(k)/qvf2
419	qvf = 1. + rvovrd*moist(i,k,j,P_QV)
420	grid%em_alt(i,k,j) = (r_d/p1000mb)(grid%em_t_1(i,k,j)+t0)qvf(((grid%em_p(i,k,j)+grid%em_pb(i,k,j))/p1000mb)*cvpm)
421	grid%em_al(i,k,j) = grid%em_alt(i,k,j) - grid%em_alb(i,k,j)
422	ENDIF
423	ENDDO
424
425	DO k = kte-2, 1, -1
426	DO i = its, min(ite,ide-1)
427	IF ( grid%imask_nostag(i,j) .EQ. 1 ) THEN
428	qvf1 = 0.5*(moist(i,k,j,P_QV)+moist(i,k+1,j,P_QV))
429	qvf2 = 1./(1.+qvf1)
430	qvf1 = qvf1*qvf2
431	grid%em_p(i,k,j) = grid%em_p(i,k+1,j) - (grid%em_mu_1(i,j) + qvf1*grid%em_mub(i,j))/qvf2/grid%em_rdn(k+1)
432	qvf = 1. + rvovrd*moist(i,k,j,P_QV)
433	grid%em_alt(i,k,j) = (r_d/p1000mb)(grid%em_t_1(i,k,j)+t0)qvf* &
434	(((grid%em_p(i,k,j)+grid%em_pb(i,k,j))/p1000mb)**cvpm)
435	grid%em_al(i,k,j) = grid%em_alt(i,k,j) - grid%em_alb(i,k,j)
436	ENDIF
437	ENDDO
438	ENDDO
439
440	ENDDO
441
442	ENDIF
443
444	IF ( ( grid%id .NE. 1 ) .AND. .NOT. ( config_flags%restart ) .AND. &
445	( ( config_flags%input_from_hires ) .OR. ( config_flags%input_from_file ) ) ) THEN
446	DO j = jts, MIN(jte,jde-1)
447	DO i = its, MIN(ite,ide-1)
448	grid%em_mu_2(i,j) = grid%em_mu_2(i,j) + grid%em_al(i,1,j) / ( grid%em_alt(i,1,j) * grid%em_alb(i,1,j) ) * &
449	g * ( grid%ht(i,j) - grid%ht_fine(i,j) )
450	END DO
451	END DO
452	DO j = jts,min(jte,jde-1)
453	DO i = its, min(ite,ide-1)
454	grid%em_mu_1(i,j)=grid%em_mu_2(i,j)
455	ENDDO
456	ENDDO
457
458	END IF
459
460	IF ( first_trip_for_this_domain ) THEN
461
462	CALL wrf_debug ( 100 , 'module_start: start_domain_rk: Before call to phy_init' )
463
464	! namelist MPDT does not exist yet, so set it here
465	! MPDT is the call frequency for microphysics in minutes (0 means every step)
466	MPDT = 0.
467
468	! set GMT outside of phy_init because phy_init may not be called on this
469	! process if, for example, it is a moving nest and if this part of the domain is not
470	! being initialized (not the leading edge).
471	CALL domain_setgmtetc( grid, start_of_simulation )
472
473	CALL set_tiles ( grid , grid%imask_nostag, ims, ime, jms, jme, ips, ipe, jps, jpe )
474
475	! Phy_init is not necessarily thread-safe; do not multi-thread this loop.
476	! The tiling is to handle the fact that we may be masking off part of the computation.
477	DO ij = 1, grid%num_tiles
478
479	CALL phy_init ( grid%id , config_flags, grid%DT, grid%RESTART, grid%em_znw, grid%em_znu, &
480	grid%p_top, grid%tsk, grid%RADT,grid%BLDT,grid%CUDT, MPDT, &
481	grid%rthcuten, grid%rqvcuten, grid%rqrcuten, &
482	grid%rqccuten, grid%rqscuten, grid%rqicuten, &
483	grid%rublten,grid%rvblten,grid%rthblten, &
484	grid%rqvblten,grid%rqcblten,grid%rqiblten, &
485	grid%rthraten,grid%rthratenlw,grid%rthratensw, &
486	grid%stepbl,grid%stepra,grid%stepcu, &
487	grid%w0avg, grid%rainnc, grid%rainc, grid%raincv, grid%rainncv, &
488	grid%nca,grid%swrad_scat, &
489	grid%cldefi,grid%lowlyr, &
490	grid%mass_flux, &
491	grid%rthften, grid%rqvften, &
492	grid%cldfra,grid%glw,grid%gsw,grid%emiss,grid%lu_index, &
493	grid%landuse_ISICE, grid%landuse_LUCATS, &
494	grid%landuse_LUSEAS, grid%landuse_ISN, &
495	grid%lu_state, &
496	grid%xlat,grid%xlong,grid%albedo,grid%albbck,grid%GMT,grid%JULYR,grid%JULDAY, &
497	grid%levsiz, num_ozmixm, num_aerosolc, grid%paerlev, &
498	grid%tmn,grid%xland,grid%znt,grid%z0,grid%ust,grid%mol,grid%pblh,grid%tke_myj, &
499	grid%exch_h,grid%thc,grid%snowc,grid%mavail,grid%hfx,grid%qfx,grid%rainbl, &
500	grid%tslb,grid%zs,grid%dzs,config_flags%num_soil_layers,grid%warm_rain, &
501	grid%adv_moist_cond, &
502	grid%apr_gr,grid%apr_w,grid%apr_mc,grid%apr_st,grid%apr_as, &
503	grid%apr_capma,grid%apr_capme,grid%apr_capmi, &
504	grid%xice,grid%vegfra,grid%snow,grid%canwat,grid%smstav, &
505	grid%smstot, grid%sfcrunoff,grid%udrunoff,grid%grdflx,grid%acsnow, &
506	grid%acsnom,grid%ivgtyp,grid%isltyp, grid%sfcevp,grid%smois, &
507	grid%sh2o, grid%snowh, grid%smfr3d, &
508	grid%DX,grid%DY,grid%f_ice_phy,grid%f_rain_phy,grid%f_rimef_phy, &
509	grid%mp_restart_state,grid%tbpvs_state,grid%tbpvs0_state,&
510	allowed_to_read, grid%moved, start_of_simulation, &
511	ids, ide, jds, jde, kds, kde, &
512	ims, ime, jms, jme, kms, kme, &
513	grid%i_start(ij), grid%i_end(ij), grid%j_start(ij), grid%j_end(ij), kts, kte, &
514	ozmixm,grid%pin, & ! Optional
515	grid%m_ps_1,grid%m_ps_2,grid%m_hybi,aerosolc_1,aerosolc_2,& ! Optional
516	grid%rundgdten,grid%rvndgdten,grid%rthndgdten, & ! Optional
517	grid%rqvndgdten,grid%rmundgdten, & ! Optional
518	grid%FGDT,grid%stepfg, & ! Optional
519	grid%DZR, grid%DZB, grid%DZG, & !Optional urban
520	grid%TR_URB2D,grid%TB_URB2D,grid%TG_URB2D,grid%TC_URB2D, & !Optional urban
521	grid%QC_URB2D, grid%XXXR_URB2D,grid%XXXB_URB2D, & !Optional urban
522	grid%XXXG_URB2D, grid%XXXC_URB2D, & !Optional urban
523	grid%TRL_URB3D, grid%TBL_URB3D, grid%TGL_URB3D, & !Optional urban
524	grid%SH_URB2D, grid%LH_URB2D, grid%G_URB2D, grid%RN_URB2D, & !Optional urban
525	grid%TS_URB2D, grid%FRC_URB2D, grid%UTYPE_URB2D, & !Optional urban
526	itimestep=grid%itimestep, &
527	fdob=grid%fdob &
528	)
529
530
531
532	ENDDO
533
534
535
536	CALL wrf_debug ( 100 , 'module_start: start_domain_rk: After call to phy_init' )
537
538	#ifdef MCELIO
539	LU_MASK = 0.
540	WHERE ( grid%lu_index .EQ. 16 ) LU_MASK = 1.
541	#endif
542
543	END IF
544
545	#if 0
546	#include "CYCLE_TEST.inc"
547	#endif
548
549	!
550	!
551
552	! set physical boundary conditions for all initialized variables
553
554	!-----------------------------------------------------------------------
555	! Stencils for patch communications (WCS, 29 June 2001)
556	! Note: the size of this halo exchange reflects the
557	! fact that we are carrying the uncoupled variables
558	! as state variables in the mass coordinate model, as
559	! opposed to the coupled variables as in the height
560	! coordinate model.
561	!
562	! * * * * *
563	! * * * * * * * * *
564	! * + * * + * * * + * *
565	! * * * * * * * * *
566	! * * * * *
567	!
568	!j grid%em_u_1 x
569	!j grid%em_u_2 x
570	!j grid%em_v_1 x
571	!j grid%em_v_2 x
572	!j grid%em_w_1 x
573	!j grid%em_w_2 x
574	!j grid%em_t_1 x
575	!j grid%em_t_2 x
576	!j grid%em_ph_1 x
577	!j grid%em_ph_2 x
578	!
579	!j grid%em_t_init x
580	!
581	!j grid%em_phb x
582	!j grid%em_ph0 x
583	!j grid%em_php x
584	!j grid%em_pb x
585	!j grid%em_al x
586	!j grid%em_alt x
587	!j grid%em_alb x
588	!
589	! the following are 2D (xy) variables
590	!
591	!j grid%em_mu_1 x
592	!j grid%em_mu_2 x
593	!j grid%em_mub x
594	!j grid%em_mu0 x
595	!j grid%ht x
596	!j grid%msft x
597	!j grid%msfu x
598	!j grid%msfv x
599	!j grid%sina x
600	!j grid%cosa x
601	!j grid%e x
602	!j grid%f x
603	!
604	! 4D variables
605	!
606	! moist x
607	! chem x
608	!scalar x
609
610	!--------------------------------------------------------------
611
612	#ifdef DM_PARALLEL
613	# include "HALO_EM_INIT_1.inc"
614	# include "HALO_EM_INIT_2.inc"
615	# include "HALO_EM_INIT_3.inc"
616	# include "HALO_EM_INIT_4.inc"
617	# include "HALO_EM_INIT_5.inc"
618	# include "PERIOD_BDY_EM_INIT.inc"
619	# include "PERIOD_BDY_EM_MOIST.inc"
620	# include "PERIOD_BDY_EM_CHEM.inc"
621	#endif
622
623
624	CALL set_physical_bc3d( grid%em_u_1 , 'U' , config_flags , &
625	ids , ide , jds , jde , kds , kde , &
626	ims , ime , jms , jme , kms , kme , &
627	its , ite , jts , jte , kts , kte , &
628	its , ite , jts , jte , kts , kte )
629	CALL set_physical_bc3d( grid%em_u_2 , 'U' , config_flags , &
630	ids , ide , jds , jde , kds , kde , &
631	ims , ime , jms , jme , kms , kme , &
632	its , ite , jts , jte , kts , kte , &
633	its , ite , jts , jte , kts , kte )
634
635	CALL set_physical_bc3d( grid%em_v_1 , 'V' , config_flags , &
636	ids , ide , jds , jde , kds , kde , &
637	ims , ime , jms , jme , kms , kme , &
638	its , ite , jts , jte , kts , kte , &
639	its , ite , jts , jte , kts , kte )
640	CALL set_physical_bc3d( grid%em_v_2 , 'V' , config_flags , &
641	ids , ide , jds , jde , kds , kde , &
642	ims , ime , jms , jme , kms , kme , &
643	its , ite , jts , jte , kts , kte , &
644	its , ite , jts , jte , kts , kte )
645
646	! set kinematic condition for w
647
648	CALL set_physical_bc2d( grid%ht , 'r' , config_flags , &
649	ids , ide , jds , jde , &
650	ims , ime , jms , jme , &
651	its , ite , jts , jte , &
652	its , ite , jts , jte )
653
654	IF ( .not. config_flags%restart ) THEN
655	CALL set_w_surface( config_flags, &
656	grid%em_w_1, grid%ht, grid%em_u_1, grid%em_v_1, grid%cf1, &
657	grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msft, &
658	ids, ide, jds, jde, kds, kde, &
659	ips, ipe, jps, jpe, kps, kpe, &
660	its, ite, jts, jte, kts, kte, &
661	ims, ime, jms, jme, kms, kme )
662	CALL set_w_surface( config_flags, &
663	grid%em_w_2, grid%ht, grid%em_u_2, grid%em_v_2, grid%cf1, &
664	grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msft, &
665	ids, ide, jds, jde, kds, kde, &
666	ips, ipe, jps, jpe, kps, kpe, &
667	its, ite, jts, jte, kts, kte, &
668	ims, ime, jms, jme, kms, kme )
669	END IF
670
671	! finished setting kinematic condition for w at the surface
672
673	CALL set_physical_bc3d( grid%em_w_1 , 'W' , config_flags , &
674	ids , ide , jds , jde , kds , kde , &
675	ims , ime , jms , jme , kms , kme , &
676	its , ite , jts , jte , kts , kte , &
677	its , ite , jts , jte , kts , kte )
678	CALL set_physical_bc3d( grid%em_w_2 , 'W' , config_flags , &
679	ids , ide , jds , jde , kds , kde , &
680	ims , ime , jms , jme , kms , kme , &
681	its , ite , jts , jte , kts , kte , &
682	its , ite , jts , jte , kts , kte )
683
684	CALL set_physical_bc3d( grid%em_ph_1 , 'W' , config_flags , &
685	ids , ide , jds , jde , kds , kde , &
686	ims , ime , jms , jme , kms , kme , &
687	its , ite , jts , jte , kts , kte , &
688	its , ite , jts , jte , kts , kte )
689
690	CALL set_physical_bc3d( grid%em_ph_2 , 'W' , config_flags , &
691	ids , ide , jds , jde , kds , kde , &
692	ims , ime , jms , jme , kms , kme , &
693	its , ite , jts , jte , kts , kte , &
694	its , ite , jts , jte , kts , kte )
695
696	CALL set_physical_bc3d( grid%em_t_1 , 't' , config_flags , &
697	ids , ide , jds , jde , kds , kde , &
698	ims , ime , jms , jme , kms , kme , &
699	its , ite , jts , jte , kts , kte , &
700	its , ite , jts , jte , kts , kte )
701
702	CALL set_physical_bc3d( grid%em_t_2 , 't' , config_flags , &
703	ids , ide , jds , jde , kds , kde , &
704	ims , ime , jms , jme , kms , kme , &
705	its , ite , jts , jte , kts , kte , &
706	its , ite , jts , jte , kts , kte )
707
708	CALL set_physical_bc2d( grid%em_mu_1, 't' , config_flags , &
709	ids , ide , jds , jde , &
710	ims , ime , jms , jme , &
711	its , ite , jts , jte , &
712	its , ite , jts , jte )
713	CALL set_physical_bc2d( grid%em_mu_2, 't' , config_flags , &
714	ids , ide , jds , jde , &
715	ims , ime , jms , jme , &
716	its , ite , jts , jte , &
717	its , ite , jts , jte )
718	CALL set_physical_bc2d( grid%em_mub , 't' , config_flags , &
719	ids , ide , jds , jde , &
720	ims , ime , jms , jme , &
721	its , ite , jts , jte , &
722	its , ite , jts , jte )
723	CALL set_physical_bc2d( grid%em_mu0 , 't' , config_flags , &
724	ids , ide , jds , jde , &
725	ims , ime , jms , jme , &
726	its , ite , jts , jte , &
727	its , ite , jts , jte )
728
729
730	CALL set_physical_bc3d( grid%em_phb , 'W' , config_flags , &
731	ids , ide , jds , jde , kds , kde , &
732	ims , ime , jms , jme , kms , kme , &
733	its , ite , jts , jte , kts , kte , &
734	its , ite , jts , jte , kts , kte )
735	CALL set_physical_bc3d( grid%em_ph0 , 'W' , config_flags , &
736	ids , ide , jds , jde , kds , kde , &
737	ims , ime , jms , jme , kms , kme , &
738	its , ite , jts , jte , kts , kte , &
739	its , ite , jts , jte , kts , kte )
740	CALL set_physical_bc3d( grid%em_php , 'W' , config_flags , &
741	ids , ide , jds , jde , kds , kde , &
742	ims , ime , jms , jme , kms , kme , &
743	its , ite , jts , jte , kts , kte , &
744	its , ite , jts , jte , kts , kte )
745
746	CALL set_physical_bc3d( grid%em_pb , 't' , config_flags , &
747	ids , ide , jds , jde , kds , kde , &
748	ims , ime , jms , jme , kms , kme , &
749	its , ite , jts , jte , kts , kte , &
750	its , ite , jts , jte , kts , kte )
751	CALL set_physical_bc3d( grid%em_al , 't' , config_flags , &
752	ids , ide , jds , jde , kds , kde , &
753	ims , ime , jms , jme , kms , kme , &
754	its , ite , jts , jte , kts , kte , &
755	its , ite , jts , jte , kts , kte )
756	CALL set_physical_bc3d( grid%em_alt , 't' , config_flags , &
757	ids , ide , jds , jde , kds , kde , &
758	ims , ime , jms , jme , kms , kme , &
759	its , ite , jts , jte , kts , kte , &
760	its , ite , jts , jte , kts , kte )
761	CALL set_physical_bc3d( grid%em_alb , 't' , config_flags , &
762	ids , ide , jds , jde , kds , kde , &
763	ims , ime , jms , jme , kms , kme , &
764	its , ite , jts , jte , kts , kte , &
765	its , ite , jts , jte , kts , kte )
766	CALL set_physical_bc3d(grid%em_t_init, 't' , config_flags , &
767	ids , ide , jds , jde , kds , kde , &
768	ims , ime , jms , jme , kms , kme , &
769	its , ite , jts , jte , kts , kte , &
770	its , ite , jts , jte , kts , kte )
771
772	IF (num_moist > 0) THEN
773
774	! use of (:,:,:,loop) not efficient on DEC, but (ims,kms,jms,loop) not portable to SGI/Cray
775
776	loop_3d_m : DO loop = 1 , num_moist
777	CALL set_physical_bc3d( moist(:,:,:,loop) , 'r' , config_flags , &
778	ids , ide , jds , jde , kds , kde , &
779	ims , ime , jms , jme , kms , kme , &
780	its , ite , jts , jte , kts , kte , &
781	its , ite , jts , jte , kts , kte )
782	END DO loop_3d_m
783
784	ENDIF
785
786	#ifdef WRF_CHEM
787	! Zero out the Mie arrays so they do not cause instabilities when not calling Fast-J.
788	if( .NOT. config_flags%restart ) then
789	do j=jts,jte
790	do k=kts,kte
791	do i=its,ite
792	grid%tauaer1(i,k,j) = 0.
793	grid%tauaer2(i,k,j) = 0.
794	grid%tauaer3(i,k,j) = 0.
795	grid%tauaer4(i,k,j) = 0.
796	grid%gaer1(i,k,j) = 0.
797	grid%gaer2(i,k,j) = 0.
798	grid%gaer3(i,k,j) = 0.
799	grid%gaer4(i,k,j) = 0.
800	grid%waer1(i,k,j) = 0.
801	grid%waer2(i,k,j) = 0.
802	grid%waer3(i,k,j) = 0.
803	grid%waer4(i,k,j) = 0.
804	end do
805	end do
806	end do
807	end if
808	#endif
809
810	IF (num_chem >= PARAM_FIRST_SCALAR ) THEN
811	#ifdef WRF_CHEM
812	!
813	! we do this here, so we only have one chem_init routine for either core....
814	!
815	do j=jts,min(jte,jde-1)
816	do i=its,min(ite,ide-1)
817	do k=kts,kte
818	z_at_w(i,k,j)=(grid%em_ph_2(i,k,j)+grid%em_phb(i,k,j))/g
819	enddo
820	do k=kts,min(kte,kde-1)
821	tempfac=(grid%em_t_1(i,k,j) + t0)((grid%em_p(i,k,j) + grid%em_pb(i,k,j))/p1000mb)*rcp
822	convfac(i,k,j) = (grid%em_p(i,k,j)+grid%em_pb(i,k,j))/rgasuniv/tempfac
823	enddo
824	enddo
825	enddo
826
827	CALL chem_init (grid%id,chem,grid%dt,grid%bioemdt,grid%photdt, &
828	grid%chemdt, &
829	grid%stepbioe,grid%stepphot,grid%stepchem, &
830	z_at_w,g,grid%aerwrf,config_flags, &
831	grid%em_alt,grid%em_t_1,grid%em_p,convfac, &
832	grid%tauaer1,grid%tauaer2,grid%tauaer3,grid%tauaer4, &
833	grid%gaer1,grid%gaer2,grid%gaer3,grid%gaer4, &
834	grid%waer1,grid%waer2,grid%waer3,grid%waer4, &
835	grid%pm2_5_dry,grid%pm2_5_water,grid%pm2_5_dry_ec, &
836	grid%chem_in_opt, &
837	ids , ide , jds , jde , kds , kde , &
838	ims , ime , jms , jme , kms , kme , &
839	its , ite , jts , jte , kts , kte )
840
841	!
842	! calculate initial pm
843	!
844	! print *,'calculating initial pm'
845	select case (config_flags%chem_opt)
846	case (RADM2SORG, RACMSORG)
847	call sum_pm_sorgam ( &
848	grid%em_alt, chem, grid%h2oaj, grid%h2oai, &
849	grid%pm2_5_dry, grid%pm2_5_water, grid%pm2_5_dry_ec, grid%pm10, &
850	ids,ide, jds,jde, kds,kde, &
851	ims,ime, jms,jme, kms,kme, &
852	its,ite, jts,jte, kts,kte )
853
854	case (CBMZ_MOSAIC_AA, CBMZ_MOSAIC_BB)
855	call sum_pm_mosaic ( &
856	grid%em_alt, chem, &
857	grid%pm2_5_dry, grid%pm2_5_water, grid%pm2_5_dry_ec, grid%pm10, &
858	ids,ide, jds,jde, kds,kde, &
859	ims,ime, jms,jme, kms,kme, &
860	its,ite, jts,jte, kts,kte )
861
862	case default
863	do j=jts,min(jte,jde-1)
864	do k=kts,min(kte,kde-1)
865	do i=its,min(ite,ide-1)
866	grid%pm2_5_dry(i,k,j) = 0.
867	grid%pm2_5_water(i,k,j) = 0.
868	grid%pm2_5_dry_ec(i,k,j) = 0.
869	grid%pm10(i,k,j) = 0.
870	enddo
871	enddo
872	enddo
873	end select
874	#endif
875
876	! use of (:,:,:,loop) not efficient on DEC, but (ims,kms,jms,loop) not portable to SGI/Cray
877
878	loop_3d_c : DO loop = PARAM_FIRST_SCALAR , num_chem
879	CALL set_physical_bc3d( chem(:,:,:,loop) , 'r' , config_flags , &
880	ids , ide , jds , jde , kds , kde , &
881	ims , ime , jms , jme , kms , kme , &
882	its , ite , jts , jte , kts , kte , &
883	its , ite , jts , jte , kts , kte )
884	END DO loop_3d_c
885
886	ENDIF
887
888	CALL set_physical_bc2d( grid%msft , 'r' , config_flags , &
889	ids , ide , jds , jde , &
890	ims , ime , jms , jme , &
891	its , ite , jts , jte , &
892	its , ite , jts , jte )
893	CALL set_physical_bc2d( grid%msfu , 'x' , config_flags , &
894	ids , ide , jds , jde , &
895	ims , ime , jms , jme , &
896	its , ite , jts , jte , &
897	its , ite , jts , jte )
898	CALL set_physical_bc2d( grid%msfv , 'y' , config_flags , &
899	ids , ide , jds , jde , &
900	ims , ime , jms , jme , &
901	its , ite , jts , jte , &
902	its , ite , jts , jte )
903	CALL set_physical_bc2d( grid%sina , 'r' , config_flags , &
904	ids , ide , jds , jde , &
905	ims , ime , jms , jme , &
906	its , ite , jts , jte , &
907	its , ite , jts , jte )
908	CALL set_physical_bc2d( grid%cosa , 'r' , config_flags , &
909	ids , ide , jds , jde , &
910	ims , ime , jms , jme , &
911	its , ite , jts , jte , &
912	its , ite , jts , jte )
913	CALL set_physical_bc2d( grid%e , 'r' , config_flags , &
914	ids , ide , jds , jde , &
915	ims , ime , jms , jme , &
916	its , ite , jts , jte , &
917	its , ite , jts , jte )
918	CALL set_physical_bc2d( grid%f , 'r' , config_flags , &
919	ids , ide , jds , jde , &
920	ims , ime , jms , jme , &
921	its , ite , jts , jte , &
922	its , ite , jts , jte )
923
924	#ifdef DM_PARALLEL
925	# include "HALO_EM_INIT_1.inc"
926	# include "HALO_EM_INIT_2.inc"
927	# include "HALO_EM_INIT_3.inc"
928	# include "HALO_EM_INIT_4.inc"
929	# include "HALO_EM_INIT_5.inc"
930	# include "PERIOD_BDY_EM_INIT.inc"
931	# include "PERIOD_BDY_EM_MOIST.inc"
932	# include "PERIOD_BDY_EM_CHEM.inc"
933	#endif
934
935	CALL wrf_debug ( 100 , 'module_start: start_domain_rk: Returning' )
936
937	RETURN
938
939	END SUBROUTINE start_domain_em
940

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

Download in other formats:

Original Format