source: codes/icosagcm/devel/src/dynamics/caldyn_hevi.f90 @ 921

Last change on this file since 921 was 921, checked in by dubos, 5 years ago

devel : unique interface for compute_caldyn_fast_X, compute_caldyn_slow_hydro_X

File size: 7.3 KB
Line 
1MODULE caldyn_hevi_mod
2  USE icosa
3  USE transfert_mod
4  USE caldyn_vars_mod
5  USE caldyn_kernels_hevi_mod
6  USE caldyn_kernels_base_mod
7  USE compute_theta_mod, ONLY : compute_theta => compute_theta_manual
8  USE compute_geopot_mod, ONLY : compute_geopot => compute_geopot_manual
9  USE compute_caldyn_kv_mod, ONLY : compute_caldyn_kv
10  USE compute_caldyn_Coriolis_mod, ONLY : compute_caldyn_Coriolis
11  USE compute_caldyn_slow_hydro_mod, ONLY : compute_caldyn_slow_hydro => compute_caldyn_slow_hydro_manual
12  USE compute_caldyn_slow_NH_mod, ONLY : compute_caldyn_slow_NH
13  USE compute_caldyn_solver_mod, ONLY : compute_caldyn_solver
14  USE compute_caldyn_fast_mod, ONLY : compute_caldyn_fast => compute_caldyn_fast_manual
15  USE compute_NH_geopot_mod, ONLY : compute_NH_geopot
16  IMPLICIT NONE
17  PRIVATE
18  PUBLIC caldyn_hevi
19
20CONTAINS
21 
22  SUBROUTINE caldyn_hevi(write_out,tau, f_phis, f_ps, f_mass, f_theta_rhodz, f_u, f_q, &
23       f_W, f_geopot, f_hflux, f_wflux, f_dps, f_dmass, f_dtheta_rhodz, &
24       f_du_slow, f_du_fast, f_dPhi_slow, f_dPhi_fast, f_dW_slow, f_dW_fast) 
25    USE icosa
26    USE observable_mod
27    USE disvert_mod, ONLY : caldyn_eta, eta_mass
28    USE vorticity_mod
29    USE kinetic_mod
30    USE theta2theta_rhodz_mod
31    USE wind_mod
32    USE mpipara
33    USE trace
34    USE omp_para
35    USE output_field_mod
36    USE checksum_mod
37    USE compute_caldyn_mod, ONLY : compute_pvort_only
38    IMPLICIT NONE
39    LOGICAL,INTENT(IN)    :: write_out
40    REAL(rstd), INTENT(IN) :: tau
41    TYPE(t_field),POINTER :: f_phis(:)
42    TYPE(t_field),POINTER :: f_ps(:)
43    TYPE(t_field),POINTER :: f_mass(:)
44    TYPE(t_field),POINTER :: f_theta_rhodz(:)
45    TYPE(t_field),POINTER :: f_u(:)
46    TYPE(t_field),POINTER :: f_q(:)
47    TYPE(t_field),POINTER :: f_W(:)
48    TYPE(t_field),POINTER :: f_geopot(:)
49    TYPE(t_field),POINTER :: f_hflux(:), f_wflux(:)
50    TYPE(t_field) :: f_dps(:)
51    TYPE(t_field) :: f_dmass(:)
52    TYPE(t_field) :: f_dtheta_rhodz(:)
53    TYPE(t_field) :: f_du_slow(:)
54    TYPE(t_field) :: f_du_fast(:)
55    TYPE(t_field) :: f_dW_slow(:)
56    TYPE(t_field) :: f_dW_fast(:)
57    TYPE(t_field) :: f_dPhi_slow(:)
58    TYPE(t_field) :: f_dPhi_fast(:)
59   
60    REAL(rstd),POINTER :: ps(:), dps(:), phis(:)
61    REAL(rstd),POINTER :: mass(:,:), theta_rhodz(:,:,:), dtheta_rhodz(:,:,:)
62    REAL(rstd),POINTER :: du(:,:), dW(:,:), dPhi(:,:), hflux(:,:), wflux(:,:)
63    REAL(rstd),POINTER :: u(:,:), w(:,:), qu(:,:), qv(:,:), Kv(:,:), hv(:,:)
64
65! temporary shared variable
66    REAL(rstd),POINTER  :: theta(:,:,:) 
67    REAL(rstd),POINTER  :: pk(:,:)
68    REAL(rstd),POINTER  :: geopot(:,:)
69    REAL(rstd),POINTER  :: convm(:,:) 
70    REAL(rstd),POINTER  :: wwuu(:,:)
71    REAL(rstd),POINTER  :: F_el(:,:), gradPhi2(:,:), w_il(:,:) , W_etadot(:,:), pres(:,:), m_il(:,:)
72    INTEGER :: ind
73    LOGICAL,SAVE :: first=.TRUE.
74!$OMP THREADPRIVATE(first)
75   
76    IF (first) THEN
77      first=.FALSE.
78      IF(caldyn_eta==eta_mass) THEN
79         CALL init_message(f_ps,req_i1,req_ps)
80      ELSE
81         CALL init_message(f_mass,req_i1,req_mass)
82      END IF
83      CALL init_message(f_theta_rhodz,req_i1,req_theta_rhodz)
84      CALL init_message(f_u,req_e1_vect,req_u)
85      CALL init_message(f_qu,req_e1_scal,req_qu)
86      IF(caldyn_kinetic==kinetic_consistent) CALL init_message(f_Kv,req_z1_scal,req_Kv)
87      IF(.NOT.hydrostatic) THEN
88         CALL init_message(f_geopot,req_i1,req_geopot)
89         CALL init_message(f_w,req_i1,req_w)
90      END IF
91    ENDIF
92   
93    CALL trace_start("caldyn")
94   
95    IF(caldyn_eta==eta_mass) THEN
96       CALL send_message(f_ps,req_ps) ! COM00
97       CALL wait_message(req_ps) ! COM00
98    ELSE
99       CALL send_message(f_mass,req_mass) ! COM00
100       CALL wait_message(req_mass) ! COM00
101    END IF
102    CALL send_message(f_theta_rhodz,req_theta_rhodz) ! COM01
103    CALL wait_message(req_theta_rhodz) ! COM01 Moved from caldyn_pvort
104
105    IF(.NOT.hydrostatic) THEN
106       CALL send_message(f_geopot,req_geopot) ! COM03
107       CALL wait_message(req_geopot) ! COM03
108       CALL send_message(f_w,req_w) ! COM04
109       CALL wait_message(req_w) ! COM04
110    END IF
111   
112    DO ind=1,ndomain
113       IF (.NOT. assigned_domain(ind)) CYCLE
114       CALL swap_dimensions(ind)
115       CALL swap_geometry(ind)
116       ps=f_ps(ind)
117       theta_rhodz=f_theta_rhodz(ind)
118       mass=f_mass(ind)
119       theta = f_theta(ind)
120       CALL compute_theta(ps,theta_rhodz, mass,theta)
121       pk = f_pk(ind)
122       geopot = f_geopot(ind)
123       du=f_du_fast(ind)
124       IF(hydrostatic) THEN
125          du(:,:)=0.
126          CALL compute_geopot(mass,theta, ps,pk,geopot)
127       ELSE
128          phis = f_phis(ind)
129          W = f_W(ind)
130          dW = f_dW_fast(ind)
131          dPhi = f_dPhi_fast(ind)
132          ! reuse buffers
133          m_il = f_wil(ind)
134          pres = f_gradPhi2(ind)
135          CALL compute_caldyn_solver(tau,phis, mass,theta,pk,geopot,W, m_il,pres, dPhi,dW,du) ! computes d(Phi,W,du)_fast and updates Phi,W
136       END IF
137       u=f_u(ind)
138       Kv=f_Kv(ind) ! buffer for partial Bernoulli function
139       CALL compute_caldyn_fast(tau,theta,geopot, pk,Kv,du,u) ! computes du_fast and updates u
140    ENDDO
141   
142    CALL send_message(f_u,req_u) ! COM02
143    CALL wait_message(req_u)   ! COM02
144   
145    DO ind=1,ndomain
146       IF (.NOT. assigned_domain(ind)) CYCLE
147       CALL swap_dimensions(ind)
148       CALL swap_geometry(ind)
149       u=f_u(ind)
150       mass=f_mass(ind)
151       qu=f_qu(ind)
152       qv=f_qv(ind)
153       hv=f_hv(ind)
154!       PRINT *, 'caldyn_hevi : u, mass, qu qv, hv', &
155!            SHAPE(u), SHAPE(mass), SHAPE(qu), SHAPE(qv), SHAPE(hv)
156       CALL compute_pvort_only(u,mass,qu,qv,hv)
157       IF(caldyn_kinetic==kinetic_consistent) THEN
158          Kv=f_Kv(ind)
159          CALL compute_caldyn_Kv(u,Kv)
160       END IF
161    ENDDO
162   
163    CALL send_message(f_qu,req_qu) ! COM03
164    CALL wait_message(req_qu) ! COM03
165
166    IF(caldyn_kinetic==kinetic_consistent) THEN
167       CALL transfert_request(f_hv,req_z1_scal)
168       CALL send_message(f_Kv,req_Kv)
169       CALL wait_message(req_Kv)
170    END IF
171
172    DO ind=1,ndomain
173       IF (.NOT. assigned_domain(ind)) CYCLE
174       CALL swap_dimensions(ind)
175       CALL swap_geometry(ind)
176       u=f_u(ind)
177       mass=f_mass(ind)
178       theta = f_theta(ind)
179       qu=f_qu(ind)
180       hflux=f_hflux(ind)
181       convm = f_dmass(ind)
182       dtheta_rhodz=f_dtheta_rhodz(ind)
183       du=f_du_slow(ind)
184
185       IF(hydrostatic) THEN
186          hv=f_hv(ind)
187          Kv=f_Kv(ind)
188          CALL compute_caldyn_slow_hydro(.TRUE., u,mass,hv,Kv, Kv, hflux,du)
189       ELSE
190          W = f_W(ind)
191          dW = f_dW_slow(ind)
192          geopot = f_geopot(ind)
193          dPhi = f_dPhi_slow(ind)
194          F_el = f_Fel(ind)
195          gradPhi2 = f_gradPhi2(ind)
196          w_il = f_wil(ind)
197          CALL compute_caldyn_slow_NH(u,mass,geopot,W, F_el,gradPhi2,w_il, hflux,du,dPhi,dW)
198       END IF
199       CALL compute_caldyn_Coriolis(hflux,theta,qu, convm,dtheta_rhodz,du)
200       IF(caldyn_eta==eta_mass) THEN
201          wflux=f_wflux(ind)
202          wwuu=f_wwuu(ind)
203          dps=f_dps(ind)
204          CALL compute_caldyn_vert(u,theta,mass,convm, wflux,wwuu, dps, dtheta_rhodz, du)
205          IF(.NOT.hydrostatic) THEN
206             W_etadot=f_Wetadot(ind)
207             CALL compute_caldyn_vert_NH(mass,geopot,W,wflux, W_etadot, du,dPhi,dW)
208          END IF
209       END IF
210    ENDDO
211   
212!$OMP BARRIER
213    !    CALL check_mass_conservation(f_ps,f_dps)
214    CALL trace_end("caldyn_hevi")
215!!$OMP BARRIER
216   
217  END SUBROUTINE caldyn_hevi
218
219END MODULE caldyn_hevi_mod
Note: See TracBrowser for help on using the repository browser.