Changeset 1363 for LMDZ4/branches/LMDZ4V5.0-dev/libf/dyn3dpar/iniacademic.F

Timestamp:

Apr 16, 2010, 11:50:10 AM (14 years ago)

Author:

Ehouarn Millour

Message:

Added possibility to run in "Shallow Water" mode. SW mode is automatically set if the number of atmospheric layers is 1.
To use SW mode, the model should be compiled without physics (makelmdz_fcm -p nophys) and/or without calls to the physics (i.e. set iflag_phys=0 in gcm.def).

-Updated some definitions and comments in run.def & gcm.def

-Fixed misplaced #ifdef CPP_EARTH in calfis.F + some write(lunout,*)

-Specific initialisation of fields for SW are done in sw_case_williamson91_6 (called by iniacademic, when read_start=.false.)

• Checked (on Vargas & Brodie) that these changes don't alter usual bench GCM outputs.

EM

File:

• LMDZ4/branches/LMDZ4V5.0-dev/libf/dyn3dpar/iniacademic.F (modified) (8 diffs)

LMDZ4/branches/LMDZ4V5.0-dev/libf/dyn3dpar/iniacademic.F

@@ -8 +8 @@
       USE filtreg_mod
       USE infotrac, ONLY : nqtot
+      USE control_mod
+ 
 
 c%W%    %G%
@@ -31 +33 @@
 c
 c=======================================================================
-      USE control_mod
       IMPLICIT NONE
 c-----------------------------------------------------------------------
@@ -46 +47 @@
 #include "temps.h"
 #include "iniprint.h"
+#include "logic.h"
 
 c   Arguments:
@@ -55 +57 @@
       REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
       REAL teta(ip1jmp1,llm)                 ! temperature potentielle
-      REAL q(ip1jmp1,llm,nqtot)              ! champs advectes
+      REAL q(ip1jmp1,llm,nqtot)               ! champs advectes
       REAL ps(ip1jmp1)                       ! pression  au sol
       REAL masse(ip1jmp1,llm)                ! masse d'air
@@ -84 +86 @@
         time_0=0.
         day_ref=0
-        annee_ref=0
+        annee_ref=0
 
         im         = iim
@@ -105 +107 @@
         ang0       = 0.
 
+        if (llm.eq.1) then
+          ! specific initializations for the shallow water case
+          kappa=1
+        endif
+        
         CALL iniconst
         CALL inigeom
         CALL inifilr
 
-        ps=0.
-        phis=0.
+        if (llm.eq.1) then
+          ! initialize fields for the shallow water case, if required
+          if (.not.read_start) then
+            phis(:)=0.
+            q(:,:,1)=1.e-10
+            q(:,:,2)=1.e-15
+            q(:,:,3:nqtot)=0.
+            CALL sw_case_williamson91_6(vcov,ucov,teta,masse,ps)
+          endif
+        endif
+
+        if (iflag_phys.eq.2) then
+          ! initializations for the academic case
+          ps(:)=1.e5
+          phis(:)=0.
 c---------------------------------------------------------------------
 
-        taurappel=10.*daysec
+          taurappel=10.*daysec
 
 c---------------------------------------------------------------------
@@ -119 +139 @@
 c   --------------------------------------
 
-        DO l=1,llm
-         zsig=ap(l)/preff+bp(l)
-         if (zsig.gt.0.3) then
-           lsup=l
-           tetarappell=1./8.*(-log(zsig)-.5)
-           DO j=1,jjp1
+          DO l=1,llm
+            zsig=ap(l)/preff+bp(l)
+            if (zsig.gt.0.3) then
+             lsup=l
+             tetarappell=1./8.*(-log(zsig)-.5)
+             DO j=1,jjp1
              ddsin=sin(rlatu(j))-sin(pi/20.)
              tetajl(j,l)=300.*(1+1./18.*(1.-3.*ddsin*ddsin)+tetarappell)
-           ENDDO
-          else
+             ENDDO
+            else
 c   Choix isotherme au-dessus de 300 mbar
-           do j=1,jjp1
-             tetajl(j,l)=tetajl(j,lsup)*(0.3/zsig)**kappa
-           enddo
-          endif ! of if (zsig.gt.0.3)
-        ENDDO ! of DO l=1,llm
-
-        do l=1,llm
-           do j=1,jjp1
+             do j=1,jjp1
+               tetajl(j,l)=tetajl(j,lsup)*(0.3/zsig)**kappa
+             enddo
+            endif ! of if (zsig.gt.0.3)
+          ENDDO ! of DO l=1,llm
+
+          do l=1,llm
+            do j=1,jjp1
               do i=1,iip1
                  ij=(j-1)*iip1+i
                  tetarappel(ij,l)=tetajl(j,l)
               enddo
-           enddo
-        enddo
+            enddo
+          enddo
 
 c       call dump2d(jjp1,llm,tetajl,'TEQ   ')
 
-        ps=1.e5
-        phis=0.
-        CALL pression ( ip1jmp1, ap, bp, ps, p       )
-        CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
-        CALL massdair(p,masse)
+          CALL pression ( ip1jmp1, ap, bp, ps, p       )
+          CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
+          CALL massdair(p,masse)
 
 c  intialisation du vent et de la temperature
-        teta(:,:)=tetarappel(:,:)
-        CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
-        call ugeostr(phi,ucov)
-        vcov=0.
-        q(:,:,1   )=1.e-10
-        q(:,:,2   )=1.e-15
-        q(:,:,3:nqtot)=0.
+          teta(:,:)=tetarappel(:,:)
+          CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
+          call ugeostr(phi,ucov)
+          vcov=0.
+          q(:,:,1   )=1.e-10
+          q(:,:,2   )=1.e-15
+          q(:,:,3:nqtot)=0.
 
 
 c   perturbation aleatoire sur la temperature
-        idum  = -1
-        zz = ran1(idum)
-        idum  = 0
-        do l=1,llm
-           do ij=iip2,ip1jm
+          idum  = -1
+          zz = ran1(idum)
+          idum  = 0
+          do l=1,llm
+            do ij=iip2,ip1jm
               teta(ij,l)=teta(ij,l)*(1.+0.005*ran1(idum))
-           enddo
-        enddo
-
-        do l=1,llm
-           do ij=1,ip1jmp1,iip1
+            enddo
+          enddo
+
+          do l=1,llm
+            do ij=1,ip1jmp1,iip1
               teta(ij+iim,l)=teta(ij,l)
-           enddo
-        enddo
+            enddo
+          enddo
 
 
@@ -187 +205 @@
 
 c   initialisation d'un traceur sur une colonne
-        j=jjp1*3/4
-        i=iip1/2
-        ij=(j-1)*iip1+i
-        q(ij,:,3)=1.
-      
+          j=jjp1*3/4
+          i=iip1/2
+          ij=(j-1)*iip1+i
+          q(ij,:,3)=1.
+        endif ! of if (iflag_phys.eq.2)
+        
       else
         write(lunout,*)"iniacademic: planet types other than earth",