Changeset 836 for codes/icosagcm/devel/src/kernels_unst/compute_geopot.k90

Timestamp:

05/03/19 13:25:21 (5 years ago)

Author:

dubos

Message:

devel : Cp(T) thermodynamics (TBC)

File:

• codes/icosagcm/devel/src/kernels_unst/compute_geopot.k90 (modified) (1 diff)

codes/icosagcm/devel/src/kernels_unst/compute_geopot.k90

@@ -70 +70 @@
       END DO
       !$OMP END DO
+   CASE(thermo_variable_Cp)
+      ! thermodynamics with variable Cp
+      ! Cp.dT = dh = Tds + vdp
+      ! pv = RT
+      ! => ds = (dh+v.dp)/T = Cp.dT/T - R dp/p
+      ! Cp(T) = Cp0 * (T/T0)^nu
+      ! => s(p,T) = Cp(T)/nu - R log(p/preff)
+      ! h = Cp(T).T/(nu+1)
+      !$OMP DO SCHEDULE(STATIC)
+      DO ij=1,primal_num
+         pk(llm,ij) = ptop + .5*g* rhodz(llm,ij)
+         DO l = llm-1,1,-1
+            pk(l,ij) = pk(l+1,ij) + (.5*g)*( rhodz(l,ij) + rhodz(l+1,ij) )
+         END DO
+         IF(caldyn_eta == eta_lag) THEN
+            ps(ij) = pk(1,ij) + .5*g* rhodz(1,ij)
+         END IF
+      END DO
+      !$OMP END DO
+      !$OMP DO SCHEDULE(STATIC)
+      DO ij=1,primal_num
+         DO l = 1,llm
+            p_ik = pk(l,ij)
+            Cp_ik = nu*( theta(l,ij,1) + Rd*log(p_ik/preff) )
+            temp_ik = Treff* (Cp_ik/cpp)**(1./nu)
+            gv = (g*Rd)*temp_ik/p_ik ! specific volume v = Rd*T/p
+            pk(l,ij) = temp_ik
+            geopot(l+1,ij) = geopot(l,ij) + gv*rhodz(l,ij)
+         END DO
+      END DO
+      !$OMP END DO
    CASE(thermo_moist)
       !$OMP DO SCHEDULE(STATIC)