Changeset 858
- Timestamp:
- 05/06/19 16:53:07 (5 years ago)
- Location:
- codes/icosagcm/devel
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/devel/Python/src/kernels_caldyn_hevi.jin
r852 r858 127 127 END_BLOCK 128 128 129 KERNEL(caldyn_fast)130 !131 SELECT CASE(caldyn_thermo)132 CASE(thermo_boussinesq)133 FORALL_CELLS()134 ON_PRIMAL135 berni(CELL) = pk(CELL)136 ! from now on pk contains the vertically-averaged geopotential137 pk(CELL) = .5*(geopot(CELL)+geopot(UP(CELL)))138 END_BLOCK139 END_BLOCK140 CASE(thermo_theta)141 FORALL_CELLS()142 ON_PRIMAL143 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL)))144 END_BLOCK145 END_BLOCK146 CASE(thermo_entropy)147 FORALL_CELLS()148 ON_PRIMAL149 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL)))150 berni(CELL) = berni(CELL) + pk(CELL)*(cpp-theta(CELL,1)) ! Gibbs = Cp.T-Ts = T(Cp-s)151 END_BLOCK152 END_BLOCK153 CASE(thermo_variable_Cp)154 ! thermodynamics with variable Cp155 ! Cp(T) = Cp0 * (T/T0)^nu156 ! => h = Cp(T).T/(nu+1)157 158 FORALL_CELLS()159 ON_PRIMAL160 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL)))161 cp_ik = cpp*(pk(CELL)/Treff)**nu162 berni(CELL) = berni(CELL) + pk(CELL)*(cp_ik/(nu+1.)-theta(CELL,1)) ! Gibbs = h-Ts = T(Cp/(nu+1)-s)163 END_BLOCK164 END_BLOCK165 CASE DEFAULT166 PRINT *, 'Unsupported value of caldyn_thermo : ',caldyn_thermo ! FIXME167 STOP168 END SELECT169 !170 FORALL_CELLS()171 ON_EDGES172 due = .5*(theta(CELL1,1)+theta(CELL2,1))*(pk(CELL2)-pk(CELL1)) + berni(CELL2)-berni(CELL1)173 du(EDGE) = du(EDGE) - SIGN*due174 u(EDGE) = u(EDGE) + tau*du(EDGE)175 END_BLOCK176 END_BLOCK177 !178 END_BLOCK179 180 129 KERNEL(coriolis) 181 130 ! -
codes/icosagcm/devel/src/dynamics/compute_caldyn_fast.F90
r855 r858 4 4 PRIVATE 5 5 6 #include "../unstructured/unstructured.h90" 7 6 8 PUBLIC :: compute_caldyn_fast 7 9 8 10 CONTAINS 11 12 #ifdef BEGIN_DYSL 13 14 KERNEL(caldyn_fast) 15 ! 16 SELECT CASE(caldyn_thermo) 17 CASE(thermo_boussinesq) 18 FORALL_CELLS() 19 ON_PRIMAL 20 berni(CELL) = pk(CELL) 21 ! from now on pk contains the vertically-averaged geopotential 22 pk(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 23 END_BLOCK 24 END_BLOCK 25 CASE(thermo_theta) 26 FORALL_CELLS() 27 ON_PRIMAL 28 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 29 END_BLOCK 30 END_BLOCK 31 CASE(thermo_entropy) 32 FORALL_CELLS() 33 ON_PRIMAL 34 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 35 berni(CELL) = berni(CELL) + pk(CELL)*(cpp-theta(CELL,1)) ! Gibbs = Cp.T-Ts = T(Cp-s) 36 END_BLOCK 37 END_BLOCK 38 CASE(thermo_variable_Cp) 39 ! thermodynamics with variable Cp 40 ! Cp(T) = Cp0 * (T/T0)^nu 41 ! => h = Cp(T).T/(nu+1) 42 43 FORALL_CELLS() 44 ON_PRIMAL 45 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 46 cp_ik = cpp*(pk(CELL)/Treff)**nu 47 berni(CELL) = berni(CELL) + pk(CELL)*(cp_ik/(nu+1.)-theta(CELL,1)) ! Gibbs = h-Ts = T(Cp/(nu+1)-s) 48 END_BLOCK 49 END_BLOCK 50 CASE DEFAULT 51 PRINT *, 'Unsupported value of caldyn_thermo : ',caldyn_thermo ! FIXME 52 STOP 53 END SELECT 54 ! 55 FORALL_CELLS() 56 ON_EDGES 57 due = .5*(theta(CELL1,1)+theta(CELL2,1))*(pk(CELL2)-pk(CELL1)) + berni(CELL2)-berni(CELL1) 58 du(EDGE) = du(EDGE) - SIGN*due 59 u(EDGE) = u(EDGE) + tau*du(EDGE) 60 END_BLOCK 61 END_BLOCK 62 ! 63 END_BLOCK 64 65 #endif END_DYSL 66 67 SUBROUTINE compute_caldyn_fast_unst(tau, pk,berni,theta,geopot, du,u) 68 USE ISO_C_BINDING, only : C_DOUBLE, C_FLOAT 69 USE earth_const 70 USE grid_param, ONLY : nqdyn 71 USE data_unstructured_mod, ONLY : enter_trace, exit_trace, & 72 id_fast, primal_num, dual_num, edge_num, & 73 dual_deg, dual_edge, dual_ne, dual_vertex, & 74 up, down, left, right, Av, fv, Riv2 75 NUM, INTENT(IN) :: tau 76 FIELD_MASS :: pk,berni ! INOUT, OUT 77 FIELD_THETA :: theta ! IN 78 FIELD_GEOPOT :: geopot ! IN 79 FIELD_U :: u,du ! INOUT,INOUT 80 DECLARE_INDICES 81 DECLARE_EDGES 82 NUM :: due, cp_ik 83 START_TRACE(id_fast, 4,0,2) ! primal, dual, edge 84 #include "../kernels_unst/caldyn_fast.k90" 85 STOP_TRACE 86 END SUBROUTINE compute_caldyn_fast_unst 9 87 10 88 SUBROUTINE compute_caldyn_fast(tau,u,rhodz,theta,pk,geopot,du)
Note: See TracChangeset
for help on using the changeset viewer.