Changeset 935 for codes/icosagcm/devel/src/dynamics
- Timestamp:
- 07/03/19 17:15:11 (5 years ago)
- Location:
- codes/icosagcm/devel/src/dynamics
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/devel/src/dynamics/caldyn_hevi.f90
r928 r935 5 5 USE compute_caldyn_vert_mod, ONLY : compute_caldyn_vert => compute_caldyn_vert_manual 6 6 USE compute_caldyn_vert_NH_mod, ONLY : compute_caldyn_vert_NH => compute_caldyn_vert_NH_manual 7 USE compute_theta_mod, ONLY : compute_theta => compute_theta_manual8 USE compute_geopot_mod, ONLY : compute_geopot => compute_geopot_manual9 7 USE compute_caldyn_kv_mod, ONLY : compute_caldyn_kv 10 8 USE compute_caldyn_Coriolis_mod, ONLY : compute_caldyn_Coriolis => compute_caldyn_Coriolis_manual … … 12 10 USE compute_caldyn_slow_NH_mod, ONLY : compute_caldyn_slow_NH 13 11 USE compute_caldyn_solver_mod, ONLY : compute_caldyn_solver 14 USE compute_caldyn_fast_mod, ONLY : compute_caldyn_fast => compute_caldyn_fast_manual15 12 USE compute_NH_geopot_mod, ONLY : compute_NH_geopot 16 13 IMPLICIT NONE … … 35 32 USE output_field_mod 36 33 USE checksum_mod 37 USE compute_caldyn_mod, ONLY : compute_pvort_only 34 USE compute_caldyn_mod, ONLY : compute_pvort_only, compute_theta, & 35 compute_geopot, compute_caldyn_fast 38 36 IMPLICIT NONE 39 37 LOGICAL,INTENT(IN) :: write_out -
codes/icosagcm/devel/src/dynamics/compute_caldyn.f90
r884 r935 1 1 MODULE compute_caldyn_mod 2 USE prec, ONLY : rstd 2 3 IMPLICIT NONE 3 4 SAVE 5 6 ! fake array dimensions, for interfaces 7 INTEGER, PARAMETER :: iim_jjm_i=1, iim_jjm_u=1, iim_jjm_v =1, llm_=1, llm1=1, nqdyn_=1 4 8 5 9 INTERFACE 10 11 SUBROUTINE comp_pvort_only(u,rhodz,qu,qv, hv) 12 IMPORT 13 REAL(rstd),INTENT(IN) :: u(iim_jjm_i, llm_) 14 REAL(rstd),INTENT(INOUT) :: rhodz(iim_jjm_i, llm_) 15 REAL(rstd),INTENT(OUT) :: qu(iim_jjm_u, llm_) 16 REAL(rstd),INTENT(OUT) :: qv(iim_jjm_v, llm_) 17 REAL(rstd),INTENT(OUT) :: hv(iim_jjm_v, llm_) 18 END SUBROUTINE comp_pvort_only 6 19 7 SUBROUTINE comp_pvort_only(u,rhodz,qu,qv, hv_) 8 USE prec, ONLY : rstd 9 REAL(rstd),INTENT(IN) :: u(1,1) 10 REAL(rstd),INTENT(INOUT) :: rhodz(1,1) 11 REAL(rstd),INTENT(OUT) :: qu(1,1) 12 REAL(rstd),INTENT(OUT) :: qv(1,1) 13 REAL(rstd),INTENT(OUT) :: hv_(1,1) 14 END SUBROUTINE comp_pvort_only 20 SUBROUTINE comp_theta(mass_col,theta_rhodz, rhodz,theta) 21 IMPORT 22 REAL(rstd),INTENT(IN) :: mass_col(iim_jjm_i) 23 REAL(rstd),INTENT(IN) :: theta_rhodz(iim_jjm_i, llm_, nqdyn_) 24 REAL(rstd),INTENT(INOUT) :: rhodz(iim_jjm_i, llm_) 25 REAL(rstd),INTENT(OUT) :: theta(iim_jjm_i, llm_, nqdyn_) 26 END SUBROUTINE comp_theta 27 28 SUBROUTINE comp_geopot(rhodz,theta, ps,pk,geopot) 29 IMPORT 30 REAL(rstd),INTENT(IN) :: rhodz(iim_jjm_i, llm_) ! rho*dz = mass per unit surface in each full model level 31 REAL(rstd),INTENT(IN) :: theta(iim_jjm_i, llm_, nqdyn_) ! active scalars : theta/entropy, moisture, ... 32 REAL(rstd),INTENT(INOUT) :: ps(iim_jjm_i) ! surface pressure 33 REAL(rstd),INTENT(OUT) :: pk(iim_jjm_i, llm_) ! Exner function (compressible) /Lagrange multiplier (Boussinesq) 34 REAL(rstd),INTENT(INOUT) :: geopot(iim_jjm_i, llm1) ! geopotential 35 END SUBROUTINE comp_geopot 36 37 SUBROUTINE comp_caldyn_fast(tau,theta,geopot, pk,berni,du,u) 38 IMPORT 39 REAL(rstd),INTENT(IN) :: tau ! "solve" u-tau*du/dt = rhs 40 REAL(rstd),INTENT(IN) :: theta(iim_jjm_i, llm_, nqdyn_) 41 REAL(rstd),INTENT(IN) :: geopot(iim_jjm_i, llm1) 42 REAL(rstd),INTENT(INOUT) :: pk(iim_jjm_i, llm_) 43 REAL(rstd),INTENT(INOUT) :: berni(iim_jjm_i, llm_) ! partial Bernoulli function 44 REAL(rstd),INTENT(INOUT) :: du(iim_jjm_u, llm_) 45 REAL(rstd),INTENT(INOUT) :: u(iim_jjm_u, llm_) ! INOUT if tau>0 46 END SUBROUTINE comp_caldyn_fast 47 48 SUBROUTINE comp_caldyn_slow_hydro(zero, u,rhodz,hv,Kv, berni, hflux,du) 49 IMPORT 50 LOGICAL, INTENT(IN) :: zero 51 REAL(rstd),INTENT(IN) :: u(iim_jjm_u, llm_) ! prognostic "velocity" 52 REAL(rstd),INTENT(IN) :: rhodz(iim_jjm_i, llm_) 53 REAL(rstd),INTENT(IN) :: hv(iim_jjm_v, llm_) ! height/mass averaged to vertices 54 REAL(rstd),INTENT(IN) :: Kv(iim_jjm_v, llm_) ! kinetic energy at vertices 55 REAL(rstd), INTENT(OUT) :: berni(iim_jjm_i, llm_) ! Bernoulli function 56 REAL(rstd),INTENT(OUT) :: hflux(iim_jjm_u, llm_) ! hflux in kg/s 57 REAL(rstd),INTENT(INOUT) :: du(iim_jjm_u, llm_) 58 END SUBROUTINE comp_caldyn_slow_hydro 15 59 16 60 END INTERFACE 17 61 18 PROCEDURE(comp_pvort_only), POINTER :: compute_pvort_only => NULL() 62 PROCEDURE(comp_pvort_only), POINTER :: compute_pvort_only => NULL() 63 PROCEDURE(comp_theta), POINTER :: compute_theta => NULL() 64 PROCEDURE(comp_geopot), POINTER :: compute_geopot => NULL() 65 PROCEDURE(comp_caldyn_fast), POINTER :: compute_caldyn_fast => NULL() 66 PROCEDURE(comp_caldyn_slow_hydro), POINTER :: compute_caldyn_slow_hydro => NULL() 19 67 20 68 END MODULE compute_caldyn_mod -
codes/icosagcm/devel/src/dynamics/compute_caldyn_fast.F90
r921 r935 17 17 #ifdef BEGIN_DYSL 18 18 19 KERNEL(caldyn_fast) 20 ! 19 {% macro case_caldyn_fast(name) %} 20 CASE({{name}}) 21 FORALL_CELLS() 22 ON_PRIMAL 23 Phi_ik = .5*(geopot(CELL)+geopot(UP(CELL))) 24 {{ caller() }} 25 END_BLOCK 26 END_BLOCK 27 {%- endmacro %} 28 29 KERNEL(caldyn_fast) 30 ! 21 31 SELECT CASE(caldyn_thermo) 22 32 CASE(thermo_boussinesq) 23 FORALL_CELLS() 24 ON_PRIMAL 25 berni(CELL) = pk(CELL) 26 ! from now on pk contains the vertically-averaged geopotential 27 pk(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 28 END_BLOCK 29 END_BLOCK 30 CASE(thermo_theta) 31 FORALL_CELLS() 32 ON_PRIMAL 33 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 34 END_BLOCK 35 END_BLOCK 36 CASE(thermo_entropy) 37 FORALL_CELLS() 38 ON_PRIMAL 39 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 40 berni(CELL) = berni(CELL) + pk(CELL)*(cpp-theta(CELL,1)) ! Gibbs = Cp.T-Ts = T(Cp-s) 41 END_BLOCK 42 END_BLOCK 43 CASE(thermo_variable_Cp) 33 FORALL_CELLS() 34 ON_PRIMAL 35 berni(CELL) = pk(CELL) 36 ! from now on pk contains the vertically-averaged geopotential 37 pk(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 38 END_BLOCK 39 END_BLOCK 40 {% call case_caldyn_fast('thermo_theta') %} 41 berni(CELL) = Phi_ik 42 {%- endcall %} 43 {% call case_caldyn_fast('thermo_entropy') %} 44 berni(CELL) = Phi_ik + pk(CELL)*(cpp-theta(CELL,1)) ! Gibbs = Cp.T-Ts = T(Cp-s) 45 {%- endcall %} 46 {% call case_caldyn_fast('thermo_variable_Cp') %} 44 47 ! thermodynamics with variable Cp 45 48 ! Cp(T) = Cp0 * (T/T0)^nu 46 49 ! => h = Cp(T).T/(nu+1) 47 48 FORALL_CELLS() 49 ON_PRIMAL 50 berni(CELL) = .5*(geopot(CELL)+geopot(UP(CELL))) 51 cp_ik = cpp*(pk(CELL)/Treff)**nu 52 berni(CELL) = berni(CELL) + pk(CELL)*(cp_ik/(nu+1.)-theta(CELL,1)) ! Gibbs = h-Ts = T(Cp/(nu+1)-s) 53 END_BLOCK 54 END_BLOCK 50 cp_ik = cpp*(pk(ij,l)/Treff)**nu 51 berni(CELL) = Phi_ik + pk(CELL)*(cp_ik/(nu+1.)-theta(CELL,1)) ! Gibbs = h-Ts = T(Cp/(nu+1)-s) 52 {%- endcall %} 55 53 CASE DEFAULT 56 54 PRINT *, 'Unsupported value of caldyn_thermo : ',caldyn_thermo ! FIXME … … 83 81 DECLARE_INDICES 84 82 DECLARE_EDGES 85 NUM :: due, cp_ik 83 NUM :: due, cp_ik, Phi_ik 86 84 START_TRACE(id_fast, 4,0,2) ! primal, dual, edge 87 85 #include "../kernels_unst/caldyn_fast.k90" … … 101 99 102 100 INTEGER :: ij,l 103 REAL(rstd) :: cp_ik, qv, temp, chi, log_p_preff, due, due_right, due_lup, due_ldown101 REAL(rstd) :: due, cp_ik, Phi_ik 104 102 105 103 CALL trace_start("compute_caldyn_fast") 106 104 #include "../kernels_hex/caldyn_fast.k90" 107 105 CALL trace_end("compute_caldyn_fast") 106 108 107 END SUBROUTINE compute_caldyn_fast_hex 109 108 -
codes/icosagcm/devel/src/dynamics/compute_caldyn_slow_hydro.F90
r921 r935 73 73 LOGICAL, INTENT(IN) :: zero 74 74 REAL(rstd),INTENT(IN) :: u(3*iim*jjm,llm) ! prognostic "velocity" 75 REAL(rstd),INTENT(IN) :: rhodz(iim*jjm,llm) 76 REAL(rstd),INTENT(IN) :: hv(2*iim*jjm,llm) ! height/mass averaged to vertices 75 77 REAL(rstd),INTENT(IN) :: Kv(2*iim*jjm,llm) ! kinetic energy at vertices 76 REAL(rstd),INTENT(IN) :: hv(2*iim*jjm,llm) ! height/mass averaged to vertices 77 REAL(rstd),INTENT(IN) :: rhodz(iim*jjm,llm) 78 REAL(rstd), INTENT(OUT) :: berni(iim*jjm,llm) ! Bernoulli function 78 79 REAL(rstd),INTENT(OUT) :: hflux(3*iim*jjm,llm) ! hflux in kg/s 79 80 REAL(rstd),INTENT(INOUT) :: du(3*iim*jjm,llm) 80 81 REAL(rstd) :: berni(iim*jjm,llm) ! Bernoulli function82 81 REAL(rstd) :: berni1(iim*jjm) ! Bernoulli function 83 82 REAL(rstd) :: uu_right, uu_lup, uu_ldown, ke, uu
Note: See TracChangeset
for help on using the changeset viewer.