source: codes/icosagcm/devel/src/unstructured/caldyn_unstructured.F90 @ 802

Last change on this file since 802 was 796, checked in by jisesh, 5 years ago

unstructured : reconstruction of velocity components (bugged)
File size: 6.9 KB
Line 
1MODULE caldyn_unstructured_mod
2  USE ISO_C_BINDING
3  USE OMP_LIB
4  USE data_unstructured_mod
5  IMPLICIT NONE
6  SAVE
7
8CONTAINS
9
10#include "unstructured.h90"
11
12#define PHI_BOT(ij) Phi_bot
13
14!----------------------------- Non-Hydrostatic -----------------------------
15
16SUBROUTINE compute_NH_geopot(tau, m_ik, m_il, theta, W_il, Phi_il)
17  FIELD_MASS   :: m_ik, theta   ! IN*2
18  FIELD_GEOPOT :: m_il, W_il, Phi_il, Phi_star_il  ! IN,INOUT*2, LOCAL
19  NUM :: tau, gamma, tau2_g, tau2_g2, g2, gm2, vreff, Rd_preff
20  INTEGER :: iter
21  LOGICAL :: debug_hevi_solver
22  DECLARE_INDICES
23  NUM :: rho_ij, X_ij, Y_ij, wil, rho_c2_mik, c2_mik, ml_g2
24#define COLUMN 0
25#if COLUMN 
26  NUM1(llm)  :: pk, Ak, Ck
27  NUM1(llm+1):: Rl, Bl, Dl, xl
28#define p_ik(l,ij) pk(l)
29#define A_ik(l,ij) Ak(l)
30#define C_ik(l,ij) Ck(l)
31#define R_il(l,ij) Rl(l)
32#define B_il(l,ij) Bl(l)
33#define D_il(l,ij) Dl(l)
34#define x_il(l,ij) xl(l)
35#else
36  FIELD_MASS :: p_ik, A_ik, C_ik
37  FIELD_GEOPOT :: R_il, B_il, D_il, x_il
38#endif
39
40  debug_hevi_solver=.FALSE.
41!$OMP MASTER
42  debug_hevi_solver = debug_hevi_solver_
43!$OMP END MASTER
44
45  START_TRACE(id_NH_geopot, 7,0,0)
46#include "../kernels_unst/compute_NH_geopot.k90"
47  STOP_TRACE 
48
49!$OMP MASTER
50  debug_hevi_solver_ = debug_hevi_solver
51!$OMP END MASTER
52
53#if COLUMN
54#undef p_ik
55#undef A_ik
56#undef C_ik
57#undef R_il
58#undef B_il
59#undef D_il
60#undef x_il
61#endif
62#undef COLUMN
63END SUBROUTINE compute_NH_geopot
64
65SUBROUTINE compute_caldyn_slow_NH(u,rhodz,Phi,W, w_il,berni,gradPhi2,DePhil,v_el,G_el,F_el, hflux,du,dPhi,dW)
66  FIELD_U      :: u, hflux, du   ! IN, OUT, OUT
67  FIELD_MASS   :: rhodz, berni   ! IN, BUF
68  FIELD_GEOPOT :: Phi,W,dPhi,dW, w_il, gradPhi2  ! IN,IN, OUT,OUT, BUF*2
69  FIELD_UL     :: DePhil, v_el, G_el, F_el ! BUF*4
70  DECLARE_INDICES
71  DECLARE_EDGES
72  NUM :: W_el, W2_el, gPhi2, dP, divG, u2, uu
73  START_TRACE(id_slow_NH, 5,0,3)
74#include "../kernels_unst/caldyn_slow_NH.k90"
75  STOP_TRACE
76END SUBROUTINE compute_caldyn_slow_NH
77
78SUBROUTINE compute_caldyn_solver(tau,rhodz,theta, berni,pres,m_il, pk,geopot,W,dPhi,dW,du)
79  NUM, INTENT(IN) :: tau
80  FIELD_MASS   :: rhodz,pk,berni,pres    ! IN, OUT, BUF*2
81  FIELD_THETA  :: theta                  ! IN
82  FIELD_GEOPOT :: geopot,W,dPhi,dW, m_il ! INOUT,INOUT, OUT,OUT, BUF
83  FIELD_U      :: du                     ! OUT
84  DECLARE_INDICES
85  NUM :: X_ij, rho_ij, T_ij, gamma, Cvd, vreff, Rd_preff
86#include "../kernels_unst/caldyn_mil.k90"
87  IF(tau>0) THEN ! solve implicit problem for geopotential
88    CALL compute_NH_geopot(tau, rhodz, m_il, theta, W, geopot)
89  END IF
90  START_TRACE(id_solver, 7,0,1)
91#include "../kernels_unst/caldyn_solver.k90"
92  STOP_TRACE
93END SUBROUTINE compute_caldyn_solver
94
95SUBROUTINE compute_caldyn_vert_NH(mass,geopot,W,wflux, eta_dot,wcov,W_etadot, du,dPhi,dW)
96  FIELD_MASS   :: mass, eta_dot, wcov, W_etadot ! IN, BUF*3
97  FIELD_GEOPOT :: geopot,W,wflux,dPhi,dW ! IN*3, INOUT*2
98  FIELD_U      :: du ! INOUT
99  DECLARE_INDICES
100  NUM :: w_ij, wflux_ij
101  START_TRACE(id_vert_NH, 6,0,1)
102#include "../kernels_unst/caldyn_vert_NH.k90"
103  STOP_TRACE
104END SUBROUTINE compute_caldyn_vert_NH
105
106!----------------------------- Hydrostatic -------------------------------
107
108SUBROUTINE compute_geopot(rhodz,theta,ps,pk,geopot)
109  FIELD_MASS  :: rhodz,pk   ! IN, OUT
110  FIELD_THETA :: theta      ! IN
111  FIELD_GEOPOT :: geopot    ! IN(l=1)/OUT(l>1)
112  FIELD_PS     :: ps        ! OUT
113  DECLARE_INDICES
114  NUM :: gdz, ke, uu, chi, gv, exner_ik, temp_ik, p_ik, qv, Rmix
115  START_TRACE(id_geopot, 3,0,3)
116#include "../kernels_unst/compute_geopot.k90"
117  STOP_TRACE
118END SUBROUTINE compute_geopot
119!
120SUBROUTINE compute_caldyn_slow_hydro(rhodz,theta,u, berni,hflux,du)
121  FIELD_MASS  :: rhodz,berni! IN
122  FIELD_THETA :: theta      ! IN
123  FIELD_U     :: u,hflux,du ! IN, OUT, OUT
124  DECLARE_INDICES
125  DECLARE_EDGES
126  LOGICAL, PARAMETER :: zero=.TRUE.
127  NUM :: ke, uu
128  START_TRACE(id_slow_hydro, 3,0,3)
129#include "../kernels_unst/caldyn_slow_hydro.k90"
130  STOP_TRACE
131END SUBROUTINE compute_caldyn_slow_hydro
132
133!---------------------------------- Generic ------------------------------
134
135SUBROUTINE caldyn_vert(convm,rhodz,theta,u, dmass_col,wflux,dtheta_rhodz,du, wwuu)
136  FIELD_PS    :: dmass_col
137  FIELD_MASS  :: convm, rhodz
138  FIELD_U     :: u,du
139  FIELD_THETA :: dtheta_rhodz, theta
140  FIELD_GEOPOT :: wflux
141  FIELD_UL     :: wwuu
142  DECLARE_INDICES
143  NUM :: dF_deta, dFu_deta
144  wwuu=0.
145  !$OMP BARRIER
146  START_TRACE(id_vert, 5,0,3)
147#include "../kernels_unst/caldyn_wflux.k90"
148#include "../kernels_unst/caldyn_dmass.k90"
149#include "../kernels_unst/caldyn_vert.k90"
150  STOP_TRACE
151END SUBROUTINE caldyn_vert
152
153SUBROUTINE compute_coriolis(hflux,theta,qu,Ftheta, convm,dtheta_rhodz,du)
154  FIELD_U     :: hflux, Ftheta, qu, du
155  FIELD_MASS  :: convm
156  FIELD_THETA :: theta, dtheta_rhodz
157  DECLARE_INDICES
158  DECLARE_EDGES
159  NUM :: divF, du_trisk
160  START_TRACE(id_coriolis, 3,4,0) ! primal, dual, edge
161#include "../kernels_unst/coriolis.k90"
162  STOP_TRACE
163END SUBROUTINE
164
165SUBROUTINE compute_theta(mass_col,rhodz,theta_rhodz, theta)
166  FIELD_PS :: mass_col
167  FIELD_MASS :: rhodz
168  FIELD_THETA :: theta, theta_rhodz
169  DECLARE_INDICES
170  NUM :: m
171  START_TRACE(id_theta, 3,0,0) ! primal, dual, edge
172#include "../kernels_unst/theta.k90"
173  STOP_TRACE
174END SUBROUTINE
175
176SUBROUTINE compute_pvort_only(rhodz,u,qv,qu)
177  FIELD_MASS :: rhodz
178  FIELD_U    :: u,qu
179  FIELD_Z    :: qv
180  DECLARE_INDICES
181  DECLARE_EDGES
182  DECLARE_VERTICES
183  NUM :: etav, hv
184  START_TRACE(id_pvort_only, 1,1,2) ! primal, dual, edge
185#include "../kernels_unst/pvort_only.k90"
186  STOP_TRACE
187END SUBROUTINE compute_pvort_only
188
189SUBROUTINE compute_caldyn_fast(tau, pk,berni,theta,geopot, du,u)
190  NUM, INTENT(IN) :: tau
191  FIELD_MASS   :: pk,berni  ! INOUT, OUT
192  FIELD_THETA  :: theta     ! IN
193  FIELD_GEOPOT :: geopot    ! IN
194  FIELD_U      :: u,du      ! INOUT,INOUT
195  DECLARE_INDICES
196  DECLARE_EDGES
197  NUM          :: due
198  START_TRACE(id_fast, 4,0,2) ! primal, dual, edge
199#include "../kernels_unst/caldyn_fast.k90"
200  STOP_TRACE
201END SUBROUTINE compute_caldyn_fast
202
203!----------------------------- Unused -----------------------------
204
205SUBROUTINE gradient(b,grad) BINDC(gradient)
206  FIELD_MASS :: b
207  FIELD_U    :: grad
208  DECLARE_INDICES
209#include "../kernels_unst/gradient.k90"
210END SUBROUTINE
211!
212SUBROUTINE div(u,divu) BINDC(div)
213  FIELD_MASS :: divu
214  FIELD_U    :: u
215  DECLARE_INDICES
216  DECLARE_EDGES
217  NUM :: div_ij
218  !$OMP PARALLEL NUM_THREADS(nb_threads)
219#include "../kernels_unst/div.k90"
220  !$OMP END PARALLEL
221END SUBROUTINE
222
223SUBROUTINE compute_scalar_laplacian(b,grad,divu)
224  FIELD_MASS   :: b, divu ! IN, OUT
225  FIELD_U      :: grad ! BUF
226  DECLARE_INDICES
227  DECLARE_EDGES
228  NUM :: div_ij
229  START_TRACE(id_scalar_laplacian, 0,0,1)
230#include "../kernels_unst/scalar_laplacian.k90"
231  STOP_TRACE
232END SUBROUTINE compute_scalar_laplacian
233
234SUBROUTINE compute_curl_laplacian(u,qv,curlcurl)
235  FIELD_Z   :: qv ! BUF
236  FIELD_U   :: u, curlcurl ! IN,OUT
237  DECLARE_INDICES
238  DECLARE_EDGES
239  DECLARE_VERTICES
240  NUM :: etav
241  START_TRACE(id_scalar_laplacian, 0,0,1)
242#include "../kernels_unst/curl_laplacian.k90"
243  STOP_TRACE
244END SUBROUTINE compute_curl_laplacian
245
246END MODULE caldyn_unstructured_mod
Note: See TracBrowser for help on using the repository browser.