Ignore:
Timestamp:
Apr 25, 2014, 12:20:14 PM (11 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.

File:
1 moved

Legend:

Unmodified
Added
Removed
  • LMDZ5/trunk/libf/dyn3dmem/exner_hyb_loc_m.F90

    r1992 r2021  
    1 c
    2 c $Id$
    3 c
    4       SUBROUTINE  exner_hyb_loc(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       USE parallel_lmdz
    32       USE mod_filtreg_p
    33       USE write_field_loc
    34       IMPLICIT NONE
    35 c
    36 #include "dimensions.h"
    37 #include "paramet.h"
    38 #include "comconst.h"
    39 #include "comgeom.h"
    40 #include "comvert.h"
    41 #include "serre.h"
     1module exner_hyb_loc_m
    422
    43       INTEGER  ngrid
    44       REAL p(ijb_u:ije_u,llmp1),pk(ijb_u:ije_u,llm)
    45       REAL pkf(ijb_u:ije_u,llm)
    46       REAL ps(ijb_u:ije_u),pks(ijb_u:ije_u)
    47       REAL alpha(ijb_u:ije_u,llm),beta(ijb_u:ije_u,llm)
     3  IMPLICIT NONE
    484
    49 c    .... variables locales   ...
     5contains
    506
    51       INTEGER l, ij
    52       REAL unpl2k,dellta
     7  SUBROUTINE  exner_hyb_loc(ngrid, ps, p, pks,pk,pkf)
    538
    54       REAL ppn(iim),pps(iim)
    55       REAL xpn, xps
    56       REAL SSUM
    57       EXTERNAL SSUM
    58       INTEGER ije,ijb,jje,jjb
    59       logical,save :: firstcall=.true.
    60 !$OMP THREADPRIVATE(firstcall)
    61       character(len=*),parameter :: modname="exner_hyb_loc"
    62 c
    63 c$OMP BARRIER           
     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    USE parallel_lmdz
     35    USE mod_filtreg_p
     36    USE write_field_loc
     37    !
     38    include "dimensions.h"
     39    include "paramet.h"
     40    include "comconst.h"
     41    include "comgeom.h"
     42    include "comvert.h"
     43    include "serre.h"
    6444
    65       ! Sanity check
    66       if (firstcall) then
    67         ! sanity checks for Shallow Water case (1 vertical layer)
    68         if (llm.eq.1) then
     45    INTEGER  ngrid
     46    REAL p(ijb_u:ije_u,llmp1),pk(ijb_u:ije_u,llm)
     47    REAL, optional:: pkf(ijb_u:ije_u,llm)
     48    REAL ps(ijb_u:ije_u),pks(ijb_u:ije_u)
     49    REAL alpha(ijb_u:ije_u,llm),beta(ijb_u:ije_u,llm)
     50
     51    !    .... variables locales   ...
     52
     53    INTEGER l, ij
     54    REAL unpl2k,dellta
     55
     56    INTEGER ije,ijb,jje,jjb
     57    logical,save :: firstcall=.true.
     58    !$OMP THREADPRIVATE(firstcall)
     59    character(len=*),parameter :: modname="exner_hyb_loc"
     60    !
     61    !$OMP BARRIER           
     62
     63    ! Sanity check
     64    if (firstcall) then
     65       ! sanity checks for Shallow Water case (1 vertical layer)
     66       if (llm.eq.1) then
    6967          if (kappa.ne.1) then
    70             call abort_gcm(modname,
    71      &      "kappa!=1 , but running in Shallow Water mode!!",42)
     68             call abort_gcm(modname, &
     69                  "kappa!=1 , but running in Shallow Water mode!!",42)
    7270          endif
    7371          if (cpp.ne.r) then
    74             call abort_gcm(modname,
    75      &      "cpp!=r , but running in Shallow Water mode!!",42)
     72             call abort_gcm(modname, &
     73                  "cpp!=r , but running in Shallow Water mode!!",42)
    7674          endif
    77         endif ! of if (llm.eq.1)
     75       endif ! of if (llm.eq.1)
    7876
    79         firstcall=.false.
    80       endif ! of if (firstcall)
     77       firstcall=.false.
     78    endif ! of if (firstcall)
    8179
    82 c$OMP BARRIER
     80    !$OMP BARRIER
    8381
    84 ! Specific behaviour for Shallow Water (1 vertical layer) case
    85       if (llm.eq.1) then
    86        
    87         ! Compute pks(:),pk(:),pkf(:)
    88         ijb=ij_begin
    89         ije=ij_end
    90 !$OMP DO SCHEDULE(STATIC)
    91         DO ij=ijb, ije
    92           pks(ij)=(cpp/preff)*ps(ij)
     82    ! Specific behaviour for Shallow Water (1 vertical layer) case:
     83    if (llm.eq.1) then
     84
     85       ! Compute pks(:),pk(:),pkf(:)
     86       ijb=ij_begin
     87       ije=ij_end
     88       !$OMP DO SCHEDULE(STATIC)
     89       DO ij=ijb, ije
     90          pks(ij) = (cpp/preff) * ps(ij)
    9391          pk(ij,1) = .5*pks(ij)
    94           pkf(ij,1)=pk(ij,1)
    95         ENDDO
    96 !$OMP ENDDO
     92          if (present(pkf)) pkf(ij,1)=pk(ij,1)
     93       ENDDO
     94       !$OMP ENDDO
    9795
    98 !$OMP MASTER
    99       if (pole_nord) then
    100         DO  ij   = 1, iim
    101           ppn(ij) = aire(   ij   ) * pks(  ij     )
    102         ENDDO
    103         xpn      = SSUM(iim,ppn,1) /apoln
    104  
    105         DO ij   = 1, iip1
    106           pks(   ij     )  =  xpn
    107           pk(ij,1) = .5*pks(ij)
    108           pkf(ij,1)=pk(ij,1)
    109         ENDDO
    110       endif
    111      
    112       if (pole_sud) then
    113         DO  ij   = 1, iim
    114           pps(ij) = aire(ij+ip1jm) * pks(ij+ip1jm )
    115         ENDDO
    116         xps      = SSUM(iim,pps,1) /apols
    117  
    118         DO ij   = 1, iip1
    119           pks( ij+ip1jm )  =  xps
    120           pk(ij+ip1jm,1)=.5*pks(ij+ip1jm)
    121           pkf(ij+ip1jm,1)=pk(ij+ip1jm,1)
    122         ENDDO
    123       endif
    124 !$OMP END MASTER
    125 !$OMP BARRIER
    126         jjb=jj_begin
    127         jje=jj_end
    128         CALL filtreg_p ( pkf,jjb_u,jje_u,jjb,jje, jmp1, llm,
    129      &                 2, 1, .TRUE., 1 )
     96       !$OMP BARRIER
     97       if (present(pkf)) then
     98          jjb=jj_begin
     99          jje=jj_end
     100          CALL filtreg_p ( pkf,jjb_u,jje_u,jjb,jje, jmp1, llm, &
     101               2, 1, .TRUE., 1 )
     102       end if
    130103
    131         ! our work is done, exit routine
    132         return
    133       endif ! of if (llm.eq.1)
     104       ! our work is done, exit routine
     105       return
     106    endif ! of if (llm.eq.1)
    134107
     108    ! General case:
    135109
    136       unpl2k    = 1.+ 2.* kappa
    137 c
    138       ijb=ij_begin
    139       ije=ij_end
     110    unpl2k    = 1.+ 2.* kappa
    140111
    141 c$OMP DO SCHEDULE(STATIC)
    142       DO   ij  = ijb, ije
    143         pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
    144       ENDDO
    145 c$OMP ENDDO
    146 c Synchro OPENMP ici
     112    !     -------------
     113    !     Calcul de pks
     114    !     -------------
    147115
    148 c$OMP MASTER
    149       if (pole_nord) then
    150         DO  ij   = 1, iim
    151           ppn(ij) = aire(   ij   ) * pks(  ij     )
    152         ENDDO
    153         xpn      = SSUM(iim,ppn,1) /apoln
    154  
    155         DO ij   = 1, iip1
    156           pks(   ij     )  =  xpn
    157         ENDDO
    158       endif
    159      
    160       if (pole_sud) then
    161         DO  ij   = 1, iim
    162           pps(ij) = aire(ij+ip1jm) * pks(ij+ip1jm )
    163         ENDDO
    164         xps      = SSUM(iim,pps,1) /apols
    165  
    166         DO ij   = 1, iip1
    167           pks( ij+ip1jm )  =  xps
    168         ENDDO
    169       endif
    170 c$OMP END MASTER
    171 c$OMP BARRIER
    172 c
    173 c
    174 c    .... Calcul des coeff. alpha et beta  pour la couche l = llm ..
    175 c
    176 c$OMP DO SCHEDULE(STATIC)
    177       DO     ij      = ijb,ije
     116    ijb=ij_begin
     117    ije=ij_end
     118
     119    !$OMP DO SCHEDULE(STATIC)
     120    DO   ij  = ijb, ije
     121       pks(ij) = cpp * ( ps(ij)/preff ) ** kappa
     122    ENDDO
     123    !$OMP ENDDO
     124    ! Synchro OPENMP ici
     125
     126    !$OMP BARRIER
     127    !
     128    !
     129    !    .... Calcul des coeff. alpha et beta  pour la couche l = llm ..
     130    !
     131    !$OMP DO SCHEDULE(STATIC)
     132    DO     ij      = ijb,ije
    178133       alpha(ij,llm) = 0.
    179134       beta (ij,llm) = 1./ unpl2k
    180       ENDDO
    181 c$OMP ENDDO NOWAIT
    182 c
    183 c     ... Calcul des coeff. alpha et beta  pour l = llm-1  a l = 2 ...
    184 c
    185       DO l = llm -1 , 2 , -1
    186 c
    187 c$OMP DO SCHEDULE(STATIC)
    188         DO ij = ijb, ije
    189         dellta = p(ij,l)* unpl2k + p(ij,l+1)* ( beta(ij,l+1)-unpl2k )
    190         alpha(ij,l)  = - p(ij,l+1) / dellta * alpha(ij,l+1)
    191         beta (ij,l)  =   p(ij,l  ) / dellta   
    192         ENDDO
    193 c$OMP ENDDO NOWAIT
    194 c
    195       ENDDO
     135    ENDDO
     136    !$OMP ENDDO NOWAIT
     137    !
     138    !     ... Calcul des coeff. alpha et beta  pour l = llm-1  a l = 2 ...
     139    !
     140    DO l = llm -1 , 2 , -1
     141       !
     142       !$OMP DO SCHEDULE(STATIC)
     143       DO ij = ijb, ije
     144          dellta = p(ij,l)* unpl2k + p(ij,l+1)* ( beta(ij,l+1)-unpl2k )
     145          alpha(ij,l)  = - p(ij,l+1) / dellta * alpha(ij,l+1)
     146          beta (ij,l)  =   p(ij,l  ) / dellta   
     147       ENDDO
     148       !$OMP ENDDO NOWAIT
     149    ENDDO
    196150
    197 c
    198 c  ***********************************************************************
    199 c     .....  Calcul de pk pour la couche 1 , pres du sol  ....
    200 c
    201 c$OMP DO SCHEDULE(STATIC)
    202       DO   ij   = ijb, ije
    203        pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )  /
    204      *    (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-unpl2k )* p(ij,2) )
    205       ENDDO
    206 c$OMP ENDDO NOWAIT
    207 c
    208 c    ..... Calcul de pk(ij,l) , pour l = 2 a l = llm  ........
    209 c
    210       DO l = 2, llm
    211 c$OMP DO SCHEDULE(STATIC)
    212         DO   ij   = ijb, ije
    213          pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1)
    214         ENDDO
    215 c$OMP ENDDO NOWAIT       
    216       ENDDO
    217 c
    218 c
    219 c      CALL SCOPY   ( ngrid * llm, pk, 1, pkf, 1 )
    220       DO l = 1, llm
    221 c$OMP DO SCHEDULE(STATIC)
    222          DO   ij   = ijb, ije
    223            pkf(ij,l)=pk(ij,l)
    224          ENDDO
    225 c$OMP ENDDO NOWAIT             
    226       ENDDO
     151    !  ***********************************************************************
     152    !     .....  Calcul de pk pour la couche 1 , pres du sol  ....
     153    !
     154    !$OMP DO SCHEDULE(STATIC)
     155    DO   ij   = ijb, ije
     156       pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )  / &
     157            (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-unpl2k )* p(ij,2) )
     158    ENDDO
     159    !$OMP ENDDO NOWAIT
     160    !
     161    !    ..... Calcul de pk(ij,l) , pour l = 2 a l = llm  ........
     162    !
     163    DO l = 2, llm
     164       !$OMP DO SCHEDULE(STATIC)
     165       DO   ij   = ijb, ije
     166          pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1)
     167       ENDDO
     168       !$OMP ENDDO NOWAIT       
     169    ENDDO
    227170
    228 c$OMP BARRIER
    229      
    230       jjb=jj_begin
    231       jje=jj_end
     171    if (present(pkf)) then
     172       !    calcul de pkf
     173
     174       DO l = 1, llm
     175          !$OMP DO SCHEDULE(STATIC)
     176          DO   ij   = ijb, ije
     177             pkf(ij,l)=pk(ij,l)
     178          ENDDO
     179          !$OMP ENDDO NOWAIT             
     180       ENDDO
     181
     182       !$OMP BARRIER
     183
     184       jjb=jj_begin
     185       jje=jj_end
    232186#ifdef DEBUG_IO   
    233       call WriteField_u('pkf',pkf)
     187       call WriteField_u('pkf',pkf)
    234188#endif
    235       CALL filtreg_p ( pkf,jjb_u,jje_u,jjb,jje, jmp1, llm,
    236      &                 2, 1, .TRUE., 1 )
     189       CALL filtreg_p ( pkf,jjb_u,jje_u,jjb,jje, jmp1, llm, &
     190            2, 1, .TRUE., 1 )
    237191#ifdef DEBUG_IO   
    238       call WriteField_u('pkf',pkf)
    239 #endif     
     192       call WriteField_u('pkf',pkf)
     193#endif     
     194    end if
    240195
    241       RETURN
    242       END
     196  END SUBROUTINE exner_hyb_loc
     197
     198end module exner_hyb_loc_m
Note: See TracChangeset for help on using the changeset viewer.