1 | MODULE caldyn_unstructured_mod |
---|
2 | USE ISO_C_BINDING |
---|
3 | USE OMP_LIB |
---|
4 | USE data_unstructured_mod |
---|
5 | IMPLICIT NONE |
---|
6 | SAVE |
---|
7 | |
---|
8 | CONTAINS |
---|
9 | |
---|
10 | #include "unstructured.h90" |
---|
11 | |
---|
12 | #define INDICES1 ij,l,iq,iedge,edge,ivertex,vertex,ij_left,ij_right |
---|
13 | #define INDICES2 ij_up,ij_down,itrisk,edge_trisk,kup,kdown |
---|
14 | #define EDGES edge1,edge2,edge3,edge4,edge5,edge6 |
---|
15 | #define VERTICES vertex1,vertex2,vertex3,vertex4,vertex5,vertex6 |
---|
16 | #define SIGNS sign1,sign2,sign3,sign4,sign5,sign6 |
---|
17 | #define EDGE_ENDS ij_up1,ij_up2,ij_up3,ij_up4,ij_up5,ij_up6,ij_down1,ij_down2,ij_down3,ij_down4,ij_down5,ij_down6 |
---|
18 | #define LENGTHS le_de1,le_de2,le_de3,le_de4,le_de5,le_de6 |
---|
19 | #define DECLARE_INDICES INTEGER INDICES1,INDICES2 |
---|
20 | #define DECLARE_EDGES NUM SIGNS,LENGTHS ; INTEGER EDGES, EDGE_ENDS |
---|
21 | #define DECLARE_VERTICES INTEGER VERTICES |
---|
22 | #define PHI_BOT(ij) Phi_bot |
---|
23 | |
---|
24 | #define HASNAN(field) (ANY(.NOT.ABS(field)<1e20)) |
---|
25 | |
---|
26 | #define START_TRACE(id,nprimal,ndual,nedge) CALL enter_trace(id, 8*llm*((nprimal)*primal_num+(ndual)*dual_num+(nedge)*edge_num) ) |
---|
27 | #define STOP_TRACE CALL exit_trace() |
---|
28 | |
---|
29 | !----------------------------- Non-Hydrostatic ----------------------------- |
---|
30 | |
---|
31 | SUBROUTINE compute_NH_geopot(tau, m_ik, m_il, theta, W_il, Phi_il) |
---|
32 | FIELD_MASS :: m_ik, theta ! IN*2 |
---|
33 | FIELD_GEOPOT :: m_il, W_il, Phi_il, Phi_star_il ! IN,INOUT*2, LOCAL |
---|
34 | NUM :: tau, gamma, tau2_g, tau2_g2, g2, gm2, vreff, Rd_preff |
---|
35 | INTEGER :: iter |
---|
36 | DECLARE_INDICES |
---|
37 | NUM :: rho_ij, X_ij, Y_ij, wil, rho_c2_mik, c2_mik, ml_g2 |
---|
38 | #define COLUMN 0 |
---|
39 | #if COLUMN |
---|
40 | NUM1(llm) :: pk, Ak, Ck |
---|
41 | NUM1(llm+1):: Rl, Bl, Dl, xl |
---|
42 | #define p_ik(l,ij) pk(l) |
---|
43 | #define A_ik(l,ij) Ak(l) |
---|
44 | #define C_ik(l,ij) Ck(l) |
---|
45 | #define R_il(l,ij) Rl(l) |
---|
46 | #define B_il(l,ij) Bl(l) |
---|
47 | #define D_il(l,ij) Dl(l) |
---|
48 | #define x_il(l,ij) xl(l) |
---|
49 | #else |
---|
50 | FIELD_MASS :: p_ik, A_ik, C_ik |
---|
51 | FIELD_GEOPOT :: R_il, B_il, D_il, x_il |
---|
52 | #endif |
---|
53 | |
---|
54 | START_TRACE(id_NH_geopot, 7,0,0) |
---|
55 | #include "../kernels_unst/compute_NH_geopot.k90" |
---|
56 | STOP_TRACE |
---|
57 | |
---|
58 | #if COLUMN |
---|
59 | #undef p_ik |
---|
60 | #undef A_ik |
---|
61 | #undef C_ik |
---|
62 | #undef R_il |
---|
63 | #undef B_il |
---|
64 | #undef D_il |
---|
65 | #undef x_il |
---|
66 | #endif |
---|
67 | #undef COLUMN |
---|
68 | END SUBROUTINE compute_NH_geopot |
---|
69 | |
---|
70 | SUBROUTINE compute_caldyn_slow_NH(u,rhodz,Phi,W, w_il,berni,gradPhi2,DePhil,v_el,G_el,F_el, hflux,du,dPhi,dW) |
---|
71 | FIELD_U :: u, hflux, du ! IN, OUT, OUT |
---|
72 | FIELD_MASS :: rhodz, berni ! IN, BUF |
---|
73 | FIELD_GEOPOT :: Phi,W,dPhi,dW, w_il, gradPhi2 ! IN,IN, OUT,OUT, BUF*2 |
---|
74 | FIELD_UL :: DePhil, v_el, G_el, F_el ! BUF*4 |
---|
75 | DECLARE_INDICES |
---|
76 | DECLARE_EDGES |
---|
77 | NUM :: W_el, W2_el, gPhi2, dP, divG, u2, uu |
---|
78 | START_TRACE(id_slow_NH, 5,0,3) |
---|
79 | #include "../kernels_unst/caldyn_slow_NH.k90" |
---|
80 | STOP_TRACE |
---|
81 | END SUBROUTINE compute_caldyn_slow_NH |
---|
82 | |
---|
83 | SUBROUTINE compute_caldyn_solver(tau,rhodz,theta, berni,pres,m_il, pk,geopot,W,dPhi,dW,du) |
---|
84 | NUM, INTENT(IN) :: tau |
---|
85 | FIELD_MASS :: rhodz,pk,berni,pres ! IN, OUT, BUF*2 |
---|
86 | FIELD_THETA :: theta ! IN |
---|
87 | FIELD_GEOPOT :: geopot,W,dPhi,dW, m_il ! INOUT,INOUT, OUT,OUT, BUF |
---|
88 | FIELD_U :: du ! OUT |
---|
89 | DECLARE_INDICES |
---|
90 | NUM :: X_ij, rho_ij, T_ij, gamma, Cvd, vreff, Rd_preff |
---|
91 | #include "../kernels_unst/caldyn_mil.k90" |
---|
92 | IF(tau>0) THEN ! solve implicit problem for geopotential |
---|
93 | CALL compute_NH_geopot(tau, rhodz, m_il, theta, W, geopot) |
---|
94 | END IF |
---|
95 | START_TRACE(id_solver, 7,0,1) |
---|
96 | #include "../kernels_unst/caldyn_solver.k90" |
---|
97 | STOP_TRACE |
---|
98 | END SUBROUTINE compute_caldyn_solver |
---|
99 | |
---|
100 | SUBROUTINE compute_caldyn_vert_NH(mass,geopot,W,wflux, eta_dot,wcov,W_etadot, du,dPhi,dW) |
---|
101 | FIELD_MASS :: mass, eta_dot, wcov, W_etadot ! IN, BUF*3 |
---|
102 | FIELD_GEOPOT :: geopot,W,wflux,dPhi,dW ! IN*3, INOUT*2 |
---|
103 | FIELD_U :: du ! INOUT |
---|
104 | DECLARE_INDICES |
---|
105 | NUM :: w_ij, wflux_ij |
---|
106 | START_TRACE(id_vert_NH, 6,0,1) |
---|
107 | #include "../kernels_unst/caldyn_vert_NH.k90" |
---|
108 | STOP_TRACE |
---|
109 | END SUBROUTINE compute_caldyn_vert_NH |
---|
110 | |
---|
111 | !----------------------------- Hydrostatic ------------------------------- |
---|
112 | |
---|
113 | SUBROUTINE compute_geopot(rhodz,theta,ps,pk,geopot) |
---|
114 | FIELD_MASS :: rhodz,pk ! IN, OUT |
---|
115 | FIELD_THETA :: theta ! IN |
---|
116 | FIELD_GEOPOT :: geopot ! IN(l=1)/OUT(l>1) |
---|
117 | FIELD_PS :: ps ! OUT |
---|
118 | DECLARE_INDICES |
---|
119 | NUM :: gdz, ke, uu, chi, gv, exner_ik, temp_ik, p_ik, qv, Rmix |
---|
120 | START_TRACE(id_geopot, 3,0,3) |
---|
121 | #include "../kernels_unst/compute_geopot.k90" |
---|
122 | STOP_TRACE |
---|
123 | END SUBROUTINE compute_geopot |
---|
124 | ! |
---|
125 | SUBROUTINE compute_caldyn_slow_hydro(rhodz,theta,u, berni,hflux,du) |
---|
126 | FIELD_MASS :: rhodz,berni! IN |
---|
127 | FIELD_THETA :: theta ! IN |
---|
128 | FIELD_U :: u,hflux,du ! IN, OUT, OUT |
---|
129 | DECLARE_INDICES |
---|
130 | DECLARE_EDGES |
---|
131 | LOGICAL, PARAMETER :: zero=.TRUE. |
---|
132 | NUM :: ke, uu |
---|
133 | START_TRACE(id_slow_hydro, 3,0,3) |
---|
134 | #include "../kernels_unst/caldyn_slow_hydro.k90" |
---|
135 | STOP_TRACE |
---|
136 | END SUBROUTINE compute_caldyn_slow_hydro |
---|
137 | |
---|
138 | !---------------------------------- Generic ------------------------------ |
---|
139 | |
---|
140 | SUBROUTINE caldyn_vert(convm,rhodz,theta,u, dmass_col,wflux,dtheta_rhodz,du, wwuu) |
---|
141 | FIELD_PS :: dmass_col |
---|
142 | FIELD_MASS :: convm, rhodz |
---|
143 | FIELD_U :: u,du |
---|
144 | FIELD_THETA :: dtheta_rhodz, theta |
---|
145 | FIELD_GEOPOT :: wflux |
---|
146 | FIELD_UL :: wwuu |
---|
147 | DECLARE_INDICES |
---|
148 | NUM :: dF_deta, dFu_deta |
---|
149 | wwuu=0. |
---|
150 | !$OMP BARRIER |
---|
151 | START_TRACE(id_vert, 5,0,3) |
---|
152 | #include "../kernels_unst/caldyn_wflux.k90" |
---|
153 | #include "../kernels_unst/caldyn_dmass.k90" |
---|
154 | #include "../kernels_unst/caldyn_vert.k90" |
---|
155 | STOP_TRACE |
---|
156 | END SUBROUTINE caldyn_vert |
---|
157 | |
---|
158 | SUBROUTINE compute_coriolis(hflux,theta,qu,Ftheta, convm,dtheta_rhodz,du) |
---|
159 | FIELD_U :: hflux, Ftheta, qu, du |
---|
160 | FIELD_MASS :: convm |
---|
161 | FIELD_THETA :: theta, dtheta_rhodz |
---|
162 | DECLARE_INDICES |
---|
163 | DECLARE_EDGES |
---|
164 | NUM :: divF, du_trisk |
---|
165 | START_TRACE(id_coriolis, 3,4,0) ! primal, dual, edge |
---|
166 | #include "../kernels_unst/coriolis.k90" |
---|
167 | STOP_TRACE |
---|
168 | END SUBROUTINE |
---|
169 | |
---|
170 | SUBROUTINE compute_theta(mass_col,rhodz,theta_rhodz, theta) |
---|
171 | FIELD_PS :: mass_col |
---|
172 | FIELD_MASS :: rhodz |
---|
173 | FIELD_THETA :: theta, theta_rhodz |
---|
174 | DECLARE_INDICES |
---|
175 | NUM :: m |
---|
176 | START_TRACE(id_theta, 3,0,0) ! primal, dual, edge |
---|
177 | #include "../kernels_unst/theta.k90" |
---|
178 | STOP_TRACE |
---|
179 | END SUBROUTINE |
---|
180 | |
---|
181 | SUBROUTINE compute_pvort_only(rhodz,u,qv,qu) |
---|
182 | FIELD_MASS :: rhodz |
---|
183 | FIELD_U :: u,qu |
---|
184 | FIELD_Z :: qv |
---|
185 | DECLARE_INDICES |
---|
186 | DECLARE_EDGES |
---|
187 | DECLARE_VERTICES |
---|
188 | NUM :: etav, hv |
---|
189 | START_TRACE(id_pvort_only, 1,1,2) ! primal, dual, edge |
---|
190 | #include "../kernels_unst/pvort_only.k90" |
---|
191 | STOP_TRACE |
---|
192 | END SUBROUTINE compute_pvort_only |
---|
193 | |
---|
194 | SUBROUTINE compute_caldyn_fast(tau, pk,berni,theta,geopot, du,u) |
---|
195 | NUM, INTENT(IN) :: tau |
---|
196 | FIELD_MASS :: pk,berni ! INOUT, OUT |
---|
197 | FIELD_THETA :: theta ! IN |
---|
198 | FIELD_GEOPOT :: geopot ! IN |
---|
199 | FIELD_U :: u,du ! INOUT,INOUT |
---|
200 | DECLARE_INDICES |
---|
201 | DECLARE_EDGES |
---|
202 | NUM :: due |
---|
203 | START_TRACE(id_fast, 4,0,2) ! primal, dual, edge |
---|
204 | #include "../kernels_unst/caldyn_fast.k90" |
---|
205 | STOP_TRACE |
---|
206 | END SUBROUTINE compute_caldyn_fast |
---|
207 | |
---|
208 | !----------------------------- Unused ----------------------------- |
---|
209 | |
---|
210 | SUBROUTINE gradient(b,grad) BINDC(gradient) |
---|
211 | FIELD_MASS :: b |
---|
212 | FIELD_U :: grad |
---|
213 | DECLARE_INDICES |
---|
214 | #include "../kernels_unst/gradient.k90" |
---|
215 | END SUBROUTINE |
---|
216 | ! |
---|
217 | SUBROUTINE div(u,divu) BINDC(div) |
---|
218 | FIELD_MASS :: divu |
---|
219 | FIELD_U :: u |
---|
220 | DECLARE_INDICES |
---|
221 | DECLARE_EDGES |
---|
222 | NUM :: div_ij |
---|
223 | !$OMP PARALLEL NUM_THREADS(nb_threads) |
---|
224 | #include "../kernels_unst/div.k90" |
---|
225 | !$OMP END PARALLEL |
---|
226 | END SUBROUTINE |
---|
227 | |
---|
228 | END MODULE caldyn_unstructured_mod |
---|