New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 2715 for trunk/NEMOGCM/NEMO/OPA_SRC/OBC/obc_oce.F90 – NEMO

Ignore:
Timestamp:
2011-03-30T17:58:35+02:00 (13 years ago)
Author:
rblod
Message:

First attempt to put dynamic allocation on the trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/OBC/obc_oce.F90

    r2528 r2715  
    44   !! Open Boundary Cond. :   define related variables 
    55   !!============================================================================== 
     6   !! history :  OPA  ! 1991-01 (CLIPPER)  Original code  
     7   !!   NEMO     1.0  ! 2002-02 (C. Talandier)  modules, F90 
     8   !!---------------------------------------------------------------------- 
     9#if defined key_obc 
    610   !!---------------------------------------------------------------------- 
    711   !!   'key_obc' :                                 Open Boundary Condition 
    812   !!---------------------------------------------------------------------- 
    9    !! history : 
    10    !!  8.0   01/91   (CLIPPER)  Original code  
    11    !!  8.5   06/02   (C. Talandier)  modules 
    12    !!        06/04   (F. Durand) ORCA2_ZIND config 
    13    !!                  
    14    !!---------------------------------------------------------------------- 
    15    !! * Modules used 
    1613   USE par_oce         ! ocean parameters 
    1714   USE obc_par         ! open boundary condition parameters 
    1815 
    19 #if defined key_obc 
    20  
    2116   IMPLICIT NONE 
    2217   PUBLIC 
     18    
     19   PUBLIC   obc_oce_alloc   ! called by obcini.F90 module 
    2320 
    2421   !!---------------------------------------------------------------------- 
     
    2623   !!---------------------------------------------------------------------- 
    2724   ! 
    28    !                                  !!* Namelist namobc: open boundary condition * 
     25   !                                            !!* Namelist namobc: open boundary condition * 
    2926   INTEGER           ::   nn_obcdta   = 0        !:  = 0 use the initial state as obc data 
    3027   !                                             !   = 1 read obc data in obcxxx.dta files 
     
    6360   !!General variables for open boundaries: 
    6461   !!-------------------------------------- 
    65    LOGICAL ::              &  !: 
    66       lfbceast, lfbcwest,  &  !: logical flag for a fixed East and West open boundaries        
    67       lfbcnorth, lfbcsouth    !: logical flag for a fixed North and South open boundaries        
    68       !                       !  These logical flags are set to 'true' if damping time  
    69       !                       !  scale are set to 0 in the namelist, for both inflow and outflow). 
     62   LOGICAL ::   lfbceast, lfbcwest      !: logical flag for a fixed East and West open boundaries        
     63   LOGICAL ::   lfbcnorth, lfbcsouth    !: logical flag for a fixed North and South open boundaries        
     64   !                                    !  These logical flags are set to 'true' if damping time  
     65   !                                    !  scale are set to 0 in the namelist, for both inflow and outflow). 
    7066 
    7167   REAL(wp), PUBLIC ::   obcsurftot       !: Total lateral surface of open boundaries 
    7268    
    73    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: &  !: 
     69   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: &  !: 
    7470      obctmsk,            &  !: mask array identical to tmask, execpt along OBC where it is set to 0 
    7571      !                      !  it used to calculate the cumulate flux E-P in the obcvol.F90 routine 
     
    8783   INTEGER ::   nje1m2, nje0m1    !: do loop index in mpp case for jpjefm1-1,jpjed 
    8884 
    89    REAL(wp), DIMENSION(jpj) ::   &  !: 
     85   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::   &  !: 
    9086      sshfoe,           & !: now climatology of the east boundary sea surface height 
    9187      ubtfoe,vbtfoe       !: now climatology of the east boundary barotropic transport 
    9288      
    93    REAL(wp), DIMENSION(jpj,jpk) ::   &  !: 
     89   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   &  !: 
    9490      ufoe, vfoe,       & !: now climatology of the east boundary velocities  
    9591      tfoe, sfoe,       & !: now climatology of the east boundary temperature and salinity 
     
    9793      !                   ! in the obcdyn.F90 routine 
    9894 
    99    REAL(wp), DIMENSION(jpi,jpj) ::   sshfoe_b      !: east boundary ssh correction averaged over the barotropic loop 
     95   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sshfoe_b   !: east boundary ssh correction averaged over the barotropic loop 
    10096      !                                            !  (if Flather's algoritm applied at open boundary) 
    10197 
     
    10399   !! Arrays for radiative East OBC:  
    104100   !!------------------------------- 
    105    REAL(wp), DIMENSION(jpj,jpk,3,3) ::   uebnd, vebnd      !: baroclinic u & v component of the velocity over 3 rows  
     101   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   uebnd, vebnd      !: baroclinic u & v component of the velocity over 3 rows  
    106102      !                                                    !  and 3 time step (now, before, and before before) 
    107    REAL(wp), DIMENSION(jpj,jpk,2,2) ::   tebnd, sebnd      !: East boundary temperature and salinity over 2 rows  
     103   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   tebnd, sebnd      !: East boundary temperature and salinity over 2 rows  
    108104      !                                                    !  and 2 time step (now and before) 
    109    REAL(wp), DIMENSION(jpj,jpk) ::   u_cxebnd, v_cxebnd    !: Zonal component of the phase speed ratio computed with  
     105   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   u_cxebnd, v_cxebnd    !: Zonal component of the phase speed ratio computed with  
    110106      !                                                    !  radiation of u and v velocity (respectively) at the  
    111107      !                                                    !  east open boundary (u_cxebnd = cx rdt ) 
    112    REAL(wp), DIMENSION(jpj,jpk) ::   uemsk, vemsk, temsk   !: 2D mask for the East OB 
     108   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   uemsk, vemsk, temsk   !: 2D mask for the East OB 
    113109 
    114110   ! Note that those arrays are optimized for mpp case  
     
    124120   INTEGER ::   njw1m2, njw0m1     !: do loop index in mpp case for jpjwfm2,jpjwd 
    125121 
    126    REAL(wp), DIMENSION(jpj) ::   &  !: 
     122   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:) ::   &  !: 
    127123      sshfow,           & !: now climatology of the west boundary sea surface height 
    128124      ubtfow,vbtfow       !: now climatology of the west boundary barotropic transport 
    129125 
    130    REAL(wp), DIMENSION(jpj,jpk) ::   &  !: 
     126   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   &  !: 
    131127      ufow, vfow,       & !: now climatology of the west velocities  
    132128      tfow, sfow,       & !: now climatology of the west temperature and salinity 
     
    134130      !                   !  in the obcdyn.F90 routine 
    135131 
    136    REAL(wp), DIMENSION(jpi,jpj) ::   sshfow_b    !: west boundary ssh correction averaged over the barotropic loop 
     132   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sshfow_b    !: west boundary ssh correction averaged over the barotropic loop 
    137133      !                                          !  (if Flather's algoritm applied at open boundary) 
    138134 
     
    140136   !! Arrays for radiative West OBC 
    141137   !!------------------------------- 
    142    REAL(wp), DIMENSION(jpj,jpk,3,3) ::   uwbnd, vwbnd     !: baroclinic u & v components of the velocity over 3 rows  
     138   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   uwbnd, vwbnd     !: baroclinic u & v components of the velocity over 3 rows  
    143139      !                                                   !  and 3 time step (now, before, and before before) 
    144    REAL(wp), DIMENSION(jpj,jpk,2,2) ::   twbnd, swbnd     !: west boundary temperature and salinity over 2 rows and  
     140   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   twbnd, swbnd     !: west boundary temperature and salinity over 2 rows and  
    145141      !                                                   !  2 time step (now and before) 
    146    REAL(wp), DIMENSION(jpj,jpk) ::   u_cxwbnd, v_cxwbnd   !: Zonal component of the phase speed ratio computed with  
     142   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   u_cxwbnd, v_cxwbnd   !: Zonal component of the phase speed ratio computed with  
    147143      !                                                   !  radiation of zonal and meridional velocity (respectively)  
    148144      !                                                   !  at the west open boundary (u_cxwbnd = cx rdt ) 
    149    REAL(wp), DIMENSION(jpj,jpk) ::   uwmsk, vwmsk, twmsk  !: 2D mask for the West OB 
     145   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   uwmsk, vwmsk, twmsk  !: 2D mask for the West OB 
    150146 
    151147   ! Note that those arrays are optimized for mpp case  
     
    162158   INTEGER ::   njn0m1, njn1m1     !: do loop index in mpp case for jpnob-1 
    163159 
    164    REAL(wp), DIMENSION(jpi) ::   &  !: 
     160   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:) ::   &  !: 
    165161      sshfon,           & !: now climatology of the north boundary sea surface height 
    166162      ubtfon,vbtfon       !: now climatology of the north boundary barotropic transport 
    167163 
    168    REAL(wp), DIMENSION(jpi,jpk) ::   &    !: 
     164   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   &    !: 
    169165      ufon, vfon,       & !: now climatology of the north boundary velocities 
    170166      tfon, sfon,       & !: now climatology of the north boundary temperature and salinity 
     
    172168      !                   !  in yhe obcdyn.F90 routine 
    173169 
    174    REAL(wp), DIMENSION(jpi,jpj) ::   sshfon_b      !: north boundary ssh correction averaged over the barotropic loop 
     170   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sshfon_b      !: north boundary ssh correction averaged over the barotropic loop 
    175171      !                                            !  (if Flather's algoritm applied at open boundary) 
    176172 
     
    178174   !! Arrays for radiative North OBC 
    179175   !!-------------------------------- 
    180    REAL(wp), DIMENSION(jpi,jpk,3,3) ::   unbnd, vnbnd      !: baroclinic u & v components of the velocity over 3 
     176   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   unbnd, vnbnd      !: baroclinic u & v components of the velocity over 3 
    181177      !                                                    !  rows and 3 time step (now, before, and before before) 
    182    REAL(wp), DIMENSION(jpi,jpk,2,2) ::   tnbnd, snbnd      !: north boundary temperature and salinity over 
     178   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   tnbnd, snbnd      !: north boundary temperature and salinity over 
    183179      !                                                    !  2 rows and 2 time step (now and before) 
    184    REAL(wp), DIMENSION(jpi,jpk) ::   u_cynbnd, v_cynbnd    !: Meridional component of the phase speed ratio compu- 
     180   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   u_cynbnd, v_cynbnd    !: Meridional component of the phase speed ratio compu- 
    185181      !                                                    !  ted with radiation of zonal and meridional velocity  
    186182      !                                                    !  (respectively) at the north OB (u_cynbnd = cx rdt ) 
    187    REAL(wp), DIMENSION(jpi,jpk) ::   unmsk, vnmsk, tnmsk   !: 2D mask for the North OB 
     183   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   unmsk, vnmsk, tnmsk   !: 2D mask for the North OB 
    188184 
    189185   ! Note that those arrays are optimized for mpp case  
     
    199195   INTEGER ::   njs0p1, njs1p1     !: do loop index in mpp case for jpsob+1 
    200196 
    201    REAL(wp), DIMENSION(jpi) ::    &   !: 
     197   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:) ::    &   !: 
    202198      sshfos,           & !: now climatology of the south boundary sea surface height 
    203199      ubtfos,vbtfos       !: now climatology of the south boundary barotropic transport 
    204200 
    205    REAL(wp), DIMENSION(jpi,jpk) ::    &   !: 
     201   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::    &   !: 
    206202      ufos, vfos,       & !: now climatology of the south boundary velocities  
    207203      tfos, sfos,       & !: now climatology of the south boundary temperature and salinity 
     
    209205      !                   !  in the obcdyn.F90 routine 
    210206 
    211    REAL(wp), DIMENSION(jpi,jpj) ::   sshfos_b     !: south boundary ssh correction averaged over the barotropic loop 
     207   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sshfos_b     !: south boundary ssh correction averaged over the barotropic loop 
    212208      !                                           !  (if Flather's algoritm applied at open boundary) 
    213209 
     
    215211   !! Arrays for radiative South OBC   (computed by the forward time step in dynspg) 
    216212   !!-------------------------------- 
    217    REAL(wp), DIMENSION(jpi,jpk,3,3) ::   usbnd, vsbnd     !: baroclinic u & v components of the velocity over 3  
     213   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   usbnd, vsbnd     !: baroclinic u & v components of the velocity over 3  
    218214      !                                                   !  rows and 3 time step (now, before, and before before) 
    219    REAL(wp), DIMENSION(jpi,jpk,2,2) ::   tsbnd, ssbnd     !: south boundary temperature and salinity over 
     215   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   tsbnd, ssbnd     !: south boundary temperature and salinity over 
    220216      !                                                   !  2 rows and 2 time step (now and before) 
    221    REAL(wp), DIMENSION(jpi,jpk) ::   u_cysbnd, v_cysbnd   !: Meridional component of the phase speed ratio 
     217   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   u_cysbnd, v_cysbnd   !: Meridional component of the phase speed ratio 
    222218      !                                                   !  computed with radiation of zonal and meridional velocity  
    223219      !                                                   !  (repsectively) at the south OB (u_cynbnd = cx rdt ) 
    224    REAL(wp), DIMENSION(jpi,jpk) ::   usmsk, vsmsk, tsmsk  !: 2D mask for the South OB 
    225  
    226 #else 
    227    !!---------------------------------------------------------------------- 
    228    !!   Default option :                                       Empty module 
    229    !!---------------------------------------------------------------------- 
    230 #endif 
     220   REAL(wp), ALLOCATABLE, SAVE,     DIMENSION(:,:) ::   usmsk, vsmsk, tsmsk  !: 2D mask for the South OB 
    231221 
    232222   !!---------------------------------------------------------------------- 
     
    234224   !! $Id$  
    235225   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     226   !!---------------------------------------------------------------------- 
     227CONTAINS 
     228 
     229   INTEGER FUNCTION obc_oce_alloc() 
     230      !!---------------------------------------------------------------------- 
     231      !!               ***  FUNCTION obc_oce_alloc  *** 
     232      !!---------------------------------------------------------------------- 
     233 
     234      ALLOCATE(                                                               & 
     235              !! East open boundary 
     236              obctmsk(jpi,jpj), obcumask(jpi,jpj), obcvmask(jpi,jpj),        & 
     237              sshfoe(jpjed:jpjef), ubtfoe(jpjed:jpjef), vbtfoe(jpjed:jpjef), & 
     238              ufoe(jpj,jpk), vfoe(jpj,jpk), tfoe(jpj,jpk), sfoe(jpj,jpk),    & 
     239              uclie(jpj,jpk), sshfoe_b(jpjed:jpjef,jpj),                     & 
     240              !! Arrays for radiative East OBC 
     241              uebnd(jpj,jpk,3,3), vebnd(jpj,jpk,3,3) ,                       & 
     242              tebnd(jpj,jpk,2,2), sebnd(jpj,jpk,2,2),                        & 
     243              u_cxebnd(jpj,jpk), v_cxebnd(jpj,jpk),                          & 
     244              uemsk(jpj,jpk), vemsk(jpj,jpk), temsk(jpj,jpk),                & 
     245              !! West open boundary 
     246              sshfow(jpjwd:jpjwf), ubtfow(jpjwd:jpjwf), vbtfow(jpjwd:jpjwf), & 
     247              ufow(jpj,jpk), vfow(jpj,jpk), tfow(jpj,jpk),                   & 
     248              sfow(jpj,jpk), ucliw(jpj,jpk), sshfow_b(jpjwd:jpjwf,jpj),      & 
     249              !! Arrays for radiative West OBC 
     250              uwbnd(jpj,jpk,3,3), vwbnd(jpj,jpk,3,3),                        & 
     251              twbnd(jpj,jpk,2,2), swbnd(jpj,jpk,2,2),                        & 
     252              u_cxwbnd(jpj,jpk), v_cxwbnd(jpj,jpk),                          & 
     253              uwmsk(jpj,jpk), vwmsk(jpj,jpk), twmsk(jpj,jpk),                & 
     254              !! North open boundary 
     255              sshfon(jpind:jpinf), ubtfon(jpind:jpinf), vbtfon(jpind:jpinf), & 
     256              ufon(jpi,jpk), vfon(jpi,jpk), tfon(jpi,jpk),                   & 
     257              sfon(jpi,jpk), vclin(jpi,jpk), sshfon_b(jpind:jpinf,jpj),      & 
     258              !! Arrays for radiative North OBC 
     259              unbnd(jpi,jpk,3,3), vnbnd(jpi,jpk,3,3),                        & 
     260              tnbnd(jpi,jpk,2,2), snbnd(jpi,jpk,2,2),                        & 
     261              u_cynbnd(jpi,jpk), v_cynbnd(jpi,jpk),                          & 
     262              unmsk(jpi,jpk), vnmsk(jpi,jpk), tnmsk (jpi,jpk),               & 
     263              !! South open boundary 
     264              sshfos(jpisd:jpisf), ubtfos(jpisd:jpisf), vbtfos(jpisd:jpisf), & 
     265              ufos(jpi,jpk), vfos(jpi,jpk), tfos(jpi,jpk),                   & 
     266              sfos(jpi,jpk), vclis(jpi,jpk),                                 & 
     267              sshfos_b(jpisd:jpisf,jpj),                                     & 
     268              !! Arrays for radiative South OBC  
     269              usbnd(jpi,jpk,3,3), vsbnd(jpi,jpk,3,3),                        & 
     270              tsbnd(jpi,jpk,2,2), ssbnd(jpi,jpk,2,2),                        & 
     271              u_cysbnd(jpi,jpk), v_cysbnd(jpi,jpk),                          & 
     272              usmsk(jpi,jpk), vsmsk(jpi,jpk), tsmsk(jpi,jpk),                & 
     273              !! 
     274              STAT=obc_oce_alloc ) 
     275      ! 
     276   END FUNCTION obc_oce_alloc 
     277    
     278#else 
     279   !!---------------------------------------------------------------------- 
     280   !!   Default option         Empty module                          No OBC 
     281   !!---------------------------------------------------------------------- 
     282#endif 
     283 
    236284   !!====================================================================== 
    237285END MODULE obc_oce 
Note: See TracChangeset for help on using the changeset viewer.