Changeset 1046 for LMDZ4/trunk/libf/dyn3dpar/guide_p.F

Timestamp:

Nov 6, 2008, 5:39:24 PM (16 years ago)

Author:

lmdzadmin

Message:

Modifs guidage pour utiliser des champs de guidage sur niveaux hybrides
Ajout cle logique ok_gradsfile (.false. par defaut) pour activer sorties grads du guidage
FC/IM

File:

• LMDZ4/trunk/libf/dyn3dpar/guide_p.F (modified) (15 diffs)

LMDZ4/trunk/libf/dyn3dpar/guide_p.F

@@ -79 +79 @@
       real alpha_T(ip1jmp1),alpha_P(ip1jmp1)
       real alpha_u(ip1jmp1),alpha_v(ip1jm)
+      real alpha_pcor(llm)
       real dday_step,toto,reste,itau_test
       INTEGER step_rea,count_no_rea
@@ -112 +113 @@
 
       save alpha_T,alpha_q,alpha_u,alpha_v,alpha_P,itau_test
+      save alpha_pcor
       save step_rea,count_no_rea
 
@@ -159 +161 @@
          
          if (mpi_rank==0) then
+         if (ok_gradsfile) then
          call inigrads(igrads,iip1
      s  ,rlonv,180./pi,-180.,180.,jjp1,rlatu,-90.,90.,180./pi
      s  ,llm,presnivs,1.
      s  ,dtgrads,file,'dyn_zon ')
+         endif !ok_gradsfile
          endif
          
@@ -208 +212 @@
 c           physic=.false.
          endif
-
+c correction de rappel dans couche limite
+c F.Codron
+         if (guide_BL) then
+           alpha_pcor(:)=1.
+         else
+           do l=1,llm
+             alpha_pcor(l)=(1.+tanh((0.85-presnivs(l)/preff)/0.05))/2.
+           enddo
+         endif
          itau_test=1001
          step_rea=1
@@ -216 +228 @@
 c lecture d'un fichier netcdf pour determiner le nombre de niveaux
          IF (mpi_rank==0) THEN
-        
+         if (guide_modele) then
+           if (ncidpl.eq.-99) ncidpl=NCOPN('apbp.nc',NCNOWRIT,rcod)
+         else   
          if (guide_u) then
            if (ncidpl.eq.-99) ncidpl=NCOPN('u.nc',NCNOWRIT,rcod)
@@ -233 +247 @@
          endif
 c
+         endif 
          if (ncep) then
           status=NF_INQ_DIMID(ncidpl,'LEVEL',rid)
@@ -239 +254 @@
          endif
           status=NF_INQ_DIMLEN(ncidpl,rid,nlev)
-         print *,'nlev', nlev 
-          call ncclos(ncidpl,rcod)
-         
-         ENDIF
-         
+         print *,'nlev guide', nlev 
+         call ncclos(ncidpl,rcod)
 c   Lecture du premier etat des reanalyses.
          call Gather_Field(ps,ip1jmp1,1,0)
@@ -316 +328 @@
       
        if (mpi_rank==0) then
+        if (ok_gradsfile) then
          call wrgrads(igrads,1,aire   ,'aire      ','aire      ' )
          call wrgrads(igrads,1,dxdys  ,'dxdy      ','dxdy      ' )
@@ -328 +341 @@
          call wrgrads(igrads,llm,q,'Q         ','Q         ' )
          call wrgrads(igrads,llm,qsat,'QSAT      ','QSAT      ' )
+       endif !(ok_gradsfile) then
       endif
       
@@ -357 +371 @@
             do ij=ijb,ije
                 a=(1.-tau)*ucovrea1(ij,l)+tau*ucovrea2(ij,l)
-                ucov(ij,l)=(1.-alpha_u(ij))*ucov(ij,l)+alpha_u(ij)*a
+                ucov(ij,l)=(1.-alpha_u(ij))*ucov(ij,l)+alpha_u(ij)
+     $                     *alpha_pcor(l)*a
                 if (first.and.ini_anal) ucov(ij,l)=a
             enddo
@@ -368 +383 @@
             do ij=ijb,ije
                 a=(1.-tau)*tetarea1(ij,l)+tau*tetarea2(ij,l)
-                teta(ij,l)=(1.-alpha_T(ij))*teta(ij,l)+alpha_T(ij)*a
+                teta(ij,l)=(1.-alpha_T(ij))*teta(ij,l)+alpha_T(ij)
+     $                     *alpha_pcor(l)*a
                 if (first.and.ini_anal) teta(ij,l)=a
             enddo
@@ -393 +409 @@
 c   hum relative en % -> hum specif
                 a=qsat(ij,l)*a*0.01
-                q(ij,l)=(1.-alpha_Q(ij))*q(ij,l)+alpha_Q(ij)*a
+                q(ij,l)=(1.-alpha_Q(ij))*q(ij,l)+alpha_Q(ij)
+     $                     *alpha_pcor(l)*a
                 if (first.and.ini_anal) q(ij,l)=a
             enddo
@@ -406 +423 @@
             do ij=ijb,ije
                 a=(1.-tau)*vcovrea1(ij,l)+tau*vcovrea2(ij,l)
-                vcov(ij,l)=(1.-alpha_v(ij))*vcov(ij,l)+alpha_v(ij)*a
+                vcov(ij,l)=(1.-alpha_v(ij))*vcov(ij,l)+alpha_v(ij)
+     $                     *alpha_pcor(l)*a
                 if (first.and.ini_anal) vcov(ij,l)=a
             enddo
@@ -493 +511 @@
          enddo
 
-         call dump2d(iip1,jjp1,dxdys,'DX2DY2 SCAL  ')
-         call dump2d(iip1,jjp1,dxdyu,'DX2DY2 U     ')
+c         call dump2d(iip1,jjp1,dxdys,'DX2DY2 SCAL  ')
+c         call dump2d(iip1,jjp1,dxdyu,'DX2DY2 U     ')
          call dump2d(iip1,jjp1,dxdyv,'DX2DY2 v     ')
 
@@ -520 +538 @@
               stop
             endif
+            gamma=log(0.5)/log(gamma)
+            if (gamma4) then 
+              gamma=min(gamma,4.)
+            endif
             print*,'gamma=',gamma
-            gamma=log(0.5)/log(gamma)
          endif
       endif