Ignore:
Timestamp:
Apr 25, 2014, 12:20:14 PM (10 years ago)
Author:
lguez
Message:

Removed unused variables pks, pk, pkf from main program unit gcm.

Encapsulated procedures exner_hyb, exner_hyb_p, exner_hyb_loc,
exner_milieu, exner_milieu_p and exner_milieu_loc into
modules. (Compulsory to allow optional arguments.)

In the procedures exner_hyb, exner_hyb_p, exner_hyb_loc, donwgraded
arguments alpha and beta to local variables. In exner_milieu,
exner_milieu_p and exner_milieu_loc, removed beta altogether. In the
six procedures exner_*, made pkf an optional argument. Made some
cosmetic modifications in order to keep the six procedures exner_* as
close as possible.

In the six procedures exner_*, removed the averaging of pks at the
poles: this is not useful because ps is already the same at all
longitudes at the poles. This modification changes the results of the
program. Motivation: had to do this for exner_hyb because we call it
from test_disvert with a few surface pressure values.

In all the procedures calling exner_*, removed the variables alpha and
beta. Also removed variables alpha and beta from module leapfrog_mod
and from module call_calfis_mod.

Removed actual argument pkf in call to exner_hyb* and exner_milieu*
from guide_interp, guide_main, iniacademic and iniacademic_loc (pkf
was not used in those procedures).

Argument workvar of startget_dyn is used only if varname is tpot or

  1. When varname is tpot or q, the actual argument associated to

workvar in etat0_netcdf is not y. So y in etat0_netcdf is a
place-holder, never used. So we remove optional argument y in the
calls to exner_hyb and exner_milieu from etat0_netcdf.

Created procedure test_disvert, called only by etat0_netcdf. This
procedure tests the order of pressure values at half-levels and full
levels.

Location:
LMDZ5/trunk/libf/dyn3d_common
Files:
2 moved

Legend:

Unmodified
Added
Removed
  • LMDZ5/trunk/libf/dyn3d_common/exner_hyb_m.F90

    r1992 r2021  
    1 !
    2 ! $Id $
    3 !
    4       SUBROUTINE  exner_hyb ( ngrid, ps, p,alpha,beta, pks, pk, pkf )
    5 c
    6 c     Auteurs :  P.Le Van  , Fr. Hourdin  .
    7 c    ..........
    8 c
    9 c    ....  ngrid, ps,p             sont des argum.d'entree  au sous-prog ...
    10 c    .... alpha,beta, pks,pk,pkf   sont des argum.de sortie au sous-prog ...
    11 c
    12 c   ************************************************************************
    13 c    Calcule la fonction d'Exner pk = Cp * p ** kappa , aux milieux des
    14 c    couches .   Pk(l) sera calcule aux milieux  des couches l ,entre les
    15 c    pressions p(l) et p(l+1) ,definis aux interfaces des llm couches .
    16 c   ************************************************************************
    17 c  .. N.B : Au sommet de l'atmosphere,  p(llm+1) = 0. , et ps et pks sont
    18 c    la pression et la fonction d'Exner  au  sol  .
    19 c
    20 c                                 -------- z                                   
    21 c    A partir des relations  ( 1 ) p*dz(pk) = kappa *pk*dz(p)      et
    22 c                            ( 2 ) pk(l) = alpha(l)+ beta(l)*pk(l-1)
    23 c    ( voir note de Fr.Hourdin )  ,
    24 c
    25 c    on determine successivement , du haut vers le bas des couches, les
    26 c    coef. alpha(llm),beta(llm) .,.,alpha(l),beta(l),,,alpha(2),beta(2),
    27 c    puis pk(ij,1). Ensuite ,on calcule,du bas vers le haut des couches, 
    28 c     pk(ij,l)  donne  par la relation (2),  pour l = 2 a l = llm .
    29 c
    30 c
    31       IMPLICIT NONE
    32 c
    33 #include "dimensions.h"
    34 #include "paramet.h"
    35 #include "comconst.h"
    36 #include "comgeom.h"
    37 #include "comvert.h"
    38 #include "serre.h"
     1module exner_hyb_m
    392
    40       INTEGER  ngrid
    41       REAL p(ngrid,llmp1),pk(ngrid,llm),pkf(ngrid,llm)
    42       REAL ps(ngrid),pks(ngrid), alpha(ngrid,llm),beta(ngrid,llm)
     3  IMPLICIT NONE
    434
    44 c    .... variables locales   ...
     5contains
    456
    46       INTEGER l, ij
    47       REAL unpl2k,dellta
     7  SUBROUTINE  exner_hyb ( ngrid, ps, p, pks, pk, pkf )
    488
    49       REAL ppn(iim),pps(iim)
    50       REAL xpn, xps
    51       REAL SSUM
    52 c
    53       logical,save :: firstcall=.true.
    54       character(len=*),parameter :: modname="exner_hyb"
    55      
    56       ! Sanity check
    57       if (firstcall) then
    58         ! sanity checks for Shallow Water case (1 vertical layer)
    59         if (llm.eq.1) then
     9    !     Auteurs :  P.Le Van  , Fr. Hourdin  .
     10    !    ..........
     11    !
     12    !    ....  ngrid, ps,p             sont des argum.d'entree  au sous-prog ...
     13    !    ....  pks,pk,pkf   sont des argum.de sortie au sous-prog ...
     14    !
     15    !   ************************************************************************
     16    !    Calcule la fonction d'Exner pk = Cp * (p/preff) ** kappa , aux milieux des
     17    !    couches .   Pk(l) sera calcule aux milieux  des couches l ,entre les
     18    !    pressions p(l) et p(l+1) ,definis aux interfaces des llm couches .
     19    !   ************************************************************************
     20    !  .. N.B : Au sommet de l'atmosphere,  p(llm+1) = 0. , et ps et pks sont
     21    !    la pression et la fonction d'Exner  au  sol  .
     22    !
     23    !                                 -------- z
     24    !    A partir des relations  ( 1 ) p*dz(pk) = kappa *pk*dz(p)      et
     25    !                            ( 2 ) pk(l) = alpha(l)+ beta(l)*pk(l-1)
     26    !    ( voir note de Fr.Hourdin )  ,
     27    !
     28    !    on determine successivement , du haut vers le bas des couches, les
     29    !    coef. alpha(llm),beta(llm) .,.,alpha(l),beta(l),,,alpha(2),beta(2),
     30    !    puis pk(ij,1). Ensuite ,on calcule,du bas vers le haut des couches, 
     31    !     pk(ij,l)  donne  par la relation (2),  pour l = 2 a l = llm .
     32    !
     33    !
     34    !
     35    include "dimensions.h"
     36    include "paramet.h"
     37    include "comconst.h"
     38    include "comgeom.h"
     39    include "comvert.h"
     40    include "serre.h"
     41
     42    INTEGER  ngrid
     43    REAL p(ngrid,llmp1),pk(ngrid,llm)
     44    real, optional:: pkf(ngrid,llm)
     45    REAL ps(ngrid),pks(ngrid), alpha(ngrid,llm),beta(ngrid,llm)
     46
     47    !    .... variables locales   ...
     48
     49    INTEGER l, ij
     50    REAL unpl2k,dellta
     51
     52    logical,save :: firstcall=.true.
     53    character(len=*),parameter :: modname="exner_hyb"
     54
     55    ! Sanity check
     56    if (firstcall) then
     57       ! sanity checks for Shallow Water case (1 vertical layer)
     58       if (llm.eq.1) then
    6059          if (kappa.ne.1) then
    61             call abort_gcm(modname,
    62      &      "kappa!=1 , but running in Shallow Water mode!!",42)
     60             call abort_gcm(modname, &
     61                  "kappa!=1 , but running in Shallow Water mode!!",42)
    6362          endif
    6463          if (cpp.ne.r) then
    65             call abort_gcm(modname,
    66      &      "cpp!=r , but running in Shallow Water mode!!",42)
     64             call abort_gcm(modname, &
     65                  "cpp!=r , but running in Shallow Water mode!!",42)
    6766          endif
    68         endif ! of if (llm.eq.1)
     67       endif ! of if (llm.eq.1)
    6968
    70         firstcall=.false.
    71       endif ! of if (firstcall)
     69       firstcall=.false.
     70    endif ! of if (firstcall)
    7271
    73       if (llm.eq.1) then
    74        
    75         ! Compute pks(:),pk(:),pkf(:)
    76        
    77         DO   ij  = 1, ngrid
    78           pks(ij) = (cpp/preff) * ps(ij)
     72    ! Specific behaviour for Shallow Water (1 vertical layer) case:
     73    if (llm.eq.1) then
     74
     75       ! Compute pks(:),pk(:),pkf(:)
     76
     77       DO   ij  = 1, ngrid
     78          pks(ij) = (cpp/preff) * ps(ij)
    7979          pk(ij,1) = .5*pks(ij)
    80         ENDDO
    81        
    82         CALL SCOPY   ( ngrid * llm, pk, 1, pkf, 1 )
    83         CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
    84        
    85         ! our work is done, exit routine
    86         return
     80       ENDDO
    8781
    88       endif ! of if (llm.eq.1)
     82       if (present(pkf)) then
     83          pkf = pk
     84          CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
     85       end if
    8986
    90 !!!! General case:
    91      
    92       unpl2k    = 1.+ 2.* kappa
    93 c
    94       DO   ij  = 1, ngrid
    95         pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
    96       ENDDO
     87       ! our work is done, exit routine
     88       return
     89    endif ! of if (llm.eq.1)
    9790
    98       DO  ij   = 1, iim
    99         ppn(ij) = aire(   ij   ) * pks(  ij     )
    100         pps(ij) = aire(ij+ip1jm) * pks(ij+ip1jm )
    101       ENDDO
    102       xpn      = SSUM(iim,ppn,1) /apoln
    103       xps      = SSUM(iim,pps,1) /apols
     91    ! General case:
    10492
    105       DO ij   = 1, iip1
    106         pks(   ij     )  =  xpn
    107         pks( ij+ip1jm )  =  xps
    108       ENDDO
    109 c
    110 c
    111 c    .... Calcul des coeff. alpha et beta  pour la couche l = llm ..
    112 c
    113       DO     ij      = 1, ngrid
     93    unpl2k    = 1.+ 2.* kappa
     94
     95    !     -------------
     96    !     Calcul de pks
     97    !     -------------
     98
     99    DO   ij  = 1, ngrid
     100       pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
     101    ENDDO
     102
     103    !    .... Calcul des coeff. alpha et beta  pour la couche l = llm ..
     104    !
     105    DO     ij      = 1, ngrid
    114106       alpha(ij,llm) = 0.
    115107       beta (ij,llm) = 1./ unpl2k
    116       ENDDO
    117 c
    118 c     ... Calcul des coeff. alpha et beta  pour l = llm-1  a l = 2 ...
    119 c
    120       DO l = llm -1 , 2 , -1
    121 c
    122         DO ij = 1, ngrid
    123         dellta = p(ij,l)* unpl2k + p(ij,l+1)* ( beta(ij,l+1)-unpl2k )
    124         alpha(ij,l)  = - p(ij,l+1) / dellta * alpha(ij,l+1)
    125         beta (ij,l)  =   p(ij,l  ) / dellta   
    126         ENDDO
    127 c
    128       ENDDO
    129 c
    130 c  ***********************************************************************
    131 c     .....  Calcul de pk pour la couche 1 , pres du sol  ....
    132 c
    133       DO   ij   = 1, ngrid
    134        pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )  /
    135      *    (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-unpl2k )* p(ij,2) )
    136       ENDDO
    137 c
    138 c    ..... Calcul de pk(ij,l) , pour l = 2 a l = llm  ........
    139 c
    140       DO l = 2, llm
    141         DO   ij   = 1, ngrid
    142          pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1)
    143         ENDDO
    144       ENDDO
    145 c
    146 c
    147       CALL SCOPY   ( ngrid * llm, pk, 1, pkf, 1 )
    148       CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
    149      
     108    ENDDO
     109    !
     110    !     ... Calcul des coeff. alpha et beta  pour l = llm-1  a l = 2 ...
     111    !
     112    DO l = llm -1 , 2 , -1
     113       !
     114       DO ij = 1, ngrid
     115          dellta = p(ij,l)* unpl2k + p(ij,l+1)* ( beta(ij,l+1)-unpl2k )
     116          alpha(ij,l)  = - p(ij,l+1) / dellta * alpha(ij,l+1)
     117          beta (ij,l)  =   p(ij,l  ) / dellta   
     118       ENDDO
     119    ENDDO
    150120
    151       RETURN
    152       END
     121    !  ***********************************************************************
     122    !     .....  Calcul de pk pour la couche 1 , pres du sol  ....
     123    !
     124    DO   ij   = 1, ngrid
     125       pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )  / &
     126            (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-unpl2k )* p(ij,2) )
     127    ENDDO
     128    !
     129    !    ..... Calcul de pk(ij,l) , pour l = 2 a l = llm  ........
     130    !
     131    DO l = 2, llm
     132       DO   ij   = 1, ngrid
     133          pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1)
     134       ENDDO
     135    ENDDO
     136
     137    if (present(pkf)) then
     138       !    calcul de pkf
     139       pkf = pk
     140       CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
     141    end if
     142
     143  END SUBROUTINE exner_hyb
     144
     145end module exner_hyb_m
  • LMDZ5/trunk/libf/dyn3d_common/exner_milieu_m.F90

    r1992 r2021  
    1 !
    2 ! $Id $
    3 !
    4       SUBROUTINE  exner_milieu ( ngrid, ps, p,beta, pks, pk, pkf )
    5 c
    6 c     Auteurs :  F. Forget , Y. Wanherdrick
    7 c P.Le Van  , Fr. Hourdin  .
    8 c    ..........
    9 c
    10 c    ....  ngrid, ps,p             sont des argum.d'entree  au sous-prog ...
    11 c    .... beta, pks,pk,pkf   sont des argum.de sortie au sous-prog ...
    12 c
    13 c   ************************************************************************
    14 c    Calcule la fonction d'Exner pk = Cp * (p/preff) ** kappa , aux milieux des
    15 c    couches .   Pk(l) sera calcule aux milieux  des couches l ,entre les
    16 c    pressions p(l) et p(l+1) ,definis aux interfaces des llm couches .
    17 c   ************************************************************************
    18 c    .. N.B : Au sommet de l'atmosphere,  p(llm+1) = 0. , et ps et pks sont
    19 c    la pression et la fonction d'Exner  au  sol  .
    20 c
    21 c     WARNING : CECI est une version speciale de exner_hyb originale
    22 c               Utilise dans la version martienne pour pouvoir
    23 c               tourner avec des coordonnees verticales complexe
    24 c              => Il ne verifie PAS la condition la proportionalite en
    25 c              energie totale/ interne / potentielle (F.Forget 2001)
    26 c    ( voir note de Fr.Hourdin )  ,
    27 c
    28       IMPLICIT NONE
    29 c
    30 #include "dimensions.h"
    31 #include "paramet.h"
    32 #include "comconst.h"
    33 #include "comgeom.h"
    34 #include "comvert.h"
    35 #include "serre.h"
     1module exner_milieu_m
    362
    37       INTEGER  ngrid
    38       REAL p(ngrid,llmp1),pk(ngrid,llm),pkf(ngrid,llm)
    39       REAL ps(ngrid),pks(ngrid), beta(ngrid,llm)
     3  IMPLICIT NONE
    404
    41 c    .... variables locales   ...
     5contains
    426
    43       INTEGER l, ij
    44       REAL dum1
     7  SUBROUTINE  exner_milieu ( ngrid, ps, p, pks, pk, pkf )
     8    !
     9    !     Auteurs :  F. Forget , Y. Wanherdrick
     10    ! P.Le Van  , Fr. Hourdin  .
     11    !    ..........
     12    !
     13    !    ....  ngrid, ps,p             sont des argum.d'entree  au sous-prog ...
     14    !    ....  pks,pk,pkf   sont des argum.de sortie au sous-prog ...
     15    !
     16    !   ************************************************************************
     17    !    Calcule la fonction d'Exner pk = Cp * (p/preff) ** kappa , aux milieux des
     18    !    couches .   Pk(l) sera calcule aux milieux  des couches l ,entre les
     19    !    pressions p(l) et p(l+1) ,definis aux interfaces des llm couches .
     20    !   ************************************************************************
     21    !  .. N.B : Au sommet de l'atmosphere,  p(llm+1) = 0. , et ps et pks sont
     22    !    la pression et la fonction d'Exner  au  sol  .
     23    !
     24    !     WARNING : CECI est une version speciale de exner_hyb originale
     25    !               Utilise dans la version martienne pour pouvoir
     26    !               tourner avec des coordonnees verticales complexe
     27    !              => Il ne verifie PAS la condition la proportionalite en
     28    !              energie totale/ interne / potentielle (F.Forget 2001)
     29    !    ( voir note de Fr.Hourdin )  ,
     30    !
     31    !
     32    include "dimensions.h"
     33    include "paramet.h"
     34    include "comconst.h"
     35    include "comgeom.h"
     36    include "comvert.h"
     37    include "serre.h"
    4538
    46       REAL ppn(iim),pps(iim)
    47       REAL xpn, xps
    48       REAL SSUM
    49       EXTERNAL SSUM
    50       logical,save :: firstcall=.true.
    51       character(len=*),parameter :: modname="exner_milieu"
     39    INTEGER  ngrid
     40    REAL p(ngrid,llmp1),pk(ngrid,llm)
     41    real, optional:: pkf(ngrid,llm)
     42    REAL ps(ngrid),pks(ngrid)
    5243
    53       ! Sanity check
    54       if (firstcall) then
    55         ! sanity checks for Shallow Water case (1 vertical layer)
    56         if (llm.eq.1) then
     44    !    .... variables locales   ...
     45
     46    INTEGER l, ij
     47    REAL dum1
     48
     49    logical,save :: firstcall=.true.
     50    character(len=*),parameter :: modname="exner_milieu"
     51
     52    ! Sanity check
     53    if (firstcall) then
     54       ! sanity checks for Shallow Water case (1 vertical layer)
     55       if (llm.eq.1) then
    5756          if (kappa.ne.1) then
    58             call abort_gcm(modname,
    59      &      "kappa!=1 , but running in Shallow Water mode!!",42)
     57             call abort_gcm(modname, &
     58                  "kappa!=1 , but running in Shallow Water mode!!",42)
    6059          endif
    6160          if (cpp.ne.r) then
    62             call abort_gcm(modname,
    63      &      "cpp!=r , but running in Shallow Water mode!!",42)
     61             call abort_gcm(modname, &
     62                  "cpp!=r , but running in Shallow Water mode!!",42)
    6463          endif
    65         endif ! of if (llm.eq.1)
     64       endif ! of if (llm.eq.1)
    6665
    67         firstcall=.false.
    68       endif ! of if (firstcall)
     66       firstcall=.false.
     67    endif ! of if (firstcall)
    6968
    70 !!!! Specific behaviour for Shallow Water (1 vertical layer) case:
    71       if (llm.eq.1) then
    72      
    73         ! Compute pks(:),pk(:),pkf(:)
    74        
    75         DO   ij  = 1, ngrid
    76           pks(ij) = (cpp/preff) * ps(ij) 
     69    ! Specific behaviour for Shallow Water (1 vertical layer) case:
     70    if (llm.eq.1) then
     71
     72       ! Compute pks(:),pk(:),pkf(:)
     73
     74       DO   ij  = 1, ngrid
     75          pks(ij) = (cpp/preff) * ps(ij)
    7776          pk(ij,1) = .5*pks(ij)
    78         ENDDO
    79        
    80         CALL SCOPY   ( ngrid * llm, pk, 1, pkf, 1 )
    81         CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
    82        
    83         ! our work is done, exit routine
    84         return
     77       ENDDO
    8578
    86       endif ! of if (llm.eq.1)
     79       if (present(pkf)) then
     80          pkf = pk
     81          CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
     82       end if
    8783
    88 !!!! General case:
     84       ! our work is done, exit routine
     85       return
     86    endif ! of if (llm.eq.1)
    8987
    90 c     -------------
    91 c     Calcul de pks
    92 c     -------------
    93    
    94       DO   ij  = 1, ngrid
    95         pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
    96       ENDDO
     88    ! General case:
    9789
    98       DO  ij   = 1, iim
    99         ppn(ij) = aire(   ij   ) * pks(  ij     )
    100         pps(ij) = aire(ij+ip1jm) * pks(ij+ip1jm )
    101       ENDDO
    102       xpn      = SSUM(iim,ppn,1) /apoln
    103       xps      = SSUM(iim,pps,1) /apols
     90    !     -------------
     91    !     Calcul de pks
     92    !     -------------
    10493
    105       DO ij   = 1, iip1
    106         pks(   ij     )  =  xpn
    107         pks( ij+ip1jm )  =  xps
    108       ENDDO
    109 c
    110 c
    111 c    .... Calcul de pk  pour la couche l
    112 c    --------------------------------------------
    113 c
    114       dum1 = cpp * (2*preff)**(-kappa)
    115       DO l = 1, llm-1
    116         DO   ij   = 1, ngrid
    117          pk(ij,l) = dum1 * (p(ij,l) + p(ij,l+1))**kappa
    118         ENDDO
    119       ENDDO
     94    DO   ij  = 1, ngrid
     95       pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
     96    ENDDO
    12097
    121 c    .... Calcul de pk  pour la couche l = llm ..
    122 c    (on met la meme distance (en log pression)  entre Pk(llm)
    123 c    et Pk(llm -1) qu'entre Pk(llm-1) et Pk(llm-2)
     98    !    .... Calcul de pk  pour la couche l
     99    !    --------------------------------------------
     100    !
     101    dum1 = cpp * (2*preff)**(-kappa)
     102    DO l = 1, llm-1
     103       DO   ij   = 1, ngrid
     104          pk(ij,l) = dum1 * (p(ij,l) + p(ij,l+1))**kappa
     105       ENDDO
     106    ENDDO
    124107
    125       DO   ij   = 1, ngrid
    126          pk(ij,llm) = pk(ij,llm-1)**2 / pk(ij,llm-2)
    127       ENDDO
     108    !    .... Calcul de pk  pour la couche l = llm ..
     109    !    (on met la meme distance (en log pression)  entre Pk(llm)
     110    !    et Pk(llm -1) qu'entre Pk(llm-1) et Pk(llm-2)
    128111
     112    DO   ij   = 1, ngrid
     113       pk(ij,llm) = pk(ij,llm-1)**2 / pk(ij,llm-2)
     114    ENDDO
    129115
    130 c    calcul de pkf
    131 c    -------------
    132       CALL SCOPY   ( ngrid * llm, pk, 1, pkf, 1 )
    133       CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
    134      
    135 c    EST-CE UTILE ?? : calcul de beta
    136 c    --------------------------------
    137       DO l = 2, llm
    138         DO   ij   = 1, ngrid
    139           beta(ij,l) = pk(ij,l) / pk(ij,l-1)   
    140         ENDDO
    141       ENDDO
     116    if (present(pkf)) then
     117       !    calcul de pkf
     118       pkf = pk
     119       CALL filtreg ( pkf, jmp1, llm, 2, 1, .TRUE., 1 )
     120    end if
    142121
    143       RETURN
    144       END
     122  END SUBROUTINE exner_milieu
     123
     124end module exner_milieu_m
Note: See TracChangeset for help on using the changeset viewer.