Changeset 366 for codes/icosagcm/trunk/src/etat0_dcmip3.f90

Timestamp:

10/30/15 15:41:06 (9 years ago)

Author:

dubos

Message:

Progress towards NH

File:

• codes/icosagcm/trunk/src/etat0_dcmip3.f90 (modified) (3 diffs)

codes/icosagcm/trunk/src/etat0_dcmip3.f90

@@ -7 +7 @@
 CONTAINS
   
-  SUBROUTINE compute_etat0(ngrid,lon,lat, phis,ps,temp,ulon,ulat,q)
+  SUBROUTINE compute_etat0(ngrid,lon,lat, phis,ps,temp,ulon,ulat,geopot,q)
     USE genmod
     USE dcmip_initial_conditions_test_1_2_3
@@ -20 +20 @@
     REAL(rstd), INTENT(OUT) :: ulat(ngrid,llm)
     REAL(rstd), INTENT(OUT) :: temp(ngrid,llm)
+    REAL(rstd), INTENT(OUT) :: geopot(ngrid,llm+1)
     REAL(rstd), INTENT(OUT) :: q(ngrid,llm,nqtot)
     REAL(rstd),PARAMETER :: Peq=1e5        ! Reference surface pressure at the equator (hPa)
-    REAL(rstd) :: dummy, pp
+    REAL(rstd) :: dummy, pp, zz
     INTEGER :: l,ij
     pp=peq
@@ -28 +29 @@
        CALL test3_gravity_wave(lon(ij),lat(ij),pp,dummy,0, &
             dummy,dummy,dummy,dummy,phis(ij),ps(ij),dummy,dummy)
+    END DO
+    DO l=ll_begin,ll_endp1
+       DO ij=1,ngrid
+          pp = ap(l) + bp(l)*ps(ij) ! half-layer pressure
+          CALL test3_gravity_wave(lon(ij),lat(ij),pp,zz,0, &
+               dummy,dummy,dummy,dummy,dummy,dummy,dummy,dummy)
+          geopot(ij,l) = g*zz ! initialize geopotential for NH
+       END DO
+       q(:,l,:)=0.
     END DO
     DO l=ll_begin,ll_end