source: ICOSA_LMDZ/src/phydev/interface_icosa_lmdz.f90 @ 5232

Last change on this file since 5232 was 4855, checked in by yann meurdesoif, 8 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
RevLine 
[4498]1MODULE 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
39CONTAINS
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
[4855]133    CALL param_gravity_wave
[4498]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) + dtphy*sum( 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
780END MODULE interface_icosa_lmdz_mod
Note: See TracBrowser for help on using the repository browser.