Context Navigation

interface_icosa_lmdz.f90 @ 4855

Last change on this file since 4855 was 4855, checked in by yann meurdesoif, 2 months ago
move dynamic etat0 initialisation into driver via new etat0 plugin functionnality. move computation of parameter for gravity wave in driver add method of computation of gravity wave parameters similar of what is done currently by lmdz (legacy method) add switch between legacy method (old) sso method in .def file ( param_gw_method = 'legacy' / 'sso' ) add corresponding XML files
File size: 21.9 KB

Line
1	MODULE interface_icosa_lmdz_mod
2
3	USE field_mod, ONLY: t_field
4	USE transfert_mod, ONLY: t_message
5
6
7	TYPE(t_message),SAVE :: req_u
8	TYPE(t_message),SAVE :: req_dps0, req_dulon0, req_dulat0, req_dTemp0, req_dq0
9
10	TYPE(t_field),POINTER,SAVE :: f_p(:)
11	TYPE(t_field),POINTER,SAVE :: f_pks(:)
12	TYPE(t_field),POINTER,SAVE :: f_pk(:)
13	TYPE(t_field),POINTER,SAVE :: f_p_layer(:)
14	TYPE(t_field),POINTER,SAVE :: f_theta(:)
15	TYPE(t_field),POINTER,SAVE :: f_phi(:)
16	TYPE(t_field),POINTER,SAVE :: f_Temp(:)
17	TYPE(t_field),POINTER,SAVE :: f_ulon(:)
18	TYPE(t_field),POINTER,SAVE :: f_ulat(:)
19	TYPE(t_field),POINTER,SAVE :: f_dulon(:)
20	TYPE(t_field),POINTER,SAVE :: f_dulat(:)
21	TYPE(t_field),POINTER,SAVE :: f_dTemp(:)
22	TYPE(t_field),POINTER,SAVE :: f_dq(:)
23	TYPE(t_field),POINTER,SAVE :: f_dps(:)
24	TYPE(t_field),POINTER,SAVE :: f_duc(:)
25	TYPE(t_field),POINTER,SAVE :: f_bounds_lon(:)
26	TYPE(t_field),POINTER,SAVE :: f_bounds_lat(:)
27
28	INTEGER :: start_clock
29	INTEGER :: stop_clock
30	INTEGER :: count_clock=0
31
32	REAL :: day_length
33
34
35	INTEGER,SAVE :: nbp_phys
36	INTEGER,SAVE :: nbp_phys_glo
37
38
39	CONTAINS
40
41	SUBROUTINE initialize_physics
42	USE distrib_icosa_lmdz_mod, ONLY : init_distrib_icosa_lmdz, transfer_icosa_to_lmdz
43	! from dynamico
44	USE domain_mod
45	USE dimensions
46	USE mpi_mod
47	USE mpipara
48	USE disvert_mod
49	USE xios_mod
50	USE transfert_mod
51
52	! from LMDZ
53	USE mod_grid_phy_lmdz, ONLY : unstructured
54	USE mod_phys_lmdz_omp_data, ONLY: klon_omp
55	USE transfert_mod
56	USE physics_distribution_mod, ONLY : init_physics_distribution
57
58
59	IMPLICIT NONE
60	INTEGER :: ind,i,j,ij,pos
61	REAL(rstd),POINTER :: bounds_lon(:,:)
62	REAL(rstd),POINTER :: bounds_lat(:,:)
63
64	REAL(rstd),ALLOCATABLE :: latfi(:)
65	REAL(rstd),ALLOCATABLE :: lonfi(:)
66	REAL(rstd),ALLOCATABLE :: airefi(:)
67	REAL(rstd),ALLOCATABLE :: bounds_latfi(:,:)
68	REAL(rstd),ALLOCATABLE :: bounds_lonfi(:,:)
69
70	INTEGER :: nbp_phys, nbp_phys_glo
71
72	!$OMP PARALLEL
73	CALL allocate_field(f_bounds_lon,field_t,type_real,6)
74	CALL allocate_field(f_bounds_lat,field_t,type_real,6)
75	CALL allocate_field(f_p,field_t,type_real,llm+1,name="p_in")
76	CALL allocate_field(f_pks,field_t,type_real)
77	CALL allocate_field(f_pk,field_t,type_real,llm)
78	CALL allocate_field(f_p_layer,field_t,type_real,llm,name="p_layer_in")
79	CALL allocate_field(f_theta,field_t,type_real,llm)
80	CALL allocate_field(f_phi,field_t,type_real,llm,name="phi_in")
81	CALL allocate_field(f_Temp,field_t,type_real,llm,name="Temp_in")
82	CALL allocate_field(f_ulon,field_t,type_real,llm,name="ulon_in")
83	CALL allocate_field(f_ulat,field_t,type_real,llm,name="ulat_in")
84	CALL allocate_field(f_dulon,field_t,type_real,llm,name="dulon_out")
85	CALL allocate_field(f_dulat,field_t,type_real,llm,name="dulat_out")
86	CALL allocate_field(f_dTemp,field_t,type_real,llm,name="dTemp_out")
87	CALL allocate_field(f_dq,field_t,type_real,llm,nqtot,name="dq_out")
88	CALL allocate_field(f_dps,field_t,type_real,name="dps_out")
89	CALL allocate_field(f_duc,field_t,type_real,3,llm)
90
91	CALL init_message(f_dps,req_i0,req_dps0)
92	CALL init_message(f_dulon,req_i0,req_dulon0)
93	CALL init_message(f_dulat,req_i0,req_dulat0)
94	CALL init_message(f_dTemp,req_i0,req_dTemp0)
95	CALL init_message(f_dq,req_i0,req_dq0)
96	!$OMP END PARALLEL
97
98	nbp_phys=0
99	DO ind=1,ndomain
100	CALL swap_dimensions(ind)
101	DO j=jj_begin,jj_end
102	DO i=ii_begin,ii_end
103	IF (domain(ind)%own(i,j)) nbp_phys=nbp_phys+1
104	ENDDO
105	ENDDO
106	ENDDO
107
108
109	!initialize LMDZ5 physic mpi decomposition
110	CALL MPI_ALLREDUCE(nbp_phys,nbp_phys_glo,1,MPI_INTEGER,MPI_SUM,comm_icosa,ierr)
111	CALL init_physics_distribution(unstructured, 6, nbp_phys, 1, nbp_phys_glo, llm, comm_icosa)
112
113	DO ind=1,ndomain
114	CALL swap_dimensions(ind)
115	CALL swap_geometry(ind)
116	bounds_lon=f_bounds_lon(ind)
117	bounds_lat=f_bounds_lat(ind)
118	DO j=jj_begin,jj_end
119	DO i=ii_begin,ii_end
120	ij=(j-1)*iim+i
121	CALL xyz2lonlat(xyz_v(ij+z_rup,:), bounds_lon(ij,1), bounds_lat(ij,1))
122	CALL xyz2lonlat(xyz_v(ij+z_up,:), bounds_lon(ij,2), bounds_lat(ij,2))
123	CALL xyz2lonlat(xyz_v(ij+z_lup,:), bounds_lon(ij,3), bounds_lat(ij,3))
124	CALL xyz2lonlat(xyz_v(ij+z_ldown,:), bounds_lon(ij,4), bounds_lat(ij,4))
125	CALL xyz2lonlat(xyz_v(ij+z_down,:), bounds_lon(ij,5), bounds_lat(ij,5))
126	CALL xyz2lonlat(xyz_v(ij+z_rdown,:), bounds_lon(ij,6), bounds_lat(ij,6))
127	ENDDO
128	ENDDO
129	ENDDO
130
131	!$OMP PARALLEL
132	CALL initialize_physics_omp
133	CALL param_gravity_wave
134	!$OMP END PARALLEL
135
136	CALL xios_set_context
137
138
139
140
141	END SUBROUTINE initialize_physics
142
143
144	SUBROUTINE initialize_physics_omp
145	USE distrib_icosa_lmdz_mod, ONLY : init_distrib_icosa_lmdz, transfer_icosa_to_lmdz
146	! from dynamico
147	USE domain_mod
148	USE dimensions
149	USE mpi_mod
150	USE mpipara
151	USE disvert_mod
152	USE earth_const, ONLY: scale_height
153	USE xios_mod
154	USE time_mod , init_time_icosa=> init_time
155	USE omp_para
156
157	! from LMDZ
158	USE mod_grid_phy_lmdz, ONLY : unstructured
159	USE mod_phys_lmdz_omp_data, ONLY: klon_omp
160	USE transfert_mod
161	USE geometry_mod, ONLY : init_geometry
162	USE vertical_layers_mod, ONLY : init_vertical_layers
163	USE infotrac_phy, ONLY : init_infotrac_phy
164	USE inifis_mod, ONLY : inifis
165	USE phyaqua_mod, ONLY : iniaqua
166
167
168	IMPLICIT NONE
169
170
171
172	INTEGER :: ind,i,j,ij,pos
173	REAL(rstd),POINTER :: bounds_lon(:,:)
174	REAL(rstd),POINTER :: bounds_lat(:,:)
175
176	REAL(rstd),ALLOCATABLE :: latfi(:)
177	REAL(rstd),ALLOCATABLE :: lonfi(:)
178	REAL(rstd),ALLOCATABLE :: airefi(:)
179	REAL(rstd),ALLOCATABLE :: bounds_latfi(:,:)
180	REAL(rstd),ALLOCATABLE :: bounds_lonfi(:,:)
181	REAL(rstd),ALLOCATABLE :: ind_cell_glo(:)
182
183	REAL(rstd) :: pseudoalt(llm)
184	REAL(rstd) :: aps(llm)
185	REAL(rstd) :: bps(llm)
186	REAL(rstd) :: scaleheight
187
188	CHARACTER(len=4) :: type_trac
189
190	TYPE(t_field),POINTER,SAVE :: f_ind_cell_glo(:)
191
192	INTEGER :: iflag_phys
193
194	CALL init_distrib_icosa_lmdz
195
196	ALLOCATE(latfi(klon_omp))
197	ALLOCATE(lonfi(klon_omp))
198	ALLOCATE(airefi(klon_omp))
199	ALLOCATE(bounds_latfi(klon_omp,6))
200	ALLOCATE(bounds_lonfi(klon_omp,6))
201	ALLOCATE(ind_cell_glo(klon_omp))
202
203	CALL transfer_icosa_to_lmdz(geom%lat_i,latfi)
204	CALL transfer_icosa_to_lmdz(geom%lon_i,lonfi)
205	CALL transfer_icosa_to_lmdz(f_bounds_lat,bounds_latfi)
206	CALL transfer_icosa_to_lmdz(f_bounds_lon,bounds_lonfi)
207	CALL transfer_icosa_to_lmdz(geom%Ai,airefi)
208
209	CALL allocate_field(f_ind_cell_glo,field_t,type_real)
210
211	DO ind=1,ndomain
212	IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master ) CYCLE
213	CALL swap_dimensions(ind)
214	CALL swap_geometry(ind)
215	DO j=jj_begin,jj_end
216	DO i=ii_begin,ii_end
217	ij=(j-1)*iim+i
218	f_ind_cell_glo(ind)%rval2d(ij)=domain(ind)%assign_cell_glo(i,j)
219	ENDDO
220	ENDDO
221	ENDDO
222
223	!$OMP BARRIER
224
225	CALL transfer_icosa_to_lmdz(f_ind_cell_glo,ind_cell_glo)
226	CALL deallocate_field(f_ind_cell_glo)
227
228
229	CALL init_geometry(klon_omp,lonfi, latfi, bounds_lonfi, bounds_latfi, airefi, INT(ind_cell_glo))
230
231	scaleheight=scale_height/1000. ! Atmospheric scale height (km)
232	aps(1:llm)=0.5*(ap(1:llm)+ap(2:llm+1))
233	bps(1:llm)=0.5*(bp(1:llm)+bp(2:llm+1))
234	pseudoalt(:)=-scaleheight*log(presnivs(:)/preff)
235	CALL init_vertical_layers(llm,preff,scaleheight,ap,bp,aps,bps,presnivs,pseudoalt)
236
237
238	! Initialize physical constant
239	CALL inifis(radius,g,kappa*cpp,cpp)
240
241	! Initialize tracer names, numbers, etc. for physics
242
243	!Config Key = type_trac
244	!Config Desc = Choix de couplage avec model de chimie INCA ou REPROBUS
245	!Config Def = lmdz
246	!Config Help =
247	!Config 'lmdz' = pas de couplage, pur LMDZ
248	!Config 'inca' = model de chime INCA
249	!Config 'repr' = model de chime REPROBUS
250	type_trac = 'lmdz'
251	CALL getin('type_trac',type_trac)
252
253	CALL init_infotrac_phy(nqtot,type_trac)
254
255	! Additional initializations for aquaplanets
256	CALL getin("iflag_phys",iflag_phys)
257	IF (iflag_phys>=100) THEN
258	CALL iniaqua(klon_omp, iflag_phys)
259	END IF
260
261
262	END SUBROUTINE initialize_physics_omp
263
264
265
266
267	SUBROUTINE physics
268	USE icosa
269	USE time_mod
270	USE disvert_mod
271	USE transfert_mod
272	USE mpipara
273	USE xios_mod
274	USE wxios
275	USE trace
276	USE distrib_icosa_lmdz_mod, ONLY : transfer_icosa_to_lmdz, transfer_lmdz_to_icosa
277	USE physics_external_mod, ONLY : it, f_phis, f_ps, f_theta_rhodz, f_u, f_wflux, f_q
278	USE write_field_mod
279	USE checksum_mod
280	! from LMDZ
281	USE mod_phys_lmdz_omp_data, ONLY: klon_omp
282	USE geometry_mod, ONLY : cell_area
283	USE physiq_mod, ONLY: physiq
284	IMPLICIT NONE
285
286	REAL(rstd),POINTER :: phis(:)
287	REAL(rstd),POINTER :: ps(:)
288	REAL(rstd),POINTER :: theta_rhodz(:,:,:)
289	REAL(rstd),POINTER :: u(:,:)
290	REAL(rstd),POINTER :: wflux(:,:)
291	REAL(rstd),POINTER :: q(:,:,:)
292	REAL(rstd),POINTER :: p(:,:)
293	REAL(rstd),POINTER :: pks(:)
294	REAL(rstd),POINTER :: pk(:,:)
295	REAL(rstd),POINTER :: p_layer(:,:)
296	REAL(rstd),POINTER :: theta(:,:)
297	REAL(rstd),POINTER :: phi(:,:)
298	REAL(rstd),POINTER :: Temp(:,:)
299	REAL(rstd),POINTER :: ulon(:,:)
300	REAL(rstd),POINTER :: ulat(:,:)
301	REAL(rstd),POINTER :: dulon(:,:)
302	REAL(rstd),POINTER :: dulat(:,:)
303	REAL(rstd),POINTER :: dTemp(:,:)
304	REAL(rstd),POINTER :: dq(:,:,:)
305	REAL(rstd),POINTER :: dps(:)
306	REAL(rstd),POINTER :: duc(:,:,:)
307
308
309	INTEGER :: ind,l
310
311	REAL(rstd),ALLOCATABLE,SAVE :: ps_phy(:)
312	!$OMP THREADPRIVATE(ps_phy)
313	REAL(rstd),ALLOCATABLE,SAVE :: p_phy(:,:)
314	!$OMP THREADPRIVATE(p_phy)
315	REAL(rstd),ALLOCATABLE,SAVE :: p_layer_phy(:,:)
316	!$OMP THREADPRIVATE(p_layer_phy)
317	REAL(rstd),ALLOCATABLE,SAVE :: Temp_phy(:,:)
318	!$OMP THREADPRIVATE(Temp_phy)
319	REAL(rstd),ALLOCATABLE,SAVE :: phis_phy(:)
320	!$OMP THREADPRIVATE(phis_phy)
321	REAL(rstd),ALLOCATABLE,SAVE :: phi_phy(:,:)
322	!$OMP THREADPRIVATE(phi_phy)
323	REAL(rstd),ALLOCATABLE,SAVE :: ulon_phy(:,:)
324	!$OMP THREADPRIVATE(ulon_phy)
325	REAL(rstd),ALLOCATABLE,SAVE :: ulat_phy(:,:)
326	!$OMP THREADPRIVATE(ulat_phy)
327	REAL(rstd),ALLOCATABLE,SAVE :: q_phy(:,:,:)
328	!$OMP THREADPRIVATE(q_phy)
329	REAL(rstd),ALLOCATABLE,SAVE :: wflux_phy(:,:)
330	!$OMP THREADPRIVATE(wflux_phy)
331	REAL(rstd),ALLOCATABLE,SAVE :: dulon_phy(:,:)
332	!$OMP THREADPRIVATE(dulon_phy)
333	REAL(rstd),ALLOCATABLE,SAVE :: dulat_phy(:,:)
334	!$OMP THREADPRIVATE(dulat_phy)
335	REAL(rstd),ALLOCATABLE,SAVE :: dTemp_phy(:,:)
336	!$OMP THREADPRIVATE(dTemp_phy)
337	REAL(rstd),ALLOCATABLE,SAVE :: dq_phy(:,:,:)
338	!$OMP THREADPRIVATE(dq_phy)
339	REAL(rstd),ALLOCATABLE,SAVE :: dps_phy(:)
340	!$OMP THREADPRIVATE(dps_phy)
341	REAL(rstd) :: dtphy
342	LOGICAL :: debut
343	LOGICAL :: lafin
344	LOGICAL,SAVE :: first=.TRUE.
345	!$OMP THREADPRIVATE(first)
346
347
348	IF(first) THEN
349	debut=.TRUE.
350	ELSE
351	debut=.FALSE.
352	ENDIF
353
354
355	IF(it-itau0>=itaumax) THEN
356	lafin=.TRUE.
357	ELSE
358	lafin=.FALSE.
359	ENDIF
360
361	IF (first) THEN
362	first=.FALSE.
363	CALL init_message(f_u,req_e1_vect,req_u)
364	ALLOCATE(ps_phy(klon_omp))
365	ALLOCATE(p_phy(klon_omp,llm+1))
366	ALLOCATE(p_layer_phy(klon_omp,llm))
367	ALLOCATE(Temp_phy(klon_omp,llm))
368	ALLOCATE(phis_phy(klon_omp))
369	ALLOCATE(phi_phy(klon_omp,llm))
370	ALLOCATE(ulon_phy(klon_omp,llm))
371	ALLOCATE(ulat_phy(klon_omp,llm))
372	ALLOCATE(q_phy(klon_omp,llm,nqtot))
373	ALLOCATE(wflux_phy(klon_omp,llm))
374	ALLOCATE(dulon_phy(klon_omp,llm))
375	ALLOCATE(dulat_phy(klon_omp,llm))
376	ALLOCATE(dTemp_phy(klon_omp,llm))
377	ALLOCATE(dq_phy(klon_omp,llm,nqtot))
378	ALLOCATE(dps_phy(klon_omp))
379	!$OMP BARRIER
380	ENDIF
381
382
383	!$OMP MASTER
384	! CALL update_calendar(it)
385	!$OMP END MASTER
386	!$OMP BARRIER
387	dtphy=itau_physics*dt
388
389
390
391	CALL transfert_message(f_u,req_u)
392
393	DO ind=1,ndomain
394	CALL swap_dimensions(ind)
395	IF (assigned_domain(ind)) THEN
396	CALL swap_geometry(ind)
397
398	phis=f_phis(ind)
399	ps=f_ps(ind)
400	theta_rhodz=f_theta_rhodz(ind)
401	u=f_u(ind)
402	q=f_q(ind)
403	wflux=f_wflux(ind)
404	p=f_p(ind)
405	pks=f_pks(ind)
406	pk=f_pk(ind)
407	p_layer=f_p_layer(ind)
408	theta=f_theta(ind)
409	phi=f_phi(ind)
410	Temp=f_Temp(ind)
411	ulon=f_ulon(ind)
412	ulat=f_ulat(ind)
413
414	CALL grid_icosa_to_physics
415
416	ENDIF
417	ENDDO
418
419	!$OMP BARRIER
420	!$OMP MASTER
421	CALL SYSTEM_CLOCK(start_clock)
422	!$OMP END MASTER
423	CALL trace_start("physic")
424	! CALL trace_off()
425
426
427	! CALL writeField("p_in",f_p)
428	! CALL writeField("p_layer_in",f_p_layer)
429	! CALL writeField("phi_in",f_phi)
430	! CALL writeField("phis_in",f_phis)
431	! CALL writeField("ulon_in",f_ulon)
432	! CALL writeField("ulat_in",f_ulat)
433	! CALL writeField("Temp_in",f_Temp)
434	! CALL writeField("q_in",f_q)
435	! CALL writeField("wflux_in",f_wflux)
436
437	! CALL checksum(f_p)
438	! CALL checksum(f_p_layer)
439	! CALL checksum(f_phi)
440	! CALL checksum(f_phis)
441	! CALL checksum(f_ulon)
442	! CALL checksum(f_ulat)
443	! CALL checksum(f_Temp)
444	! CALL checksum(f_q)
445	! CALL checksum(f_wflux)
446
447	CALL transfer_icosa_to_lmdz(f_p , p_phy)
448	CALL transfer_icosa_to_lmdz(f_p_layer, p_layer_phy)
449	CALL transfer_icosa_to_lmdz(f_phi , phi_phy)
450	CALL transfer_icosa_to_lmdz(f_phis , phis_phy )
451	CALL transfer_icosa_to_lmdz(f_ulon , ulon_phy )
452	CALL transfer_icosa_to_lmdz(f_ulat , ulat_phy)
453	CALL transfer_icosa_to_lmdz(f_Temp , Temp_phy)
454	CALL transfer_icosa_to_lmdz(f_q , q_phy)
455	CALL transfer_icosa_to_lmdz(f_wflux , wflux_phy)
456
457	DO l=1,llm
458	wflux_phy(:,l)=wflux_phy(:,l)*cell_area(:)
459	phi_phy(:,l)=phi_phy(:,l)-phis_phy(:)
460	ENDDO
461
462	CALL wxios_set_context()
463
464	CALL physiq(klon_omp, llm, debut, lafin, dtphy, &
465	p_phy, p_layer_phy, phi_phy, phis_phy, presnivs, &
466	ulon_phy, ulat_phy, Temp_phy, q_phy, wflux_phy, &
467	dulon_phy, dulat_phy, dTemp_phy, dq_phy, dps_phy)
468
469	CALL transfer_lmdz_to_icosa(dulon_phy, f_dulon )
470	CALL transfer_lmdz_to_icosa(dulat_phy, f_dulat )
471	CALL transfer_lmdz_to_icosa(dTemp_phy, f_dTemp )
472	CALL transfer_lmdz_to_icosa(dq_phy , f_dq )
473	CALL transfer_lmdz_to_icosa(dps_phy , f_dps )
474
475	! CALL writeField("dulon_out",f_dulon)
476	! CALL writeField("dulat_out",f_dulat)
477	! CALL writeField("dTemp_out",f_dTemp)
478	! CALL writeField("dq_out",f_dq)
479	! CALL writeField("dps_out",f_dps)
480
481	! CALL checksum(f_dulon)
482	! CALL checksum(f_dulat)
483	! CALL checksum(f_dTemp)
484	! CALL checksum(f_dq)
485	! CALL checksum(f_dps)
486
487	CALL send_message(f_dps,req_dps0)
488	CALL send_message(f_dulon,req_dulon0)
489	CALL send_message(f_dulat,req_dulat0)
490	CALL send_message(f_dTemp,req_dTemp0)
491	CALL send_message(f_dq,req_dq0)
492
493	CALL wait_message(req_dps0)
494	CALL wait_message(req_dulon0)
495	CALL wait_message(req_dulat0)
496	CALL wait_message(req_dTemp0)
497	CALL wait_message(req_dq0)
498
499
500	! CALL trace_on()
501	CALL trace_end("physic")
502	!$OMP MASTER
503	CALL SYSTEM_CLOCK(stop_clock)
504	count_clock=count_clock+stop_clock-start_clock
505	!$OMP END MASTER
506
507	!$OMP BARRIER
508
509	DO ind=1,ndomain
510	CALL swap_dimensions(ind)
511	IF (assigned_domain(ind)) THEN
512	CALL swap_geometry(ind)
513
514	theta_rhodz=f_theta_rhodz(ind)
515	u=f_u(ind)
516	q=f_q(ind)
517	ps=f_ps(ind)
518	dulon=f_dulon(ind)
519	dulat=f_dulat(ind)
520	Temp=f_temp(ind)
521	dTemp=f_dTemp(ind)
522	dq=f_dq(ind)
523	dps=f_dps(ind)
524	duc=f_duc(ind)
525	p=f_p(ind)
526	pks=f_pks(ind)
527	pk=f_pk(ind)
528
529	CALL grid_physics_to_icosa
530	ENDIF
531	ENDDO
532
533	!$OMP BARRIER
534	CALL xios_set_context
535
536
537	CONTAINS
538
539	SUBROUTINE grid_icosa_to_physics
540	USE pression_mod
541	USE exner_mod
542	USE theta2theta_rhodz_mod
543	USE geopotential_mod
544	USE wind_mod
545	USE omp_para
546	IMPLICIT NONE
547
548	REAL(rstd) :: uc(3)
549	INTEGER :: i,j,ij,l
550
551
552	! compute pression
553
554	DO l = ll_begin,ll_endp1
555	DO j=jj_begin,jj_end
556	DO i=ii_begin,ii_end
557	ij=(j-1)*iim+i
558	p(ij,l) = ap(l) + bp(l) * ps(ij)
559	ENDDO
560	ENDDO
561	ENDDO
562
563	!$OMP BARRIER
564
565	! compute exner
566
567	IF (is_omp_first_level) THEN
568	DO j=jj_begin,jj_end
569	DO i=ii_begin,ii_end
570	ij=(j-1)*iim+i
571	pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
572	ENDDO
573	ENDDO
574	ENDIF
575
576	! 3D : pk
577	DO l = ll_begin,ll_end
578	DO j=jj_begin,jj_end
579	DO i=ii_begin,ii_end
580	ij=(j-1)*iim+i
581	pk(ij,l) = cpp * ((.5/preff)(p(ij,l)+p(ij,l+1))) * kappa
582	ENDDO
583	ENDDO
584	ENDDO
585
586	!$OMP BARRIER
587
588	! compute theta, temperature and pression at layer
589	DO l = ll_begin, ll_end
590	DO j=jj_begin,jj_end
591	DO i=ii_begin,ii_end
592	ij=(j-1)*iim+i
593	theta(ij,l) = theta_rhodz(ij,l,1) / ((p(ij,l)-p(ij,l+1))/g)
594	Temp(ij,l) = theta(ij,l) * pk(ij,l) / cpp
595	p_layer(ij,l)=preff(pk(ij,l)/cpp)*(1./kappa)
596	ENDDO
597	ENDDO
598	ENDDO
599
600
601	!!! Compute geopotential
602
603	! for first layer
604	IF (is_omp_first_level) THEN
605	DO j=jj_begin,jj_end
606	DO i=ii_begin,ii_end
607	ij=(j-1)*iim+i
608	phi( ij,1 ) = phis( ij ) + theta(ij,1) * ( pks(ij) - pk(ij,1) )
609	ENDDO
610	ENDDO
611	ENDIF
612	!!-> implicit flush on phi(:,1)
613
614	! for other layers
615	DO l = ll_beginp1, ll_end
616	DO j=jj_begin,jj_end
617	DO i=ii_begin,ii_end
618	ij=(j-1)*iim+i
619	phi(ij,l) = 0.5 * ( theta(ij,l) + theta(ij,l-1) ) &
620	* ( pk(ij,l-1) - pk(ij,l) )
621	ENDDO
622	ENDDO
623	ENDDO
624
625	!$OMP BARRIER
626
627
628	IF (is_omp_first_level) THEN
629	DO l = 2, llm
630	DO j=jj_begin,jj_end
631	! ---> Bug compilo intel ici en openmp
632	! ---> Couper la boucle
633	IF (j==jj_end+1) PRINT*,"this message must not be printed"
634	DO i=ii_begin,ii_end
635	ij=(j-1)*iim+i
636	phi(ij,l) = phi(ij,l)+ phi(ij,l-1)
637	ENDDO
638	ENDDO
639	ENDDO
640	! --> IMPLICIT FLUSH on phi --> non
641	ENDIF
642
643	! compute wind centered lon lat compound
644	DO l=ll_begin,ll_end
645	DO j=jj_begin,jj_end
646	DO i=ii_begin,ii_end
647	ij=(j-1)*iim+i
648	uc(:)=1/Ai(ij)* &
649	( ne(ij,right)u(ij+u_right,l)le(ij+u_right)*((xyz_v(ij+z_rdown,:)+xyz_v(ij+z_rup,:))/2-centroid(ij,:)) &
650	+ ne(ij,rup)u(ij+u_rup,l)le(ij+u_rup)*((xyz_v(ij+z_rup,:)+xyz_v(ij+z_up,:))/2-centroid(ij,:)) &
651	+ ne(ij,lup)u(ij+u_lup,l)le(ij+u_lup)*((xyz_v(ij+z_up,:)+xyz_v(ij+z_lup,:))/2-centroid(ij,:)) &
652	+ ne(ij,left)u(ij+u_left,l)le(ij+u_left)*((xyz_v(ij+z_lup,:)+xyz_v(ij+z_ldown,:))/2-centroid(ij,:)) &
653	+ ne(ij,ldown)u(ij+u_ldown,l)le(ij+u_ldown)*((xyz_v(ij+z_ldown,:)+xyz_v(ij+z_down,:))/2-centroid(ij,:))&
654	+ ne(ij,rdown)u(ij+u_rdown,l)le(ij+u_rdown)*((xyz_v(ij+z_down,:)+xyz_v(ij+z_rdown,:))/2-centroid(ij,:)))
655	ulon(ij,l)=sum(uc(:)*elon_i(ij,:))
656	ulat(ij,l)=sum(uc(:)*elat_i(ij,:))
657	ENDDO
658	ENDDO
659	ENDDO
660
661	!$OMP BARRIER
662	END SUBROUTINE grid_icosa_to_physics
663
664
665	SUBROUTINE grid_physics_to_icosa
666	USE theta2theta_rhodz_mod
667	USE omp_para
668	IMPLICIT NONE
669	INTEGER :: i,j,ij,l,iq
670
671	DO l=ll_begin,ll_end
672	DO j=jj_begin,jj_end
673	DO i=ii_begin,ii_end
674	ij=(j-1)*iim+i
675	duc(ij,:,l)=dulon(ij,l)elon_i(ij,:)+dulat(ij,l)elat_i(ij,:)
676	ENDDO
677	ENDDO
678	ENDDO
679
680	DO l=ll_begin,ll_end
681	DO j=jj_begin,jj_end
682	DO i=ii_begin,ii_end
683	ij=(j-1)*iim+i
684	u(ij+u_right,l) = u(ij+u_right,l) + dtphy * sum( 0.5(duc(ij,:,l) + duc(ij+t_right,:,l))ep_e(ij+u_right,:) )
685	u(ij+u_lup,l) = u(ij+u_lup,l) + dtphy * sum( 0.5(duc(ij,:,l) + duc(ij+t_lup,:,l))ep_e(ij+u_lup,:) )
686	u(ij+u_ldown,l) = u(ij+u_ldown,l) + dtphysum( 0.5(duc(ij,:,l) + duc(ij+t_ldown,:,l))*ep_e(ij+u_ldown,:) )
687	ENDDO
688	ENDDO
689	ENDDO
690
691	DO l=ll_begin,ll_end
692	DO j=jj_begin,jj_end
693	DO i=ii_begin,ii_end
694	ij=(j-1)*iim+i
695	Temp(ij,l)=Temp(ij,l)+ dtphy * dTemp(ij,l)
696	ENDDO
697	ENDDO
698	ENDDO
699
700	DO iq=1,nqtot
701	DO l=ll_begin,ll_end
702	DO j=jj_begin,jj_end
703	DO i=ii_begin,ii_end
704	ij=(j-1)*iim+i
705	q(ij,l,iq)=q(ij,l,iq)+ dtphy * dq(ij,l,iq)
706	ENDDO
707	ENDDO
708	ENDDO
709	ENDDO
710
711	!$OMP BARRIER
712
713	IF (is_omp_first_level) THEN
714	DO j=jj_begin,jj_end
715	DO i=ii_begin,ii_end
716	ij=(j-1)*iim+i
717	ps(ij)=ps(ij)+ dtphy * dps(ij)
718	ENDDO
719	ENDDO
720	ENDIF
721
722	! CALL compute_temperature2theta_rhodz(ps,Temp,theta_rhodz,0)
723
724	! compute pression
725	!$OMP BARRIER
726	DO l = ll_begin,ll_endp1
727	DO j=jj_begin,jj_end
728	DO i=ii_begin,ii_end
729	ij=(j-1)*iim+i
730	p(ij,l) = ap(l) + bp(l) * ps(ij)
731	ENDDO
732	ENDDO
733	ENDDO
734
735	!$OMP BARRIER
736
737	! compute exner
738
739	IF (is_omp_first_level) THEN
740	DO j=jj_begin,jj_end
741	DO i=ii_begin,ii_end
742	ij=(j-1)*iim+i
743	pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
744	ENDDO
745	ENDDO
746	ENDIF
747
748	! 3D : pk
749	DO l = ll_begin,ll_end
750	DO j=jj_begin,jj_end
751	DO i=ii_begin,ii_end
752	ij=(j-1)*iim+i
753	pk(ij,l) = cpp * ((.5/preff)(p(ij,l)+p(ij,l+1))) * kappa
754	ENDDO
755	ENDDO
756	ENDDO
757
758	!$OMP BARRIER
759
760	! compute theta, temperature and pression at layer
761	DO l = ll_begin, ll_end
762	DO j=jj_begin,jj_end
763	DO i=ii_begin,ii_end
764	ij=(j-1)*iim+i
765	theta_rhodz(ij,l,1) = temp(ij,l) * ((p(ij,l)-p(ij,l+1))/g) / (pk(ij,l) / cpp )
766	ENDDO
767	ENDDO
768	ENDDO
769
770	END SUBROUTINE grid_physics_to_icosa
771
772
773
774	END SUBROUTINE physics
775
776
777
778
779
780	END MODULE interface_icosa_lmdz_mod

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

Download in other formats:

Original Format