1 | MODULE compute_caldyn_fast_mod |
---|
2 | USE grid_param, ONLY : llm |
---|
3 | IMPLICIT NONE |
---|
4 | PRIVATE |
---|
5 | |
---|
6 | PUBLIC :: compute_caldyn_fast |
---|
7 | |
---|
8 | CONTAINS |
---|
9 | |
---|
10 | SUBROUTINE compute_caldyn_fast(tau,u,rhodz,theta,pk,geopot,du) |
---|
11 | USE icosa |
---|
12 | USE trace |
---|
13 | USE caldyn_vars_mod |
---|
14 | USE omp_para, ONLY : ll_begin, ll_end |
---|
15 | REAL(rstd),INTENT(IN) :: tau ! "solve" u-tau*du/dt = rhs |
---|
16 | REAL(rstd),INTENT(INOUT) :: u(iim*3*jjm,llm) ! OUT if tau>0 |
---|
17 | REAL(rstd),INTENT(IN) :: rhodz(iim*jjm,llm) |
---|
18 | REAL(rstd),INTENT(IN) :: theta(iim*jjm,llm,nqdyn) |
---|
19 | REAL(rstd),INTENT(INOUT) :: pk(iim*jjm,llm) |
---|
20 | REAL(rstd),INTENT(INOUT) :: geopot(iim*jjm,llm+1) |
---|
21 | REAL(rstd),INTENT(INOUT) :: du(iim*3*jjm,llm) |
---|
22 | REAL(rstd) :: berni(iim*jjm,llm) ! Bernoulli function |
---|
23 | REAL(rstd) :: berniv(iim*jjm,llm) ! moist Bernoulli function |
---|
24 | |
---|
25 | INTEGER :: i,j,ij,l |
---|
26 | !REAL(rstd) :: cp_ik, qv, temp, chi, nu, due, due_right, due_lup, due_ldown |
---|
27 | REAL(rstd) :: cp_ik, qv, temp, chi, due, due_right, due_lup, due_ldown |
---|
28 | |
---|
29 | CALL trace_start("compute_caldyn_fast") |
---|
30 | |
---|
31 | IF(dysl_caldyn_fast) THEN |
---|
32 | #include "../kernels_hex/caldyn_fast.k90" |
---|
33 | ELSE |
---|
34 | |
---|
35 | ! Compute Bernoulli term |
---|
36 | IF(boussinesq) THEN |
---|
37 | DO l=ll_begin,ll_end |
---|
38 | !DIR$ SIMD |
---|
39 | DO ij=ij_begin,ij_end |
---|
40 | berni(ij,l) = pk(ij,l) |
---|
41 | ! from now on pk contains the vertically-averaged geopotential |
---|
42 | pk(ij,l) = .5*(geopot(ij,l)+geopot(ij,l+1)) |
---|
43 | END DO |
---|
44 | END DO |
---|
45 | ELSE ! compressible |
---|
46 | |
---|
47 | DO l=ll_begin,ll_end |
---|
48 | SELECT CASE(caldyn_thermo) |
---|
49 | CASE(thermo_theta) ! vdp = theta.dpi => B = Phi |
---|
50 | !DIR$ SIMD |
---|
51 | DO ij=ij_begin,ij_end |
---|
52 | berni(ij,l) = .5*(geopot(ij,l)+geopot(ij,l+1)) |
---|
53 | END DO |
---|
54 | CASE(thermo_entropy) ! vdp = dG + sdT => B = Phi + G, G=h-Ts=T*(cpp-s) |
---|
55 | !DIR$ SIMD |
---|
56 | DO ij=ij_begin,ij_end |
---|
57 | berni(ij,l) = .5*(geopot(ij,l)+geopot(ij,l+1)) & |
---|
58 | + pk(ij,l)*(cpp-theta(ij,l,1)) ! pk=temperature, theta=entropy |
---|
59 | END DO |
---|
60 | CASE(thermo_moist) |
---|
61 | !DIR$ SIMD |
---|
62 | DO ij=ij_begin,ij_end |
---|
63 | ! du/dt = grad(Bd)+rv.grad(Bv)+s.grad(T) |
---|
64 | ! Bd = Phi + gibbs_d |
---|
65 | ! Bv = Phi + gibbs_v |
---|
66 | ! pk=temperature, theta=entropy |
---|
67 | qv = theta(ij,l,2) |
---|
68 | temp = pk(ij,l) |
---|
69 | chi = log(temp/Treff) |
---|
70 | nu = (chi*(cpp+qv*cppv)-theta(ij,l,1))/(Rd+qv*Rv) ! log(p/preff) |
---|
71 | berni(ij,l) = .5*(geopot(ij,l)+geopot(ij,l+1)) & |
---|
72 | + temp*(cpp*(1.-chi)+Rd*nu) |
---|
73 | berniv(ij,l) = .5*(geopot(ij,l)+geopot(ij,l+1)) & |
---|
74 | + temp*(cppv*(1.-chi)+Rv*nu) |
---|
75 | END DO |
---|
76 | END SELECT |
---|
77 | END DO |
---|
78 | |
---|
79 | END IF ! Boussinesq/compressible |
---|
80 | |
---|
81 | !!! u:=u+tau*du, du = -grad(B)-theta.grad(pi) |
---|
82 | DO l=ll_begin,ll_end |
---|
83 | IF(caldyn_thermo == thermo_moist) THEN |
---|
84 | !DIR$ SIMD |
---|
85 | DO ij=ij_begin,ij_end |
---|
86 | due_right = berni(ij+t_right,l)-berni(ij,l) & |
---|
87 | + 0.5*(theta(ij,l,1)+theta(ij+t_right,l,1)) & |
---|
88 | *(pk(ij+t_right,l)-pk(ij,l)) & |
---|
89 | + 0.5*(theta(ij,l,2)+theta(ij+t_right,l,2)) & |
---|
90 | *(berniv(ij+t_right,l)-berniv(ij,l)) |
---|
91 | |
---|
92 | due_lup = berni(ij+t_lup,l)-berni(ij,l) & |
---|
93 | + 0.5*(theta(ij,l,1)+theta(ij+t_lup,l,1)) & |
---|
94 | *(pk(ij+t_lup,l)-pk(ij,l)) & |
---|
95 | + 0.5*(theta(ij,l,2)+theta(ij+t_lup,l,2)) & |
---|
96 | *(berniv(ij+t_lup,l)-berniv(ij,l)) |
---|
97 | |
---|
98 | due_ldown = berni(ij+t_ldown,l)-berni(ij,l) & |
---|
99 | + 0.5*(theta(ij,l,1)+theta(ij+t_ldown,l,1)) & |
---|
100 | *(pk(ij+t_ldown,l)-pk(ij,l)) & |
---|
101 | + 0.5*(theta(ij,l,2)+theta(ij+t_ldown,l,2)) & |
---|
102 | *(berniv(ij+t_ldown,l)-berniv(ij,l)) |
---|
103 | |
---|
104 | du(ij+u_right,l) = du(ij+u_right,l) - ne_right*due_right |
---|
105 | du(ij+u_lup,l) = du(ij+u_lup,l) - ne_lup*due_lup |
---|
106 | du(ij+u_ldown,l) = du(ij+u_ldown,l) - ne_ldown*due_ldown |
---|
107 | u(ij+u_right,l) = u(ij+u_right,l) + tau*du(ij+u_right,l) |
---|
108 | u(ij+u_lup,l) = u(ij+u_lup,l) + tau*du(ij+u_lup,l) |
---|
109 | u(ij+u_ldown,l) = u(ij+u_ldown,l) + tau*du(ij+u_ldown,l) |
---|
110 | END DO |
---|
111 | ELSE |
---|
112 | !DIR$ SIMD |
---|
113 | DO ij=ij_begin,ij_end |
---|
114 | due_right = 0.5*(theta(ij,l,1)+theta(ij+t_right,l,1)) & |
---|
115 | *(pk(ij+t_right,l)-pk(ij,l)) & |
---|
116 | + berni(ij+t_right,l)-berni(ij,l) |
---|
117 | due_lup = 0.5*(theta(ij,l,1)+theta(ij+t_lup,l,1)) & |
---|
118 | *(pk(ij+t_lup,l)-pk(ij,l)) & |
---|
119 | + berni(ij+t_lup,l)-berni(ij,l) |
---|
120 | due_ldown = 0.5*(theta(ij,l,1)+theta(ij+t_ldown,l,1)) & |
---|
121 | *(pk(ij+t_ldown,l)-pk(ij,l)) & |
---|
122 | + berni(ij+t_ldown,l)-berni(ij,l) |
---|
123 | du(ij+u_right,l) = du(ij+u_right,l) - ne_right*due_right |
---|
124 | du(ij+u_lup,l) = du(ij+u_lup,l) - ne_lup*due_lup |
---|
125 | du(ij+u_ldown,l) = du(ij+u_ldown,l) - ne_ldown*due_ldown |
---|
126 | u(ij+u_right,l) = u(ij+u_right,l) + tau*du(ij+u_right,l) |
---|
127 | u(ij+u_lup,l) = u(ij+u_lup,l) + tau*du(ij+u_lup,l) |
---|
128 | u(ij+u_ldown,l) = u(ij+u_ldown,l) + tau*du(ij+u_ldown,l) |
---|
129 | END DO |
---|
130 | END IF |
---|
131 | END DO |
---|
132 | |
---|
133 | END IF ! dysl |
---|
134 | CALL trace_end("compute_caldyn_fast") |
---|
135 | |
---|
136 | END SUBROUTINE compute_caldyn_fast |
---|
137 | |
---|
138 | END MODULE compute_caldyn_fast_mod |
---|