Ignore:
Timestamp:
Oct 21, 2024, 2:58:45 PM (23 hours ago)
Author:
abarral
Message:

Convert fixed-form to free-form sources .F -> .{f,F}90
(WIP: some .F remain, will be handled in subsequent commits)

File:
1 moved

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/dyn3d/integrd.f90

    r5245 r5246  
    22! $Id$
    33!
    4       SUBROUTINE integrd
    5      $  (  nq,vcovm1,ucovm1,tetam1,psm1,massem1,
    6      $     dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis !,finvmaold
    7      &  )
    8 
    9       use control_mod, only : planet_type
    10       use comconst_mod, only: pi
    11       USE logic_mod, ONLY: leapf
    12       use comvert_mod, only: ap, bp
    13       USE temps_mod, ONLY: dt
    14 
    15       IMPLICIT NONE
    16 
    17 
    18 c=======================================================================
    19 c
    20 c   Auteur:  P. Le Van
    21 c   -------
    22 c
    23 c   objet:
    24 c   ------
    25 c
    26 c   Incrementation des tendances dynamiques
    27 c
    28 c=======================================================================
    29 c-----------------------------------------------------------------------
    30 c   Declarations:
    31 c   -------------
    32 
    33       include "dimensions.h"
    34       include "paramet.h"
    35       include "comgeom.h"
    36       include "iniprint.h"
    37 
    38 c   Arguments:
    39 c   ----------
    40 
    41       integer,intent(in) :: nq ! number of tracers to handle in this routine
    42       real,intent(inout) :: vcov(ip1jm,llm) ! covariant meridional wind
    43       real,intent(inout) :: ucov(ip1jmp1,llm) ! covariant zonal wind
    44       real,intent(inout) :: teta(ip1jmp1,llm) ! potential temperature
    45       real,intent(inout) :: q(ip1jmp1,llm,nq) ! advected tracers
    46       real,intent(inout) :: ps(ip1jmp1) ! surface pressure
    47       real,intent(inout) :: masse(ip1jmp1,llm) ! atmospheric mass
    48       real,intent(in) :: phis(ip1jmp1) ! ground geopotential !!! unused
    49       ! values at previous time step
    50       real,intent(inout) :: vcovm1(ip1jm,llm)
    51       real,intent(inout) :: ucovm1(ip1jmp1,llm)
    52       real,intent(inout) :: tetam1(ip1jmp1,llm)
    53       real,intent(inout) :: psm1(ip1jmp1)
    54       real,intent(inout) :: massem1(ip1jmp1,llm)
    55       ! the tendencies to add
    56       real,intent(in) :: dv(ip1jm,llm)
    57       real,intent(in) :: du(ip1jmp1,llm)
    58       real,intent(in) :: dteta(ip1jmp1,llm)
    59       real,intent(in) :: dp(ip1jmp1)
    60       real,intent(in) :: dq(ip1jmp1,llm,nq) !!! unused
    61 !      real,intent(out) :: finvmaold(ip1jmp1,llm) !!! unused
    62 
    63 c   Local:
    64 c   ------
    65 
    66       REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
    67       REAL massescr( ip1jmp1,llm )
    68 !      REAL finvmasse(ip1jmp1,llm)
    69       REAL p(ip1jmp1,llmp1)
    70       REAL tpn,tps,tppn(iim),tpps(iim)
    71       REAL qpn,qps,qppn(iim),qpps(iim)
    72       REAL deltap( ip1jmp1,llm )
    73 
    74       INTEGER  l,ij,iq,i,j
    75 
    76       REAL SSUM
    77 
    78 c-----------------------------------------------------------------------
    79 
    80       DO  l = 1,llm
    81         DO  ij = 1,iip1
    82          ucov(    ij    , l) = 0.
    83          ucov( ij +ip1jm, l) = 0.
    84          uscr(     ij      ) = 0.
    85          uscr( ij +ip1jm   ) = 0.
    86         ENDDO
     4SUBROUTINE integrd &
     5        (  nq,vcovm1,ucovm1,tetam1,psm1,massem1, &
     6        dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis & !,finvmaold
     7        )
     8
     9  use control_mod, only : planet_type
     10  use comconst_mod, only: pi
     11  USE logic_mod, ONLY: leapf
     12  use comvert_mod, only: ap, bp
     13  USE temps_mod, ONLY: dt
     14
     15  IMPLICIT NONE
     16
     17
     18  !=======================================================================
     19  !
     20  !   Auteur:  P. Le Van
     21  !   -------
     22  !
     23  !   objet:
     24  !   ------
     25  !
     26  !   Incrementation des tendances dynamiques
     27  !
     28  !=======================================================================
     29  !-----------------------------------------------------------------------
     30  !   Declarations:
     31  !   -------------
     32
     33  include "dimensions.h"
     34  include "paramet.h"
     35  include "comgeom.h"
     36  include "iniprint.h"
     37
     38  !   Arguments:
     39  !   ----------
     40
     41  integer,intent(in) :: nq ! number of tracers to handle in this routine
     42  real,intent(inout) :: vcov(ip1jm,llm) ! covariant meridional wind
     43  real,intent(inout) :: ucov(ip1jmp1,llm) ! covariant zonal wind
     44  real,intent(inout) :: teta(ip1jmp1,llm) ! potential temperature
     45  real,intent(inout) :: q(ip1jmp1,llm,nq) ! advected tracers
     46  real,intent(inout) :: ps(ip1jmp1) ! surface pressure
     47  real,intent(inout) :: masse(ip1jmp1,llm) ! atmospheric mass
     48  real,intent(in) :: phis(ip1jmp1) ! ground geopotential !!! unused
     49  ! ! values at previous time step
     50  real,intent(inout) :: vcovm1(ip1jm,llm)
     51  real,intent(inout) :: ucovm1(ip1jmp1,llm)
     52  real,intent(inout) :: tetam1(ip1jmp1,llm)
     53  real,intent(inout) :: psm1(ip1jmp1)
     54  real,intent(inout) :: massem1(ip1jmp1,llm)
     55  ! ! the tendencies to add
     56  real,intent(in) :: dv(ip1jm,llm)
     57  real,intent(in) :: du(ip1jmp1,llm)
     58  real,intent(in) :: dteta(ip1jmp1,llm)
     59  real,intent(in) :: dp(ip1jmp1)
     60  real,intent(in) :: dq(ip1jmp1,llm,nq) !!! unused
     61   ! real,intent(out) :: finvmaold(ip1jmp1,llm) !!! unused
     62
     63  !   Local:
     64  !   ------
     65
     66  REAL :: vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
     67  REAL :: massescr( ip1jmp1,llm )
     68   ! REAL finvmasse(ip1jmp1,llm)
     69  REAL :: p(ip1jmp1,llmp1)
     70  REAL :: tpn,tps,tppn(iim),tpps(iim)
     71  REAL :: qpn,qps,qppn(iim),qpps(iim)
     72  REAL :: deltap( ip1jmp1,llm )
     73
     74  INTEGER :: l,ij,iq,i,j
     75
     76  REAL :: SSUM
     77
     78  !-----------------------------------------------------------------------
     79
     80  DO  l = 1,llm
     81    DO  ij = 1,iip1
     82     ucov(    ij    , l) = 0.
     83     ucov( ij +ip1jm, l) = 0.
     84     uscr(     ij      ) = 0.
     85     uscr( ij +ip1jm   ) = 0.
     86    ENDDO
     87  ENDDO
     88
     89
     90  !    ............    integration  de       ps         ..............
     91
     92  CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
     93
     94  DO ij = 1,ip1jmp1
     95   pscr (ij)    = ps(ij)
     96   ps (ij)      = psm1(ij) + dt * dp(ij)
     97  ENDDO
     98  !
     99  DO ij = 1,ip1jmp1
     100    IF( ps(ij).LT.0. ) THEN
     101     write(lunout,*) "integrd: negative surface pressure ",ps(ij)
     102     write(lunout,*) " at node ij =", ij
     103     ! ! since ij=j+(i-1)*jjp1 , we have
     104     j=modulo(ij,jjp1)
     105     i=1+(ij-j)/jjp1
     106     write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg", &
     107           " lat = ",rlatu(j)*180./pi, " deg"
     108     call abort_gcm("integrd", "", 1)
     109    ENDIF
     110  ENDDO
     111  !
     112  DO  ij    = 1, iim
     113   tppn(ij) = aire(   ij   ) * ps(  ij    )
     114   tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
     115  ENDDO
     116   tpn      = SSUM(iim,tppn,1)/apoln
     117   tps      = SSUM(iim,tpps,1)/apols
     118  DO ij   = 1, iip1
     119   ps(   ij   )  = tpn
     120   ps(ij+ip1jm)  = tps
     121  ENDDO
     122  !
     123  !  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...
     124  !
     125  CALL pression ( ip1jmp1, ap, bp, ps, p )
     126  CALL massdair (     p  , masse         )
     127
     128  ! Ehouarn : we don't use/need finvmaold and finvmasse,
     129        ! so might as well not compute them
     130   ! CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
     131   ! CALL filtreg( finvmasse, jjp1, llm, -2, 2, .TRUE., 1  )
     132  !
     133
     134  !    ............   integration  de  ucov, vcov,  h     ..............
     135
     136  DO l = 1,llm
     137
     138   DO ij = iip2,ip1jm
     139    uscr( ij )   =  ucov( ij,l )
     140    ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
     141   ENDDO
     142
     143   DO ij = 1,ip1jm
     144    vscr( ij )   =  vcov( ij,l )
     145    vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
     146   ENDDO
     147
     148   DO ij = 1,ip1jmp1
     149    hscr( ij )    =  teta(ij,l)
     150    teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l) &
     151          + dt * dteta(ij,l) / masse(ij,l)
     152   ENDDO
     153
     154  !   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......
     155  !
     156  !
     157   DO  ij   = 1, iim
     158    tppn(ij) = aire(   ij   ) * teta(  ij    ,l)
     159    tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
     160   ENDDO
     161    tpn      = SSUM(iim,tppn,1)/apoln
     162    tps      = SSUM(iim,tpps,1)/apols
     163
     164   DO ij   = 1, iip1
     165    teta(   ij   ,l)  = tpn
     166    teta(ij+ip1jm,l)  = tps
     167   ENDDO
     168  !
     169
     170   IF(leapf)  THEN
     171     CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
     172     CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
     173     CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
     174   END IF
     175
     176  ENDDO ! of DO l = 1,llm
     177
     178
     179  !
     180  !   .......  integration de   q   ......
     181  !
     182  !$$$      IF( iadv(1).NE.3.AND.iadv(2).NE.3 )    THEN
     183  !$$$c
     184  !$$$       IF( forward.OR. leapf )  THEN
     185  !$$$        DO iq = 1,2
     186  !$$$        DO  l = 1,llm
     187  !$$$        DO ij = 1,ip1jmp1
     188  !$$$        q(ij,l,iq) = ( q(ij,l,iq)*finvmaold(ij,l) + dtvr *dq(ij,l,iq) )/
     189  !$$$     $                            finvmasse(ij,l)
     190  !$$$        ENDDO
     191  !$$$        ENDDO
     192  !$$$        ENDDO
     193  !$$$       ELSE
     194  !$$$         DO iq = 1,2
     195  !$$$         DO  l = 1,llm
     196  !$$$         DO ij = 1,ip1jmp1
     197  !$$$         q( ij,l,iq ) = q( ij,l,iq ) * finvmaold(ij,l) / finvmasse(ij,l)
     198  !$$$         ENDDO
     199  !$$$         ENDDO
     200  !$$$         ENDDO
     201  !$$$
     202  !$$$       END IF
     203  !$$$c
     204  !$$$      ENDIF
     205
     206  if (planet_type.eq."earth") then
     207  ! Earth-specific treatment of first 2 tracers (water)
     208    DO l = 1, llm
     209      DO ij = 1, ip1jmp1
     210        deltap(ij,l) =  p(ij,l) - p(ij,l+1)
    87211      ENDDO
    88 
    89 
    90 c    ............    integration  de       ps         ..............
    91 
    92       CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
    93 
    94       DO ij = 1,ip1jmp1
    95        pscr (ij)    = ps(ij)
    96        ps (ij)      = psm1(ij) + dt * dp(ij)
    97       ENDDO
    98 c
    99       DO ij = 1,ip1jmp1
    100         IF( ps(ij).LT.0. ) THEN
    101          write(lunout,*) "integrd: negative surface pressure ",ps(ij)
    102          write(lunout,*) " at node ij =", ij
    103          ! since ij=j+(i-1)*jjp1 , we have
    104          j=modulo(ij,jjp1)
    105          i=1+(ij-j)/jjp1
    106          write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg",
    107      &                   " lat = ",rlatu(j)*180./pi, " deg"
    108          call abort_gcm("integrd", "", 1)
    109         ENDIF
    110       ENDDO
    111 c
    112       DO  ij    = 1, iim
    113        tppn(ij) = aire(   ij   ) * ps(  ij    )
    114        tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
    115       ENDDO
    116        tpn      = SSUM(iim,tppn,1)/apoln
    117        tps      = SSUM(iim,tpps,1)/apols
    118       DO ij   = 1, iip1
    119        ps(   ij   )  = tpn
    120        ps(ij+ip1jm)  = tps
    121       ENDDO
    122 c
    123 c  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...
    124 c
    125       CALL pression ( ip1jmp1, ap, bp, ps, p )
    126       CALL massdair (     p  , masse         )
    127 
    128 ! Ehouarn : we don't use/need finvmaold and finvmasse,
    129 !           so might as well not compute them
    130 !      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
    131 !      CALL filtreg( finvmasse, jjp1, llm, -2, 2, .TRUE., 1  )
    132 c
    133 
    134 c    ............   integration  de  ucov, vcov,  h     ..............
    135 
    136       DO l = 1,llm
    137 
    138        DO ij = iip2,ip1jm
    139         uscr( ij )   =  ucov( ij,l )
    140         ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
     212    ENDDO
     213
     214    CALL qminimum( q, nq, deltap )
     215
     216  !
     217  !    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
     218  !
     219
     220   DO iq = 1, nq
     221    DO l = 1, llm
     222
     223       DO ij = 1, iim
     224         qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)
     225         qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
    141226       ENDDO
    142 
    143        DO ij = 1,ip1jm
    144         vscr( ij )   =  vcov( ij,l )
    145         vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
     227         qpn  =  SSUM(iim,qppn,1)/apoln
     228         qps  =  SSUM(iim,qpps,1)/apols
     229
     230       DO ij = 1, iip1
     231         q(   ij   ,l,iq)  = qpn
     232         q(ij+ip1jm,l,iq)  = qps
    146233       ENDDO
    147234
    148        DO ij = 1,ip1jmp1
    149         hscr( ij )    =  teta(ij,l)
    150         teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l)
    151      &                + dt * dteta(ij,l) / masse(ij,l)
    152        ENDDO
    153 
    154 c   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......
    155 c
    156 c
    157        DO  ij   = 1, iim
    158         tppn(ij) = aire(   ij   ) * teta(  ij    ,l)
    159         tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
    160        ENDDO
    161         tpn      = SSUM(iim,tppn,1)/apoln
    162         tps      = SSUM(iim,tpps,1)/apols
    163 
    164        DO ij   = 1, iip1
    165         teta(   ij   ,l)  = tpn
    166         teta(ij+ip1jm,l)  = tps
    167        ENDDO
    168 c
    169 
    170        IF(leapf)  THEN
    171          CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
    172          CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
    173          CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
    174        END IF
    175 
    176       ENDDO ! of DO l = 1,llm
    177 
    178 
    179 c
    180 c   .......  integration de   q   ......
    181 c
    182 c$$$      IF( iadv(1).NE.3.AND.iadv(2).NE.3 )    THEN
    183 c$$$c
    184 c$$$       IF( forward. OR . leapf )  THEN
    185 c$$$        DO iq = 1,2
    186 c$$$        DO  l = 1,llm
    187 c$$$        DO ij = 1,ip1jmp1
    188 c$$$        q(ij,l,iq) = ( q(ij,l,iq)*finvmaold(ij,l) + dtvr *dq(ij,l,iq) )/
    189 c$$$     $                            finvmasse(ij,l)
    190 c$$$        ENDDO
    191 c$$$        ENDDO
    192 c$$$        ENDDO
    193 c$$$       ELSE
    194 c$$$         DO iq = 1,2
    195 c$$$         DO  l = 1,llm
    196 c$$$         DO ij = 1,ip1jmp1
    197 c$$$         q( ij,l,iq ) = q( ij,l,iq ) * finvmaold(ij,l) / finvmasse(ij,l)
    198 c$$$         ENDDO
    199 c$$$         ENDDO
    200 c$$$         ENDDO
    201 c$$$
    202 c$$$       END IF
    203 c$$$c
    204 c$$$      ENDIF
    205 
    206       if (planet_type.eq."earth") then
    207 ! Earth-specific treatment of first 2 tracers (water)
    208         DO l = 1, llm
    209           DO ij = 1, ip1jmp1
    210             deltap(ij,l) =  p(ij,l) - p(ij,l+1)
    211           ENDDO
    212         ENDDO
    213 
    214         CALL qminimum( q, nq, deltap )
    215 
    216 c
    217 c    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
    218 c
    219 
    220        DO iq = 1, nq
    221         DO l = 1, llm
    222 
    223            DO ij = 1, iim
    224              qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)
    225              qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
    226            ENDDO
    227              qpn  =  SSUM(iim,qppn,1)/apoln
    228              qps  =  SSUM(iim,qpps,1)/apols
    229 
    230            DO ij = 1, iip1
    231              q(   ij   ,l,iq)  = qpn
    232              q(ij+ip1jm,l,iq)  = qps
    233            ENDDO
    234 
    235         ENDDO
    236        ENDDO
    237 
    238 ! Ehouarn: forget about finvmaold
    239 !      CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
    240 
    241       endif ! of if (planet_type.eq."earth")
    242 c
    243 c
    244 c     .....   FIN  de l'integration  de   q    .......
    245 
    246 c    .................................................................
    247 
    248 
    249       IF( leapf )  THEN
    250          CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )
    251          CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )
    252       END IF
    253 
    254       RETURN
    255       END
     235    ENDDO
     236   ENDDO
     237
     238  ! Ehouarn: forget about finvmaold
     239   ! CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
     240
     241  endif ! of if (planet_type.eq."earth")
     242  !
     243  !
     244  ! .....   FIN  de l'integration  de   q    .......
     245
     246  !    .................................................................
     247
     248
     249  IF( leapf )  THEN
     250     CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )
     251     CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )
     252  END IF
     253
     254  RETURN
     255END SUBROUTINE integrd
Note: See TracChangeset for help on using the changeset viewer.