Ignore:
Timestamp:
Jul 23, 2024, 7:14:34 PM (8 weeks ago)
Author:
abarral
Message:

Replace 1DUTILS.h by module lmdz_1dutils.f90
Replace 1DConv.h by module lmdz_old_1dconv.f90 (it's only used by old_* files)
Convert *.F to *.f90
Fix gradsdef.h formatting
Remove unnecessary "RETURN" at the end of functions/subroutines

File:
1 moved

Legend:

Unmodified
Added
Removed
  • LMDZ6/branches/Amaury_dev/libf/dyn3dmem/integrd_loc.f90

    r5104 r5105  
    22! $Id: integrd_p.F 1299 2010-01-20 14:27:21Z fairhead $
    33
    4       SUBROUTINE integrd_loc
    5      $  (  nq,vcovm1,ucovm1,tetam1,psm1,massem1,
    6      $     dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps0,masse,phis) !,finvmaold)
    7       USE parallel_lmdz
    8       USE control_mod
    9       USE mod_filtreg_p
    10       USE write_field_loc
    11       USE write_field
    12       USE integrd_mod
    13       USE comconst_mod, ONLY: pi
    14       USE logic_mod, ONLY: leapf
    15       USE comvert_mod, ONLY: ap, bp
    16       USE temps_mod, ONLY: dt
    17       USE strings_mod, ONLY: int2str
    18      
    19       IMPLICIT NONE
    20 
    21 
    22 c=======================================================================
    23 c
    24 c   Auteur:  P. Le Van
    25 c   -------
    26 c
    27 c   objet:
    28 c   ------
    29 c
    30 c   Incrementation des tendances dynamiques
    31 c
    32 c=======================================================================
    33 c-----------------------------------------------------------------------
    34 c   Declarations:
    35 c   -------------
    36 
    37       include "dimensions.h"
    38       include "paramet.h"
    39       include "comgeom.h"
    40       include "iniprint.h"
    41 
    42 c   Arguments:
    43 c   ----------
    44 
    45       INTEGER,intent(in) :: nq ! number of tracers to handle in this routine
    46 
    47       REAL,INTENT(INOUT) :: vcov(ijb_v:ije_v,llm) ! covariant meridional wind
    48       REAL,INTENT(INOUT) :: ucov(ijb_u:ije_u,llm) ! covariant zonal wind
    49       REAL,INTENT(INOUT) :: teta(ijb_u:ije_u,llm) ! potential temperature
    50       REAL,INTENT(INOUT) :: q(ijb_u:ije_u,llm,nq) ! advected tracers
    51       REAL,INTENT(INOUT) :: ps0(ijb_u:ije_u) ! surface pressure
    52       REAL,INTENT(INOUT) :: masse(ijb_u:ije_u,llm) ! atmospheric mass
    53       REAL,INTENT(INOUT) :: phis(ijb_u:ije_u) ! ground geopotential !!! unused
    54       ! values at previous time step
    55       REAL,INTENT(INOUT) :: vcovm1(ijb_v:ije_v,llm)
    56       REAL,INTENT(INOUT) :: ucovm1(ijb_u:ije_u,llm)
    57       REAL,INTENT(INOUT) :: tetam1(ijb_u:ije_u,llm)
    58       REAL,INTENT(INOUT) :: psm1(ijb_u:ije_u)
    59       REAL,INTENT(INOUT) :: massem1(ijb_u:ije_u,llm)
    60       ! the tendencies to add
    61       REAL,INTENT(INOUT) :: dv(ijb_v:ije_v,llm)
    62       REAL,INTENT(INOUT) :: du(ijb_u:ije_u,llm)
    63       REAL,INTENT(INOUT) :: dteta(ijb_u:ije_u,llm)
    64       REAL,INTENT(INOUT) :: dp(ijb_u:ije_u)
    65       REAL,INTENT(INOUT) :: dq(ijb_u:ije_u,llm,nq) !!! unused
    66 !      REAL,INTENT(INOUT) ::finvmaold(ijb_u:ije_u,llm) !!! unused
    67 
    68 c   Local:
    69 c   ------
    70 
    71       REAL vscr( ijb_v:ije_v ),uscr( ijb_u:ije_u )
    72       REAL hscr( ijb_u:ije_u ),pscr(ijb_u:ije_u)
    73       REAL massescr( ijb_u:ije_u,llm )
    74 !      REAL finvmasse(ijb_u:ije_u,llm)
    75       REAL tpn,tps,tppn(iim),tpps(iim)
    76       REAL qpn,qps,qppn(iim),qpps(iim)
    77 
    78       INTEGER  l,ij,iq,i,j
    79 
    80       REAL SSUM
    81       EXTERNAL SSUM
    82       INTEGER ijb,ije,jjb,jje
    83       LOGICAL :: checksum
    84       LOGICAL,SAVE :: checksum_all=.TRUE.
    85       INTEGER :: stop_it
    86       INTEGER :: ierr
    87 
    88       !write(*,*) 'integrd 88: entree, nq=',nq
    89 c-----------------------------------------------------------------------
    90 
    91 c$OMP BARRIER     
    92       if (pole_nord) THEN
    93 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    94         DO  l = 1,llm
    95           DO  ij = 1,iip1
    96            ucov(    ij    , l) = 0.
    97            uscr(     ij      ) = 0.
    98            ENDDO
    99         ENDDO
    100 c$OMP END DO NOWAIT       
    101       ENDIF
    102 
    103       if (pole_sud) THEN
    104 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    105         DO  l = 1,llm
    106           DO  ij = 1,iip1
    107            ucov( ij +ip1jm, l) = 0.
    108            uscr( ij +ip1jm   ) = 0.
    109           ENDDO
    110         ENDDO
    111 c$OMP END DO NOWAIT     
    112       ENDIF
    113 
    114 c    ............    integration  de       ps         ..............
    115 
    116 c      CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
    117 
    118       ijb=ij_begin
    119       ije=ij_end
    120 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    121       DO  l = 1,llm
    122         massescr(ijb:ije,l)=masse(ijb:ije,l)
     4SUBROUTINE integrd_loc &
     5        (  nq,vcovm1,ucovm1,tetam1,psm1,massem1, &
     6        dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps0,masse,phis) !,finvmaold)
     7  USE parallel_lmdz
     8  USE control_mod
     9  USE mod_filtreg_p
     10  USE write_field_loc
     11  USE write_field
     12  USE integrd_mod
     13  USE comconst_mod, ONLY: pi
     14  USE logic_mod, ONLY: leapf
     15  USE comvert_mod, ONLY: ap, bp
     16  USE temps_mod, ONLY: dt
     17  USE strings_mod, ONLY: int2str
     18
     19  IMPLICIT NONE
     20
     21
     22  !=======================================================================
     23  !
     24  !   Auteur:  P. Le Van
     25  !   -------
     26  !
     27  !   objet:
     28  !   ------
     29  !
     30  !   Incrementation des tendances dynamiques
     31  !
     32  !=======================================================================
     33  !-----------------------------------------------------------------------
     34  !   Declarations:
     35  !   -------------
     36
     37  include "dimensions.h"
     38  include "paramet.h"
     39  include "comgeom.h"
     40  include "iniprint.h"
     41
     42  !   Arguments:
     43  !   ----------
     44
     45  INTEGER,intent(in) :: nq ! number of tracers to handle in this routine
     46
     47  REAL,INTENT(INOUT) :: vcov(ijb_v:ije_v,llm) ! covariant meridional wind
     48  REAL,INTENT(INOUT) :: ucov(ijb_u:ije_u,llm) ! covariant zonal wind
     49  REAL,INTENT(INOUT) :: teta(ijb_u:ije_u,llm) ! potential temperature
     50  REAL,INTENT(INOUT) :: q(ijb_u:ije_u,llm,nq) ! advected tracers
     51  REAL,INTENT(INOUT) :: ps0(ijb_u:ije_u) ! surface pressure
     52  REAL,INTENT(INOUT) :: masse(ijb_u:ije_u,llm) ! atmospheric mass
     53  REAL,INTENT(INOUT) :: phis(ijb_u:ije_u) ! ground geopotential !!! unused
     54  ! ! values at previous time step
     55  REAL,INTENT(INOUT) :: vcovm1(ijb_v:ije_v,llm)
     56  REAL,INTENT(INOUT) :: ucovm1(ijb_u:ije_u,llm)
     57  REAL,INTENT(INOUT) :: tetam1(ijb_u:ije_u,llm)
     58  REAL,INTENT(INOUT) :: psm1(ijb_u:ije_u)
     59  REAL,INTENT(INOUT) :: massem1(ijb_u:ije_u,llm)
     60  ! ! the tendencies to add
     61  REAL,INTENT(INOUT) :: dv(ijb_v:ije_v,llm)
     62  REAL,INTENT(INOUT) :: du(ijb_u:ije_u,llm)
     63  REAL,INTENT(INOUT) :: dteta(ijb_u:ije_u,llm)
     64  REAL,INTENT(INOUT) :: dp(ijb_u:ije_u)
     65  REAL,INTENT(INOUT) :: dq(ijb_u:ije_u,llm,nq) !!! unused
     66   ! REAL,INTENT(INOUT) ::finvmaold(ijb_u:ije_u,llm) !!! unused
     67
     68  !   Local:
     69  !   ------
     70
     71  REAL :: vscr( ijb_v:ije_v ),uscr( ijb_u:ije_u )
     72  REAL :: hscr( ijb_u:ije_u ),pscr(ijb_u:ije_u)
     73  REAL :: massescr( ijb_u:ije_u,llm )
     74   ! REAL finvmasse(ijb_u:ije_u,llm)
     75  REAL :: tpn,tps,tppn(iim),tpps(iim)
     76  REAL :: qpn,qps,qppn(iim),qpps(iim)
     77
     78  INTEGER :: l,ij,iq,i,j
     79
     80  REAL :: SSUM
     81  EXTERNAL SSUM
     82  INTEGER :: ijb,ije,jjb,jje
     83  LOGICAL :: checksum
     84  LOGICAL,SAVE :: checksum_all=.TRUE.
     85  INTEGER :: stop_it
     86  INTEGER :: ierr
     87
     88  ! !write(*,*) 'integrd 88: entree, nq=',nq
     89  !-----------------------------------------------------------------------
     90
     91!$OMP BARRIER
     92  if (pole_nord) THEN
     93!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     94    DO  l = 1,llm
     95      DO  ij = 1,iip1
     96       ucov(    ij    , l) = 0.
     97       uscr(     ij      ) = 0.
     98       ENDDO
     99    ENDDO
     100!$OMP END DO NOWAIT
     101  ENDIF
     102
     103  if (pole_sud) THEN
     104!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     105    DO  l = 1,llm
     106      DO  ij = 1,iip1
     107       ucov( ij +ip1jm, l) = 0.
     108       uscr( ij +ip1jm   ) = 0.
    123109      ENDDO
    124 c$OMP END DO NOWAIT
    125 
    126 c$OMP DO SCHEDULE(STATIC)
    127       DO ij = ijb,ije
    128        pscr (ij)    = ps0(ij)
    129        ps (ij)      = psm1(ij) + dt * dp(ij)     
    130 
    131       END DO
    132 
    133 c$OMP END DO 
    134 c$OMP BARRIER
    135 c --> ici synchro OPENMP pour ps
    136        
    137       checksum=.TRUE.
    138       stop_it=0
    139 
    140 c$OMP MASTER
    141 !c$OMP DO SCHEDULE(STATIC)
    142       DO ij = ijb,ije
    143          IF( ps(ij)<0. ) THEN
    144            IF (checksum) stop_it=ij
    145            checksum=.FALSE.
    146          ENDIF
     110    ENDDO
     111!$OMP END DO NOWAIT
     112  ENDIF
     113
     114  !    ............    integration  de       ps         ..............
     115
     116   ! CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
     117
     118  ijb=ij_begin
     119  ije=ij_end
     120!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     121  DO  l = 1,llm
     122    massescr(ijb:ije,l)=masse(ijb:ije,l)
     123  ENDDO
     124!$OMP END DO NOWAIT
     125
     126!$OMP DO SCHEDULE(STATIC)
     127  DO ij = ijb,ije
     128   pscr (ij)    = ps0(ij)
     129   ps (ij)      = psm1(ij) + dt * dp(ij)
     130
     131  END DO
     132
     133!$OMP END DO
     134!$OMP BARRIER
     135  ! --> ici synchro OPENMP pour ps
     136
     137  checksum=.TRUE.
     138  stop_it=0
     139
     140!$OMP MASTER
     141  !c$OMP DO SCHEDULE(STATIC)
     142  DO ij = ijb,ije
     143     IF( ps(ij)<0. ) THEN
     144       IF (checksum) stop_it=ij
     145       checksum=.FALSE.
     146     ENDIF
     147   ENDDO
     148  !c$OMP END DO NOWAIT
     149
     150   ! CALL MPI_ALLREDUCE(checksum,checksum_all,1,
     151  ! &                   MPI_LOGICAL,MPI_LOR,COMM_LMDZ,ierr)
     152  IF( .NOT. checksum ) THEN
     153     write(lunout,*) "integrd: ps = ", ps(stop_it)
     154     write(lunout,*) " at node ij =", stop_it
     155     ! ! since ij=j+(i-1)*jjp1 , we have
     156      j=modulo(stop_it,jjp1)
     157      i=1+(stop_it-j)/jjp1
     158      write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg", &
     159            " lat = ",rlatu(j)*180./pi, " deg"
     160     CALL abort_gcm("integrd_loc", "negative surface pressure", 1)
     161  ENDIF
     162
     163!$OMP END MASTER
     164!$OMP BARRIER
     165    ! !write(*,*) 'integrd 170'
     166  IF (.NOT. Checksum_all) THEN
     167    CALL WriteField_v('int_vcov',vcov)
     168    CALL WriteField_u('int_ucov',ucov)
     169    CALL WriteField_u('int_teta',teta)
     170    CALL WriteField_u('int_ps0',ps0)
     171    CALL WriteField_u('int_masse',masse)
     172    CALL WriteField_u('int_phis',phis)
     173    CALL WriteField_v('int_vcovm1',vcovm1)
     174    CALL WriteField_u('int_ucovm1',ucovm1)
     175    CALL WriteField_u('int_tetam1',tetam1)
     176    CALL WriteField_u('int_psm1',psm1)
     177    CALL WriteField_u('int_massem1',massem1)
     178
     179    CALL WriteField_v('int_dv',dv)
     180    CALL WriteField_u('int_du',du)
     181    CALL WriteField_u('int_dteta',dteta)
     182    CALL WriteField_u('int_dp',dp)
     183     ! CALL WriteField_u('int_finvmaold',finvmaold)
     184    do j=1,nq
     185      CALL WriteField_u('int_q'//trim(int2str(j)), &
     186            q(:,:,j))
     187      CALL WriteField_u('int_dq'//trim(int2str(j)), &
     188            dq(:,:,j))
     189    enddo
     190    CALL abort_gcm("integrd_loc", "", 1)
     191  ENDIF
     192
     193
     194  !
     195  !   !write(*,*) 'integrd 200'
     196!$OMP MASTER
     197  if (pole_nord) THEN
     198
     199    DO  ij    = 1, iim
     200     tppn(ij) = aire(   ij   ) * ps(  ij    )
     201    ENDDO
     202     tpn      = SSUM(iim,tppn,1)/apoln
     203    DO ij   = 1, iip1
     204     ps(   ij   )  = tpn
     205    ENDDO
     206
     207  ENDIF
     208
     209  if (pole_sud) THEN
     210
     211    DO  ij    = 1, iim
     212     tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
     213    ENDDO
     214     tps      = SSUM(iim,tpps,1)/apols
     215    DO ij   = 1, iip1
     216     ps(ij+ip1jm)  = tps
     217    ENDDO
     218
     219  ENDIF
     220!$OMP END MASTER
     221!$OMP BARRIER
     222  ! !write(*,*) 'integrd 217'
     223  !
     224  !  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...
     225  !
     226
     227  CALL pression_loc ( ip1jmp1, ap, bp, ps, p )
     228
     229!$OMP BARRIER
     230  CALL massdair_loc (     p  , masse         )
     231
     232  ! Ehouarn : we don't use/need finvmaold and finvmasse,
     233        ! so might as well not compute them
     234  !c      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
     235   ! ijb=ij_begin
     236   ! ije=ij_end
     237
     238  !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     239   ! DO  l = 1,llm
     240   !   finvmasse(ijb:ije,l)=masse(ijb:ije,l)
     241   ! ENDDO
     242  !c$OMP END DO NOWAIT
     243
     244   ! jjb=jj_begin
     245   ! jje=jj_end
     246   ! CALL filtreg_p( finvmasse,jjb_u,jje_u,jjb,jje, jjp1, llm,
     247  ! &                -2, 2, .TRUE., 1  )
     248  !
     249
     250  !    ............   integration  de  ucov, vcov,  h     ..............
     251
     252!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     253  DO l = 1,llm
     254
     255  ijb=ij_begin
     256  ije=ij_end
     257  if (pole_nord) ijb=ij_begin+iip1
     258  if (pole_sud)  ije=ij_end-iip1
     259
     260  DO ij = ijb,ije
     261  uscr( ij )   =  ucov( ij,l )
     262  ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
     263  END DO
     264
     265  ijb=ij_begin
     266  ije=ij_end
     267  if (pole_sud)  ije=ij_end-iip1
     268
     269  DO ij = ijb,ije
     270  vscr( ij )   =  vcov( ij,l )
     271  vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
     272  END DO
     273
     274  ijb=ij_begin
     275  ije=ij_end
     276
     277  DO ij = ijb,ije
     278  hscr( ij )    =  teta(ij,l)
     279  teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l) &
     280        + dt * dteta(ij,l) / masse(ij,l)
     281  END DO
     282
     283  !   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......
     284  !
     285  !
     286  !   !write(*,*) 'integrd 291'
     287  IF (pole_nord) THEN
     288
     289    DO  ij   = 1, iim
     290      tppn(ij) = aire(   ij   ) * teta(  ij    ,l)
     291    ENDDO
     292      tpn      = SSUM(iim,tppn,1)/apoln
     293
     294    DO ij   = 1, iip1
     295      teta(   ij   ,l)  = tpn
     296    ENDDO
     297
     298  ENDIF
     299
     300  IF (pole_sud) THEN
     301
     302    DO  ij   = 1, iim
     303      tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
     304    ENDDO
     305      tps      = SSUM(iim,tpps,1)/apols
     306
     307    DO ij   = 1, iip1
     308      teta(ij+ip1jm,l)  = tps
     309    ENDDO
     310
     311  ENDIF
     312  !
     313
     314  IF(leapf)  THEN
     315      ! CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
     316      ! CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
     317      ! CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
     318    ijb=ij_begin
     319    ije=ij_end
     320    ucovm1(ijb:ije,l)=uscr(ijb:ije)
     321    tetam1(ijb:ije,l)=hscr(ijb:ije)
     322    if (pole_sud) ije=ij_end-iip1
     323    vcovm1(ijb:ije,l)=vscr(ijb:ije)
     324
     325  END IF
     326
     327  END DO
     328!$OMP END DO NOWAIT
     329
     330  !
     331  !   .......  integration de   q   ......
     332  !
     333  ijb=ij_begin
     334  ije=ij_end
     335
     336     if (planet_type=="earth") then
     337  ! Earth-specific treatment of first 2 tracers (water)
     338!$OMP BARRIER
     339!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     340      DO l = 1, llm
     341       DO ij = ijb, ije
     342        deltap(ij,l) =  p(ij,l) - p(ij,l+1)
    147343       ENDDO
    148 !c$OMP END DO NOWAIT
    149        
    150 !      CALL MPI_ALLREDUCE(checksum,checksum_all,1,
    151 !     &                   MPI_LOGICAL,MPI_LOR,COMM_LMDZ,ierr)
    152       IF( .NOT. checksum ) THEN
    153          write(lunout,*) "integrd: ps = ", ps(stop_it)
    154          write(lunout,*) " at node ij =", stop_it
    155          ! since ij=j+(i-1)*jjp1 , we have
    156           j=modulo(stop_it,jjp1)
    157           i=1+(stop_it-j)/jjp1
    158           write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg",
    159      &                    " lat = ",rlatu(j)*180./pi, " deg"
    160          CALL abort_gcm("integrd_loc", "negative surface pressure", 1)
    161       ENDIF
    162 
    163 c$OMP END MASTER
    164 c$OMP BARRIER
    165         !write(*,*) 'integrd 170'
    166       IF (.NOT. Checksum_all) THEN
    167         CALL WriteField_v('int_vcov',vcov)
    168         CALL WriteField_u('int_ucov',ucov)
    169         CALL WriteField_u('int_teta',teta)
    170         CALL WriteField_u('int_ps0',ps0)
    171         CALL WriteField_u('int_masse',masse)
    172         CALL WriteField_u('int_phis',phis)
    173         CALL WriteField_v('int_vcovm1',vcovm1)
    174         CALL WriteField_u('int_ucovm1',ucovm1)
    175         CALL WriteField_u('int_tetam1',tetam1)
    176         CALL WriteField_u('int_psm1',psm1)
    177         CALL WriteField_u('int_massem1',massem1)
    178 
    179         CALL WriteField_v('int_dv',dv)
    180         CALL WriteField_u('int_du',du)
    181         CALL WriteField_u('int_dteta',dteta)
    182         CALL WriteField_u('int_dp',dp)
    183 !        CALL WriteField_u('int_finvmaold',finvmaold)
    184         do j=1,nq
    185           CALL WriteField_u('int_q'//trim(int2str(j)),
    186      .                q(:,:,j))
    187           CALL WriteField_u('int_dq'//trim(int2str(j)),
    188      .                dq(:,:,j))
    189         enddo
    190         CALL abort_gcm("integrd_loc", "", 1)
    191       ENDIF
    192    
    193        
    194 c
    195         !write(*,*) 'integrd 200'
    196 C$OMP MASTER
    197       if (pole_nord) THEN
    198      
    199         DO  ij    = 1, iim
    200          tppn(ij) = aire(   ij   ) * ps(  ij    )
    201         ENDDO
    202          tpn      = SSUM(iim,tppn,1)/apoln
    203         DO ij   = 1, iip1
    204          ps(   ij   )  = tpn
    205         ENDDO
    206      
    207       ENDIF
    208      
    209       if (pole_sud) THEN
    210      
    211         DO  ij    = 1, iim
    212          tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
    213         ENDDO
    214          tps      = SSUM(iim,tpps,1)/apols
    215         DO ij   = 1, iip1
    216          ps(ij+ip1jm)  = tps
    217         ENDDO
    218      
    219       ENDIF
    220 c$OMP END MASTER
    221 c$OMP BARRIER
    222       !write(*,*) 'integrd 217' 
    223 c
    224 c  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...
    225 c
    226 
    227       CALL pression_loc ( ip1jmp1, ap, bp, ps, p )
    228 
    229 c$OMP BARRIER
    230       CALL massdair_loc (     p  , masse         )
    231 
    232 ! Ehouarn : we don't use/need finvmaold and finvmasse,
    233 !           so might as well not compute them
    234 !c      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
    235 !      ijb=ij_begin
    236 !      ije=ij_end
    237 
    238 !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    239 !      DO  l = 1,llm
    240 !        finvmasse(ijb:ije,l)=masse(ijb:ije,l)
    241 !      ENDDO
    242 !c$OMP END DO NOWAIT
    243 
    244 !      jjb=jj_begin
    245 !      jje=jj_end
    246 !      CALL filtreg_p( finvmasse,jjb_u,jje_u,jjb,jje, jjp1, llm,
    247 !     &                -2, 2, .TRUE., 1  )
    248 c
    249 
    250 c    ............   integration  de  ucov, vcov,  h     ..............
    251 
    252 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    253       DO l = 1,llm
    254      
    255       ijb=ij_begin
    256       ije=ij_end
    257       if (pole_nord) ijb=ij_begin+iip1
    258       if (pole_sud)  ije=ij_end-iip1
    259      
    260       DO ij = ijb,ije
    261       uscr( ij )   =  ucov( ij,l )
    262       ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
    263       END DO
    264 
    265       ijb=ij_begin
    266       ije=ij_end
    267       if (pole_sud)  ije=ij_end-iip1
    268      
    269       DO ij = ijb,ije
    270       vscr( ij )   =  vcov( ij,l )
    271       vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
    272       END DO
    273      
    274       ijb=ij_begin
    275       ije=ij_end
    276      
    277       DO ij = ijb,ije
    278       hscr( ij )    =  teta(ij,l)
    279       teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l)
    280      $                + dt * dteta(ij,l) / masse(ij,l)
    281       END DO
    282 
    283 c   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......
    284 c
    285 c
    286         !write(*,*) 'integrd 291'
    287       IF (pole_nord) THEN
    288        
    289         DO  ij   = 1, iim
    290           tppn(ij) = aire(   ij   ) * teta(  ij    ,l)
    291         ENDDO
    292           tpn      = SSUM(iim,tppn,1)/apoln
    293 
    294         DO ij   = 1, iip1
    295           teta(   ij   ,l)  = tpn
    296         ENDDO
    297      
    298       ENDIF
    299      
    300       IF (pole_sud) THEN
    301        
    302         DO  ij   = 1, iim
    303           tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
    304         ENDDO
    305           tps      = SSUM(iim,tpps,1)/apols
    306 
    307         DO ij   = 1, iip1
    308           teta(ij+ip1jm,l)  = tps
    309         ENDDO
    310      
    311       ENDIF
    312 c
    313 
    314       IF(leapf)  THEN
    315 c         CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
    316 c         CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
    317 c         CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
    318         ijb=ij_begin
    319         ije=ij_end
    320         ucovm1(ijb:ije,l)=uscr(ijb:ije)
    321         tetam1(ijb:ije,l)=hscr(ijb:ije)
    322         if (pole_sud) ije=ij_end-iip1
    323         vcovm1(ijb:ije,l)=vscr(ijb:ije)
    324      
    325       END IF
    326 
    327       END DO
    328 c$OMP END DO NOWAIT
    329 
    330 c
    331 c   .......  integration de   q   ......
    332 c
    333       ijb=ij_begin
    334       ije=ij_end
    335 
    336          if (planet_type=="earth") then
    337 ! Earth-specific treatment of first 2 tracers (water)
    338 c$OMP BARRIER
    339 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    340           DO l = 1, llm
    341            DO ij = ijb, ije
    342             deltap(ij,l) =  p(ij,l) - p(ij,l+1)
    343            ENDDO
    344           ENDDO
    345          
    346 c$OMP END DO NOWAIT
    347 c$OMP BARRIER
    348 
    349         CALL check_isotopes(q,ijb,ije,'integrd 342')
    350 
    351         !write(*,*) 'integrd 341'
    352         CALL qminimum_loc( q, nq, deltap )
    353         !write(*,*) 'integrd 343'
    354 
    355         CALL check_isotopes(q,ijb,ije,'integrd 346')
    356 c
    357 c    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
    358 c
    359 c$OMP BARRIER
    360       IF (pole_nord) THEN
    361      
    362         DO iq = 1, nq
    363        
    364 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    365           DO l = 1, llm
    366  
    367              DO ij = 1, iim
    368                qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)
    369              ENDDO
    370                qpn  =  SSUM(iim,qppn,1)/apoln
    371      
    372              DO ij = 1, iip1
    373                q(   ij   ,l,iq)  = qpn
    374              ENDDO   
    375  
    376           ENDDO
    377 c$OMP END DO NOWAIT
    378 
    379         ENDDO
    380      
    381       ENDIF
    382 
    383       IF (pole_sud) THEN
    384      
    385         DO iq = 1, nq
    386 
    387 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    388           DO l = 1, llm
    389  
    390              DO ij = 1, iim
    391                qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
    392              ENDDO
    393                qps  =  SSUM(iim,qpps,1)/apols
    394  
    395              DO ij = 1, iip1
    396                q(ij+ip1jm,l,iq)  = qps
    397              ENDDO   
    398  
    399           ENDDO
    400 c$OMP END DO NOWAIT
    401 
    402         ENDDO
    403      
    404       ENDIF
    405 
    406       CALL check_isotopes(q,ijb,ije,'integrd 409')
    407      
    408 ! Ehouarn: forget about finvmaold
    409 !c         CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
    410 
    411 !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    412 !      DO l = 1, llm     
    413 !        finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)       
    414 !      ENDDO
    415 !c$OMP END DO NOWAIT
    416 
    417       endif ! of if (planet_type.eq."earth")
    418 
    419 c
    420 c
    421 c     .....   FIN  de l'integration  de   q    .......
    422 
    423 15    continue
    424           !write(*,*) 'integrd 410'
    425 
    426 c$OMP DO SCHEDULE(STATIC)
    427       DO ij=ijb,ije 
    428         ps0(ij)=ps(ij)
    429344      ENDDO
    430 c$OMP END DO NOWAIT
    431 
    432 c    .................................................................
    433 
    434 
    435       IF( leapf )  THEN
    436 c       CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )
    437 c       CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )
    438 c$OMP DO SCHEDULE(STATIC)
    439       DO ij=ijb,ije 
    440         psm1(ij)=pscr(ij)
     345
     346!$OMP END DO NOWAIT
     347!$OMP BARRIER
     348
     349    CALL check_isotopes(q,ijb,ije,'integrd 342')
     350
     351    ! !write(*,*) 'integrd 341'
     352    CALL qminimum_loc( q, nq, deltap )
     353    ! !write(*,*) 'integrd 343'
     354
     355    CALL check_isotopes(q,ijb,ije,'integrd 346')
     356  !
     357  !    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
     358  !
     359!$OMP BARRIER
     360  IF (pole_nord) THEN
     361
     362    DO iq = 1, nq
     363
     364!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     365      DO l = 1, llm
     366
     367         DO ij = 1, iim
     368           qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)
     369         ENDDO
     370           qpn  =  SSUM(iim,qppn,1)/apoln
     371
     372         DO ij = 1, iip1
     373           q(   ij   ,l,iq)  = qpn
     374         ENDDO
     375
    441376      ENDDO
    442 c$OMP END DO NOWAIT
    443 
    444 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    445           DO l = 1, llm
    446             massem1(ijb:ije,l)=massescr(ijb:ije,l)
    447           ENDDO
    448 c$OMP END DO NOWAIT         
    449       END IF
    450 c$OMP BARRIER
    451       RETURN
    452       END
     377!$OMP END DO NOWAIT
     378
     379    ENDDO
     380
     381  ENDIF
     382
     383  IF (pole_sud) THEN
     384
     385    DO iq = 1, nq
     386
     387!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     388      DO l = 1, llm
     389
     390         DO ij = 1, iim
     391           qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
     392         ENDDO
     393           qps  =  SSUM(iim,qpps,1)/apols
     394
     395         DO ij = 1, iip1
     396           q(ij+ip1jm,l,iq)  = qps
     397         ENDDO
     398
     399      ENDDO
     400!$OMP END DO NOWAIT
     401
     402    ENDDO
     403
     404  ENDIF
     405
     406  CALL check_isotopes(q,ijb,ije,'integrd 409')
     407
     408  ! Ehouarn: forget about finvmaold
     409  !c         CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
     410
     411  !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     412   ! DO l = 1, llm
     413   !   finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)
     414   ! ENDDO
     415  !c$OMP END DO NOWAIT
     416
     417  endif ! of if (planet_type.eq."earth")
     418
     419  !
     420  !
     421  ! .....   FIN  de l'integration  de   q    .......
     422
     423      ! !write(*,*) 'integrd 410'
     424
     425!$OMP DO SCHEDULE(STATIC)
     426  DO ij=ijb,ije
     427    ps0(ij)=ps(ij)
     428  ENDDO
     429!$OMP END DO NOWAIT
     430
     431  !    .................................................................
     432
     433
     434  IF( leapf )  THEN
     435    ! CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )
     436    ! CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )
     437!$OMP DO SCHEDULE(STATIC)
     438  DO ij=ijb,ije
     439    psm1(ij)=pscr(ij)
     440  ENDDO
     441!$OMP END DO NOWAIT
     442
     443!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     444      DO l = 1, llm
     445        massem1(ijb:ije,l)=massescr(ijb:ije,l)
     446      ENDDO
     447!$OMP END DO NOWAIT
     448  END IF
     449!$OMP BARRIER
     450
     451END SUBROUTINE integrd_loc
Note: See TracChangeset for help on using the changeset viewer.