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 7412 – NEMO

Changeset 7412


Ignore:
Timestamp:
2016-12-01T11:30:29+01:00 (8 years ago)
Author:
lovato
Message:

Merge dev_NOC_CMCC_merge_2016 into branch

Location:
branches/2016/dev_merge_2016
Files:
15 added
1 deleted
80 edited
2 copied

Legend:

Unmodified
Added
Removed

  • branches/2016/dev_merge_2016/DOC/Namelists/nam_tide

    r6997 r7412  
    11!----------------------------------------------------------------------- 
    2 &nam_tide      !   tide parameters                                      ("key_tide") 
     2&nam_tide      !   tide parameters 
    33!----------------------------------------------------------------------- 
    4    ln_tide_pot = .true.    !  use tidal potential forcing 
    5    ln_tide_ramp= .false.   ! 
    6    rdttideramp =    0.     ! 
    7    clname(1)   = 'DUMMY'   !  name of constituent - all tidal components must be set in namelist_cfg 
     4   ln_tide       = .true.   !  Activate tide module 
     5   ln_tide_pot   = .true.   !  use tidal potential forcing 
     6   ln_tide_ramp  = .false.  ! 
     7   rdttideramp   =    0.    ! 
     8   clname(1)     = 'DUMMY'  !  name of constituent - all tidal components must be set in namelist_cfg 
    89/ 

  • branches/2016/dev_merge_2016/DOC/Namelists/nambdy

    r6140 r7412  
    11!----------------------------------------------------------------------- 
    2 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
     2&nambdy        !  unstructured open boundaries 
    33!----------------------------------------------------------------------- 
     4    ln_bdy         = .true.               !  Activate BDY module 
    45    nb_bdy         = 0                    !  number of open boundary sets 
    56    ln_coords_file = .true.               !  =T : read bdy coordinates from file 

  • branches/2016/dev_merge_2016/DOC/Namelists/nambdy_dta

    r6997 r7412  
    11!----------------------------------------------------------------------- 
    2 &nambdy_dta    !  open boundaries - external data                       ("key_bdy") 
     2&nambdy_dta    !  open boundaries - external data 
    33!----------------------------------------------------------------------- 
    44!              !  file name      ! frequency (hours) ! variable  ! time interp.!  clim   ! 'yearly'/ ! weights  ! rotation ! land/sea mask ! 

  • branches/2016/dev_merge_2016/DOC/TexFiles/Chapters/Chap_CFG.tex

    r6997 r7412  
    299299In particular, the AMM uses $S$-coordinates in the vertical rather than 
    300300$z$-coordinates and is forced with tidal lateral boundary conditions 
    301 using a flather boundary condition from the BDY module (key\_bdy). 
     301using a flather boundary condition from the BDY module. 
    302302The AMM configuration  uses the GLS (key\_zdfgls) turbulence scheme, the 
    303303VVL non-linear free surface(key\_vvl) and time-splitting 
    … …  
    306306In addition to the tidal boundary condition the model may also take 
    307307open boundary conditions from a North Atlantic model. Boundaries may be 
    308 completely ommited by removing the BDY key (key\_bdy). 
     308completely omitted by setting \np{ln\_bdy} to false. 
    309309Sample surface fluxes, river forcing and a sample initial restart file 
    310310are included to test a realistic model run. The Baltic boundary is 

  • branches/2016/dev_merge_2016/DOC/TexFiles/Chapters/Chap_DYN.tex

    r6997 r7412  
    12091209into account when computing the surface pressure gradient. 
    12101210 
    1211 (2) When \np{ln\_tide\_pot}~=~true and \key{tide} is defined (see \S\ref{SBC_tide}),  
     1211(2) When \np{ln\_tide\_pot}~=~true and \np{ln\_tide}~=~true (see \S\ref{SBC_tide}),  
    12121212the tidal potential is taken into account when computing the surface pressure gradient. 
    12131213 

  • branches/2016/dev_merge_2016/DOC/TexFiles/Chapters/Chap_LBC.tex

    r6997 r7412  
    350350% Unstructured open boundaries BDY  
    351351% ==================================================================== 
    352 \section{Unstructured Open Boundary Conditions (\key{bdy}) (BDY)} 
     352\section{Unstructured Open Boundary Conditions (BDY)} 
    353353\label{LBC_bdy} 
    354354 
    … …  
    368368Options are defined through the \ngn{nambdy} \ngn{nambdy\_index}  
    369369\ngn{nambdy\_dta} \ngn{nambdy\_dta2} namelist variables. 
    370 The BDY module is an alternative implementation of open boundary 
     370The BDY module is the core implementation of open boundary 
    371371conditions for regional configurations. It implements the Flow 
    372372Relaxation Scheme algorithm for temperature, salinity, velocities and 
    … …  
    376376an isobath or other irregular contour.  
    377377 
    378 The BDY module was modelled on the OBC module and shares many features 
    379 and a similar coding structure \citep{Chanut2005}. 
    380  
    381 The BDY module is completely rewritten at NEMO 3.4 and there is a new 
    382 set of namelists. Boundary data files used with earlier versions of 
    383 NEMO may need to be re-ordered to work with this version. See the 
     378The BDY module was modelled on the OBC module (see NEMO 3.4) and shares many 
     379features and a similar coding structure \citep{Chanut2005}. 
     380 
     381Boundary data files used with earlier versions of NEMO may need 
     382to be re-ordered to work with this version. See the 
    384383section on the Input Boundary Data Files for details. 
    385384 
    … …  
    388387\label{BDY_namelist} 
    389388 
     389The BDY module is activated by setting \np{ln\_bdy} to true. 
    390390It is possible to define more than one boundary ``set'' and apply 
    391391different boundary conditions to each set. The number of boundary 

  • branches/2016/dev_merge_2016/DOC/TexFiles/Chapters/Chap_MISC.tex

    r6997 r7412  
    243243   b \qquad \ &= sum_2 \\ 
    244244\end{align*} 
    245 This feature can be found in \mdl{lib\_fortran} module and is effective when \key{mpp\_rep}. 
    246 In that case, all calls to glob\_sum function (summation over the entire basin excluding  
     245An example of this feature can be found in \mdl{lib\_fortran} module. 
     246It is systematicallt used in glob\_sum function (summation over the entire basin excluding  
    247247duplicated rows and columns due to cyclic or north fold boundary condition as well as  
    248 overlap MPP areas).  
    249 Note this implementation may be sensitive to the optimization level.  
     248overlap MPP areas). The self-compensated summation method should be used in all summation 
     249in i- and/or j-direction. See closea.F90 module for an example. 
     250Note also that this implementation may be sensitive to the optimization level.  
    250251 
    251252\subsection{MPP scalability} 

  • branches/2016/dev_merge_2016/DOC/TexFiles/Chapters/Chap_SBC.tex

    r6997 r7412  
    776776 
    777777A module is available to compute the tidal potential and use it in the momentum equation. 
    778 This option is activated when \key{tide} is defined. 
     778This option is activated when \np{ln\_tide} is set to true in \ngn{nam\_tide}. 
    779779 
    780780Some parameters are available in namelist \ngn{nam\_tide}: 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/AMM12/EXP00/namelist_cfg

    r6489 r7412  
    189189/ 
    190190!----------------------------------------------------------------------- 
    191 &nam_tide      !   tide parameters (#ifdef key_tide) 
    192 !----------------------------------------------------------------------- 
     191&nam_tide      !   tide parameters 
     192!----------------------------------------------------------------------- 
     193   ln_tide     = .true. 
    193194   clname(1)     =   'Q1'   !  name of constituent 
    194195   clname(2)     =   'O1' 
    … …  
    208209/ 
    209210!----------------------------------------------------------------------- 
    210 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
     211&nambdy        !  unstructured open boundaries 
     212    ln_bdy         = .true. 
    211213    nb_bdy         =  1 
    212214    cn_dyn2d       = 'flather' 
    … …  
    216218/ 
    217219!----------------------------------------------------------------------- 
    218 &nambdy_dta      !  open boundaries - external data           ("key_bdy") 
     220&nambdy_dta      !  open boundaries - external data 
    219221!----------------------------------------------------------------------- 
    220222!          !  file name  ! frequency (hours) ! variable  ! time interp. !  clim  ! 'yearly'/ ! weights  ! rotation ! land/sea mask ! 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/AMM12/cpp_AMM12.fcm

    r6140 r7412  
    1  bld::tool::fppkeys  key_bdy key_tide key_zdfgls key_diainstant key_mpp_mpi key_iomput 
     1 bld::tool::fppkeys key_zdfgls key_diainstant key_mpp_mpi key_iomput 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/C1D_PAPA/EXP00/namelist_cfg

    r6489 r7412  
    174174/ 
    175175!----------------------------------------------------------------------- 
    176 &nam_tide      !    tide parameters (#ifdef key_tide) 
    177 !----------------------------------------------------------------------- 
    178 / 
    179 !----------------------------------------------------------------------- 
    180 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
    181 !----------------------------------------------------------------------- 
    182 / 
    183 !----------------------------------------------------------------------- 
    184 &nambdy_dta      !  open boundaries - external data           ("key_bdy") 
     176&nam_tide      !    tide parameters 
     177!----------------------------------------------------------------------- 
     178/ 
     179!----------------------------------------------------------------------- 
     180&nambdy        !  unstructured open boundaries                           
     181!----------------------------------------------------------------------- 
     182/ 
     183!----------------------------------------------------------------------- 
     184&nambdy_dta      !  open boundaries - external data            
    185185!----------------------------------------------------------------------- 
    186186/ 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/GYRE/EXP00/namelist_cfg

    r6489 r7412  
    154154/ 
    155155!----------------------------------------------------------------------- 
    156 &nam_tide      !    tide parameters (#ifdef key_tide) 
    157 !----------------------------------------------------------------------- 
    158 / 
    159 !----------------------------------------------------------------------- 
    160 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
    161 !----------------------------------------------------------------------- 
    162 / 
    163 !----------------------------------------------------------------------- 
    164 &nambdy_dta      !  open boundaries - external data           ("key_bdy") 
     156&nam_tide      !    tide parameters 
     157!----------------------------------------------------------------------- 
     158/ 
     159!----------------------------------------------------------------------- 
     160&nambdy        !  unstructured open boundaries                           
     161!----------------------------------------------------------------------- 
     162/ 
     163!----------------------------------------------------------------------- 
     164&nambdy_dta      !  open boundaries - external data            
    165165!----------------------------------------------------------------------- 
    166166/ 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/GYRE_BFM/EXP00/namelist_cfg

    r6489 r7412  
    159159/ 
    160160!----------------------------------------------------------------------- 
    161 &nam_tide      !    tide parameters (#ifdef key_tide) 
    162 !----------------------------------------------------------------------- 
    163 / 
    164 !----------------------------------------------------------------------- 
    165 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
    166 !----------------------------------------------------------------------- 
    167 / 
    168 !----------------------------------------------------------------------- 
    169 &nambdy_dta      !  open boundaries - external data           ("key_bdy") 
     161&nam_tide      !    tide parameters 
     162!----------------------------------------------------------------------- 
     163/ 
     164!----------------------------------------------------------------------- 
     165&nambdy        !  unstructured open boundaries                           
     166!----------------------------------------------------------------------- 
     167/ 
     168!----------------------------------------------------------------------- 
     169&nambdy_dta      !  open boundaries - external data            
    170170!----------------------------------------------------------------------- 
    171171/ 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/GYRE_XIOS/EXP00/namelist_cfg

    r6489 r7412  
    152152/ 
    153153!----------------------------------------------------------------------- 
    154 &nam_tide      !    tide parameters (#ifdef key_tide) 
    155 !----------------------------------------------------------------------- 
    156 / 
    157 !----------------------------------------------------------------------- 
    158 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
    159 !----------------------------------------------------------------------- 
    160 / 
    161 !----------------------------------------------------------------------- 
    162 &nambdy_dta      !  open boundaries - external data           ("key_bdy") 
     154&nam_tide      !    tide parameters 
     155!----------------------------------------------------------------------- 
     156/ 
     157!----------------------------------------------------------------------- 
     158&nambdy        !  unstructured open boundaries                           
     159!----------------------------------------------------------------------- 
     160/ 
     161!----------------------------------------------------------------------- 
     162&nambdy_dta      !  open boundaries - external data            
    163163!----------------------------------------------------------------------- 
    164164/ 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/SHARED/namelist_ref

    r7403 r7412  
    607607!!   namagrif      agrif nested grid ( read by child model only )       ("key_agrif") 
    608608!!   nam_tide      Tidal forcing  
    609 !!   nambdy        Unstructured open boundaries                         ("key_bdy") 
    610 !!   nambdy_dta    Unstructured open boundaries - external data         ("key_bdy") 
    611 !!   nambdy_tide   tidal forcing at open boundaries                     ("key_bdy_tides") 
     609!!   nambdy        Unstructured open boundaries                          
     610!!   nambdy_dta    Unstructured open boundaries - external data          
     611!!   nambdy_tide   tidal forcing at open boundaries                      
    612612!!====================================================================== 
    613613! 
    … …  
    629629/ 
    630630!----------------------------------------------------------------------- 
    631 &nam_tide      !   tide parameters                                      ("key_tide") 
    632 !----------------------------------------------------------------------- 
     631&nam_tide      !   tide parameters 
     632!----------------------------------------------------------------------- 
     633   ln_tide     = .false. 
    633634   ln_tide_pot = .true.    !  use tidal potential forcing 
    634635   ln_tide_ramp= .false.   ! 
    … …  
    637638/ 
    638639!----------------------------------------------------------------------- 
    639 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
    640 !----------------------------------------------------------------------- 
     640&nambdy        !  unstructured open boundaries                           
     641!----------------------------------------------------------------------- 
     642    ln_bdy         = .false.              !  Use unstructured open boundaries 
    641643    nb_bdy         = 0                    !  number of open boundary sets 
    642644    ln_coords_file = .true.               !  =T : read bdy coordinates from file 
    … …  
    669671    ln_vol        = .false.               !  total volume correction (see nn_volctl parameter) 
    670672    nn_volctl     = 1                     !  = 0, the total water flux across open boundaries is zero 
    671 / 
    672 !----------------------------------------------------------------------- 
    673 &nambdy_dta    !  open boundaries - external data                       ("key_bdy") 
     673    nb_jpk_bdy    = -1                    ! number of levels in the bdy data (set < 0 if consistent with planned run) 
     674/ 
     675!----------------------------------------------------------------------- 
     676&nambdy_dta    !  open boundaries - external data                        
    674677!----------------------------------------------------------------------- 
    675678!              !  file name      ! frequency (hours) ! variable  ! time interp.!  clim   ! 'yearly'/ ! weights  ! rotation ! land/sea mask ! 
    … …  
    958961   !                                !  = 30  F(i,j,k)=c2d*c1d 
    959962   !                                !  = 31  F(i,j,k)=F(grid spacing and local velocity) 
     963   !                                !  = 32  F(i,j,k)=F(local gridscale and deformation rate) 
     964   ! Caution in 20 and 30 cases the coefficient have to be given for a 1 degree grid (~111km) 
    960965   rn_ahm_0      =  40000.     !  horizontal laplacian eddy viscosity   [m2/s] 
    961966   rn_ahm_b      =      0.     !  background eddy viscosity for ldf_iso [m2/s] 
    962967   rn_bhm_0      = 1.e+12      !  horizontal bilaplacian eddy viscosity [m4/s] 
    963    ! 
    964    ! Caution in 20 and 30 cases the coefficient have to be given for a 1 degree grid (~111km) 
     968   !                       !  Smagorinsky settings (nn_ahm_ijk_t  = 32) : 
     969   rn_csmc       = 3.5         !  Smagorinsky constant of proportionality 
     970   rn_minfac     = 1.0         !  multiplier of theorectical lower limit 
     971   rn_maxfac     = 1.0         !  multiplier of theorectical upper limit 
    965972/ 
    966973 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/EXP00/namelist_cfg

    r7403 r7412  
    55&namrun        !   parameters of the run 
    66!----------------------------------------------------------------------- 
    7    cn_exp      =  "GYRE"   !  experience name 
     7   cn_exp      =  "WAD"    !  experience name 
    88   nn_it000    =       1   !  first time step 
    9    nn_itend    =    4320   !  last  time step 
     9   nn_itend    =      5760  !  last  time step 
    1010   nn_leapy    =      30   !  Leap year calendar (1) or not (0) 
    11    nn_stock    =    4320   !  frequency of creation of a restart file (modulo referenced to 1) 
    12    nn_write    =      60   !  frequency of write in the output file   (modulo referenced to nn_it000) 
     11   nn_stock    =    48000  !  frequency of creation of a restart file (modulo referenced to 1) 
    1312 
    1413   ln_clobber  = .true.    !  clobber (overwrite) an existing file 
     14   nn_istate   =       0   !  output the initial state (1) or not (0) 
    1515 
    1616/ 
    … …  
    1818&namcfg     !   parameters of the configuration    
    1919!----------------------------------------------------------------------- 
    20    cp_cfg      =  "gyre"                 !  name of the configuration 
     20   cp_cfg      =  "wad"                  !  name of the configuration 
    2121   jp_cfg      =       1                 !  resolution of the configuration 
    22    jpidta      =      32                 !  1st lateral dimension ( >= jpi ) = 30*jp_cfg+2 
    23    jpjdta      =      22                 !  2nd    "         "    ( >= jpj ) = 20*jp_cfg+2  
    24    jpkdta      =      31                 !  number of levels      ( >= jpk ) 
    25    jpiglo      =      32                 !  1st dimension of global domain --> i  = jpidta 
    26    jpjglo      =      22                 !  2nd    -                  -    --> j  = jpjdta 
     22   jpidta      =      51                 !  1st lateral dimension ( >= jpi ) = 30*jp_cfg+2 
     23   jpjdta      =      34                 !  2nd    "         "    ( >= jpj ) = 20*jp_cfg+2  
     24   jpkdta      =      10                 !  number of levels      ( >= jpk ) 
     25   jpiglo      =      51                 !  1st dimension of global domain --> i  = jpidta 
     26   jpjglo      =      34                 !  2nd    -                  -    --> j  = jpjdta 
    2727   jpizoom     =       1                 !  left bottom (i,j) indices of the zoom 
    2828   jpjzoom     =       1                 !  in data domain indices 
    … …  
    3232&namzgr        !   vertical coordinate 
    3333!----------------------------------------------------------------------- 
    34    ln_zco      = .true.    !  z-coordinate - full    steps 
    35    ln_linssh   = .true.    !  linear free surface 
     34   ln_sco      = .true.    !  s- or hybrid z-s-coordinate 
     35   ln_linssh   = .false.   !  linear free surface 
    3636/ 
    3737!----------------------------------------------------------------------- 
    3838&namzgr_sco    !   s-coordinate or hybrid z-s-coordinate 
    3939!----------------------------------------------------------------------- 
     40   ln_s_sh94   = .false.   !  Song & Haidvogel 1994 hybrid S-sigma   (T)| 
     41   ln_s_sf12   = .true.    !  Siddorn & Furner 2012 hybrid S-z-sigma (T)| if both are false the NEMO tanh stretching is applied 
     42   ln_sigcrit  = .true.    !  use sigma coordinates below critical depth (T) or Z coordinates (F) for Siddorn & Furner stretch 
     43                           !  stretching coefficients for all functions 
     44   rn_sbot_min =   0.01     !  minimum depth of s-bottom surface (>0) (m) 
     45   rn_sbot_max =   15.0    !  maximum depth of s-bottom surface (= ocean depth) (>0) (m) 
     46   rn_hc       =    3.0    !  critical depth for transition to stretched coordinates 
    4047/ 
    4148!----------------------------------------------------------------------- 
    4249&namdom        !   space and time domain (bathymetry, mesh, timestep) 
    4350!----------------------------------------------------------------------- 
     51   nn_msh      =    1      !  create (=1) a mesh file or not (=0) 
    4452   nn_bathy    =    0      !  compute (=0) or read (=1) the bathymetry file 
    45    rn_rdt      = 7200.     !  time step for the dynamics  
    46    jphgr_msh   =       5                 !  type of horizontal mesh 
     53   rn_bathy    =    10.    !  value of the bathymetry. if (=0) bottom flat at jpkm1 
     54   rn_rdt      =    12.    !  time step for the dynamics  
     55   jphgr_msh   =       1                 !  type of horizontal mesh 
    4756   ppglam0     =       0.0               !  longitude of first raw and column T-point (jphgr_msh = 1) 
    48    ppgphi0     =      29.0               ! latitude  of first raw and column T-point (jphgr_msh = 1) 
    49    ppe1_deg    =  999999.0               !  zonal      grid-spacing (degrees) 
    50    ppe2_deg    =  999999.0               !  meridional grid-spacing (degrees) 
     57   ppgphi0     =       10.0               ! latitude  of first raw and column T-point (jphgr_msh = 1) 
     58   ppe1_deg    =       0.01              !  zonal      grid-spacing (degrees) 
     59   ppe2_deg    =       0.01              !  meridional grid-spacing (degrees) 
    5160   ppe1_m      =  999999.0               !  zonal      grid-spacing (degrees) 
    5261   ppe2_m      =  999999.0               !  meridional grid-spacing (degrees) 
    53    ppsur       =   -2033.194295283385    !  ORCA r4, r2 and r05 coefficients 
    54    ppa0        =     155.8325369664153   ! (default coefficients) 
    55    ppa1        =     146.3615918601890   ! 
    56    ppkth       =      17.28520372419791  ! 
    57    ppacr       =       5.0               ! 
    58    ppdzmin     =  999999.0               !  Minimum vertical spacing 
    59    pphmax      =  999999.0               !  Maximum depth 
     62   ppsur       =  999999.0               !  ORCA r4, r2 and r05 coefficients 
     63   ppa0        =  999999.0               ! (default coefficients) 
     64   ppa1        =  999999.0               ! 
     65   ppkth       =  999999.0               ! 
     66   ppacr       =  999999.0               ! 
     67   ppdzmin     =  0.2                    !  Minimum vertical spacing 
     68   pphmax      =  10.0                   !  Maximum depth 
    6069   ldbletanh   =  .FALSE.                !  Use/do not use double tanf function for vertical coordinates 
    6170   ppa2        =  999999.0               !  Double tanh function parameters 
    … …  
    91100!----------------------------------------------------------------------- 
    92101   nn_tau000   =   100     !  gently increase the stress over the first ntau_rst time-steps 
    93    rn_utau0    =   0.1e0   !  uniform value for the i-stress 
     102   rn_utau0    =   0.0e0   !  uniform value for the i-stress 
    94103/ 
    95104!----------------------------------------------------------------------- 
    … …  
    158167/ 
    159168!----------------------------------------------------------------------- 
    160 &nambdy        !  unstructured open boundaries                          ("key_bdy") 
    161 !----------------------------------------------------------------------- 
    162 / 
    163 !----------------------------------------------------------------------- 
    164 &nambdy_dta      !  open boundaries - external data           ("key_bdy") 
    165 !----------------------------------------------------------------------- 
     169&nambdy        !  unstructured open boundaries 
     170!----------------------------------------------------------------------- 
     171    ln_bdy         = .true.              
     172    nb_bdy         = 0                    !  number of open boundary sets 
     173    ln_coords_file = .false.              !  =T : read bdy coordinates from file 
     174    cn_coords_file = 'coordinates.bdy.nc' !  bdy coordinates files 
     175    ln_mask_file   = .false.              !  =T : read mask from file 
     176    cn_mask_file   = ''                   !  name of mask file (if ln_mask_file=.TRUE.) 
     177    cn_dyn2d       = 'flather'            ! 
     178    nn_dyn2d_dta   =  1                   !  = 0, bdy data are equal to the initial state 
     179                                          !  = 1, bdy data are read in 'bdydata   .nc' files 
     180                                          !  = 2, use tidal harmonic forcing data from files 
     181                                          !  = 3, use external data AND tidal harmonic forcing 
     182    cn_dyn3d      =  'none'               ! 
     183    nn_dyn3d_dta  =  0                    !  = 0, bdy data are equal to the initial state 
     184                                          !  = 1, bdy data are read in 'bdydata   .nc' files 
     185    cn_tra        =  'frs'                ! 
     186    nn_tra_dta    =  0                    !  = 0, bdy data are equal to the initial state 
     187                                          !  = 1, bdy data are read in 'bdydata   .nc' files 
     188    cn_ice_lim      =  'none'             ! 
     189    nn_ice_lim_dta  =  0                  !  = 0, bdy data are equal to the initial state 
     190                                          !  = 1, bdy data are read in 'bdydata   .nc' files 
     191    rn_ice_tem      = 270.                !  lim3 only: arbitrary temperature of incoming sea ice 
     192    rn_ice_sal      = 10.                 !  lim3 only:      --   salinity           -- 
     193    rn_ice_age      = 30.                 !  lim3 only:      --   age                -- 
     194 
     195    ln_tra_dmp    =.false.                !  open boudaries conditions for tracers 
     196    ln_dyn3d_dmp  =.false.                !  open boundary condition for baroclinic velocities 
     197    rn_time_dmp   =  1.                   ! Damping time scale in days 
     198    rn_time_dmp_out =  1.                 ! Outflow damping time scale 
     199    nn_rimwidth   = 10                    !  width of the relaxation zone 
     200    ln_vol        = .false.               !  total volume correction (see nn_volctl parameter) 
     201    nn_volctl     = 1                     !  = 0, the total water flux across open boundaries is zero 
     202/ 
     203!----------------------------------------------------------------------- 
     204&nambdy_index 
     205!----------------------------------------------------------------------- 
     206    ctypebdy = 'E' 
     207    nbdyind  = 49 
     208    nbdybeg  = 1 
     209    nbdyend  = 34 
     210    !ctypebdy = 'W' 
     211    !nbdyind  = 2 
     212    !nbdybeg  = 1 
     213    !nbdyend  = 34 
     214/ 
     215!----------------------------------------------------------------------- 
     216&nambdy_dta      !  open boundaries - external data 
     217!----------------------------------------------------------------------- 
     218!              !  file name      ! frequency (hours) ! variable  ! time interp. !  clim   ! 'yearly'/ ! weights  ! rotation ! land/sea mask ! 
     219!              !                 !  (if <0  months)  !   name    !  (logical)   !  (T/F ) ! 'monthly' ! filename ! pairing  ! filename      ! 
     220   bn_ssh =     'bdyssh_2.5_slow_stop' ,         1        , 'sshbdy',       .true.   , .true.  ,  'daily'  ,    ''    ,   ''     , '' 
     221   bn_u2d =     'bdyuv_2.5_slow'  ,         1        , 'ubdy'  ,     .true.     , .true.  ,  'daily'  ,    ''    ,   ''     , '' 
     222   bn_v2d =     'bdyuv_2.5_slow'  ,         1        , 'vbdy'  ,     .true.     , .true.  ,  'daily'  ,    ''    ,   ''     , '' 
     223!   bn_u3d  =    'amm12_bdyU_u3d' ,         24        , 'vozocrtx',     .true.   , .false. ,  'daily'  ,    ''    ,   ''     , '' 
     224!   bn_v3d  =    'amm12_bdyV_u3d' ,         24        , 'vomecrty',     .true.   , .false. ,  'daily'  ,    ''    ,   ''     , '' 
     225!   bn_tem  =    'amm12_bdyT_tra' ,         24        , 'votemper',     .true.   , .false. ,  'daily'  ,    ''    ,   ''     , '' 
     226!   bn_sal  =    'amm12_bdyT_tra' ,         24        , 'vosaline',     .true.   , .false. ,  'daily'  ,    ''    ,   ''     , '' 
     227   cn_dir      =    './'   !  root directory for the location of the bulk files 
     228   ln_full_vel = .false.        ! 
    166229/ 
    167230!----------------------------------------------------------------------- 
    … …  
    173236!----------------------------------------------------------------------- 
    174237   nn_bfr      =    2      !  type of bottom friction :   = 0 : free slip,  = 1 : linear friction 
     238   !rn_bfri2    =    1.e-5  !  bottom drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
     239   !rn_bfri2_max =   1.e-4  !  max. bottom drag coefficient (non linear case and ln_loglayer=T) 
     240   rn_bfri2    =    1.e-5  !  bottom drag coefficient (non linear case). Minimum coeft if ln_loglayer=T 
     241   rn_bfri2_max =   1.e-4  !  max. bottom drag coefficient (non linear case and ln_loglayer=T) 
     242   !rn_bfeb2    =    2.5e-3 !  bottom turbulent kinetic energy background  (m2/s2) 
     243   !rn_bfrz0    =    3.e-3  !  bottom roughness [m] if ln_loglayer=T 
     244   ln_loglayer = .true.    !  logarithmic formulation (non linear case) 
    175245/ 
    176246!----------------------------------------------------------------------- 
    … …  
    187257&nameos        !   ocean physical parameters 
    188258!----------------------------------------------------------------------- 
    189    ln_eos80    = .true.         !  = Use EOS80 equation of state 
     259   nn_eos      =  0       !  type of equation of state and Brunt-Vaisala frequency 
     260                                 !  =-1, TEOS-10  
     261                                 !  = 0, EOS-80  
     262                                 !  = 1, S-EOS   (simplified eos) 
     263   ln_useCT    = .false.  ! use of Conservative Temp. ==> surface CT converted in Pot. Temp. in sbcssm 
    190264   !                             ! 
    191265   !                      ! S-EOS coefficients : 
    … …  
    205279&namtra_adv    !   advection scheme for tracer 
    206280!----------------------------------------------------------------------- 
     281   ln_traadv_cen =  .false.  !  2nd order centered scheme 
     282   ln_traadv_mus =  .false.  !  MUSCL scheme 
    207283   ln_traadv_fct =  .true.   !  FCT scheme 
    208284      nn_fct_h   =  2               !  =2/4, horizontal 2nd / 4th order  
    … …  
    272348&namdyn_hpg    !   Hydrostatic pressure gradient option 
    273349!----------------------------------------------------------------------- 
    274    ln_hpg_zco  = .true.    !  z-coordinate - full steps 
     350   ln_hpg_zco  = .false.   !  z-coordinate - full steps 
    275351   ln_hpg_zps  = .false.   !  z-coordinate - partial steps (interpolation) 
     352   ln_hpg_sco  = .true.    !  s-coordinate 
    276353/ 
    277354!----------------------------------------------------------------------- 
    … …  
    300377   !                                !  = 30  F(i,j,k)=c2d*c1d 
    301378   !                                !  = 31  F(i,j,k)=F(grid spacing and local velocity) 
    302    rn_ahm_0      = 100000.     !  horizontal laplacian eddy viscosity   [m2/s] 
     379   rn_ahm_0      = 1000.        !  horizontal laplacian eddy viscosity   [m2/s] 
    303380   rn_ahm_b      =      0.     !  background eddy viscosity for ldf_iso [m2/s] 
    304381   rn_bhm_0      =      0.     !  horizontal bilaplacian eddy viscosity [m4/s] 
    … …  
    395472!----------------------------------------------------------------------- 
    396473/ 
     474!----------------------------------------------------------------------- 
     475&namwad  !   Wetting and drying 
     476!----------------------------------------------------------------------- 
     477   ln_wd             = .true.   ! T/F activation of wetting and drying 
     478   !rn_wdmin1         =  0.25    ! Minimum wet depth on dried cells 
     479   rn_wdmin1         =  0.4     ! Minimum wet depth on dried cells 
     480   rn_wdmin2         =  0.00001    ! Tolerance of min wet depth on dried cells 
     481   rn_wdld           =  10.0    ! Land elevation below which wetting/drying is allowed 
     482   nn_wdit           =  50      ! Max iterations for W/D limiter 
     483/ 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/MY_SRC/iom.F90

    r7403 r7412  
    7878   !!---------------------------------------------------------------------- 
    7979   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    80    !! $Id$ 
     80   !! $Id: iom.F90 6140 2015-12-21 11:35:23Z timgraham $ 
    8181   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    8282   !!---------------------------------------------------------------------- 
    … …  
    114114      CASE (30)   ;   CALL xios_set_context_attr(TRIM(clname), calendar_type= "D360") 
    115115      END SELECT 
    116       WRITE(cldate,"(i4.4,'-',i2.2,'-',i2.2,' ',i2.2,':',i2.2,':00')") nyear,nmonth,nday,nhour,nminute 
     116      WRITE(cldate,"(i4.4,'-',i2.2,'-',i2.2,' 00:00:00')") nyear,nmonth,nday  
    117117      CALL xios_set_context_attr(TRIM(clname), start_date=cldate ) 
    118118 
    … …  
    172172      z_bnds(1:jpkm1,2) = gdepw_1d(2:jpk) 
    173173      z_bnds(jpk:   ,2) = gdepw_1d(jpk) + e3t_1d(jpk) 
    174       CALL iom_set_axis_attr( "deptht", bounds=z_bnds ) 
    175       CALL iom_set_axis_attr( "depthu", bounds=z_bnds ) 
    176       CALL iom_set_axis_attr( "depthv", bounds=z_bnds ) 
     174      !CALL iom_set_axis_attr( "deptht", bounds=z_bnds ) 
     175      !CALL iom_set_axis_attr( "depthu", bounds=z_bnds ) 
     176      !CALL iom_set_axis_attr( "depthv", bounds=z_bnds ) 
    177177      z_bnds(:    ,2) = gdept_1d(:) 
    178178      z_bnds(2:jpk,1) = gdept_1d(1:jpkm1) 
    179179      z_bnds(1    ,1) = gdept_1d(1) - e3w_1d(1) 
    180       CALL iom_set_axis_attr( "depthw", bounds=z_bnds ) 
     180      !CALL iom_set_axis_attr( "depthw", bounds=z_bnds ) 
    181181 
    182182# if defined key_floats 

  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/cfg.txt

    r7403 r7412  
    1212ORCA2_LIM_PISCES OPA_SRC LIM_SRC_2 NST_SRC TOP_SRC 
    1313ORCA2_LIM3_TRC OPA_SRC LIM_SRC_3 NST_SRC TOP_SRC 
     14WAD_TEST_CASES OPA_SRC 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_2/limdyn_2.F90

    r5123 r7412  
    8787         ! --------------------------------------------------- 
    8888          
    89          IF( lk_mpp .OR. lk_mpp_rep ) THEN                    ! mpp: compute over the whole domain 
     89         IF( lk_mpp ) THEN                    ! mpp: compute over the whole domain 
    9090            i_j1 = 1    
    9191            i_jpj = jpj 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90

    r6416 r7412  
    286286      REAL(wp), PARAMETER             :: zconv = 1.e-9 ! convert W to GW and kg to Mt 
    287287 
    288 #if ! defined key_bdy 
    289288      ! heat flux 
    290289      zhfx  = glob_sum( ( hfx_in - hfx_out - diag_heat - diag_trp_ei - diag_trp_es - SUM( qevap_ice * a_i_b, dim=3 ) )  & 
    … …  
    304303      IF( ABS( zsfx ) > zs_sill ) WRITE(numout,*) 'violation sfx    [psu*Mt/day]   (',cd_routine,')  = ',(zsfx) 
    305304      IF( ABS( zhfx ) > zh_sill ) WRITE(numout,*) 'violation hfx    [GW]           (',cd_routine,')  = ',(zhfx) 
    306 #endif 
    307305 
    308306   END SUBROUTINE lim_cons_final 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90

    r5836 r7412  
    9494         ! --------------------------------------------------- 
    9595 
    96          IF( lk_mpp .OR. lk_mpp_rep ) THEN                    ! mpp: compute over the whole domain 
     96         IF( lk_mpp ) THEN                    ! mpp: compute over the whole domain 
    9797            i_j1 = 1 
    9898            i_jpj = jpj 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

    r6416 r7412  
    4141   USE agrif_lim2_interp 
    4242#endif 
    43 #if defined key_bdy 
     43   USE bdy_oce   , ONLY: ln_bdy 
    4444   USE bdyice_lim 
    45 #endif 
    4645 
    4746   IMPLICIT NONE 
    … …  
    460459            CALL agrif_rhg_lim2( jter, nn_nevp, 'U' ) 
    461460#endif 
    462 #if defined key_bdy 
    463          CALL bdy_ice_lim_dyn( 'U' ) 
    464 #endif          
     461         IF( ln_bdy ) CALL bdy_ice_lim_dyn( 'U' ) 
    465462 
    466463            DO jj = k_j1+1, k_jpj-1 
    … …  
    486483            CALL agrif_rhg_lim2( jter, nn_nevp, 'V' ) 
    487484#endif 
    488 #if defined key_bdy 
    489          CALL bdy_ice_lim_dyn( 'V' ) 
    490 #endif          
     485         IF( ln_bdy ) CALL bdy_ice_lim_dyn( 'V' ) 
    491486 
    492487         ELSE  
    … …  
    513508            CALL agrif_rhg_lim2( jter, nn_nevp, 'V' ) 
    514509#endif 
    515 #if defined key_bdy 
    516          CALL bdy_ice_lim_dyn( 'V' ) 
    517 #endif          
     510         IF( ln_bdy ) CALL bdy_ice_lim_dyn( 'V' ) 
    518511 
    519512            DO jj = k_j1+1, k_jpj-1 
    … …  
    538531            CALL agrif_rhg_lim2( jter, nn_nevp, 'U' ) 
    539532#endif 
    540 #if defined key_bdy 
    541          CALL bdy_ice_lim_dyn( 'U' ) 
    542 #endif          
     533         IF( ln_bdy ) CALL bdy_ice_lim_dyn( 'U' ) 
    543534 
    544535         ENDIF 
    … …  
    577568      CALL agrif_rhg_lim2( nn_nevp , nn_nevp, 'V' ) 
    578569#endif 
    579 #if defined key_bdy 
    580       CALL bdy_ice_lim_dyn( 'U' ) 
    581       CALL bdy_ice_lim_dyn( 'V' ) 
    582 #endif          
     570      IF( ln_bdy ) THEN 
     571         CALL bdy_ice_lim_dyn( 'U' ) ; CALL bdy_ice_lim_dyn( 'V' ) 
     572      ENDIF 
    583573 
    584574      DO jj = k_j1+1, k_jpj-1  

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

    r6416 r7412  
    3636   USE limctl         !  
    3737   USE limcons        !  
     38   USE bdy_oce  , ONLY: ln_bdy 
    3839   ! 
    3940   USE in_out_manager ! I/O manager 
    … …  
    221222 
    222223      ! conservation test 
    223       IF( ln_limdiahsb )   CALL lim_cons_final( 'limsbc' ) 
     224      IF( ln_limdiahsb .AND. .NOT. ln_bdy)   CALL lim_cons_final( 'limsbc' ) 
    224225 
    225226      ! control prints 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/NST_SRC/agrif_user.F90

    r6140 r7412  
    6161   USE nemogcm 
    6262   USE tradmp 
    63    USE bdy_par 
     63   USE bdy_oce   , ONLY: ln_bdy 
    6464 
    6565   IMPLICIT NONE 
    … …  
    7878   ln_tradmp = .FALSE. 
    7979   ! no open boundary on fine grids 
    80    lk_bdy = .FALSE. 
     80   ln_bdy = .FALSE. 
    8181 
    8282 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdy_oce.F90

    r6140 r7412  
    1010   !!            3.6  !  2014-01  (C. Rousset) add ice boundary conditions for lim3 
    1111   !!---------------------------------------------------------------------- 
    12 #if defined key_bdy  
    13    !!---------------------------------------------------------------------- 
    14    !!   'key_bdy'                      Unstructured Open Boundary Condition 
    15    !!---------------------------------------------------------------------- 
    1612   USE par_oce         ! ocean parameters 
    17    USE bdy_par         ! Unstructured boundary parameters 
    1813   USE lib_mpp         ! distributed memory computing 
    1914 
    2015   IMPLICIT NONE 
    2116   PUBLIC 
     17 
     18   INTEGER, PUBLIC, PARAMETER ::   jp_bdy  = 10       !: Maximum number of bdy sets 
     19   INTEGER, PUBLIC, PARAMETER ::   jpbgrd  = 3        !: Number of horizontal grid types used  (T, U, V) 
    2220 
    2321   TYPE, PUBLIC ::   OBC_INDEX    !: Indices and weights which define the open boundary 
    … …  
    4947      LOGICAL                           ::  ll_tem 
    5048      LOGICAL                           ::  ll_sal 
     49      LOGICAL                           ::  ll_fvl 
    5150      REAL(wp), POINTER, DIMENSION(:)   ::  ssh 
    5251      REAL(wp), POINTER, DIMENSION(:)   ::  u2d 
    … …  
    8281   !! Namelist variables 
    8382   !!---------------------------------------------------------------------- 
     83   LOGICAL, PUBLIC            ::   ln_bdy                   !: Unstructured Ocean Boundary Condition 
     84 
    8485   CHARACTER(len=80), DIMENSION(jp_bdy) ::   cn_coords_file !: Name of bdy coordinates file 
    8586   CHARACTER(len=80)                    ::   cn_mask_file   !: Name of bdy mask file 
    … …  
    9192   ! 
    9293   INTEGER                    ::   nb_bdy                   !: number of open boundary sets 
     94   INTEGER                    ::   nb_jpk_bdy               !: number of levels in the bdy data (set < 0 if consistent with planned run) 
    9395   INTEGER, DIMENSION(jp_bdy) ::   nn_rimwidth              !: boundary rim width for Flow Relaxation Scheme 
    9496   INTEGER                    ::   nn_volctl                !: = 0 the total volume will have the variability of the surface Flux E-P  
    … …  
    134136                                                                          !: =1 => some data to be read in from data files 
    135137   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), TARGET ::   dta_global        !: workspace for reading in global data arrays (unstr.  bdy) 
     138   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), TARGET ::   dta_global_z      !: workspace for reading in global depth arrays (unstr.  bdy) 
     139   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), TARGET ::   dta_global_dz     !: workspace for reading in global depth arrays (unstr.  bdy) 
    136140   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), TARGET ::   dta_global2       !: workspace for reading in global data arrays (struct. bdy) 
     141   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), TARGET ::   dta_global2_z     !: workspace for reading in global depth arrays (struct. bdy) 
     142   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), TARGET ::   dta_global2_dz    !: workspace for reading in global depth arrays (struct. bdy) 
    137143!$AGRIF_DO_NOT_TREAT 
    138144   TYPE(OBC_INDEX), DIMENSION(jp_bdy), TARGET      ::   idx_bdy           !: bdy indices (local process) 
    … …  
    166172   END FUNCTION bdy_oce_alloc 
    167173 
    168 #else 
    169    !!---------------------------------------------------------------------- 
    170    !!   Dummy module                NO Unstructured Open Boundary Condition 
    171    !!---------------------------------------------------------------------- 
    172    LOGICAL ::   ln_tides = .false.  !: =T apply tidal harmonic forcing along open boundaries 
    173 #endif 
    174  
    175174   !!====================================================================== 
    176175END MODULE bdy_oce 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90

    r6140 r7412  
    1212   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge 
    1313   !!            3.6  !  2012-01  (C. Rousset) add ice boundary conditions for lim3 
    14    !!---------------------------------------------------------------------- 
    15 #if defined key_bdy 
    16    !!---------------------------------------------------------------------- 
    17    !!   'key_bdy'                     Open Boundary Conditions 
    1814   !!---------------------------------------------------------------------- 
    1915   !!    bdy_dta        : read external data along open boundaries from file 
    … …  
    3632#endif 
    3733   USE sbcapr 
     34   USE sbctide         ! Tidal forcing or not 
    3835 
    3936   IMPLICIT NONE 
    … …  
    267264 
    268265                        jend = jstart + dta%nread(2) - 1 
    269                         CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  & 
    270                                      & kit=jit, kt_offset=time_offset ) 
     266                        IF( ln_full_vel_array(ib_bdy) ) THEN 
     267                           CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  & 
     268                                     & kit=jit, kt_offset=time_offset , jpk_bdy=nb_jpk_bdy, fvl=ln_full_vel_array(ib_bdy)  ) 
     269                        ELSE 
     270                           CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  & 
     271                                     & kit=jit, kt_offset=time_offset  ) 
     272                        ENDIF 
    271273 
    272274                        ! If full velocities in boundary data then extract barotropic velocities from 3D fields 
    … …  
    333335                     jend = jstart + dta%nread(1) - 1 
    334336                     CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 
    335                                   & map=nbmap_ptr(jstart:jend), kt_offset=time_offset ) 
     337                                  & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy, fvl=ln_full_vel_array(ib_bdy) ) 
    336338                  ENDIF 
    337339                  ! If full velocities in boundary data then split into barotropic and baroclinic data 
    … …  
    381383      END DO  ! ib_bdy 
    382384 
    383 #if defined key_tide 
    384       IF (ln_dynspg_ts) THEN      ! Fill temporary arrays with slow-varying bdy data                            
    385          DO ib_bdy = 1, nb_bdy    ! Tidal component added in ts loop 
    386             IF ( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN 
    387                nblen => idx_bdy(ib_bdy)%nblen 
    388                nblenrim => idx_bdy(ib_bdy)%nblenrim 
    389                IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN; ilen1(:)=nblen(:) ; ELSE ; ilen1(:)=nblenrim(:) ; ENDIF  
    390                IF ( dta_bdy(ib_bdy)%ll_ssh ) dta_bdy_s(ib_bdy)%ssh(1:ilen1(1)) = dta_bdy(ib_bdy)%ssh(1:ilen1(1)) 
    391                IF ( dta_bdy(ib_bdy)%ll_u2d ) dta_bdy_s(ib_bdy)%u2d(1:ilen1(2)) = dta_bdy(ib_bdy)%u2d(1:ilen1(2)) 
    392                IF ( dta_bdy(ib_bdy)%ll_v2d ) dta_bdy_s(ib_bdy)%v2d(1:ilen1(3)) = dta_bdy(ib_bdy)%v2d(1:ilen1(3)) 
    393             ENDIF 
    394          END DO 
    395       ELSE ! Add tides if not split-explicit free surface else this is done in ts loop 
    396          ! 
    397          CALL bdy_dta_tides( kt=kt, time_offset=time_offset ) 
     385      IF ( ln_tide ) THEN 
     386         IF (ln_dynspg_ts) THEN      ! Fill temporary arrays with slow-varying bdy data                            
     387            DO ib_bdy = 1, nb_bdy    ! Tidal component added in ts loop 
     388               IF ( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN 
     389                  nblen => idx_bdy(ib_bdy)%nblen 
     390                  nblenrim => idx_bdy(ib_bdy)%nblenrim 
     391                  IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN; ilen1(:)=nblen(:) ; ELSE ; ilen1(:)=nblenrim(:) ; ENDIF  
     392                  IF ( dta_bdy(ib_bdy)%ll_ssh ) dta_bdy_s(ib_bdy)%ssh(1:ilen1(1)) = dta_bdy(ib_bdy)%ssh(1:ilen1(1)) 
     393                  IF ( dta_bdy(ib_bdy)%ll_u2d ) dta_bdy_s(ib_bdy)%u2d(1:ilen1(2)) = dta_bdy(ib_bdy)%u2d(1:ilen1(2)) 
     394                  IF ( dta_bdy(ib_bdy)%ll_v2d ) dta_bdy_s(ib_bdy)%v2d(1:ilen1(3)) = dta_bdy(ib_bdy)%v2d(1:ilen1(3)) 
     395               ENDIF 
     396            END DO 
     397         ELSE ! Add tides if not split-explicit free surface else this is done in ts loop 
     398            ! 
     399            CALL bdy_dta_tides( kt=kt, time_offset=time_offset ) 
     400         ENDIF 
    398401      ENDIF 
    399 #endif 
    400402 
    401403      IF ( ln_apr_obc ) THEN 
    … …  
    459461      NAMELIST/nambdy_dta/ bn_a_i, bn_ht_i, bn_ht_s 
    460462#endif 
    461       NAMELIST/nambdy_dta/ ln_full_vel 
     463      NAMELIST/nambdy_dta/ ln_full_vel, nb_jpk_bdy 
    462464      !!--------------------------------------------------------------------------- 
    463465      ! 
    … …  
    899901   END SUBROUTINE bdy_dta_init 
    900902 
    901 #else 
    902    !!---------------------------------------------------------------------- 
    903    !!   Dummy module                   NO Open Boundary Conditions 
    904    !!---------------------------------------------------------------------- 
    905 CONTAINS 
    906    SUBROUTINE bdy_dta( kt, jit, time_offset ) ! Empty routine 
    907       INTEGER, INTENT( in )           ::   kt     
    908       INTEGER, INTENT( in ), OPTIONAL ::   jit    
    909       INTEGER, INTENT( in ), OPTIONAL ::   time_offset 
    910       WRITE(*,*) 'bdy_dta: You should not have seen this print! error?', kt 
    911    END SUBROUTINE bdy_dta 
    912    SUBROUTINE bdy_dta_init()                  ! Empty routine 
    913       WRITE(*,*) 'bdy_dta_init: You should not have seen this print! error?' 
    914    END SUBROUTINE bdy_dta_init 
    915 #endif 
    916  
    917903   !!============================================================================== 
    918904END MODULE bdydta 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90

    r6140 r7412  
    1111   !!            3.3  !  2010-09  (D.Storkey) add ice boundary conditions 
    1212   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge 
    13    !!---------------------------------------------------------------------- 
    14 #if defined key_bdy  
    15    !!---------------------------------------------------------------------- 
    16    !!   'key_bdy' :                    Unstructured Open Boundary Condition 
    1713   !!---------------------------------------------------------------------- 
    1814   !!   bdy_dyn        : split velocities into barotropic and baroclinic parts 
    … …  
    137133   END SUBROUTINE bdy_dyn 
    138134 
    139 #else 
    140    !!---------------------------------------------------------------------- 
    141    !!   Dummy module                   NO Unstruct Open Boundary Conditions 
    142    !!---------------------------------------------------------------------- 
    143 CONTAINS 
    144    SUBROUTINE bdy_dyn( kt )      ! Empty routine 
    145       WRITE(*,*) 'bdy_dyn: You should not have seen this print! error?', kt 
    146    END SUBROUTINE bdy_dyn 
    147 #endif 
    148  
    149135   !!====================================================================== 
    150136END MODULE bdydyn 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn2d.F90

    r5930 r7412  
    77   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) Optimization of BDY communications 
    88   !!            3.5  !  2013-07  (J. Chanut) Compliant with time splitting changes 
    9    !!---------------------------------------------------------------------- 
    10 #if defined key_bdy  
    11    !!---------------------------------------------------------------------- 
    12    !!   'key_bdy' :                    Unstructured Open Boundary Condition 
    139   !!---------------------------------------------------------------------- 
    1410   !!   bdy_dyn2d          : Apply open boundary conditions to barotropic variables. 
    … …  
    310306   END SUBROUTINE bdy_ssh 
    311307 
    312 #else 
    313    !!---------------------------------------------------------------------- 
    314    !!   Dummy module                   NO Unstruct Open Boundary Conditions 
    315    !!---------------------------------------------------------------------- 
    316 CONTAINS 
    317    SUBROUTINE bdy_dyn2d( kt )      ! Empty routine 
    318       INTEGER, intent(in) :: kt 
    319       WRITE(*,*) 'bdy_dyn2d: You should not have seen this print! error?', kt 
    320    END SUBROUTINE bdy_dyn2d 
    321  
    322 #endif 
    323  
    324308   !!====================================================================== 
    325309END MODULE bdydyn2d 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn3d.F90

    r6140 r7412  
    66   !! History :  3.4  !  2011     (D. Storkey) new module as part of BDY rewrite  
    77   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) Optimization of BDY communications 
    8    !!---------------------------------------------------------------------- 
    9 #if defined key_bdy  
    10    !!---------------------------------------------------------------------- 
    11    !!   'key_bdy' :                    Unstructured Open Boundary Condition 
    128   !!---------------------------------------------------------------------- 
    139   !!   bdy_dyn3d        : apply open boundary conditions to baroclinic velocities 
    … …  
    5753         CASE('orlanski' )   ;   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.false. ) 
    5854         CASE('orlanski_npo');   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.true. ) 
     55         CASE('zerograd')    ;   CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
     56         CASE('neumann')     ;   CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy ) 
    5957         CASE DEFAULT        ;   CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 
    6058         END SELECT 
    … …  
    110108   END SUBROUTINE bdy_dyn3d_spe 
    111109 
     110   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt , ib_bdy ) 
     111      !!---------------------------------------------------------------------- 
     112      !!                  ***  SUBROUTINE bdy_dyn3d_zgrad  *** 
     113      !! 
     114      !! ** Purpose : - Enforce a zero gradient of normal velocity 
     115      !! 
     116      !!---------------------------------------------------------------------- 
     117      INTEGER                     ::   kt 
     118      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices 
     119      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data 
     120      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index 
     121      !! 
     122      INTEGER  ::   jb, jk         ! dummy loop indices 
     123      INTEGER  ::   ii, ij, igrd   ! local integers 
     124      REAL(wp) ::   zwgt           ! boundary weight 
     125      INTEGER  ::   fu, fv 
     126      !!---------------------------------------------------------------------- 
     127      ! 
     128      IF( nn_timing == 1 ) CALL timing_start('bdy_dyn3d_zgrad') 
     129      ! 
     130      igrd = 2                      ! Copying tangential velocity into bdy points 
     131      DO jb = 1, idx%nblenrim(igrd) 
     132         DO jk = 1, jpkm1 
     133            ii   = idx%nbi(jb,igrd) 
     134            ij   = idx%nbj(jb,igrd) 
     135            fu   = ABS( ABS (NINT( idx%flagu(jb,igrd) ) ) - 1 ) 
     136            ua(ii,ij,jk) = ua(ii,ij,jk) * REAL( 1 - fu) + ( ua(ii,ij+fu,jk) * umask(ii,ij+fu,jk) & 
     137                        &+ ua(ii,ij-fu,jk) * umask(ii,ij-fu,jk) ) * umask(ii,ij,jk) * REAL( fu ) 
     138         END DO 
     139      END DO 
     140      ! 
     141      igrd = 3                      ! Copying tangential velocity into bdy points 
     142      DO jb = 1, idx%nblenrim(igrd) 
     143         DO jk = 1, jpkm1 
     144            ii   = idx%nbi(jb,igrd) 
     145            ij   = idx%nbj(jb,igrd) 
     146            fv   = ABS( ABS (NINT( idx%flagv(jb,igrd) ) ) - 1 ) 
     147            va(ii,ij,jk) = va(ii,ij,jk) * REAL( 1 - fv ) + ( va(ii+fv,ij,jk) * vmask(ii+fv,ij,jk) & 
     148                        &+ va(ii-fv,ij,jk) * vmask(ii-fv,ij,jk) ) * vmask(ii,ij,jk) * REAL( fv ) 
     149         END DO 
     150      END DO 
     151      CALL lbc_bdy_lnk( ua, 'U', -1., ib_bdy )   ! Boundary points should be updated   
     152      CALL lbc_bdy_lnk( va, 'V', -1., ib_bdy )    
     153      ! 
     154      IF( kt .eq. nit000 ) CLOSE( unit = 102 ) 
     155 
     156      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn3d_zgrad') 
     157 
     158   END SUBROUTINE bdy_dyn3d_zgrad 
    112159 
    113160   SUBROUTINE bdy_dyn3d_zro( idx, dta, kt, ib_bdy ) 
    … …  
    296343   END SUBROUTINE bdy_dyn3d_dmp 
    297344 
    298 #else 
    299    !!---------------------------------------------------------------------- 
    300    !!   Dummy module                   NO Unstruct Open Boundary Conditions 
    301    !!---------------------------------------------------------------------- 
    302 CONTAINS 
    303    SUBROUTINE bdy_dyn3d( kt )      ! Empty routine 
    304       WRITE(*,*) 'bdy_dyn3d: You should not have seen this print! error?', kt 
    305    END SUBROUTINE bdy_dyn3d 
    306    SUBROUTINE bdy_dyn3d_dmp( kt )      ! Empty routine 
    307       WRITE(*,*) 'bdy_dyn3d_dmp: You should not have seen this print! error?', kt 
    308    END SUBROUTINE bdy_dyn3d_dmp 
    309 #endif 
     345   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy ) 
     346      !!---------------------------------------------------------------------- 
     347      !!                 ***  SUBROUTINE bdy_dyn3d_nmn  *** 
     348      !!              
     349      !!              - Apply Neumann condition to baroclinic velocities.  
     350      !!              - Wrapper routine for bdy_nmn 
     351      !!  
     352      !! 
     353      !!---------------------------------------------------------------------- 
     354      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices 
     355      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index 
     356 
     357      INTEGER  ::   jb, igrd                               ! dummy loop indices 
     358      !!---------------------------------------------------------------------- 
     359 
     360      IF( nn_timing == 1 ) CALL timing_start('bdy_dyn3d_nmn') 
     361      ! 
     362      !! Note that at this stage the ub and ua arrays contain the baroclinic velocities.  
     363      ! 
     364      igrd = 2      ! Neumann bc on u-velocity;  
     365      !             
     366      CALL bdy_nmn( idx, igrd, ua ) 
     367 
     368      igrd = 3      ! Neumann bc on v-velocity 
     369      !   
     370      CALL bdy_nmn( idx, igrd, va ) 
     371      ! 
     372      CALL lbc_bdy_lnk( ua, 'U', -1., ib_bdy )    ! Boundary points should be updated 
     373      CALL lbc_bdy_lnk( va, 'V', -1., ib_bdy ) 
     374      ! 
     375      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn3d_nmn') 
     376      ! 
     377   END SUBROUTINE bdy_dyn3d_nmn 
    310378 
    311379   !!====================================================================== 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90

    r5836 r7412  
    88   !!              -   !  2012-01 (C. Rousset)  add lim3 and remove useless jk loop  
    99   !!---------------------------------------------------------------------- 
    10 #if defined   key_bdy   &&  ( defined key_lim2 || defined key_lim3 ) 
    11    !!---------------------------------------------------------------------- 
    12    !!   'key_bdy'            and                 Unstructured Open Boundary Conditions 
     10#if defined key_lim2 || defined key_lim3 
     11   !!---------------------------------------------------------------------- 
    1312   !!   'key_lim2'                                                 LIM-2 sea ice model 
    1413   !!   'key_lim3'                                                 LIM-3 sea ice model 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90

    r6140 r7412  
    1313   !!            3.4  !  2012     (J. Chanut) straight open boundary case update 
    1414   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) optimization of BDY communications 
    15    !!---------------------------------------------------------------------- 
    16 #if defined key_bdy 
    17    !!---------------------------------------------------------------------- 
    18    !!   'key_bdy'                     Unstructured Open Boundary Conditions 
     15   !!            3.7  !  2016     (T. Lovato) Remove bdy macro, call here init for dta and tides 
    1916   !!---------------------------------------------------------------------- 
    2017   !!   bdy_init      : Initialization of unstructured open boundaries 
    … …  
    2320   USE dom_oce        ! ocean space and time domain 
    2421   USE bdy_oce        ! unstructured open boundary conditions 
    25    USE sbctide  , ONLY: lk_tide ! Tidal forcing or not 
     22   USE bdydta         ! open boundary cond. setting   (bdy_dta_init routine) 
     23   USE bdytides       ! open boundary cond. setting   (bdytide_init routine) 
     24   USE sbctide        ! Tidal forcing or not 
    2625   USE phycst   , ONLY: rday 
    2726   ! 
    … …  
    5352   !!---------------------------------------------------------------------- 
    5453CONTAINS 
    55     
     54 
    5655   SUBROUTINE bdy_init 
    5756      !!---------------------------------------------------------------------- 
    5857      !!                 ***  ROUTINE bdy_init  *** 
     58      !! 
     59      !! ** Purpose :   Initialization of the dynamics and tracer fields with 
     60      !!              unstructured open boundaries. 
     61      !! 
     62      !! ** Method  :   Read initialization arrays (mask, indices) to identify 
     63      !!              an unstructured open boundary 
     64      !! 
     65      !! ** Input   :  bdy_init.nc, input file for unstructured open boundaries 
     66      !!---------------------------------------------------------------------- 
     67      NAMELIST/nambdy/ ln_bdy, nb_bdy, ln_coords_file, cn_coords_file,         & 
     68         &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,     & 
     69         &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             & 
     70         &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, & 
     71         &             cn_ice_lim, nn_ice_lim_dta,                             & 
     72         &             rn_ice_tem, rn_ice_sal, rn_ice_age,                     & 
     73         &             ln_vol, nn_volctl, nn_rimwidth, nb_jpk_bdy 
     74         ! 
     75      INTEGER  ::   ios                 ! Local integer output status for namelist read 
     76      !!---------------------------------------------------------------------- 
     77      ! 
     78      IF( nn_timing == 1 )   CALL timing_start('bdy_init') 
     79 
     80      ! ------------------------ 
     81      ! Read namelist parameters 
     82      ! ------------------------ 
     83      REWIND( numnam_ref )              ! Namelist nambdy in reference namelist :Unstructured open boundaries 
     84      READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 901) 
     85901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in reference namelist', lwp ) 
     86      ! 
     87      REWIND( numnam_cfg )              ! Namelist nambdy in configuration namelist :Unstructured open boundaries 
     88      READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 902 ) 
     89902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in configuration namelist', lwp ) 
     90      IF(lwm) WRITE ( numond, nambdy ) 
     91 
     92      ! ----------------------------------------- 
     93      ! unstructured open boundaries use control 
     94      ! ----------------------------------------- 
     95      IF ( ln_bdy ) THEN 
     96         IF(lwp) WRITE(numout,*) 
     97         IF(lwp) WRITE(numout,*) 'bdy_init : initialization of open boundaries' 
     98         IF(lwp) WRITE(numout,*) '~~~~~~~~' 
     99         ! 
     100         ! Open boundaries definition (arrays and masks) 
     101         CALL bdy_segs 
     102         ! 
     103         ! Open boundaries initialisation of external data arrays 
     104         CALL bdy_dta_init 
     105         ! 
     106         ! Open boundaries initialisation of tidal harmonic forcing 
     107         IF( ln_tide ) CALL bdytide_init 
     108         ! 
     109      ELSE 
     110         IF(lwp) WRITE(numout,*) 
     111         IF(lwp) WRITE(numout,*) 'bdy_init : open boundaries not used (ln_bdy = F)' 
     112         IF(lwp) WRITE(numout,*) '~~~~~~~~' 
     113         ! 
     114      ENDIF 
     115      ! 
     116      IF( nn_timing == 1 )   CALL timing_stop('bdy_init') 
     117      ! 
     118   END SUBROUTINE bdy_init 
     119    
     120   SUBROUTINE bdy_segs 
     121      !!---------------------------------------------------------------------- 
     122      !!                 ***  ROUTINE bdy_init  *** 
    59123      !!          
    60       !! ** Purpose :   Initialization of the dynamics and tracer fields with  
    61       !!              unstructured open boundaries. 
     124      !! ** Purpose :   Definition of unstructured open boundaries. 
    62125      !! 
    63126      !! ** Method  :   Read initialization arrays (mask, indices) to identify  
    … …  
    90153      REAL(wp), POINTER, DIMENSION(:,:)       ::   zfmask  ! temporary fmask array excluding coastal boundary condition (shlat) 
    91154      !! 
    92       CHARACTER(LEN=80),DIMENSION(jpbgrd)  ::   clfile     ! Namelist variables 
    93155      CHARACTER(LEN=1)                     ::   ctypebdy   !     -        -  
    94156      INTEGER                              ::   nbdyind, nbdybeg, nbdyend 
    95157      !! 
    96       NAMELIST/nambdy/ nb_bdy, ln_coords_file, cn_coords_file,                 & 
    97          &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,     & 
    98          &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             &   
    99          &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, & 
    100          &             cn_ice_lim, nn_ice_lim_dta,                           & 
    101          &             rn_ice_tem, rn_ice_sal, rn_ice_age,                 & 
    102          &             ln_vol, nn_volctl, nn_rimwidth 
    103          ! 
    104158      NAMELIST/nambdy_index/ ctypebdy, nbdyind, nbdybeg, nbdyend 
    105159      INTEGER  ::   ios                 ! Local integer output status for namelist read 
    106160      !!---------------------------------------------------------------------- 
    107161      ! 
    108       IF( nn_timing == 1 )   CALL timing_start('bdy_init') 
    109       ! 
    110       IF(lwp) WRITE(numout,*) 
    111       IF(lwp) WRITE(numout,*) 'bdy_init : initialization of open boundaries' 
    112       IF(lwp) WRITE(numout,*) '~~~~~~~~' 
    113       ! 
    114       IF( jperio /= 0 )   CALL ctl_stop( 'Cyclic or symmetric,',   & 
    115          &                               ' and general open boundary condition are not compatible' ) 
    116  
     162      IF( nn_timing == 1 )   CALL timing_start('bdy_segs') 
     163      ! 
    117164      cgrid = (/'t','u','v'/) 
    118        
    119       ! ------------------------ 
    120       ! Read namelist parameters 
    121       ! ------------------------ 
    122       REWIND( numnam_ref )              ! Namelist nambdy in reference namelist :Unstructured open boundaries   
    123       READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 901) 
    124 901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in reference namelist', lwp ) 
    125       ! 
    126       REWIND( numnam_cfg )              ! Namelist nambdy in configuration namelist :Unstructured open boundaries 
    127       READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 902 ) 
    128 902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in configuration namelist', lwp ) 
    129       IF(lwm) WRITE ( numond, nambdy ) 
    130165 
    131166      ! ----------------------------------------- 
    132167      ! Check and write out namelist parameters 
    133168      ! ----------------------------------------- 
    134       !                                   ! control prints 
    135       IF(lwp) WRITE(numout,*) '   nambdy' 
     169      IF( jperio /= 0 )   CALL ctl_stop( 'bdy_segs: Cyclic or symmetric,',   & 
     170         &                               ' and general open boundary condition are not compatible' ) 
    136171 
    137172      IF( nb_bdy == 0 ) THEN  
    … …  
    189224              CASE DEFAULT   ;   CALL ctl_stop( 'nn_dyn2d_dta must be between 0 and 3' ) 
    190225           END SELECT 
    191            IF (( nn_dyn2d_dta(ib_bdy) .ge. 2 ).AND.(.NOT.lk_tide)) THEN 
    192              CALL ctl_stop( 'You must activate key_tide to add tidal forcing at open boundaries' ) 
     226           IF (( nn_dyn2d_dta(ib_bdy) .ge. 2 ).AND.(.NOT.ln_tide)) THEN 
     227             CALL ctl_stop( 'You must activate with ln_tide to add tidal forcing at open boundaries' ) 
    193228           ENDIF 
    194229        ENDIF 
    … …  
    209244             dta_bdy(ib_bdy)%ll_u3d = .true. 
    210245             dta_bdy(ib_bdy)%ll_v3d = .true. 
     246          CASE('neumann') 
     247             IF(lwp) WRITE(numout,*) '      Neumann conditions' 
     248             dta_bdy(ib_bdy)%ll_u3d = .false. 
     249             dta_bdy(ib_bdy)%ll_v3d = .false. 
     250          CASE('zerograd') 
     251             IF(lwp) WRITE(numout,*) '      Zero gradient for baroclinic velocities' 
     252             dta_bdy(ib_bdy)%ll_u3d = .false. 
     253             dta_bdy(ib_bdy)%ll_v3d = .false. 
    211254          CASE('zero') 
    212255             IF(lwp) WRITE(numout,*) '      Zero baroclinic velocities (runoff case)' 
    … …  
    377420          IF(lwp) WRITE(numout,*) 'No volume correction applied at open boundaries' 
    378421          IF(lwp) WRITE(numout,*) 
     422        ENDIF 
     423        IF( nb_jpk_bdy > 0 ) THEN 
     424           IF(lwp) WRITE(numout,*) '*** open boundary will be interpolate in the vertical onto the native grid ***' 
     425        ELSE 
     426           IF(lwp) WRITE(numout,*) '*** open boundary will be read straight onto the native grid without vertical interpolation ***' 
    379427        ENDIF 
    380428     ENDIF 
    … …  
    499547            &      nbrdta(jpbdta, jpbgrd, nb_bdy) ) 
    500548 
    501          ALLOCATE( dta_global(jpbdtau, 1, jpk) ) 
    502          IF ( icount>0 ) ALLOCATE( dta_global2(jpbdtas, nrimmax, jpk) ) 
     549         IF( nb_jpk_bdy>0 ) THEN 
     550            ALLOCATE( dta_global(jpbdtau, 1, nb_jpk_bdy) ) 
     551            ALLOCATE( dta_global_z(jpbdtau, 1, nb_jpk_bdy) ) 
     552            ALLOCATE( dta_global_dz(jpbdtau, 1, nb_jpk_bdy) ) 
     553         ELSE 
     554            ALLOCATE( dta_global(jpbdtau, 1, jpk) ) 
     555            ALLOCATE( dta_global_z(jpbdtau, 1, jpk) ) ! needed ?? TODO 
     556            ALLOCATE( dta_global_dz(jpbdtau, 1, jpk) )! needed ?? TODO 
     557         ENDIF 
     558 
     559         IF ( icount>0 ) THEN 
     560            IF( nb_jpk_bdy>0 ) THEN 
     561               ALLOCATE( dta_global2(jpbdtas, nrimmax, nb_jpk_bdy) ) 
     562               ALLOCATE( dta_global2_z(jpbdtas, nrimmax, nb_jpk_bdy) ) 
     563               ALLOCATE( dta_global2_dz(jpbdtas, nrimmax, nb_jpk_bdy) ) 
     564            ELSE 
     565               ALLOCATE( dta_global2(jpbdtas, nrimmax, jpk) ) 
     566               ALLOCATE( dta_global2_z(jpbdtas, nrimmax, jpk) ) ! needed ?? TODO 
     567               ALLOCATE( dta_global2_dz(jpbdtas, nrimmax, jpk) )! needed ?? TODO   
     568            ENDIF 
     569         ENDIF 
    503570         !  
    504571      ENDIF 
    … …  
    839906               IF(lwp) THEN         ! Since all procs read global data only need to do this check on one proc... 
    840907                  IF( nbrdta(ib,igrd,ib_bdy) < nbrdta(ibm1,igrd,ib_bdy) ) THEN 
    841                      CALL ctl_stop('bdy_init : ERROR : boundary data in file must be defined ', & 
     908                     CALL ctl_stop('bdy_segs : ERROR : boundary data in file must be defined ', & 
    842909                          &        ' in order of distance from edge nbr A utility for re-ordering ', & 
    843910                          &        ' boundary coordinates and data files exists in the TOOLS/OBC directory') 
    … …  
    10921159      !          = 0  elsewhere    
    10931160  
     1161      bdytmask(:,:) = ssmask(:,:) 
     1162 
    10941163      IF( ln_mask_file ) THEN 
    10951164         CALL iom_open( cn_mask_file, inum ) 
    … …  
    11081177         CALL lbc_lnk( bdyumask(:,:), 'U', 1. )   ;   CALL lbc_lnk( bdyvmask(:,:), 'V', 1. )      ! Lateral boundary cond. 
    11091178 
    1110  
    1111          ! Mask corrections 
    1112          ! ---------------- 
    1113          DO ik = 1, jpkm1 
    1114             DO ij = 1, jpj 
    1115                DO ii = 1, jpi 
    1116                   tmask(ii,ij,ik) = tmask(ii,ij,ik) * bdytmask(ii,ij) 
    1117                   umask(ii,ij,ik) = umask(ii,ij,ik) * bdyumask(ii,ij) 
    1118                   vmask(ii,ij,ik) = vmask(ii,ij,ik) * bdyvmask(ii,ij) 
    1119                END DO       
    1120             END DO 
    1121             DO ij = 2, jpjm1 
    1122                DO ii = 2, jpim1 
    1123                   fmask(ii,ij,ik) = fmask(ii,ij,ik) * bdytmask(ii,ij  ) * bdytmask(ii+1,ij  )   & 
    1124                      &                              * bdytmask(ii,ij+1) * bdytmask(ii+1,ij+1) 
    1125                END DO       
    1126             END DO 
    1127          END DO 
    1128          tmask_i (:,:) = ssmask(:,:) * tmask_i(:,:) 
    1129          ! 
    11301179      ENDIF ! ln_mask_file=.TRUE. 
    11311180       
    1132       bdytmask(:,:) = ssmask(:,:) 
    11331181      IF( .NOT.ln_mask_file ) THEN 
    11341182         ! If .not. ln_mask_file then we need to derive mask on U and V grid from mask on T grid here. 
    … …  
    13001348      CALL wrk_dealloc(jpi,jpj,   zfmask )  
    13011349      ! 
    1302       IF( nn_timing == 1 )   CALL timing_stop('bdy_init') 
    1303       ! 
    1304    END SUBROUTINE bdy_init 
    1305  
     1350      IF( nn_timing == 1 )   CALL timing_stop('bdy_segs') 
     1351      ! 
     1352   END SUBROUTINE bdy_segs 
    13061353 
    13071354   SUBROUTINE bdy_ctl_seg 
    … …  
    17131760   END SUBROUTINE bdy_ctl_corn 
    17141761 
    1715 #else 
    1716    !!--------------------------------------------------------------------------------- 
    1717    !!   Dummy module                                   NO open boundaries 
    1718    !!--------------------------------------------------------------------------------- 
    1719 CONTAINS 
    1720    SUBROUTINE bdy_init      ! Dummy routine 
    1721    END SUBROUTINE bdy_init 
    1722 #endif 
    1723  
    17241762   !!================================================================================= 
    17251763END MODULE bdyini 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdylib.F90

    r6140 r7412  
    55   !!====================================================================== 
    66   !! History :  3.6  !  2013     (D. Storkey) original code 
     7   !!            4.0  !  2014     (T. Lovato) Generalize OBC structure 
    78   !!---------------------------------------------------------------------- 
    8 #if defined key_bdy  
    9    !!---------------------------------------------------------------------- 
    10    !!   'key_bdy' :                    Unstructured Open Boundary Condition 
    119   !!---------------------------------------------------------------------- 
    1210   !!   bdy_orlanski_2d 
    … …  
    2523   PRIVATE 
    2624 
    27    PUBLIC   bdy_orlanski_2d     ! routine called where? 
    28    PUBLIC   bdy_orlanski_3d     ! routine called where? 
     25   PUBLIC   bdy_frs, bdy_spe, bdy_nmn, bdy_orl 
     26   PUBLIC   bdy_orlanski_2d 
     27   PUBLIC   bdy_orlanski_3d 
    2928 
    3029   !!---------------------------------------------------------------------- 
    31    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     30   !! NEMO/OPA 4.0 , NEMO Consortium (2016) 
    3231   !! $Id$  
    3332   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    3433   !!---------------------------------------------------------------------- 
    3534CONTAINS 
     35 
     36   SUBROUTINE bdy_frs( idx, pta, dta ) 
     37      !!---------------------------------------------------------------------- 
     38      !!                 ***  SUBROUTINE bdy_frs  *** 
     39      !! 
     40      !! ** Purpose : Apply the Flow Relaxation Scheme for tracers at open boundaries. 
     41      !! 
     42      !! Reference : Engedahl H., 1995, Tellus, 365-382. 
     43      !!---------------------------------------------------------------------- 
     44      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices 
     45      REAL(wp), DIMENSION(:,:),            INTENT(in) ::   dta  ! OBC external data 
     46      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
     47      !! 
     48      REAL(wp) ::   zwgt           ! boundary weight 
     49      INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
     50      INTEGER  ::   ii, ij         ! 2D addresses 
     51      !!---------------------------------------------------------------------- 
     52      ! 
     53      IF( nn_timing == 1 ) CALL timing_start('bdy_frs') 
     54      !  
     55      igrd = 1                       ! Everything is at T-points here 
     56      DO ib = 1, idx%nblen(igrd) 
     57         DO ik = 1, jpkm1 
     58            ii = idx%nbi(ib,igrd)  
     59            ij = idx%nbj(ib,igrd) 
     60            zwgt = idx%nbw(ib,igrd) 
     61            pta(ii,ij,ik) = ( pta(ii,ij,ik) + zwgt * (dta(ib,ik) - pta(ii,ij,ik) ) ) * tmask(ii,ij,ik) 
     62         END DO 
     63      END DO 
     64      ! 
     65      IF( nn_timing == 1 ) CALL timing_stop('bdy_frs') 
     66      ! 
     67   END SUBROUTINE bdy_frs 
     68 
     69   SUBROUTINE bdy_spe( idx, pta, dta ) 
     70      !!---------------------------------------------------------------------- 
     71      !!                 ***  SUBROUTINE bdy_spe  *** 
     72      !! 
     73      !! ** Purpose : Apply a specified value for tracers at open boundaries. 
     74      !! 
     75      !!---------------------------------------------------------------------- 
     76      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices 
     77      REAL(wp), DIMENSION(:,:),            INTENT(in) ::   dta  ! OBC external data 
     78      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
     79      !! 
     80      REAL(wp) ::   zwgt           ! boundary weight 
     81      INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
     82      INTEGER  ::   ii, ij         ! 2D addresses 
     83      !!---------------------------------------------------------------------- 
     84      ! 
     85      IF( nn_timing == 1 ) CALL timing_start('bdy_spe') 
     86      ! 
     87      igrd = 1                       ! Everything is at T-points here 
     88      DO ib = 1, idx%nblenrim(igrd) 
     89         ii = idx%nbi(ib,igrd) 
     90         ij = idx%nbj(ib,igrd) 
     91         DO ik = 1, jpkm1 
     92            pta(ii,ij,ik) = dta(ib,ik) * tmask(ii,ij,ik) 
     93         END DO 
     94      END DO 
     95      ! 
     96      IF( nn_timing == 1 ) CALL timing_stop('bdy_spe') 
     97      ! 
     98   END SUBROUTINE bdy_spe 
     99 
     100   SUBROUTINE bdy_orl( idx, ptb, pta, dta, ll_npo ) 
     101      !!---------------------------------------------------------------------- 
     102      !!                 ***  SUBROUTINE bdy_orl  *** 
     103      !! 
     104      !! ** Purpose : Apply Orlanski radiation for tracers at open boundaries. 
     105      !!              This is a wrapper routine for bdy_orlanski_3d below 
     106      !! 
     107      !!---------------------------------------------------------------------- 
     108      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices 
     109      REAL(wp), DIMENSION(:,:),            INTENT(in) ::   dta  ! OBC external data 
     110      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptb  ! before tracer field 
     111      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
     112      LOGICAL,                             INTENT(in) ::   ll_npo  ! switch for NPO version 
     113      !! 
     114      INTEGER  ::   igrd                                    ! grid index 
     115      !!---------------------------------------------------------------------- 
     116      ! 
     117      IF( nn_timing == 1 ) CALL timing_start('bdy_orl') 
     118      ! 
     119      igrd = 1                       ! Everything is at T-points here 
     120      ! 
     121      CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, ll_npo ) 
     122      ! 
     123      IF( nn_timing == 1 ) CALL timing_stop('bdy_orl') 
     124      ! 
     125   END SUBROUTINE bdy_orl 
    36126 
    37127   SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, ll_npo ) 
    … …  
    355445   END SUBROUTINE bdy_orlanski_3d 
    356446 
    357  
    358 #else 
    359    !!---------------------------------------------------------------------- 
    360    !!   Dummy module                   NO Unstruct Open Boundary Conditions 
    361    !!---------------------------------------------------------------------- 
    362 CONTAINS 
    363    SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext  )      ! Empty routine 
    364       WRITE(*,*) 'bdy_orlanski_2d: You should not have seen this print! error?', kt 
    365    END SUBROUTINE bdy_orlanski_2d 
    366    SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext  )      ! Empty routine 
    367       WRITE(*,*) 'bdy_orlanski_3d: You should not have seen this print! error?', kt 
    368    END SUBROUTINE bdy_orlanski_3d 
    369 #endif 
     447   SUBROUTINE bdy_nmn( idx, igrd, phia ) 
     448      !!---------------------------------------------------------------------- 
     449      !!                 ***  SUBROUTINE bdy_nmn  *** 
     450      !!                     
     451      !! ** Purpose : Duplicate the value at open boundaries, zero gradient. 
     452      !!  
     453      !!---------------------------------------------------------------------- 
     454      INTEGER,                    INTENT(in)     ::   igrd     ! grid index 
     455      REAL(wp), DIMENSION(:,:,:), INTENT(inout)  ::   phia     ! model after 3D field (to be updated) 
     456      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices 
     457      !!  
     458      REAL(wp) ::   zcoef, zcoef1, zcoef2 
     459      REAL(wp), POINTER, DIMENSION(:,:,:)        :: pmask      ! land/sea mask for field 
     460      REAL(wp), POINTER, DIMENSION(:,:)        :: bdypmask      ! land/sea mask for field 
     461      INTEGER  ::   ib, ik   ! dummy loop indices 
     462      INTEGER  ::   ii, ij, ip, jp   ! 2D addresses 
     463      !!---------------------------------------------------------------------- 
     464      !! 
     465      IF( nn_timing == 1 ) CALL timing_start('bdy_nmn') 
     466      ! 
     467      SELECT CASE(igrd) 
     468         CASE(1) 
     469            pmask => tmask(:,:,:) 
     470            bdypmask => bdytmask(:,:) 
     471         CASE(2) 
     472            pmask => umask(:,:,:) 
     473            bdypmask => bdyumask(:,:) 
     474         CASE(3) 
     475            pmask => vmask(:,:,:) 
     476            bdypmask => bdyvmask(:,:) 
     477         CASE DEFAULT ;   CALL ctl_stop( 'unrecognised value for igrd in bdy_nmn' ) 
     478      END SELECT 
     479      DO ib = 1, idx%nblenrim(igrd) 
     480         ii = idx%nbi(ib,igrd) 
     481         ij = idx%nbj(ib,igrd) 
     482         DO ik = 1, jpkm1 
     483            ! search the sense of the gradient 
     484            zcoef1 = bdypmask(ii-1,ij  )*pmask(ii-1,ij,ik) +  bdypmask(ii+1,ij  )*pmask(ii+1,ij,ik) 
     485            zcoef2 = bdypmask(ii  ,ij-1)*pmask(ii,ij-1,ik) +  bdypmask(ii  ,ij+1)*pmask(ii,ij+1,ik) 
     486            IF ( nint(zcoef1+zcoef2) == 0) THEN 
     487               ! corner **** we probably only want to set the tangentail component for the dynamics here 
     488               zcoef = pmask(ii-1,ij,ik) + pmask(ii+1,ij,ik) +  pmask(ii,ij-1,ik) +  pmask(ii,ij+1,ik) 
     489               IF (zcoef > .5_wp) THEN ! Only set none isolated points. 
     490                 phia(ii,ij,ik) = phia(ii-1,ij  ,ik) * pmask(ii-1,ij  ,ik) + & 
     491                   &              phia(ii+1,ij  ,ik) * pmask(ii+1,ij  ,ik) + & 
     492                   &              phia(ii  ,ij-1,ik) * pmask(ii  ,ij-1,ik) + & 
     493                   &              phia(ii  ,ij+1,ik) * pmask(ii  ,ij+1,ik) 
     494                 phia(ii,ij,ik) = ( phia(ii,ij,ik) / zcoef ) * pmask(ii,ij,ik) 
     495               ELSE 
     496                 phia(ii,ij,ik) = phia(ii,ij  ,ik) * pmask(ii,ij  ,ik) 
     497               ENDIF 
     498            ELSEIF ( nint(zcoef1+zcoef2) == 2) THEN 
     499               ! oblique corner **** we probably only want to set the normal component for the dynamics here 
     500               zcoef = pmask(ii-1,ij,ik)*bdypmask(ii-1,ij  ) + pmask(ii+1,ij,ik)*bdypmask(ii+1,ij  ) + & 
     501                   &   pmask(ii,ij-1,ik)*bdypmask(ii,ij -1 ) +  pmask(ii,ij+1,ik)*bdypmask(ii,ij+1  ) 
     502               phia(ii,ij,ik) = phia(ii-1,ij  ,ik) * pmask(ii-1,ij  ,ik)*bdypmask(ii-1,ij  ) + & 
     503                   &            phia(ii+1,ij  ,ik) * pmask(ii+1,ij  ,ik)*bdypmask(ii+1,ij  )  + & 
     504                   &            phia(ii  ,ij-1,ik) * pmask(ii  ,ij-1,ik)*bdypmask(ii,ij -1 ) + & 
     505                   &            phia(ii  ,ij+1,ik) * pmask(ii  ,ij+1,ik)*bdypmask(ii,ij+1  ) 
     506  
     507               phia(ii,ij,ik) = ( phia(ii,ij,ik) / MAX(1._wp, zcoef) ) * pmask(ii,ij,ik) 
     508            ELSE 
     509               ip = nint(bdypmask(ii+1,ij  )*pmask(ii+1,ij,ik) - bdypmask(ii-1,ij  )*pmask(ii-1,ij,ik)) 
     510               jp = nint(bdypmask(ii  ,ij+1)*pmask(ii,ij+1,ik) - bdypmask(ii  ,ij-1)*pmask(ii,ij-1,ik)) 
     511               phia(ii,ij,ik) = phia(ii+ip,ij+jp,ik) * pmask(ii+ip,ij+jp,ik) 
     512            ENDIF 
     513         END DO 
     514      END DO 
     515      ! 
     516      IF( nn_timing == 1 ) CALL timing_stop('bdy_nmn') 
     517      ! 
     518   END SUBROUTINE bdy_nmn 
    370519 
    371520   !!====================================================================== 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90

    r6140 r7412  
    1111   !!            3.5  !  2013-07  (J. Chanut) Compliant with time splitting changes 
    1212   !!---------------------------------------------------------------------- 
    13 #if defined key_bdy 
    14    !!---------------------------------------------------------------------- 
    15    !!   'key_bdy'     Open Boundary Condition 
    16    !!---------------------------------------------------------------------- 
    1713   !!   bdytide_init  : read of namelist and initialisation of tidal harmonics data 
    1814   !!   tide_update   : calculation of tidal forcing at each timestep 
    … …  
    2117   USE dom_oce        ! ocean space and time domain 
    2218   USE phycst         ! physical constants 
    23    USE bdy_par        ! Unstructured boundary parameters 
    2419   USE bdy_oce        ! ocean open boundary conditions 
    2520   USE tideini        !  
    … …  
    598593  END SUBROUTINE tide_init_velocities 
    599594 
    600 #else 
    601    !!---------------------------------------------------------------------- 
    602    !!   Dummy module         NO Unstruct Open Boundary Conditions for tides 
    603    !!---------------------------------------------------------------------- 
    604 CONTAINS 
    605    SUBROUTINE bdytide_init             ! Empty routine 
    606       WRITE(*,*) 'bdytide_init: You should not have seen this print! error?' 
    607    END SUBROUTINE bdytide_init 
    608    SUBROUTINE bdytide_update( kt, jit )   ! Empty routine 
    609       WRITE(*,*) 'bdytide_update: You should not have seen this print! error?', kt, jit 
    610    END SUBROUTINE bdytide_update 
    611    SUBROUTINE bdy_dta_tides( kt, kit, time_offset )     ! Empty routine 
    612       INTEGER, INTENT( in )            ::   kt          ! Dummy argument empty routine       
    613       INTEGER, INTENT( in ),OPTIONAL   ::   kit         ! Dummy argument empty routine 
    614       INTEGER, INTENT( in ),OPTIONAL   ::   time_offset ! Dummy argument empty routine 
    615       WRITE(*,*) 'bdy_dta_tides: You should not have seen this print! error?', kt, jit 
    616    END SUBROUTINE bdy_dta_tides 
    617 #endif 
    618  
    619595   !!====================================================================== 
    620596END MODULE bdytides 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdytra.F90

    r6140 r7412  
    88   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge 
    99   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) Optimization of BDY communications 
     10   !!            4.0  !  2016     (T. Lovato) Generalize OBC structure 
    1011   !!---------------------------------------------------------------------- 
    11 #if defined key_bdy 
    12    !!---------------------------------------------------------------------- 
    13    !!   'key_bdy'                     Unstructured Open Boundary Conditions 
    14    !!---------------------------------------------------------------------- 
    15    !!   bdy_tra            : Apply open boundary conditions to T and S 
    16    !!   bdy_tra_frs        : Apply Flow Relaxation Scheme 
     12   !!   bdy_tra       : Apply open boundary conditions & damping to T and S 
    1713   !!---------------------------------------------------------------------- 
    1814   USE oce            ! ocean dynamics and tracers variables 
    … …  
    2016   USE bdy_oce        ! ocean open boundary conditions 
    2117   USE bdylib         ! for orlanski library routines 
    22    USE bdydta   , ONLY:   bf   !  
    2318   ! 
    2419   USE in_out_manager ! I/O manager 
    … …  
    2924   PRIVATE 
    3025 
     26   ! Local structure to rearrange tracers data 
     27   TYPE, PUBLIC ::   ztrabdy 
     28      REAL(wp), POINTER, DIMENSION(:,:) ::  tra 
     29   END TYPE 
     30 
    3131   PUBLIC   bdy_tra      ! called in tranxt.F90  
    3232   PUBLIC   bdy_tra_dmp  ! called in step.F90  
    3333 
    3434   !!---------------------------------------------------------------------- 
    35    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     35   !! NEMO/OPA 4.0, NEMO Consortium (2016) 
    3636   !! $Id$  
    3737   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    … …  
    4848      INTEGER, INTENT(in) ::   kt   ! Main time step counter 
    4949      ! 
    50       INTEGER ::   ib_bdy   ! Loop index 
     50      INTEGER                        :: ib_bdy, jn, igrd   ! Loop indeces 
     51      TYPE(ztrabdy), DIMENSION(jpts) :: zdta               ! Temporary data structure 
    5152      !!---------------------------------------------------------------------- 
     53      igrd = 1  
    5254 
    5355      DO ib_bdy=1, nb_bdy 
    5456         ! 
    55          SELECT CASE( cn_tra(ib_bdy) ) 
    56          CASE('none'        )   ;   CYCLE 
    57          CASE('frs'         )   ;   CALL bdy_tra_frs     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    58          CASE('specified'   )   ;   CALL bdy_tra_spe     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    59          CASE('neumann'     )   ;   CALL bdy_tra_nmn     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    60          CASE('orlanski'    )   ;   CALL bdy_tra_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ll_npo=.false. ) 
    61          CASE('orlanski_npo')   ;   CALL bdy_tra_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ll_npo=.true. ) 
    62          CASE('runoff'      )   ;   CALL bdy_tra_rnf     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    63          CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 
    64          END SELECT 
    65          ! Boundary points should be updated 
    66          CALL lbc_bdy_lnk( tsa(:,:,:,jp_tem), 'T', 1., ib_bdy ) 
    67          CALL lbc_bdy_lnk( tsa(:,:,:,jp_sal), 'T', 1., ib_bdy ) 
     57         zdta(1)%tra => dta_bdy(ib_bdy)%tem 
     58         zdta(2)%tra => dta_bdy(ib_bdy)%sal 
     59         ! 
     60         DO jn = 1, jpts 
     61            ! 
     62            SELECT CASE( TRIM(cn_tra(ib_bdy)) ) 
     63            CASE('none'        )   ;   CYCLE 
     64            CASE('frs'         )   ;   CALL bdy_frs ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
     65            CASE('specified'   )   ;   CALL bdy_spe ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
     66            CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , tsa(:,:,:,jn)               ) 
     67            CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.false. ) 
     68            CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.true. ) 
     69            CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                tsa(:,:,:,jn),               jn ) 
     70            CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 
     71            END SELECT 
     72            ! Boundary points should be updated 
     73            CALL lbc_bdy_lnk( tsa(:,:,:,jn), 'T', 1., ib_bdy ) 
     74            !  
     75         END DO 
    6876      END DO 
    6977      ! 
    7078   END SUBROUTINE bdy_tra 
    7179 
    72  
    73    SUBROUTINE bdy_tra_frs( idx, dta, kt ) 
     80   SUBROUTINE bdy_rnf( idx, pta, jpa ) 
    7481      !!---------------------------------------------------------------------- 
    75       !!                 ***  SUBROUTINE bdy_tra_frs  *** 
     82      !!                 ***  SUBROUTINE bdy_rnf  *** 
    7683      !!                     
    77       !! ** Purpose : Apply the Flow Relaxation Scheme for tracers at open boundaries. 
    78       !!  
    79       !! Reference : Engedahl H., 1995, Tellus, 365-382. 
    80       !!---------------------------------------------------------------------- 
    81       INTEGER,         INTENT(in) ::   kt    ! 
    82       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    83       TYPE(OBC_DATA),  INTENT(in) ::   dta   ! OBC external data 
    84       ! 
    85       REAL(wp) ::   zwgt           ! boundary weight 
    86       INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
    87       INTEGER  ::   ii, ij         ! 2D addresses 
    88       !!---------------------------------------------------------------------- 
    89       ! 
    90       IF( nn_timing == 1 )   CALL timing_start('bdy_tra_frs') 
    91       ! 
    92       igrd = 1                       ! Everything is at T-points here 
    93       DO ib = 1, idx%nblen(igrd) 
    94          DO ik = 1, jpkm1 
    95             ii = idx%nbi(ib,igrd) 
    96             ij = idx%nbj(ib,igrd) 
    97             zwgt = idx%nbw(ib,igrd) 
    98             tsa(ii,ij,ik,jp_tem) = ( tsa(ii,ij,ik,jp_tem) + zwgt * ( dta%tem(ib,ik) - tsa(ii,ij,ik,jp_tem) ) ) * tmask(ii,ij,ik)          
    99             tsa(ii,ij,ik,jp_sal) = ( tsa(ii,ij,ik,jp_sal) + zwgt * ( dta%sal(ib,ik) - tsa(ii,ij,ik,jp_sal) ) ) * tmask(ii,ij,ik) 
    100          END DO 
    101       END DO  
    102       ! 
    103       IF( kt .eq. nit000 )   CLOSE( unit = 102 ) 
    104       ! 
    105       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_frs') 
    106       ! 
    107    END SUBROUTINE bdy_tra_frs 
    108  
    109  
    110    SUBROUTINE bdy_tra_spe( idx, dta, kt ) 
    111       !!---------------------------------------------------------------------- 
    112       !!                 ***  SUBROUTINE bdy_tra_frs  *** 
    113       !!                     
    114       !! ** Purpose : Apply a specified value for tracers at open boundaries. 
     84      !! ** Purpose : Specialized routine to apply TRA runoff values at OBs: 
     85      !!                  - duplicate the neighbour value for the temperature 
     86      !!                  - specified to 0.1 PSU for the salinity 
    11587      !!  
    11688      !!---------------------------------------------------------------------- 
    117       INTEGER,         INTENT(in) ::   kt    ! 
    118       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    119       TYPE(OBC_DATA),  INTENT(in) ::   dta   ! OBC external data 
    120       ! 
    121       REAL(wp) ::   zwgt           ! boundary weight 
    122       INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
    123       INTEGER  ::   ii, ij         ! 2D addresses 
    124       !!---------------------------------------------------------------------- 
    125       ! 
    126       IF( nn_timing == 1 ) CALL timing_start('bdy_tra_spe') 
    127       ! 
    128       igrd = 1                       ! Everything is at T-points here 
    129       DO ib = 1, idx%nblenrim(igrd) 
    130          ii = idx%nbi(ib,igrd) 
    131          ij = idx%nbj(ib,igrd) 
    132          DO ik = 1, jpkm1 
    133             tsa(ii,ij,ik,jp_tem) = dta%tem(ib,ik) * tmask(ii,ij,ik) 
    134             tsa(ii,ij,ik,jp_sal) = dta%sal(ib,ik) * tmask(ii,ij,ik) 
    135          END DO 
    136       END DO 
    137       ! 
    138       IF( kt == nit000 )   CLOSE( unit = 102 ) 
    139       ! 
    140       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_spe') 
    141       ! 
    142    END SUBROUTINE bdy_tra_spe 
    143  
    144  
    145    SUBROUTINE bdy_tra_nmn( idx, dta, kt ) 
    146       !!---------------------------------------------------------------------- 
    147       !!                 ***  SUBROUTINE bdy_tra_nmn  *** 
    148       !!                     
    149       !! ** Purpose : Duplicate the value for tracers at open boundaries. 
    150       !!  
    151       !!---------------------------------------------------------------------- 
    152       INTEGER,         INTENT(in) ::   kt    !  
    153       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    154       TYPE(OBC_DATA) , INTENT(in) ::   dta   ! OBC external data 
    155       ! 
    156       REAL(wp) ::   zwgt           ! boundary weight 
    157       INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
    158       INTEGER  ::   ii, ij,zcoef, zcoef1,zcoef2, ip, jp   ! 2D addresses 
    159       !!---------------------------------------------------------------------- 
    160       ! 
    161       IF( nn_timing == 1 )   CALL timing_start('bdy_tra_nmn') 
    162       ! 
    163       igrd = 1                       ! Everything is at T-points here 
    164       DO ib = 1, idx%nblenrim(igrd) 
    165          ii = idx%nbi(ib,igrd) 
    166          ij = idx%nbj(ib,igrd) 
    167          DO ik = 1, jpkm1 
    168             ! search the sense of the gradient 
    169             zcoef1 = bdytmask(ii-1,ij  ) +  bdytmask(ii+1,ij  ) 
    170             zcoef2 = bdytmask(ii  ,ij-1) +  bdytmask(ii  ,ij+1) 
    171             IF ( zcoef1+zcoef2 == 0) THEN 
    172                ! corner 
    173                zcoef = tmask(ii-1,ij,ik) + tmask(ii+1,ij,ik) +  tmask(ii,ij-1,ik) +  tmask(ii,ij+1,ik) 
    174                tsa(ii,ij,ik,jp_tem) = tsa(ii-1,ij  ,ik,jp_tem) * tmask(ii-1,ij  ,ik) + & 
    175                  &                    tsa(ii+1,ij  ,ik,jp_tem) * tmask(ii+1,ij  ,ik) + & 
    176                  &                    tsa(ii  ,ij-1,ik,jp_tem) * tmask(ii  ,ij-1,ik) + & 
    177                  &                    tsa(ii  ,ij+1,ik,jp_tem) * tmask(ii  ,ij+1,ik) 
    178                tsa(ii,ij,ik,jp_tem) = ( tsa(ii,ij,ik,jp_tem) / MAX( 1, zcoef) ) * tmask(ii,ij,ik) 
    179                tsa(ii,ij,ik,jp_sal) = tsa(ii-1,ij  ,ik,jp_sal) * tmask(ii-1,ij  ,ik) + & 
    180                  &                    tsa(ii+1,ij  ,ik,jp_sal) * tmask(ii+1,ij  ,ik) + & 
    181                  &                    tsa(ii  ,ij-1,ik,jp_sal) * tmask(ii  ,ij-1,ik) + & 
    182                  &                    tsa(ii  ,ij+1,ik,jp_sal) * tmask(ii  ,ij+1,ik) 
    183                tsa(ii,ij,ik,jp_sal) = ( tsa(ii,ij,ik,jp_sal) / MAX( 1, zcoef) ) * tmask(ii,ij,ik) 
    184             ELSE 
    185                ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  ) 
    186                jp = bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1) 
    187                tsa(ii,ij,ik,jp_tem) = tsa(ii+ip,ij+jp,ik,jp_tem) * tmask(ii+ip,ij+jp,ik) 
    188                tsa(ii,ij,ik,jp_sal) = tsa(ii+ip,ij+jp,ik,jp_sal) * tmask(ii+ip,ij+jp,ik) 
    189             ENDIF 
    190          END DO 
    191       END DO 
    192       ! 
    193       IF( kt == nit000 )   CLOSE( unit = 102 ) 
    194       ! 
    195       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_nmn') 
    196       ! 
    197    END SUBROUTINE bdy_tra_nmn 
    198   
    199  
    200    SUBROUTINE bdy_tra_orlanski( idx, dta, ll_npo ) 
    201       !!---------------------------------------------------------------------- 
    202       !!                 ***  SUBROUTINE bdy_tra_orlanski  *** 
    203       !!              
    204       !!              - Apply Orlanski radiation to temperature and salinity.  
    205       !!              - Wrapper routine for bdy_orlanski_3d 
    206       !!  
    207       !! 
    208       !! References:  Marchesiello, McWilliams and Shchepetkin, Ocean Modelling vol. 3 (2001)     
    209       !!---------------------------------------------------------------------- 
    210       TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices 
    211       TYPE(OBC_DATA) , INTENT(in) ::   dta     ! OBC external data 
    212       LOGICAL        , INTENT(in) ::   ll_npo  ! switch for NPO version 
    213       ! 
    214       INTEGER  ::   igrd                                    ! grid index 
    215       !!---------------------------------------------------------------------- 
    216       ! 
    217       IF( nn_timing == 1 ) CALL timing_start('bdy_tra_orlanski') 
    218       ! 
    219       igrd = 1      ! Orlanski bc on temperature;  
    220       !             
    221       CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_tem), tsa(:,:,:,jp_tem), dta%tem, ll_npo ) 
    222  
    223       igrd = 1      ! Orlanski bc on salinity; 
    224       !   
    225       CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_sal), tsa(:,:,:,jp_sal), dta%sal, ll_npo ) 
    226       ! 
    227       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_orlanski') 
    228       ! 
    229    END SUBROUTINE bdy_tra_orlanski 
    230  
    231  
    232    SUBROUTINE bdy_tra_rnf( idx, dta, kt ) 
    233       !!---------------------------------------------------------------------- 
    234       !!                 ***  SUBROUTINE bdy_tra_rnf  *** 
    235       !!                     
    236       !! ** Purpose : Apply the runoff values for tracers at open boundaries: 
    237       !!                  - specified to 0.1 PSU for the salinity 
    238       !!                  - duplicate the value for the temperature 
    239       !!  
    240       !!---------------------------------------------------------------------- 
    241       INTEGER        , INTENT(in) ::   kt    !  
    242       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    243       TYPE(OBC_DATA) , INTENT(in) ::   dta   ! OBC external data 
     89      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices 
     90      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
     91      INTEGER,                             INTENT(in) ::   jpa  ! TRA index 
    24492      ! 
    24593      REAL(wp) ::   zwgt           ! boundary weight 
    … …  
    24896      !!---------------------------------------------------------------------- 
    24997      ! 
    250       IF( nn_timing == 1 )   CALL timing_start('bdy_tra_rnf') 
     98      IF( nn_timing == 1 )   CALL timing_start('bdy_rnf') 
    25199      ! 
    252100      igrd = 1                       ! Everything is at T-points here 
    … …  
    257105            ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  ) 
    258106            jp = bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1) 
    259             tsa(ii,ij,ik,jp_tem) = tsa(ii+ip,ij+jp,ik,jp_tem) * tmask(ii,ij,ik) 
    260             tsa(ii,ij,ik,jp_sal) =                        0.1 * tmask(ii,ij,ik) 
     107            if (jpa == jp_tem) pta(ii,ij,ik) = pta(ii+ip,ij+jp,ik) * tmask(ii,ij,ik) 
     108            if (jpa == jp_sal) pta(ii,ij,ik) =                 0.1 * tmask(ii,ij,ik) 
    261109         END DO 
    262110      END DO 
    263111      ! 
    264       IF( kt == nit000 )   CLOSE( unit = 102 ) 
     112      IF( nn_timing == 1 )   CALL timing_stop('bdy_rnf') 
    265113      ! 
    266       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_rnf') 
    267       ! 
    268    END SUBROUTINE bdy_tra_rnf 
    269  
     114   END SUBROUTINE bdy_rnf 
    270115 
    271116   SUBROUTINE bdy_tra_dmp( kt ) 
    … …  
    308153   END SUBROUTINE bdy_tra_dmp 
    309154  
    310 #else 
    311    !!---------------------------------------------------------------------- 
    312    !!   Dummy module                   NO Unstruct Open Boundary Conditions 
    313    !!---------------------------------------------------------------------- 
    314 CONTAINS 
    315    SUBROUTINE bdy_tra(kt)      ! Empty routine 
    316       WRITE(*,*) 'bdy_tra: You should not have seen this print! error?', kt 
    317    END SUBROUTINE bdy_tra 
    318  
    319    SUBROUTINE bdy_tra_dmp(kt)      ! Empty routine 
    320       WRITE(*,*) 'bdy_tra_dmp: You should not have seen this print! error?', kt 
    321    END SUBROUTINE bdy_tra_dmp 
    322 #endif 
    323  
    324155   !!====================================================================== 
    325156END MODULE bdytra 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90

    r6140 r7412  
    99   !!            3.0  !  2008-04  (NEMO team)  add in the reference version 
    1010   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge 
    11    !!---------------------------------------------------------------------- 
    12 #if defined key_bdy 
    13    !!---------------------------------------------------------------------- 
    14    !!   'key_bdy'                     unstructured open boundary conditions 
    1511   !!---------------------------------------------------------------------- 
    1612   USE oce            ! ocean dynamics and tracers  
    … …  
    175171   END SUBROUTINE bdy_vol 
    176172 
    177 #else 
    178    !!---------------------------------------------------------------------- 
    179    !!   Dummy module                   NO Unstruct Open Boundary Conditions 
    180    !!---------------------------------------------------------------------- 
    181 CONTAINS 
    182    SUBROUTINE bdy_vol( kt )        ! Empty routine 
    183       WRITE(*,*) 'bdy_vol: You should not have seen this print! error?', kt 
    184    END SUBROUTINE bdy_vol 
    185 #endif 
    186  
    187173   !!====================================================================== 
    188174END MODULE bdyvol 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DIA/diaharm.F90

    r6140 r7412  
    66   !! History :  3.1  !  2007  (O. Le Galloudec, J. Chanut)  Original code 
    77   !!---------------------------------------------------------------------- 
    8 #if defined key_diaharm && defined key_tide 
     8#if defined key_diaharm 
    99   !!---------------------------------------------------------------------- 
    1010   !!   'key_diaharm' 
    11    !!   'key_tide' 
    1211   !!---------------------------------------------------------------------- 
    1312   USE oce             ! ocean dynamics and tracers variables 
    … …  
    1615   USE daymod 
    1716   USE tide_mod 
     17   USE sbctide         ! Tidal forcing or not 
    1818   ! 
    1919   USE in_out_manager  ! I/O units 
    … …  
    8282         WRITE(numout,*) '~~~~~~~ ' 
    8383      ENDIF 
     84      ! 
     85      IF( .NOT. ln_tide )   CALL ctl_stop( 'dia_harm_init : ln_tide must be true for harmonic analysis') 
    8486      ! 
    8587      CALL tide_init_Wave 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r6140 r7412  
    2323   USE trabbc          ! bottom boundary condition  
    2424   USE trabbc          ! bottom boundary condition 
    25    USE bdy_par         ! (for lk_bdy) 
    2625   USE restart         ! ocean restart 
     26   USE bdy_oce   , ONLY: ln_bdy 
    2727   ! 
    2828   USE iom             ! I/O manager 
    … …  
    372372 
    373373      IF( .NOT. ln_diahsb )   RETURN 
    374          !      IF( .NOT. lk_mpp_rep ) & 
    375          !        CALL ctl_stop (' Your global mpp_sum if performed in single precision - 64 bits -', & 
    376          !             &         ' whereas the global sum to be precise must be done in double precision ',& 
    377          !             &         ' please add key_mpp_rep') 
    378374 
    379375      ! ------------------- ! 
    … …  
    399395      surf_tot  = glob_sum( surf(:,:) )                                       ! total ocean surface area 
    400396 
    401       IF( lk_bdy ) CALL ctl_warn( 'dia_hsb does not take open boundary fluxes into account' )          
     397      IF( ln_bdy ) CALL ctl_warn( 'dia_hsb does not take open boundary fluxes into account' )          
    402398      ! 
    403399      ! ---------------------------------- ! 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90

    r5836 r7412  
    220220         ! 
    221221         !                                        ! surface of closed seas  
    222          IF( lk_mpp_rep ) THEN                         ! MPP reproductible calculation 
    223             DO jc = 1, jpncs 
    224                ctmp = CMPLX( 0.e0, 0.e0, wp ) 
    225                DO jj = ncsj1(jc), ncsj2(jc) 
    226                   DO ji = ncsi1(jc), ncsi2(jc) 
    227                      ztmp = e1e2t(ji,jj) * tmask_i(ji,jj) 
    228                      CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
    229                   END DO  
    230                END DO  
    231                IF( lk_mpp )   CALL mpp_sum( ctmp ) 
    232                surf(jc) = REAL(ctmp,wp) 
    233             END DO 
    234          ELSE                                          ! Standard calculation            
    235             DO jc = 1, jpncs 
    236                DO jj = ncsj1(jc), ncsj2(jc) 
    237                   DO ji = ncsi1(jc), ncsi2(jc) 
    238                      surf(jc) = surf(jc) + e1e2t(ji,jj) * tmask_i(ji,jj)      ! surface of closed seas 
    239                   END DO  
     222         DO jc = 1, jpncs 
     223            ctmp = CMPLX( 0.e0, 0.e0, wp ) 
     224            DO jj = ncsj1(jc), ncsj2(jc) 
     225               DO ji = ncsi1(jc), ncsi2(jc) 
     226                  ztmp = e1e2t(ji,jj) * tmask_i(ji,jj) 
     227                  CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
    240228               END DO  
    241229            END DO  
    242             IF( lk_mpp )   CALL mpp_sum ( surf, jpncs )       ! mpp: sum over all the global domain 
    243          ENDIF 
     230            IF( lk_mpp )   CALL mpp_sum( ctmp ) 
     231            surf(jc) = REAL(ctmp,wp) 
     232         END DO 
    244233 
    245234         IF(lwp) WRITE(numout,*)'     Closed sea surfaces' 
    … …  
    257246      !                                                   !  update emp        ! 
    258247      zfwf = 0.e0_wp                                      !--------------------! 
    259       IF( lk_mpp_rep ) THEN                         ! MPP reproductible calculation 
    260          DO jc = 1, jpncs 
    261             ctmp = CMPLX( 0.e0, 0.e0, wp ) 
    262             DO jj = ncsj1(jc), ncsj2(jc) 
    263                DO ji = ncsi1(jc), ncsi2(jc) 
    264                   ztmp = e1e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj) 
    265                   CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
    266                END DO   
    267             END DO  
    268             IF( lk_mpp )   CALL mpp_sum( ctmp ) 
    269             zfwf(jc) = REAL(ctmp,wp) 
    270          END DO 
    271       ELSE                                          ! Standard calculation            
    272          DO jc = 1, jpncs 
    273             DO jj = ncsj1(jc), ncsj2(jc) 
    274                DO ji = ncsi1(jc), ncsi2(jc) 
    275                   zfwf(jc) = zfwf(jc) + e1e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj)  
    276                END DO   
    277             END DO  
    278          END DO 
    279          IF( lk_mpp )   CALL mpp_sum ( zfwf(:) , jpncs )       ! mpp: sum over all the global domain 
    280       ENDIF 
     248      DO jc = 1, jpncs 
     249         ctmp = CMPLX( 0.e0, 0.e0, wp ) 
     250         DO jj = ncsj1(jc), ncsj2(jc) 
     251            DO ji = ncsi1(jc), ncsi2(jc) 
     252               ztmp = e1e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj) 
     253               CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
     254            END DO   
     255         END DO  
     256         IF( lk_mpp )   CALL mpp_sum( ctmp ) 
     257         zfwf(jc) = REAL(ctmp,wp) 
     258      END DO 
    281259 
    282260      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN      ! Black Sea case for ORCA_R2 configuration 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90

    r6140 r7412  
    270270 
    271271   !!---------------------------------------------------------------------- 
    272    !! mpp reproducibility 
    273    !!---------------------------------------------------------------------- 
    274 #if defined key_mpp_rep 
    275    LOGICAL, PUBLIC, PARAMETER ::   lk_mpp_rep = .TRUE.    !: agrif flag 
    276 #else 
    277    LOGICAL, PUBLIC, PARAMETER ::   lk_mpp_rep = .FALSE.   !: agrif flag 
    278 #endif 
    279  
    280    !!---------------------------------------------------------------------- 
    281272   !! agrif domain 
    282273   !!---------------------------------------------------------------------- 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90

    r6140 r7412  
    2424   USE oce             ! ocean dynamics and tracers 
    2525   USE dom_oce         ! ocean space and time domain 
    26    ! 
     26   USE bdy_oce       
    2727   USE in_out_manager  ! I/O manager 
    2828   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2929   USE lib_mpp         ! 
     30   USE iom 
    3031   USE wrk_nemo        ! Memory allocation 
    3132   USE timing          ! Timing 
    … …  
    8889      !!      are defined with the proper value at lateral domain boundaries. 
    8990      !! 
    90       !!      In case of open boundaries (lk_bdy=T): 
     91      !!      In case of open boundaries (ln_bdy=T): 
    9192      !!        - tmask is set to 1 on the points to be computed bay the open 
    9293      !!          boundaries routines. 
    … …  
    102103      INTEGER  ::   iif, iil, ii0, ii1, ii   ! local integers 
    103104      INTEGER  ::   ijf, ijl, ij0, ij1       !   -       - 
    104       INTEGER  ::   ios 
     105      INTEGER  ::   ios, inum 
    105106      INTEGER  ::   isrow                    ! index for ORCA1 starting row 
    106107      INTEGER , POINTER, DIMENSION(:,:) ::  imsk 
    … …  
    108109      !! 
    109110      NAMELIST/namlbc/ rn_shlat, ln_vorlat 
     111      NAMELIST/nambdy/ ln_bdy ,nb_bdy, ln_coords_file, cn_coords_file,         & 
     112         &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,     & 
     113         &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             & 
     114         &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, & 
     115         &             cn_ice_lim, nn_ice_lim_dta,                           & 
     116         &             rn_ice_tem, rn_ice_sal, rn_ice_age,                 & 
     117         &             ln_vol, nn_volctl, nn_rimwidth, nb_jpk_bdy 
    110118      !!--------------------------------------------------------------------- 
    111119      ! 
    … …  
    155163      END DO   
    156164       
     165      REWIND( numnam_ref )              ! Namelist nambdy in reference namelist :Unstructured open boundaries   
     166      READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 903) 
     167903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy in reference namelist', lwp ) 
     168 
     169      REWIND( numnam_cfg )              ! Namelist nambdy in configuration namelist :Unstructured open boundaries 
     170      READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 904 ) 
     171904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy in configuration namelist', lwp ) 
     172      IF(lwm) WRITE ( numond, nambdy ) 
     173 
     174     IF( ln_bdy .AND. ln_mask_file ) THEN ! correct for bdy mask 
     175         CALL iom_open( cn_mask_file, inum ) 
     176         CALL iom_get ( inum, jpdom_data, 'bdy_msk', bdytmask(:,:) ) 
     177         CALL iom_close( inum ) 
     178 
     179         ! Mask corrections 
     180         ! ---------------- 
     181         DO jk = 1, jpkm1 
     182            DO jj = 1, jpj 
     183               DO ji = 1, jpi 
     184                  tmask(ji,jj,jk) = tmask(ji,jj,jk) * bdytmask(ji,jj) 
     185               END DO 
     186            END DO 
     187         END DO 
     188      ENDIF 
     189 
    157190      ! (ISF) define barotropic mask and mask the ice shelf point 
    158191      ssmask(:,:)=tmask(:,:,1) ! at this stage ice shelf is not masked 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90

    r6351 r7412  
    874874            ! 
    875875         ELSE                                   !* Initialize at "rest" 
    876             e3t_b(:,:,:) = e3t_0(:,:,:) 
    877             e3t_n(:,:,:) = e3t_0(:,:,:) 
    878             sshn(:,:) = 0.0_wp 
    879  
    880             IF( ln_wd ) THEN 
     876            ! 
     877            IF( ln_wd .AND. ( cp_cfg == 'wad' ) ) THEN 
     878               ! 
     879               CALL wad_istate                  ! WAD test configuration : start from  
     880                                                ! uniform T-S fields and initial ssh slope 
     881               ! needs to be called here and in istate which is called later. 
     882               ! Adjust vertical metrics 
     883               DO jk = 1, jpk 
     884                  e3t_n(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshn(:,:) ) & 
     885                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
     886                    &            + e3t_0(:,:,jk)                               * (1._wp -tmask(:,:,jk)) 
     887               END DO 
     888               e3t_b(:,:,:) = e3t_n(:,:,:) 
     889               ! 
     890            ELSEIF( ln_wd ) THEN 
     891               ! 
    881892              DO jj = 1, jpj 
    882893                DO ji = 1, jpi 
    883894                  IF( e3t_0(ji,jj,1) <= 0.5_wp * rn_wdmin1 ) THEN 
    884                      e3t_b(ji,jj,:) = 0.5_wp * rn_wdmin1  
    885                      e3t_n(ji,jj,:) = 0.5_wp * rn_wdmin1  
    886                      e3t_a(ji,jj,:) = 0.5_wp * rn_wdmin1  
     895                     e3t_b(ji,jj,:) = 0.5_wp * rn_wdmin1 
     896                     e3t_n(ji,jj,:) = 0.5_wp * rn_wdmin1 
     897                     e3t_a(ji,jj,:) = 0.5_wp * rn_wdmin1 
    887898                     sshb(ji,jj) = rn_wdmin1 - bathy(ji,jj) 
    888899                     sshn(ji,jj) = rn_wdmin1 - bathy(ji,jj) 
    … …  
    891902                ENDDO 
    892903              ENDDO 
     904               ! 
     905            ELSE 
     906               ! 
     907               e3t_b(:,:,:) = e3t_0(:,:,:) 
     908               e3t_n(:,:,:) = e3t_0(:,:,:) 
     909               sshn(:,:) = 0.0_wp 
     910               ! 
    893911            END IF 
    894912 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r6492 r7412  
    421421               IF(lwp) WRITE(numout,*) '         Depth = rn_bathy read in namelist' 
    422422               zdta(:,:) = rn_bathy 
     423! 
     424               IF( cp_cfg == 'wad' ) THEN 
     425                  SELECT CASE ( jp_cfg ) 
     426                     !                                        ! ==================== 
     427                     CASE ( 1 )                               ! WAD 1 configuration 
     428                        !                                     ! ==================== 
     429                        ! 
     430                        IF(lwp) WRITE(numout,*) 
     431                        IF(lwp) WRITE(numout,*) 'zgr_bat : Closed box with EW linear bottom slope' 
     432                        IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     433                        ! 
     434                        zdta = 1.5_wp 
     435                        DO ji = 10, jpidta 
     436                          zi = MIN(FLOAT(ji - 10)/FLOAT(jpidta - 10), 1.0 ) 
     437                          zdta(ji,:) = MAX(rn_bathy*zi, 1.5)  
     438                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     439                        END DO 
     440                        !!DO ji = 1, jpidta  
     441                        !!  zi = 1.0-EXP(-0.045*(ji-25.0)**2) 
     442                        !!  zdta(ji,:) = MAX(rn_bathy*zi, 1.5) 
     443                        !!  IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     444                        !!END DO 
     445                        zdta(1:2,:) = -2._wp 
     446                        zdta(jpidta-1:jpidta,:) = -2._wp 
     447                        zdta(:,1) = -2._wp 
     448                        zdta(:,jpjdta) = -2._wp 
     449                        zdta(:,1:3) = -2._wp 
     450                        zdta(:,jpjdta-2:jpjdta) = -2._wp 
     451                        !                                     ! ==================== 
     452                     CASE ( 2, 3 )                            ! WAD 2 or 3  configuration 
     453                        !                                     ! ==================== 
     454                        ! 
     455                        IF(lwp) WRITE(numout,*) 
     456                        IF(lwp) WRITE(numout,*) 'zgr_bat : Parobolic EW channel' 
     457                        IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     458                        ! 
     459                        DO ji = 1, jpidta 
     460                          zi = MAX(1.0-FLOAT((ji-25)**2)/484.0, 0.0 ) 
     461                          zi = MAX(1.0-FLOAT((ji-25)**2)/484.0, -2.0 ) 
     462                          zdta(ji,:) = MAX(rn_bathy*zi, -20.0) 
     463                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     464                        END DO 
     465                        zdta(1:2,:) = -2._wp 
     466                        zdta(jpidta-1:jpidta,:) = -2._wp 
     467                        zdta(:,1) = -2._wp 
     468                        zdta(:,jpjdta) = -2._wp 
     469                        zdta(:,1:3) = -2._wp 
     470                        zdta(:,jpjdta-2:jpjdta) = -2._wp 
     471                        !                                     ! ==================== 
     472                     CASE ( 4 )                               ! WAD 4 configuration 
     473                        !                                     ! ==================== 
     474                        ! 
     475                        IF(lwp) WRITE(numout,*) 
     476                        IF(lwp) WRITE(numout,*) 'zgr_bat : Parobolic bowl' 
     477                        IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     478                        ! 
     479                        DO ji = 1, jpidta 
     480                          zi = MAX(1.0-FLOAT((ji-25)**2)/484.0, -2.0 ) 
     481                        DO jj = 1, jpjdta 
     482                          zj = MAX(1.0-FLOAT((jj-17)**2)/196.0, -2.0 ) 
     483                          zdta(ji,jj) = MAX(rn_bathy*zi*zj, -2.0) 
     484                        END DO 
     485                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     486                        END DO 
     487                        zdta(1:2,:) = -2._wp 
     488                        zdta(jpidta-1:jpidta,:) = -2._wp 
     489                        zdta(:,1) = -2._wp 
     490                        zdta(:,jpjdta) = -2._wp 
     491                        zdta(:,1:3) = -2._wp 
     492                        zdta(:,jpjdta-2:jpjdta) = -2._wp 
     493                        !                                    ! =========================== 
     494                     CASE ( 5 )                              ! WAD 5 configuration 
     495                        !                                    ! ==================== 
     496                        ! 
     497                        IF(lwp) WRITE(numout,*) 
     498                        IF(lwp) WRITE(numout,*) 'zgr_bat : Double slope with shelf' 
     499                        IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     500                        ! 
     501                        DO ji = 1, jpidta 
     502                          zi = MIN(FLOAT(ji)/FLOAT(jpidta - 5), 1.0 ) 
     503                          zdta(ji,:) = MAX(rn_bathy*zi, 0.5) 
     504                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     505                        END DO 
     506                        DO ji = jpidta,46,-1 
     507                          zdta(ji,:) = 10.0 
     508                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     509                        END DO 
     510                        DO ji = 46,20,-1 
     511                          zi = 7.5/25. 
     512                          zdta(ji,:) = MAX(10. - zi*(47.-ji),2.5) 
     513                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     514                        END DO 
     515                        DO ji = 19,15,-1 
     516                          zdta(ji,:) = 2.5 
     517                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     518                        END DO 
     519                        DO ji = 15,4,-1 
     520                          zi = 2.0/11.0 
     521                          zdta(ji,:) = MAX(2.5 - zi*(16-ji), 0.5) 
     522                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     523                        END DO 
     524                        DO ji = 4,1,-1 
     525                          zdta(ji,:) = 0.5 
     526                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     527                        END DO 
     528                        !                                    ! =========================== 
     529                        zdta(1:2,:) = -4._wp 
     530                        zdta(jpidta-1:jpidta,:) = -4._wp 
     531                        zdta(:,1) = -4._wp 
     532                        zdta(:,jpjdta) = -4._wp 
     533                        zdta(:,1:3) = -4._wp 
     534                        zdta(:,jpjdta-2:jpjdta) = -4._wp 
     535                        !                                    ! =========================== 
     536                     CASE ( 6 )                              ! WAD 6 configuration 
     537                        !                                    ! ==================== 
     538                        ! 
     539                        IF(lwp) WRITE(numout,*) 
     540                        IF(lwp) WRITE(numout,*) 'zgr_bat : Parabolic channel with gaussian ridge' 
     541                        IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     542                        ! 
     543                        DO ji = 1, jpidta 
     544                          zi = MAX(1.0-FLOAT((ji-25)**2)/484.0, -2.0 ) 
     545                          zj = 0.95*MAX(EXP(-1.0*FLOAT((ji-25)**2)/32.0) , 0.0 ) 
     546                          zdta(ji,:) = MAX(rn_bathy*(zi-zj), -2.0) 
     547                          IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zdta(ji,1) 
     548                        END DO 
     549                        zdta(1:2,:) = -4._wp 
     550                        zdta(jpidta-1:jpidta,:) = -4._wp 
     551                        zdta(:,1) = -4._wp 
     552                        zdta(:,jpjdta) = -4._wp 
     553                        zdta(:,1:3) = -4._wp 
     554                        zdta(:,jpjdta-2:jpjdta) = -4._wp 
     555                        !                                    ! =========================== 
     556                     CASE DEFAULT 
     557                        !                                    ! =========================== 
     558                        WRITE(ctmp1,*) 'WAD test with a ', jp_cfg,' option is not coded' 
     559                        ! 
     560                        CALL ctl_stop( ctmp1 ) 
     561                        ! 
     562                  END SELECT 
     563               END IF 
     564               ! 
    423565               IF( ln_sco ) THEN                                   ! s-coordinate (zsc       ): idta()=jpk 
    424566                  idta(:,:) = jpkm1 
    … …  
    21932335      CALL lbc_lnk( e3vw_0, 'V', 1._wp ) 
    21942336      ! 
    2195       IF( .NOT.ln_wd ) THEN 
    2196         WHERE( e3t_0 (:,:,:) == 0._wp )   e3t_0 (:,:,:) = 1._wp 
    2197         WHERE( e3u_0 (:,:,:) == 0._wp )   e3u_0 (:,:,:) = 1._wp 
    2198         WHERE( e3v_0 (:,:,:) == 0._wp )   e3v_0 (:,:,:) = 1._wp 
    2199         WHERE( e3f_0 (:,:,:) == 0._wp )   e3f_0 (:,:,:) = 1._wp 
    2200         WHERE( e3w_0 (:,:,:) == 0._wp )   e3w_0 (:,:,:) = 1._wp 
    2201         WHERE( e3uw_0(:,:,:) == 0._wp )   e3uw_0(:,:,:) = 1._wp 
    2202         WHERE( e3vw_0(:,:,:) == 0._wp )   e3vw_0(:,:,:) = 1._wp 
    2203       END IF 
     2337      WHERE( e3t_0 (:,:,:) == 0._wp )   e3t_0 (:,:,:) = 1._wp 
     2338      WHERE( e3u_0 (:,:,:) == 0._wp )   e3u_0 (:,:,:) = 1._wp 
     2339      WHERE( e3v_0 (:,:,:) == 0._wp )   e3v_0 (:,:,:) = 1._wp 
     2340      WHERE( e3f_0 (:,:,:) == 0._wp )   e3f_0 (:,:,:) = 1._wp 
     2341      WHERE( e3w_0 (:,:,:) == 0._wp )   e3w_0 (:,:,:) = 1._wp 
     2342      WHERE( e3uw_0(:,:,:) == 0._wp )   e3uw_0(:,:,:) = 1._wp 
     2343      WHERE( e3vw_0(:,:,:) == 0._wp )   e3vw_0(:,:,:) = 1._wp 
    22042344 
    22052345#if defined key_agrif 
    … …  
    23032443               DO jk = 1, mbathy(ji,jj) 
    23042444                 ! check coordinate is monotonically increasing 
    2305                  IF (e3w_n(ji,jj,jk) <= 0._wp .OR. e3t_n(ji,jj,jk) <= 0._wp ) THEN 
     2445                 IF (e3w_0(ji,jj,jk) <= 0._wp .OR. e3t_0(ji,jj,jk) <= 0._wp ) THEN 
    23062446                    WRITE(ctmp1,*) 'ERROR zgr_sco :   e3w   or e3t   =< 0  at point (i,j,k)= ', ji, jj, jk 
    23072447                    WRITE(numout,*) 'ERROR zgr_sco :   e3w   or e3t   =< 0  at point (i,j,k)= ', ji, jj, jk 
    2308                     WRITE(numout,*) 'e3w',e3w_n(ji,jj,:) 
    2309                     WRITE(numout,*) 'e3t',e3t_n(ji,jj,:) 
     2448                    WRITE(numout,*) 'e3w',e3w_0(ji,jj,:) 
     2449                    WRITE(numout,*) 'e3t',e3t_0(ji,jj,:) 
    23102450                    CALL ctl_stop( ctmp1 ) 
    23112451                 ENDIF 
    23122452                 ! and check it has never gone negative 
    2313                  IF( gdepw_n(ji,jj,jk) < 0._wp .OR. gdept_n(ji,jj,jk) < 0._wp ) THEN 
     2453                 IF( gdepw_0(ji,jj,jk) < 0._wp .OR. gdept_0(ji,jj,jk) < 0._wp ) THEN 
    23142454                    WRITE(ctmp1,*) 'ERROR zgr_sco :   gdepw or gdept =< 0  at point (i,j,k)= ', ji, jj, jk 
    23152455                    WRITE(numout,*) 'ERROR zgr_sco :   gdepw   or gdept   =< 0  at point (i,j,k)= ', ji, jj, jk 
    2316                     WRITE(numout,*) 'gdepw',gdepw_n(ji,jj,:) 
    2317                     WRITE(numout,*) 'gdept',gdept_n(ji,jj,:) 
     2456                    WRITE(numout,*) 'gdepw',gdepw_0(ji,jj,:) 
     2457                    WRITE(numout,*) 'gdept',gdept_0(ji,jj,:) 
    23182458                    CALL ctl_stop( ctmp1 ) 
    23192459                 ENDIF 
    23202460                 ! and check it never exceeds the total depth 
    2321                  IF( gdepw_n(ji,jj,jk) > hbatt(ji,jj) ) THEN 
     2461                 IF( gdepw_0(ji,jj,jk) > hbatt(ji,jj) ) THEN 
    23222462                    WRITE(ctmp1,*) 'ERROR zgr_sco :   gdepw > hbatt  at point (i,j,k)= ', ji, jj, jk 
    23232463                    WRITE(numout,*) 'ERROR zgr_sco :   gdepw > hbatt  at point (i,j,k)= ', ji, jj, jk 
    2324                     WRITE(numout,*) 'gdepw',gdepw_n(ji,jj,:) 
     2464                    WRITE(numout,*) 'gdepw',gdepw_0(ji,jj,:) 
    23252465                    CALL ctl_stop( ctmp1 ) 
    23262466                 ENDIF 
    … …  
    23292469               DO jk = 1, mbathy(ji,jj)-1 
    23302470                 ! and check it never exceeds the total depth 
    2331                 IF( gdept_n(ji,jj,jk) > hbatt(ji,jj) ) THEN 
     2471                IF( gdept_0(ji,jj,jk) > hbatt(ji,jj) ) THEN 
    23322472                    WRITE(ctmp1,*) 'ERROR zgr_sco :   gdept > hbatt  at point (i,j,k)= ', ji, jj, jk 
    23332473                    WRITE(numout,*) 'ERROR zgr_sco :   gdept > hbatt  at point (i,j,k)= ', ji, jj, jk 
    2334                     WRITE(numout,*) 'gdept',gdept_n(ji,jj,:) 
     2474                    WRITE(numout,*) 'gdept',gdept_0(ji,jj,:) 
    23352475                    CALL ctl_stop( ctmp1 ) 
    23362476                 ENDIF 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r6140 r7412  
    3636   USE domvvl          ! varying vertical mesh 
    3737   USE iscplrst        ! ice sheet coupling 
     38   USE wet_dry         ! wetting and drying (needed for wad_istate) 
    3839   ! 
    3940   USE in_out_manager  ! I/O manager 
    … …  
    105106         ELSEIF( cp_cfg == 'gyre' ) THEN          
    106107            CALL istate_gyre                     ! GYRE  configuration : start from pre-defined T-S fields 
     108         ELSEIF( cp_cfg == 'wad' ) THEN          
     109            CALL wad_istate                      ! WAD test configuration : start from pre-defined T-S fields and initial ssh slope 
    107110         ELSE                                    ! Initial T-S, U-V fields read in files 
    108111            IF ( ln_tsd_init ) THEN              ! read 3D T and S data at nit000 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90

    r6152 r7412  
    432432      INTEGER  ::   ji, jj, jk, jii, jjj                 ! dummy loop indices 
    433433      REAL(wp) ::   zcoef0, zuap, zvap, znad, ztmp       ! temporary scalars 
    434       LOGICAL  ::   ll_tmp1, ll_tmp2, ll_tmp3            ! local logical variables 
     434      LOGICAL  ::   ll_tmp1, ll_tmp2                     ! local logical variables 
    435435      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zhpi, zhpj 
    436436      REAL(wp), POINTER, DIMENSION(:,:)   ::  zcpx, zcpy !W/D pressure filter 
    … …  
    438438      ! 
    439439      CALL wrk_alloc( jpi,jpj,jpk, zhpi, zhpj ) 
    440       IF(ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
     440      IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
    441441      ! 
    442442      IF( kt == nit000 ) THEN 
    … …  
    451451      ENDIF 
    452452      ! 
    453       IF(ln_wd) THEN 
     453      IF( ln_wd ) THEN 
    454454        DO jj = 2, jpjm1 
    455455           DO ji = 2, jpim1  
    456              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj))  
    457              ll_tmp2 = MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji+1,jj) + bathy(ji+1,jj)) > rn_wdmin1 + rn_wdmin2 
    458              ll_tmp3 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj)) + & 
    459                                                        & rn_wdmin1 + rn_wdmin2 
    460  
    461              IF(ll_tmp1.AND.ll_tmp2) THEN 
     456             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     457                  &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) .AND.            & 
     458                  &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji+1,jj) + bathy(ji+1,jj) ) & 
     459                  &                                                         > rn_wdmin1 + rn_wdmin2 
     460             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     461                  &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
     462 
     463             IF(ll_tmp1) THEN 
    462464               zcpx(ji,jj) = 1.0_wp 
    463                wduflt(ji,jj) = 1.0_wp 
    464              ELSE IF(ll_tmp3) THEN 
    465                ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here 
    466                zcpx(ji,jj) = ABS((sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) / & 
    467                            &     (sshn(ji+1,jj) - sshn(ji,jj))) 
    468                wduflt(ji,jj) = 1.0_wp 
     465             ELSE IF(ll_tmp2) THEN 
     466               ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     467               zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) & 
     468                           &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    469469             ELSE 
    470470               zcpx(ji,jj) = 0._wp 
    471                wduflt(ji,jj) = 0.0_wp 
    472471             END IF 
    473472       
    474              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1))  
    475              ll_tmp2 = MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji,jj+1) + bathy(ji,jj+1)) > rn_wdmin1 + rn_wdmin2 
    476              ll_tmp3 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) + & 
    477                                                        & rn_wdmin1 + rn_wdmin2 
    478  
    479              IF(ll_tmp1.AND.ll_tmp2) THEN 
     473             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     474                  &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) .AND.            & 
     475                  &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji,jj+1) + bathy(ji,jj+1) ) & 
     476                  &                                                         > rn_wdmin1 + rn_wdmin2 
     477             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     478                  &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     479 
     480             IF(ll_tmp1) THEN 
    480481               zcpy(ji,jj) = 1.0_wp 
    481                wdvflt(ji,jj) = 1.0_wp 
    482              ELSE IF(ll_tmp3) THEN 
    483                ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here 
    484                zcpy(ji,jj) = ABS((sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) / & 
    485                            &     (sshn(ji,jj+1) - sshn(ji,jj))) 
    486                wdvflt(ji,jj) = 1.0_wp 
     482             ELSE IF(ll_tmp2) THEN 
     483               ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     484               zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) & 
     485                           &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    487486             ELSE 
    488487               zcpy(ji,jj) = 0._wp 
    489                wdvflt(ji,jj) = 0.0_wp 
    490488             END IF 
    491489           END DO 
    492490        END DO 
    493491        CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    494       ENDIF 
    495  
     492      END IF 
    496493 
    497494      ! Surface value 
    … …  
    510507 
    511508 
    512             IF(ln_wd) THEN 
     509            IF( ln_wd ) THEN 
    513510 
    514511              zhpi(ji,jj,1) = zhpi(ji,jj,1) * zcpx(ji,jj) 
    … …  
    541538                  &           * ( gde3w_n(ji  ,jj+1,jk) - gde3w_n(ji,jj,jk) ) * r1_e2v(ji,jj) 
    542539 
    543                IF(ln_wd) THEN 
     540               IF( ln_wd ) THEN 
    544541                 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) * zcpx(ji,jj) 
    545542                 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) * zcpy(ji,jj)  
    … …  
    556553      ! 
    557554      CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj ) 
    558       IF(ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
     555      IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
    559556      ! 
    560557   END SUBROUTINE hpg_sco 
    … …  
    701698      CALL wrk_alloc( jpi, jpj, jpk, drhox, drhoy, drhoz, drhou, drhov, drhow ) 
    702699      CALL wrk_alloc( jpi, jpj, jpk, rho_i, rho_j, rho_k,  zhpi,  zhpj        ) 
    703       IF(ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
    704       ! 
    705       ! 
    706       IF(ln_wd) THEN 
     700      IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
     701      ! 
     702      ! 
     703      IF( ln_wd ) THEN 
    707704        DO jj = 2, jpjm1 
    708705           DO ji = 2, jpim1  
    709              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj)) & 
    710                      & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji+1,jj) + bathy(ji+1,jj)) & 
    711                      &  > rn_wdmin1 + rn_wdmin2 
    712              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj)) +& 
    713                                                        & rn_wdmin1 + rn_wdmin2 
     706             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     707                  &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) .AND.            & 
     708                  &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji+1,jj) + bathy(ji+1,jj) ) & 
     709                  &                                                         > rn_wdmin1 + rn_wdmin2 
     710             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     711                  &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
    714712 
    715713             IF(ll_tmp1) THEN 
    716714               zcpx(ji,jj) = 1.0_wp 
    717715             ELSE IF(ll_tmp2) THEN 
    718                ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here 
    719                zcpx(ji,jj) = ABS((sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) /& 
    720                            &     (sshn(ji+1,jj) - sshn(ji,jj))) 
     716               ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     717               zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) & 
     718                           &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    721719             ELSE 
    722720               zcpx(ji,jj) = 0._wp 
    723721             END IF 
    724722       
    725              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) & 
    726                      & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji,jj+1) + bathy(ji,jj+1)) & 
    727                      &  > rn_wdmin1 + rn_wdmin2 
    728              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) +& 
    729                                                        & rn_wdmin1 + rn_wdmin2 
     723             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     724                  &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) .AND.            & 
     725                  &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji,jj+1) + bathy(ji,jj+1) ) & 
     726                  &                                                         > rn_wdmin1 + rn_wdmin2 
     727             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     728                  &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
    730729 
    731730             IF(ll_tmp1) THEN 
    732731               zcpy(ji,jj) = 1.0_wp 
    733732             ELSE IF(ll_tmp2) THEN 
    734                ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here 
    735                zcpy(ji,jj) = ABS((sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) /& 
    736                            &     (sshn(ji,jj+1) - sshn(ji,jj))) 
     733               ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     734               zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) & 
     735                           &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    737736             ELSE 
    738737               zcpy(ji,jj) = 0._wp 
    … …  
    741740        END DO 
    742741        CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    743       ENDIF 
    744  
     742      END IF 
    745743 
    746744      IF( kt == nit000 ) THEN 
    … …  
    913911            zhpi(ji,jj,1) = ( rho_k(ji+1,jj  ,1) - rho_k(ji,jj,1) - rho_i(ji,jj,1) ) * r1_e1u(ji,jj) 
    914912            zhpj(ji,jj,1) = ( rho_k(ji  ,jj+1,1) - rho_k(ji,jj,1) - rho_j(ji,jj,1) ) * r1_e2v(ji,jj) 
    915             IF(ln_wd) THEN 
     913            IF( ln_wd ) THEN 
    916914              zhpi(ji,jj,1) = zhpi(ji,jj,1) * zcpx(ji,jj) 
    917915              zhpj(ji,jj,1) = zhpj(ji,jj,1) * zcpy(ji,jj)  
    … …  
    936934                  &           + (  ( rho_k(ji,jj+1,jk) - rho_k(ji,jj,jk  ) )    & 
    937935                  &               -( rho_j(ji,jj  ,jk) - rho_j(ji,jj,jk-1) )  ) * r1_e2v(ji,jj) 
    938                IF(ln_wd) THEN 
     936               IF( ln_wd ) THEN 
    939937                 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) * zcpx(ji,jj) 
    940938                 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) * zcpy(ji,jj)  
    … …  
    950948      CALL wrk_dealloc( jpi, jpj, jpk, drhox, drhoy, drhoz, drhou, drhov, drhow ) 
    951949      CALL wrk_dealloc( jpi, jpj, jpk, rho_i, rho_j, rho_k,  zhpi,  zhpj        ) 
    952       IF(ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
     950      IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
    953951      ! 
    954952   END SUBROUTINE hpg_djc 
    … …  
    987985      CALL wrk_alloc( jpi,jpj,jpk,   zdept, zrhh ) 
    988986      CALL wrk_alloc( jpi,jpj,       zsshu_n, zsshv_n ) 
    989       IF(ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
     987      IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
    990988      ! 
    991989      IF( kt == nit000 ) THEN 
    … …  
    1000998      IF( ln_linssh )   znad = 0._wp 
    1001999 
    1002       IF(ln_wd) THEN 
     1000      IF( ln_wd ) THEN 
    10031001        DO jj = 2, jpjm1 
    10041002           DO ji = 2, jpim1  
    1005              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj)) & 
    1006                      & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji+1,jj) + bathy(ji+1,jj)) & 
    1007                      &  > rn_wdmin1 + rn_wdmin2 
    1008              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj)) +& 
    1009                                                        & rn_wdmin1 + rn_wdmin2 
     1003             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     1004                  &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) .AND.            & 
     1005                  &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji+1,jj) + bathy(ji+1,jj) ) & 
     1006                  &                                                         > rn_wdmin1 + rn_wdmin2 
     1007             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     1008                  &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
    10101009 
    10111010             IF(ll_tmp1) THEN 
    10121011               zcpx(ji,jj) = 1.0_wp 
    10131012             ELSE IF(ll_tmp2) THEN 
    1014                ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here 
    1015                zcpx(ji,jj) = ABS((sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) /& 
    1016                            &     (sshn(ji+1,jj) - sshn(ji,jj))) 
     1013               ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     1014               zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) & 
     1015                           &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    10171016             ELSE 
    10181017               zcpx(ji,jj) = 0._wp 
    10191018             END IF 
    10201019       
    1021              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) & 
    1022                      & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji,jj+1) + bathy(ji,jj+1)) & 
    1023                      &  > rn_wdmin1 + rn_wdmin2 
    1024              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1)) +& 
    1025                                                        & rn_wdmin1 + rn_wdmin2 
    1026  
    1027              IF(ll_tmp1.OR.ll_tmp2) THEN 
     1020             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     1021                  &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) .AND.            & 
     1022                  &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji,jj+1) + bathy(ji,jj+1) ) & 
     1023                  &                                                         > rn_wdmin1 + rn_wdmin2 
     1024             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     1025                  &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     1026 
     1027             IF(ll_tmp1) THEN 
    10281028               zcpy(ji,jj) = 1.0_wp 
    10291029             ELSE IF(ll_tmp2) THEN 
    1030                ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here 
    1031                zcpy(ji,jj) = ABS((sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) /& 
    1032                            &     (sshn(ji,jj+1) - sshn(ji,jj))) 
     1030               ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     1031               zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) & 
     1032                           &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    10331033             ELSE 
    10341034               zcpy(ji,jj) = 0._wp 
    … …  
    10371037        END DO 
    10381038        CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    1039       ENDIF 
     1039      END IF 
    10401040 
    10411041      ! Clean 3-D work arrays 
    … …  
    12211221                 zdpdx2 = zcoef0 * r1_e1u(ji,jj) * REAL(jis-jid, wp) * (zpwes + zpwed) 
    12221222               ENDIF 
    1223                IF(ln_wd) THEN 
     1223               IF( ln_wd ) THEN 
    12241224                  zdpdx1 = zdpdx1 * zcpx(ji,jj) 
    12251225                  zdpdx2 = zdpdx2 * zcpx(ji,jj) 
    … …  
    12801280                  zdpdy2 = zcoef0 * r1_e2v(ji,jj) * REAL(jjs-jjd, wp) * (zpnss + zpnsd ) 
    12811281               ENDIF 
    1282                IF(ln_wd) THEN 
     1282               IF( ln_wd ) THEN 
    12831283                  zdpdy1 = zdpdy1 * zcpy(ji,jj) 
    12841284                  zdpdy2 = zdpdy2 * zcpy(ji,jj) 
    … …  
    12951295      CALL wrk_dealloc( jpi,jpj,jpk,   zdept, zrhh ) 
    12961296      CALL wrk_dealloc( jpi,jpj,       zsshu_n, zsshv_n ) 
    1297       IF(ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
     1297      IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
    12981298      ! 
    12991299   END SUBROUTINE hpg_prj 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynkeg.F90

    r5328 r7412  
    2424   USE wrk_nemo        ! Memory Allocation 
    2525   USE timing          ! Timing 
     26   USE bdy_oce         ! ocean open boundary conditions 
    2627 
    2728   IMPLICIT NONE 
    … …  
    7879      REAL(wp), POINTER, DIMENSION(:,:,:) :: zhke 
    7980      REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdu, ztrdv  
     81      INTEGER  ::   jb                 ! dummy loop indices 
     82      INTEGER  ::   ii, ij, igrd, ib_bdy   ! local integers 
     83      INTEGER  ::   fu, fv 
    8084      !!---------------------------------------------------------------------- 
    8185      ! 
    … …  
    98102      zhke(:,:,jpk) = 0._wp 
    99103       
     104      IF (ln_bdy) THEN 
     105         ! Maria Luneva & Fred Wobus: July-2016 
     106         ! compensate for lack of turbulent kinetic energy on liquid bdy points 
     107         DO ib_bdy = 1, nb_bdy 
     108            IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 
     109               igrd = 2           ! Copying normal velocity into points outside bdy 
     110               DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     111                  DO jk = 1, jpkm1 
     112                     ii   = idx_bdy(ib_bdy)%nbi(jb,igrd) 
     113                     ij   = idx_bdy(ib_bdy)%nbj(jb,igrd) 
     114                     fu   = NINT( idx_bdy(ib_bdy)%flagu(jb,igrd) ) 
     115                     un(ii-fu,ij,jk) = un(ii,ij,jk) * umask(ii,ij,jk) 
     116                  END DO 
     117               END DO 
     118               ! 
     119               igrd = 3           ! Copying normal velocity into points outside bdy 
     120               DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     121                  DO jk = 1, jpkm1 
     122                     ii   = idx_bdy(ib_bdy)%nbi(jb,igrd) 
     123                     ij   = idx_bdy(ib_bdy)%nbj(jb,igrd) 
     124                     fv   = NINT( idx_bdy(ib_bdy)%flagv(jb,igrd) ) 
     125                     vn(ii,ij-fv,jk) = vn(ii,ij,jk) * vmask(ii,ij,jk) 
     126                  END DO 
     127               END DO 
     128            ENDIF 
     129         ENDDO   
     130      ENDIF  
     131 
    100132      SELECT CASE ( kscheme )             !== Horizontal kinetic energy at T-point  ==! 
    101133      ! 
    … …  
    133165         ! 
    134166      END SELECT 
     167 
     168      IF (ln_bdy) THEN 
     169         ! restore velocity masks at points outside boundary 
     170         un(:,:,:) = un(:,:,:) * umask(:,:,:) 
     171         vn(:,:,:) = vn(:,:,:) * vmask(:,:,:) 
     172      ENDIF       
     173 
     174 
    135175      ! 
    136176      DO jk = 1, jpkm1                    !==  grad( KE ) added to the general momentum trends  ==! 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90

    r6140 r7412  
    3232   USE dynspg_ts      ! surface pressure gradient: split-explicit scheme 
    3333   USE domvvl         ! variable volume 
    34    USE bdy_oce        ! ocean open boundary conditions 
     34   USE bdy_oce   , ONLY: ln_bdy 
    3535   USE bdydta         ! ocean open boundary conditions 
    3636   USE bdydyn         ! ocean open boundary conditions 
    … …  
    7777      !!              * Apply lateral boundary conditions on after velocity  
    7878      !!             at the local domain boundaries through lbc_lnk call, 
    79       !!             at the one-way open boundaries (lk_bdy=T), 
     79      !!             at the one-way open boundaries (ln_bdy=T), 
    8080      !!             at the AGRIF zoom   boundaries (lk_agrif=T) 
    8181      !! 
    … …  
    147147      CALL lbc_lnk( va, 'V', -1. )  
    148148      ! 
    149 # if defined key_bdy 
    150149      !                                !* BDY open boundaries 
    151       IF( lk_bdy .AND. ln_dynspg_exp )   CALL bdy_dyn( kt ) 
    152       IF( lk_bdy .AND. ln_dynspg_ts  )   CALL bdy_dyn( kt, dyn3d_only=.true. ) 
     150      IF( ln_bdy .AND. ln_dynspg_exp )   CALL bdy_dyn( kt ) 
     151      IF( ln_bdy .AND. ln_dynspg_ts  )   CALL bdy_dyn( kt, dyn3d_only=.true. ) 
    153152 
    154153!!$   Do we need a call to bdy_vol here?? 
    155       ! 
    156 # endif 
    157154      ! 
    158155      IF( l_trddyn ) THEN             ! prepare the atf trend computation + some diagnostics 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90

    r6981 r7412  
    8888      ! 
    8989      IF(      ln_apr_dyn                                                &   ! atmos. pressure 
    90          .OR.  ( .NOT.ln_dynspg_ts .AND. (ln_tide_pot .AND. lk_tide) )   &   ! tide potential (no time slitting) 
     90         .OR.  ( .NOT.ln_dynspg_ts .AND. (ln_tide_pot .AND. ln_tide) )   &   ! tide potential (no time slitting) 
    9191         .OR.  nn_ice_embd == 2  ) THEN                                      ! embedded sea-ice 
    9292         ! 
    … …  
    111111         ! 
    112112         !                                    !==  tide potential forcing term  ==! 
    113          IF( .NOT.ln_dynspg_ts .AND. ( ln_tide_pot .AND. lk_tide )  ) THEN   ! N.B. added directly at sub-time-step in ts-case 
     113         IF( .NOT.ln_dynspg_ts .AND. ( ln_tide_pot .AND. ln_tide )  ) THEN   ! N.B. added directly at sub-time-step in ts-case 
    114114            ! 
    115115            CALL upd_tide( kt )                      ! update tide potential 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r6152 r7412  
    3333   USE dynvor          ! vorticity term 
    3434   USE wet_dry         ! wetting/drying flux limter 
    35    USE bdy_par         ! for lk_bdy 
     35   USE bdy_oce   , ONLY: ln_bdy 
    3636   USE bdytides        ! open boundary condition data 
    3737   USE bdydyn2d        ! open boundary conditions on barotropic variables 
    … …  
    156156      REAL(wp), POINTER, DIMENSION(:,:) :: zhf 
    157157      REAL(wp), POINTER, DIMENSION(:,:) :: zcpx, zcpy                 ! Wetting/Dying gravity filter coef. 
    158       REAL(wp), POINTER, DIMENSION(:,:) :: wduflt1, wdvflt1           ! Wetting/Dying velocity filter coef. 
    159158      !!---------------------------------------------------------------------- 
    160159      ! 
    … …  
    168167      CALL wrk_alloc( jpi,jpj,   zsshu_a, zsshv_a                  ) 
    169168      CALL wrk_alloc( jpi,jpj,   zhf ) 
    170       IF( ln_wd ) CALL wrk_alloc( jpi, jpj, zcpx, zcpy, wduflt1, wdvflt1 ) 
     169      IF( ln_wd ) CALL wrk_alloc( jpi, jpj, zcpx, zcpy ) 
    171170      ! 
    172171      zmdi=1.e+20                               !  missing data indicator for masking 
    … …  
    374373      IF( .NOT.ln_linssh ) THEN                 ! Variable volume : remove surface pressure gradient 
    375374        IF( ln_wd ) THEN                        ! Calculating and applying W/D gravity filters 
    376           wduflt1(:,:) = 1.0_wp 
    377           wdvflt1(:,:) = 1.0_wp 
    378           DO jj = 2, jpjm1 
    379              DO ji = 2, jpim1 
    380                 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj))  & 
    381                         & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji+1,jj) + bathy(ji+1,jj))   & 
    382                         &  > rn_wdmin1 + rn_wdmin2 
    383                 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-bathy(ji,jj), -bathy(ji+1,jj))   & 
    384                         &                                   + rn_wdmin1 + rn_wdmin2 
     375           DO jj = 2, jpjm1 
     376              DO ji = 2, jpim1  
     377                ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     378                     &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) .AND.            & 
     379                     &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji+1,jj) + bathy(ji+1,jj) ) & 
     380                     &                                                         > rn_wdmin1 + rn_wdmin2 
     381                ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     382                     &    MAX( -bathy(ji,jj)               , -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
     383    
    385384                IF(ll_tmp1) THEN 
    386                   zcpx(ji,jj)    = 1.0_wp 
    387                 ELSEIF(ll_tmp2) THEN 
    388                    ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen here 
    389                   zcpx(ji,jj) = ABS((sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) & 
    390                         &          /(sshn(ji+1,jj) - sshn(ji,jj))) 
     385                  zcpx(ji,jj) = 1.0_wp 
     386                ELSE IF(ll_tmp2) THEN 
     387                  ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     388                  zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + bathy(ji+1,jj) - sshn(ji,jj) - bathy(ji,jj)) & 
     389                              &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    391390                ELSE 
    392                   zcpx(ji,jj)    = 0._wp 
    393                   wduflt1(ji,jj) = 0.0_wp 
     391                  zcpx(ji,jj) = 0._wp 
    394392                END IF 
    395  
    396                 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1))   & 
    397                         & .and. MAX(sshn(ji,jj) + bathy(ji,jj), sshn(ji,jj+1) + bathy(ji,jj+1))   & 
    398                         &  > rn_wdmin1 + rn_wdmin2 
    399                 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-bathy(ji,jj), -bathy(ji,jj+1))   & 
    400                         &                                   + rn_wdmin1 + rn_wdmin2 
     393          
     394                ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     395                     &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) .AND.            & 
     396                     &    MAX(   sshn(ji,jj) + bathy(ji,jj),   sshn(ji,jj+1) + bathy(ji,jj+1) ) & 
     397                     &                                                         > rn_wdmin1 + rn_wdmin2 
     398                ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     399                     &    MAX( -bathy(ji,jj)               , -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     400    
    401401                IF(ll_tmp1) THEN 
    402                    zcpy(ji,jj)    = 1.0_wp 
    403                 ELSEIF(ll_tmp2) THEN 
    404                    ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen here 
    405                   zcpy(ji,jj) = ABS((sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) & 
    406                         &          /(sshn(ji,jj+1) - sshn(ji,jj))) 
     402                  zcpy(ji,jj) = 1.0_wp 
     403                ELSE IF(ll_tmp2) THEN 
     404                  ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     405                  zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + bathy(ji,jj+1) - sshn(ji,jj) - bathy(ji,jj)) & 
     406                              &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    407407                ELSE 
    408                   zcpy(ji,jj)    = 0._wp 
    409                   wdvflt1(ji,jj) = 0.0_wp 
    410                 ENDIF 
    411  
    412              END DO 
     408                  zcpy(ji,jj) = 0._wp 
     409                END IF 
     410              END DO 
    413411           END DO 
    414  
    415            CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    416  
     412  
    417413           DO jj = 2, jpjm1 
    418414              DO ji = 2, jpim1 
    419                  zu_trd(ji,jj) = ( zu_trd(ji,jj) - grav * ( sshn(ji+1,jj  ) - sshn(ji  ,jj ) )   & 
    420                         &                        * r1_e1u(ji,jj) ) * zcpx(ji,jj) * wduflt1(ji,jj) 
    421                  zv_trd(ji,jj) = ( zv_trd(ji,jj) - grav * ( sshn(ji  ,jj+1) - sshn(ji  ,jj ) )   & 
    422                         &                        * r1_e2v(ji,jj) ) * zcpy(ji,jj) * wdvflt1(ji,jj) 
     415                 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj  ) - sshn(ji  ,jj ) )   & 
     416                        &                        * r1_e1u(ji,jj) * zcpx(ji,jj) 
     417                 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji  ,jj+1) - sshn(ji  ,jj ) )   & 
     418                        &                        * r1_e2v(ji,jj) * zcpy(ji,jj) 
    423419              END DO 
    424420           END DO 
    … …  
    567563      ENDIF 
    568564 
    569       IF( ln_wd ) THEN      !preserve the positivity of water depth 
    570                           !ssh[b,n,a] should have already been processed for this 
    571          sshbb_e(:,:) = MAX(sshbb_e(:,:), rn_wdmin1 - bathy(:,:)) 
    572          sshb_e(:,:)  = MAX(sshb_e(:,:) , rn_wdmin1 - bathy(:,:)) 
    573       ENDIF 
    574565      ! 
    575566      IF (ln_bt_fw) THEN                  ! FORWARD integration: start from NOW fields                     
    … …  
    607598         !                                                !  ------------------ 
    608599         ! Update only tidal forcing at open boundaries 
    609 #if defined key_tide 
    610          IF( lk_bdy      .AND. lk_tide )   CALL bdy_dta_tides( kt, kit=jn, time_offset= noffset+1 ) 
    611          IF( ln_tide_pot .AND. lk_tide )   CALL upd_tide     ( kt, kit=jn, time_offset= noffset   ) 
    612 #endif 
     600         IF( ln_bdy      .AND. ln_tide )   CALL bdy_dta_tides( kt, kit=jn, time_offset= noffset+1 ) 
     601         IF( ln_tide_pot .AND. ln_tide )   CALL upd_tide     ( kt, kit=jn, time_offset= noffset   ) 
    613602         ! 
    614603         ! Set extrapolation coefficients for predictor step: 
    … …  
    646635            zhup2_e (:,:) = hu_0(:,:) + zwx(:,:)                ! Ocean depth at U- and V-points 
    647636            zhvp2_e (:,:) = hv_0(:,:) + zwy(:,:) 
    648             IF( ln_wd ) THEN 
    649               zhup2_e(:,:) = MAX(zhup2_e (:,:), rn_wdmin1) 
    650               zhvp2_e(:,:) = MAX(zhvp2_e (:,:), rn_wdmin1) 
    651             END IF 
    652637         ELSE 
    653638            zhup2_e (:,:) = hu_n(:,:) 
    … …  
    701686            END DO 
    702687         END DO 
     688 
    703689         ssha_e(:,:) = (  sshn_e(:,:) - rdtbt * ( zssh_frc(:,:) + zhdiv(:,:) )  ) * ssmask(:,:) 
    704          IF( ln_wd ) ssha_e(:,:) = MAX(ssha_e(:,:), rn_wdmin1 - bathy(:,:))  
     690          
    705691         CALL lbc_lnk( ssha_e, 'T',  1._wp ) 
    706692 
    707 #if defined key_bdy 
    708693         ! Duplicate sea level across open boundaries (this is only cosmetic if linssh=T) 
    709          IF( lk_bdy )   CALL bdy_ssh( ssha_e ) 
    710 #endif 
     694         IF( ln_bdy )   CALL bdy_ssh( ssha_e ) 
     695 
    711696#if defined key_agrif 
    712697         IF( .NOT.Agrif_Root() )   CALL agrif_ssh_ts( jn ) 
    … …  
    749734         zsshp2_e(:,:) = za0 *  ssha_e(:,:) + za1 *  sshn_e (:,:) & 
    750735          &            + za2 *  sshb_e(:,:) + za3 *  sshbb_e(:,:) 
     736 
    751737         IF( ln_wd ) THEN                   ! Calculating and applying W/D gravity filters 
    752            wduflt1(:,:) = 1._wp 
    753            wdvflt1(:,:) = 1._wp 
    754738           DO jj = 2, jpjm1 
    755               DO ji = 2, jpim1 
    756                  ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -bathy(ji,jj), -bathy(ji+1,jj) ) & 
    757                         & .AND. MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji+1,jj) + bathy(ji+1,jj) )    & 
    758                         &                                  > rn_wdmin1 + rn_wdmin2 
    759                  ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -bathy(ji,jj), -bathy(ji+1,jj) ) & 
    760                         &                                  + rn_wdmin1 + rn_wdmin2 
    761                  IF(ll_tmp1) THEN 
    762                     zcpx(ji,jj) = 1._wp 
    763                  ELSE IF(ll_tmp2) THEN 
    764                     ! no worries about zsshp2_e(ji+1,jj)-zsshp2_e(ji,jj) = 0, it won't happen here 
    765                     zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) + bathy(ji+1,jj) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 
    766                         &             / (zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj)) ) 
    767                  ELSE 
    768                     zcpx(ji,jj)    = 0._wp 
    769                     wduflt1(ji,jj) = 0._wp 
    770                  END IF 
    771  
    772                  ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -bathy(ji,jj), -bathy(ji,jj+1) ) & 
    773                         & .AND. MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji,jj+1) + bathy(ji,jj+1) )    & 
    774                         &                                  > rn_wdmin1 + rn_wdmin2 
    775                  ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -bathy(ji,jj), -bathy(ji,jj+1) ) & 
    776                         &                                  + rn_wdmin1 + rn_wdmin2 
    777                  IF(ll_tmp1) THEN 
    778                     zcpy(ji,jj) = 1._wp 
    779                  ELSE IF(ll_tmp2) THEN 
    780                     ! no worries about zsshp2_e(ji,jj+1)-zsshp2_e(ji,jj) = 0, it won't happen here 
    781                     zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) + bathy(ji,jj+1) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 
    782                         &             / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj)) ) 
    783                  ELSE 
    784                     zcpy(ji,jj)    = 0._wp 
    785                     wdvflt1(ji,jj) = 0._wp 
    786                  END IF 
     739              DO ji = 2, jpim1  
     740                ll_tmp1 = MIN( zsshp2_e(ji,jj)               , zsshp2_e(ji+1,jj) ) >                & 
     741                     &    MAX(   -bathy(ji,jj)               ,   -bathy(ji+1,jj) ) .AND.            & 
     742                     &    MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji+1,jj) + bathy(ji+1,jj) ) & 
     743                     &                                                             > rn_wdmin1 + rn_wdmin2 
     744                ll_tmp2 = MAX( zsshp2_e(ji,jj)               , zsshp2_e(ji+1,jj) ) >                & 
     745                     &    MAX(   -bathy(ji,jj)               ,   -bathy(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
     746    
     747                IF(ll_tmp1) THEN 
     748                  zcpx(ji,jj) = 1.0_wp 
     749                ELSE IF(ll_tmp2) THEN 
     750                  ! no worries about  zsshp2_e(ji+1,jj) - zsshp2_e(ji  ,jj) = 0, it won't happen ! here 
     751                  zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) +    bathy(ji+1,jj) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 
     752                              &    / (zsshp2_e(ji+1,jj) - zsshp2_e(ji  ,jj)) ) 
     753                ELSE 
     754                  zcpx(ji,jj) = 0._wp 
     755                END IF 
     756          
     757                ll_tmp1 = MIN( zsshp2_e(ji,jj)               , zsshp2_e(ji,jj+1) ) >                & 
     758                     &    MAX(   -bathy(ji,jj)               ,   -bathy(ji,jj+1) ) .AND.            & 
     759                     &    MAX( zsshp2_e(ji,jj) + bathy(ji,jj), zsshp2_e(ji,jj+1) + bathy(ji,jj+1) ) & 
     760                     &                                                             > rn_wdmin1 + rn_wdmin2 
     761                ll_tmp2 = MAX( zsshp2_e(ji,jj)               , zsshp2_e(ji,jj+1) ) >                & 
     762                     &    MAX(   -bathy(ji,jj)               ,   -bathy(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     763    
     764                IF(ll_tmp1) THEN 
     765                  zcpy(ji,jj) = 1.0_wp 
     766                ELSE IF(ll_tmp2) THEN 
     767                  ! no worries about  zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj  ) = 0, it won't happen ! here 
     768                  zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) +    bathy(ji,jj+1) - zsshp2_e(ji,jj) - bathy(ji,jj)) & 
     769                              &    / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj  )) ) 
     770                ELSE 
     771                  zcpy(ji,jj) = 0._wp 
     772                END IF 
    787773              END DO 
    788             END DO 
    789             CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    790          ENDIF 
     774           END DO 
     775         END IF 
    791776         ! 
    792777         ! Compute associated depths at U and V points: 
    … …  
    806791            END DO 
    807792 
    808             IF( ln_wd ) THEN 
    809               zhust_e(:,:) = MAX(zhust_e (:,:), rn_wdmin1 ) 
    810               zhvst_e(:,:) = MAX(zhvst_e (:,:), rn_wdmin1 ) 
    811             END IF 
    812  
    813793         ENDIF 
    814794         ! 
    … …  
    861841         ! 
    862842         ! Add tidal astronomical forcing if defined 
    863          IF ( lk_tide.AND.ln_tide_pot ) THEN 
     843         IF ( ln_tide.AND.ln_tide_pot ) THEN 
    864844            DO jj = 2, jpjm1 
    865845               DO ji = fs_2, fs_jpim1   ! vector opt. 
    … …  
    888868                 zu_spg = - grav * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) 
    889869                 zv_spg = - grav * ( zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj) ) * r1_e2v(ji,jj) 
    890                  zwx(ji,jj) = zu_spg * zcpx(ji,jj)  
    891                  zwy(ji,jj) = zv_spg * zcpy(ji,jj) 
     870                 zwx(ji,jj) = zu_spg * zcpx(ji,jj) * wdmask(ji,jj) * wdmask(ji+1, jj)  
     871                 zwy(ji,jj) = zv_spg * zcpy(ji,jj) * wdmask(ji,jj) * wdmask(ji, jj+1) 
    892872              END DO 
    893873           END DO 
    … …  
    927907               DO ji = fs_2, fs_jpim1   ! vector opt. 
    928908 
    929                   IF( ln_wd ) THEN 
    930                     zhura = MAX(hu_0(ji,jj) + zsshu_a(ji,jj), rn_wdmin1) 
    931                     zhvra = MAX(hv_0(ji,jj) + zsshv_a(ji,jj), rn_wdmin1) 
    932                   ELSE 
    933                     zhura = hu_0(ji,jj) + zsshu_a(ji,jj) 
    934                     zhvra = hv_0(ji,jj) + zsshv_a(ji,jj) 
    935                   END IF 
     909                  zhura = hu_0(ji,jj) + zsshu_a(ji,jj) 
     910                  zhvra = hv_0(ji,jj) + zsshv_a(ji,jj) 
    936911                  zhura = ssumask(ji,jj)/(zhura + 1._wp - ssumask(ji,jj)) 
    937912                  zhvra = ssvmask(ji,jj)/(zhvra + 1._wp - ssvmask(ji,jj)) 
    … …  
    953928         ! 
    954929         IF( .NOT.ln_linssh ) THEN                     !* Update ocean depth (variable volume case only) 
    955             IF( ln_wd ) THEN 
    956               hu_e (:,:) = MAX(hu_0(:,:) + zsshu_a(:,:), rn_wdmin1) 
    957               hv_e (:,:) = MAX(hv_0(:,:) + zsshv_a(:,:), rn_wdmin1) 
    958             ELSE 
    959               hu_e (:,:) = hu_0(:,:) + zsshu_a(:,:) 
    960               hv_e (:,:) = hv_0(:,:) + zsshv_a(:,:) 
    961             END IF 
     930            hu_e (:,:) = hu_0(:,:) + zsshu_a(:,:) 
     931            hv_e (:,:) = hv_0(:,:) + zsshv_a(:,:) 
    962932            hur_e(:,:) = ssumask(:,:) / ( hu_e(:,:) + 1._wp - ssumask(:,:) ) 
    963933            hvr_e(:,:) = ssvmask(:,:) / ( hv_e(:,:) + 1._wp - ssvmask(:,:) ) 
    … …  
    967937         CALL lbc_lnk_multi( ua_e, 'U', -1._wp, va_e , 'V', -1._wp ) 
    968938         ! 
    969 #if defined key_bdy   
    970939         !                                                 ! open boundaries 
    971          IF( lk_bdy )   CALL bdy_dyn2d( jn, ua_e, va_e, un_e, vn_e, hur_e, hvr_e, ssha_e ) 
    972 #endif 
     940         IF( ln_bdy )   CALL bdy_dyn2d( jn, ua_e, va_e, un_e, vn_e, hur_e, hvr_e, ssha_e ) 
     941 
    973942#if defined key_agrif                                                            
    974943         IF( .NOT.Agrif_Root() )  CALL agrif_dyn_ts( jn )  ! Agrif 
    … …  
    1024993      ! 
    1025994      ! Update barotropic trend: 
    1026       IF( ln_dynadv_vec .OR. ln_linssh ) THEN 
    1027          DO jk=1,jpkm1 
    1028             ua(:,:,jk) = ua(:,:,jk) + ( ua_b(:,:) - ub_b(:,:) ) * z1_2dt_b 
    1029             va(:,:,jk) = va(:,:,jk) + ( va_b(:,:) - vb_b(:,:) ) * z1_2dt_b 
    1030          END DO 
     995      IF(ln_wd) THEN 
     996        IF( ln_dynadv_vec .OR. ln_linssh ) THEN 
     997           DO jk=1,jpkm1 
     998              ua(:,:,jk) = ua(:,:,jk) + ( ua_b(:,:) - ub_b(:,:) ) * z1_2dt_b * wdmask(:,:) 
     999              va(:,:,jk) = va(:,:,jk) + ( va_b(:,:) - vb_b(:,:) ) * z1_2dt_b * wdmask(:,:) 
     1000           END DO 
     1001        ELSE 
     1002           ! At this stage, ssha has been corrected: compute new depths at velocity points 
     1003           DO jj = 1, jpjm1 
     1004              DO ji = 1, jpim1      ! NO Vector Opt. 
     1005                 zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1)  * r1_e1e2u(ji,jj) & 
     1006                    &              * ( e1e2t(ji  ,jj) * ssha(ji  ,jj)    & 
     1007                    &              +   e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 
     1008                 zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1)  * r1_e1e2v(ji,jj) & 
     1009                    &              * ( e1e2t(ji,jj  ) * ssha(ji,jj  )    & 
     1010                    &              +   e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 
     1011              END DO 
     1012           END DO 
     1013           CALL lbc_lnk_multi( zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 
     1014           ! 
     1015           DO jk=1,jpkm1 
     1016              ua(:,:,jk) = ua(:,:,jk) + r1_hu_n(:,:) * ( ua_b(:,:) - ub_b(:,:) * hu_b(:,:) ) * z1_2dt_b * wdmask(:,:) 
     1017              va(:,:,jk) = va(:,:,jk) + r1_hv_n(:,:) * ( va_b(:,:) - vb_b(:,:) * hv_b(:,:) ) * z1_2dt_b * wdmask(:,:) 
     1018           END DO 
     1019           ! Save barotropic velocities not transport: 
     1020           ua_b(:,:) =  ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 
     1021           va_b(:,:) =  va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     1022        ENDIF 
    10311023      ELSE 
    1032          ! At this stage, ssha has been corrected: compute new depths at velocity points 
    1033          DO jj = 1, jpjm1 
    1034             DO ji = 1, jpim1      ! NO Vector Opt. 
    1035                zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1)  * r1_e1e2u(ji,jj) & 
    1036                   &              * ( e1e2t(ji  ,jj) * ssha(ji  ,jj)    & 
    1037                   &              +   e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 
    1038                zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1)  * r1_e1e2v(ji,jj) & 
    1039                   &              * ( e1e2t(ji,jj  ) * ssha(ji,jj  )    & 
    1040                   &              +   e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 
    1041             END DO 
    1042          END DO 
    1043          CALL lbc_lnk_multi( zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 
    1044          ! 
    1045          DO jk=1,jpkm1 
    1046             ua(:,:,jk) = ua(:,:,jk) + r1_hu_n(:,:) * ( ua_b(:,:) - ub_b(:,:) * hu_b(:,:) ) * z1_2dt_b 
    1047             va(:,:,jk) = va(:,:,jk) + r1_hv_n(:,:) * ( va_b(:,:) - vb_b(:,:) * hv_b(:,:) ) * z1_2dt_b 
    1048          END DO 
    1049          ! Save barotropic velocities not transport: 
    1050          ua_b(:,:) =  ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 
    1051          va_b(:,:) =  va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 
    1052       ENDIF 
     1024        IF( ln_dynadv_vec .OR. ln_linssh ) THEN 
     1025           DO jk=1,jpkm1 
     1026              ua(:,:,jk) = ua(:,:,jk) + ( ua_b(:,:) - ub_b(:,:) ) * z1_2dt_b 
     1027              va(:,:,jk) = va(:,:,jk) + ( va_b(:,:) - vb_b(:,:) ) * z1_2dt_b 
     1028           END DO 
     1029        ELSE 
     1030           ! At this stage, ssha has been corrected: compute new depths at velocity points 
     1031           DO jj = 1, jpjm1 
     1032              DO ji = 1, jpim1      ! NO Vector Opt. 
     1033                 zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1)  * r1_e1e2u(ji,jj) & 
     1034                    &              * ( e1e2t(ji  ,jj) * ssha(ji  ,jj)    & 
     1035                    &              +   e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 
     1036                 zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1)  * r1_e1e2v(ji,jj) & 
     1037                    &              * ( e1e2t(ji,jj  ) * ssha(ji,jj  )    & 
     1038                    &              +   e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 
     1039              END DO 
     1040           END DO 
     1041           CALL lbc_lnk_multi( zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 
     1042           ! 
     1043           DO jk=1,jpkm1 
     1044              ua(:,:,jk) = ua(:,:,jk) + r1_hu_n(:,:) * ( ua_b(:,:) - ub_b(:,:) * hu_b(:,:) ) * z1_2dt_b 
     1045              va(:,:,jk) = va(:,:,jk) + r1_hv_n(:,:) * ( va_b(:,:) - vb_b(:,:) * hv_b(:,:) ) * z1_2dt_b 
     1046           END DO 
     1047           ! Save barotropic velocities not transport: 
     1048           ua_b(:,:) =  ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 
     1049           va_b(:,:) =  va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     1050        ENDIF 
     1051 
     1052      END IF 
    10531053      ! 
    10541054      DO jk = 1, jpkm1 
    … …  
    10861086      CALL wrk_dealloc( jpi,jpj,   zsshu_a, zsshv_a                                   ) 
    10871087      CALL wrk_dealloc( jpi,jpj,   zhf ) 
    1088       IF( ln_wd ) CALL wrk_dealloc( jpi, jpj, zcpx, zcpy, wduflt1, wdvflt1 ) 
     1088      IF( ln_wd ) CALL wrk_dealloc( jpi, jpj, zcpx, zcpy ) 
    10891089      ! 
    10901090      IF ( ln_diatmb ) THEN 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

    r6152 r7412  
    2222   USE divhor         ! horizontal divergence 
    2323   USE phycst         ! physical constants 
    24    USE bdy_oce        !  
    25    USE bdy_par        ! 
     24   USE bdy_oce   , ONLY: ln_bdy, bdytmask 
    2625   USE bdydyn2d       ! bdy_ssh routine 
    2726#if defined key_agrif 
    … …  
    8887      ENDIF 
    8988      ! 
    90       CALL div_hor( kt )                              ! Horizontal divergence 
    91       ! 
    92       z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog) 
     89      z2dt = 2._wp * rdt                          ! set time step size (Euler/Leapfrog) 
    9390      IF( neuler == 0 .AND. kt == nit000 )   z2dt = rdt 
     91      zcoef = 0.5_wp * r1_rau0 
    9492 
    9593      !                                           !------------------------------! 
    9694      !                                           !   After Sea Surface Height   ! 
    9795      !                                           !------------------------------! 
     96      IF(ln_wd) THEN 
     97         CALL wad_lmt(sshb, zcoef * (emp_b(:,:) + emp(:,:)), z2dt) 
     98      ENDIF 
     99 
     100      CALL div_hor( kt )                               ! Horizontal divergence 
     101      ! 
    98102      zhdiv(:,:) = 0._wp 
    99103      DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports 
    … …  
    104108      ! compute the vertical velocity which can be used to compute the non-linear terms of the momentum equations. 
    105109      !  
    106       zcoef = 0.5_wp * r1_rau0 
    107  
    108       IF(ln_wd) CALL wad_lmt(sshb, zcoef * (emp_b(:,:) + emp(:,:)), z2dt) 
    109  
    110110      ssha(:,:) = (  sshb(:,:) - z2dt * ( zcoef * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * ssmask(:,:) 
    111111 
    … …  
    116116         CALL agrif_ssh( kt ) 
    117117# endif 
    118 # if defined key_bdy 
    119          IF( lk_bdy ) THEN 
     118         IF( ln_bdy ) THEN 
    120119            CALL lbc_lnk( ssha, 'T', 1. )    ! Not sure that's necessary 
    121120            CALL bdy_ssh( ssha )             ! Duplicate sea level across open boundaries 
    122121         ENDIF 
    123 # endif 
    124122      ENDIF 
    125123 
    … …  
    211209      ENDIF 
    212210 
    213 #if defined key_bdy 
    214       IF( lk_bdy ) THEN 
     211      IF( ln_bdy ) THEN 
    215212         DO jk = 1, jpkm1 
    216213            wn(:,:,jk) = wn(:,:,jk) * bdytmask(:,:) 
    217214         END DO 
    218215      ENDIF 
    219 #endif 
    220216      ! 
    221217      IF( nn_timing == 1 )  CALL timing_stop('wzv') 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/wet_dry.F90

    r6152 r7412  
    3333   !! --------------------------------------------------------------------- 
    3434 
    35    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wduflt, wdvflt !: u- and v- filter 
    3635   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wdmask         !: u- and v- limiter  
    3736 
    … …  
    4645   PUBLIC   wad_lmt                   ! routine called by sshwzv.F90 
    4746   PUBLIC   wad_lmt_bt                ! routine called by dynspg_ts.F90 
     47   PUBLIC   wad_istate                ! routine called by istate.F90 and domvvl.F90 
    4848 
    4949   !! * Substitutions 
    … …  
    8787 
    8888      IF(ln_wd) THEN 
    89          ALLOCATE( wduflt(jpi,jpj), wdvflt(jpi,jpj), wdmask(jpi,jpj), STAT=ierr ) 
     89         ALLOCATE( wdmask(jpi,jpj), STAT=ierr ) 
    9090         IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error') 
    9191      ENDIF 
    … …  
    145145        ! Horizontal Flux in u and v direction 
    146146        DO jk = 1, jpkm1   
    147            DO jj = 1, jpjm1 
    148               DO ji = 1, jpim1 
     147           DO jj = 1, jpj 
     148              DO ji = 1, jpi 
    149149                 zflxu(ji,jj) = zflxu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk) 
    150150                 zflxv(ji,jj) = zflxv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk) 
    … …  
    156156        zflxv(:,:) = zflxv(:,:) * e1v(:,:) 
    157157        
    158         DO jj = 2, jpjm1 
    159            DO ji = 2, jpim1  
     158        wdmask(:,:) = 1 
     159        DO jj = 2, jpj 
     160           DO ji = 2, jpi  
    160161 
    161162             IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE   ! we don't care about land cells 
    … …  
    168169 
    169170              zdep2 = bathy(ji,jj) + sshb1(ji,jj) - rn_wdmin1 
    170               IF(zdep2 < 0._wp) THEN  !add more safty, but not necessary 
     171              IF(zdep2 .le. 0._wp) THEN  !add more safty, but not necessary 
    171172                !zdep2 = 0._wp 
    172173                sshb1(ji,jj) = rn_wdmin1 - bathy(ji,jj) 
     174                wdmask(ji,jj) = 0._wp 
    173175              END IF 
    174176           ENDDO 
    … …  
    183185           zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:) 
    184186           
    185            DO jj = 2, jpjm1 
    186               DO ji = 2, jpim1  
     187           DO jj = 2, jpj 
     188              DO ji = 2, jpi  
    187189         
    188                  wdmask(ji,jj) = 0 
    189190                 IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE  
    190191                 IF(bathy(ji,jj) > zdepwd) CYCLE 
    … …  
    202203                 IF(zdep1 > zdep2) THEN 
    203204                   zflag = 1 
    204                    wdmask(ji, jj) = 1 
     205                   wdmask(ji, jj) = 0 
    205206                   zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
    206207                   zcoef = max(zcoef, 0._wp) 
    … …  
    209210                   IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef 
    210211                   IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef 
    211                    IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji-1,jj) = zcoef 
     212                   IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef 
    212213                 END IF 
    213214              END DO ! ji loop 
    … …  
    231232        CALL lbc_lnk( un, 'U', -1. ) 
    232233        CALL lbc_lnk( vn, 'V', -1. ) 
     234      ! 
     235        un_b(:,:) = un_b(:,:) * zwdlmtu(:, :) 
     236        vn_b(:,:) = vn_b(:,:) * zwdlmtv(:, :) 
     237        CALL lbc_lnk( un_b, 'U', -1. ) 
     238        CALL lbc_lnk( vn_b, 'V', -1. ) 
    233239        
    234240        IF(zflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt!!!' 
    … …  
    291297        zflxp(:,:)   = 0._wp 
    292298        zflxn(:,:)   = 0._wp 
    293         !zflxu(:,:)   = 0._wp 
    294         !zflxv(:,:)   = 0._wp 
    295299 
    296300        zwdlmtu(:,:)  = 1._wp 
    … …  
    299303        ! Horizontal Flux in u and v direction 
    300304        
    301         !zflxu(:,:) = zflxu(:,:) * e2u(:,:) 
    302         !zflxv(:,:) = zflxv(:,:) * e1v(:,:) 
    303         
    304         DO jj = 2, jpjm1 
    305            DO ji = 2, jpim1  
     305        DO jj = 2, jpj 
     306           DO ji = 2, jpi  
    306307 
    307308             IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE   ! we don't care about land cells 
    … …  
    314315 
    315316              zdep2 = bathy(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 
    316               IF(zdep2 < 0._wp) THEN  !add more safty, but not necessary 
    317                 !zdep2 = 0._wp 
    318                sshn_e(ji,jj) = rn_wdmin1 - bathy(ji,jj) 
    319               END IF 
    320317           ENDDO 
    321318        END DO 
    … …  
    329326           zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:) 
    330327           
    331            DO jj = 2, jpjm1 
    332               DO ji = 2, jpim1  
     328           DO jj = 2, jpj 
     329              DO ji = 2, jpi  
    333330         
    334                  wdmask(ji,jj) = 0 
    335331                 IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE  
    336332                 IF(bathy(ji,jj) > zdepwd) CYCLE 
    … …  
    349345                 IF(zdep1 > zdep2) THEN 
    350346                   zflag = 1 
    351                    !wdmask(ji, jj) = 1 
    352347                   zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
    353348                   zcoef = max(zcoef, 0._wp) 
    … …  
    356351                   IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef 
    357352                   IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef 
    358                    IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji-1,jj) = zcoef 
     353                   IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef 
    359354                 END IF 
    360355              END DO ! ji loop 
    … …  
    379374        IF(zflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt_bt!!!' 
    380375        
    381         !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field) 
    382         !IF( nn_cla == 1 )   CALL cla_div    ( kt )             ! Cross Land Advection (update hdivn field) 
    383376        ! 
    384377        ! 
    … …  
    390383      IF( nn_timing == 1 )  CALL timing_stop('wad_lmt') 
    391384   END SUBROUTINE wad_lmt_bt 
     385 
     386   SUBROUTINE wad_istate 
     387      !!---------------------------------------------------------------------- 
     388      !!                   ***  ROUTINE wad_istate  *** 
     389      !!  
     390      !! ** Purpose :   Initialization of the dynamics and tracers for WAD test 
     391      !!      configurations (channels or bowls with initial ssh gradients) 
     392      !! 
     393      !! ** Method  : - set temperature field 
     394      !!              - set salinity field 
     395      !!              - set ssh slope (needs to be repeated in domvvl_rst_init to 
     396      !!                               set vertical metrics ) 
     397      !!---------------------------------------------------------------------- 
     398      ! 
     399      INTEGER  ::   ji, jj            ! dummy loop indices 
     400      REAL(wp) ::   zi, zj 
     401      !!---------------------------------------------------------------------- 
     402      ! 
     403      ! Uniform T & S in all test cases 
     404      tsn(:,:,:,jp_tem) = 10._wp 
     405      tsb(:,:,:,jp_tem) = 10._wp 
     406      tsn(:,:,:,jp_sal) = 35._wp 
     407      tsb(:,:,:,jp_sal) = 35._wp 
     408      SELECT CASE ( jp_cfg )  
     409         !                                        ! ==================== 
     410         CASE ( 1 )                               ! WAD 1 configuration 
     411            !                                     ! ==================== 
     412            ! 
     413            IF(lwp) WRITE(numout,*) 
     414            IF(lwp) WRITE(numout,*) 'istate_wad : Closed box with EW linear bottom slope' 
     415            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     416            ! 
     417            do ji = 1,jpi 
     418             sshn(ji,:) = ( -5.5_wp + 5.5_wp*FLOAT(mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     419            end do 
     420            !                                     ! ==================== 
     421         CASE ( 2 )                               ! WAD 2 configuration 
     422            !                                     ! ==================== 
     423            ! 
     424            IF(lwp) WRITE(numout,*) 
     425            IF(lwp) WRITE(numout,*) 'istate_wad : Parobolic EW channel, mid-range initial ssh slope' 
     426            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     427            ! 
     428            do ji = 1,jpi 
     429             sshn(ji,:) = ( -5.5_wp + 3.9_wp*FLOAT(jpidta - mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     430            end do 
     431            !                                     ! ==================== 
     432         CASE ( 3 )                               ! WAD 3 configuration 
     433            !                                     ! ==================== 
     434            ! 
     435            IF(lwp) WRITE(numout,*) 
     436            IF(lwp) WRITE(numout,*) 'istate_wad : Parobolic EW channel, extreme initial ssh slope'  
     437            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     438            ! 
     439            do ji = 1,jpi 
     440             sshn(ji,:) = ( -7.5_wp + 6.9_wp*FLOAT(jpidta - mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     441            end do 
     442 
     443            ! 
     444            !                                     ! ==================== 
     445         CASE ( 4 )                               ! WAD 4 configuration 
     446            !                                     ! ==================== 
     447            ! 
     448            IF(lwp) WRITE(numout,*) 
     449            IF(lwp) WRITE(numout,*) 'istate_wad : Parobolic bowl, mid-range initial ssh slope'  
     450            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     451            ! 
     452            DO ji = 1, jpi 
     453               zi = MAX(1.0-FLOAT((mig(ji)-25)**2)/400.0, 0.0 ) 
     454               DO jj = 1, jpj 
     455                  zj = MAX(1.0-FLOAT((mjg(jj)-17)**2)/144.0, 0.0 ) 
     456                  sshn(ji,jj) = -8.5_wp + 8.5_wp*zi*zj 
     457               END DO 
     458            END DO 
     459 
     460            ! 
     461            !                                    ! =========================== 
     462         CASE ( 5 )                              ! WAD 5 configuration 
     463            !                                    ! ==================== 
     464            ! 
     465            IF(lwp) WRITE(numout,*) 
     466            IF(lwp) WRITE(numout,*) 'istate_wad : Double slope with shelf' 
     467            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     468            ! 
     469            ! Needed rn_wdmin2 increased to 0.01 for this case? 
     470            do ji = 1,jpi 
     471             sshn(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     472            end do 
     473 
     474            ! 
     475            !                                     ! =========================== 
     476         CASE ( 6 )                               ! WAD 6 configuration 
     477            !                                     ! ==================== 
     478            ! 
     479            IF(lwp) WRITE(numout,*) 
     480            IF(lwp) WRITE(numout,*) 'istate_wad : Parobolic EW channel with gaussian ridge'  
     481            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     482            ! 
     483            do ji = 1,jpi 
     484             !6a 
     485             sshn(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(jpidta - mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     486             !Some variations in initial slope that have been tested 
     487             !6b 
     488             !sshn(ji,:) = ( -5.5_wp + 6.5_wp*FLOAT(jpidta - mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     489             !6c 
     490             !sshn(ji,:) = ( -5.5_wp + 7.5_wp*FLOAT(jpidta - mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     491             !6d 
     492             !sshn(ji,:) = ( -4.5_wp + 8.0_wp*FLOAT(jpidta - mig(ji))/FLOAT(jpidta-1))*tmask(ji,:,1) 
     493            end do 
     494 
     495            ! 
     496            !                                    ! =========================== 
     497         CASE DEFAULT                            ! NONE existing configuration 
     498            !                                    ! =========================== 
     499            WRITE(ctmp1,*) 'WAD test with a ', jp_cfg,' option is not coded' 
     500            ! 
     501            CALL ctl_stop( ctmp1 ) 
     502            ! 
     503      END SELECT 
     504      ! 
     505      ! Apply minimum wetdepth criterion 
     506      ! 
     507      do jj = 1,jpj 
     508         do ji = 1,jpi 
     509            IF( bathy(ji,jj) + sshn(ji,jj) < rn_wdmin1 ) THEN 
     510               sshn(ji,jj) = tmask(ji,jj,1)*( rn_wdmin1 - bathy(ji,jj) ) 
     511            ENDIF 
     512         end do 
     513      end do 
     514      sshb = sshn 
     515      ssha = sshn 
     516      ! 
     517   END SUBROUTINE wad_istate 
     518 
     519   !!===================================================================== 
    392520END MODULE wet_dry 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/IOM/iom.F90

    r6519 r7412  
    99   !!            3.4  ! 2012-12  (R. Bourdalle-Badie and G. Reffray)  add C1D case   
    1010   !!            3.6  ! 2014-15  DIMG format removed 
     11   !!            3.6  ! 2015-15  (J. Harle) Added procedure to read REAL attributes 
    1112   !!-------------------------------------------------------------------- 
    1213 
    … …  
    6768   END INTERFACE 
    6869   INTERFACE iom_getatt 
    69       MODULE PROCEDURE iom_g0d_intatt 
     70      MODULE PROCEDURE iom_g0d_intatt, iom_g0d_ratt 
    7071   END INTERFACE 
    7172   INTERFACE iom_rstput 
    … …  
    979980         IF( iom_file(kiomid)%nfid > 0 ) THEN 
    980981            SELECT CASE (iom_file(kiomid)%iolib) 
    981             CASE (jpnf90   )   ;   CALL iom_nf90_getatt( kiomid, cdatt, pvar ) 
     982            CASE (jpnf90   )   ;   CALL iom_nf90_getatt( kiomid, cdatt, pv_i0d=pvar ) 
    982983            CASE DEFAULT 
    983984               CALL ctl_stop( 'iom_g0d_att: accepted IO library is only jpnf90' ) 
    … …  
    987988   END SUBROUTINE iom_g0d_intatt 
    988989 
     990   SUBROUTINE iom_g0d_ratt( kiomid, cdatt, pvar, cdvar ) 
     991      INTEGER         , INTENT(in   )                 ::   kiomid    ! Identifier of the file 
     992      CHARACTER(len=*), INTENT(in   )                 ::   cdatt     ! Name of the attribute 
     993      REAL(wp)        , INTENT(  out)                 ::   pvar      ! written field 
     994      CHARACTER(len=*), INTENT(in   ), OPTIONAL       ::   cdvar     ! Name of the variable 
     995      ! 
     996      IF( kiomid > 0 ) THEN 
     997         IF( iom_file(kiomid)%nfid > 0 ) THEN 
     998            SELECT CASE (iom_file(kiomid)%iolib) 
     999            CASE (jpnf90   )   ;   IF( PRESENT(cdvar) ) THEN 
     1000                                      CALL iom_nf90_getatt( kiomid, cdatt, pv_r0d=pvar, cdvar=cdvar ) 
     1001                                   ELSE 
     1002                                      CALL iom_nf90_getatt( kiomid, cdatt, pv_r0d=pvar ) 
     1003                                   ENDIF 
     1004            CASE DEFAULT     
     1005               CALL ctl_stop( 'iom_g0d_att: accepted IO library is only jpnf90' ) 
     1006            END SELECT 
     1007         ENDIF 
     1008      ENDIF 
     1009   END SUBROUTINE iom_g0d_ratt 
    9891010 
    9901011   !!---------------------------------------------------------------------- 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/IOM/iom_nf90.F90

    r6140 r7412  
    77   !!            9.0  ! 06 02  (S. Masson) Adaptation to NEMO 
    88   !!             "   ! 07 07  (D. Storkey) Changes to iom_nf90_gettime 
     9   !!            3.6  ! 2015-15  (J. Harle) Added procedure to read REAL attributes 
    910   !!-------------------------------------------------------------------- 
    1011   !!gm  caution add !DIR nec: improved performance to be checked as well as no result changes 
    … …  
    3536   END INTERFACE 
    3637   INTERFACE iom_nf90_getatt 
    37       MODULE PROCEDURE iom_nf90_intatt 
     38      MODULE PROCEDURE iom_nf90_att 
    3839   END INTERFACE 
    3940   INTERFACE iom_nf90_rstput 
    … …  
    313314 
    314315 
    315    SUBROUTINE iom_nf90_intatt( kiomid, cdatt, pvar ) 
    316       !!----------------------------------------------------------------------- 
    317       !!                  ***  ROUTINE  iom_nf90_intatt  *** 
     316   SUBROUTINE iom_nf90_att( kiomid, cdatt, pv_i0d, pv_r0d, cdvar) 
     317      !!----------------------------------------------------------------------- 
     318      !!                  ***  ROUTINE  iom_nf90_att  *** 
    318319      !! 
    319320      !! ** Purpose : read an integer attribute with NF90 
    … …  
    321322      INTEGER         , INTENT(in   ) ::   kiomid   ! Identifier of the file 
    322323      CHARACTER(len=*), INTENT(in   ) ::   cdatt    ! attribute name 
    323       INTEGER         , INTENT(  out) ::   pvar     ! read field 
     324      INTEGER         , INTENT(  out), OPTIONAL       ::   pv_i0d    ! read field 
     325      REAL(wp),         INTENT(  out), OPTIONAL       ::   pv_r0d    ! read field  
     326      CHARACTER(len=*), INTENT(in   ), OPTIONAL       ::   cdvar     ! name of the variable 
    324327      ! 
    325328      INTEGER                         ::   if90id   ! temporary integer 
     329      INTEGER                         ::   ivarid           ! NetCDF  variable Id 
    326330      LOGICAL                         ::   llok     ! temporary logical 
    327331      CHARACTER(LEN=100)              ::   clinfo   ! info character 
    … …  
    329333      !  
    330334      if90id = iom_file(kiomid)%nfid 
    331       llok = NF90_Inquire_attribute(if90id, NF90_GLOBAL, cdatt) == nf90_noerr 
     335      IF( PRESENT(cdvar) ) THEN 
     336         llok = NF90_INQ_VARID( if90id, TRIM(cdvar), ivarid ) == nf90_noerr   ! does the variable exist in the file 
     337         IF( llok ) THEN 
     338            llok = NF90_Inquire_attribute(if90id, ivarid, cdatt) == nf90_noerr 
     339         ELSE 
     340            CALL ctl_warn('iom_nf90_getatt: no variable '//cdvar//' found') 
     341         ENDIF 
     342      ELSE 
     343         llok = NF90_Inquire_attribute(if90id, NF90_GLOBAL, cdatt) == nf90_noerr 
     344      ENDIF  
     345! 
    332346      IF( llok) THEN 
    333347         clinfo = 'iom_nf90_getatt, file: '//TRIM(iom_file(kiomid)%name)//', att: '//TRIM(cdatt) 
    334          CALL iom_nf90_check(NF90_GET_ATT(if90id, NF90_GLOBAL, cdatt, values=pvar), clinfo) 
     348         IF(     PRESENT(pv_r0d) ) THEN 
     349            IF( PRESENT(cdvar) ) THEN 
     350               CALL iom_nf90_check(NF90_GET_ATT(if90id, ivarid, cdatt, values=pv_r0d), clinfo) 
     351            ELSE 
     352               CALL iom_nf90_check(NF90_GET_ATT(if90id, NF90_GLOBAL, cdatt, values=pv_r0d), clinfo) 
     353            ENDIF 
     354         ELSE 
     355            IF( PRESENT(cdvar) ) THEN 
     356               CALL iom_nf90_check(NF90_GET_ATT(if90id, ivarid, cdatt, values=pv_i0d), clinfo) 
     357            ELSE 
     358               CALL iom_nf90_check(NF90_GET_ATT(if90id, NF90_GLOBAL, cdatt, values=pv_i0d), clinfo) 
     359            ENDIF 
     360         ENDIF 
    335361      ELSE 
    336362         CALL ctl_warn('iom_nf90_getatt: no attribute '//cdatt//' found') 
    337          pvar = -999 
     363         IF(     PRESENT(pv_r0d) ) THEN 
     364            pv_r0d = -999._wp 
     365         ELSE 
     366            pv_i0d = -999 
     367         ENDIF 
    338368      ENDIF 
    339369      !  
    340    END SUBROUTINE iom_nf90_intatt 
     370   END SUBROUTINE iom_nf90_att 
    341371 
    342372 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/LBC/mppini_2.h90

    r6412 r7412  
    4141      USE in_out_manager  ! I/O Manager 
    4242      USE iom 
     43      USE bdy_oce 
    4344      !!  
    4445      INTEGER :: ji, jj, jn, jproc, jarea     ! dummy loop indices 
    … …  
    7374      ! read namelist for ln_zco 
    7475      NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco, ln_isfcav, ln_linssh 
    75  
     76      NAMELIST/nambdy/ ln_bdy, nb_bdy, ln_coords_file, cn_coords_file,         & 
     77         &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,     & 
     78         &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             &   
     79         &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, & 
     80         &             cn_ice_lim, nn_ice_lim_dta,                           & 
     81         &             rn_ice_tem, rn_ice_sal, rn_ice_age,                 & 
     82         &             ln_vol, nn_volctl, nn_rimwidth, nb_jpk_bdy 
    7683      !!---------------------------------------------------------------------- 
    7784      !!  OPA 9.0 , LOCEAN-IPSL (2005)  
    … …  
    137144      imask(:,:)=1 
    138145      WHERE ( zdta(:,:) - zdtaisf(:,:) <= rn_isfhmin ) imask = 0 
     146 
     147      ! Adjust imask with bdy_msk if exists 
     148 
     149      REWIND( numnam_ref )              ! Namelist nambdy in reference namelist : BDY 
     150      READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 903) 
     151903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy in reference namelist (mppini_2)', lwp ) 
     152 
     153      REWIND( numnam_cfg )              ! Namelist nambdy in configuration namelist : BDY 
     154      READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 904 ) 
     155904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy in configuration namelist (mppini_2)', lwp ) 
     156      IF(lwm) WRITE ( numond, namzgr ) 
     157 
     158      IF( ln_bdy .AND. ln_mask_file ) THEN 
     159         CALL iom_open( cn_mask_file, inum ) 
     160         CALL iom_get ( inum, jpdom_unknown, 'bdy_msk', zdta(:,:), kstart=(/jpizoom,jpjzoom/), kcount=(/jpiglo,jpjglo/) ) 
     161         CALL iom_close( inum ) 
     162         WHERE ( zdta(:,:) <= 0. ) imask = 0 
     163      ENDIF 
    139164 
    140165      !  1. Dimension arrays for subdomains 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/LDF/ldfdyn.F90

    r6140 r7412  
    4242   REAL(wp), PUBLIC ::   rn_ahm_b        !: lateral laplacian background eddy viscosity [m2/s] 
    4343   REAL(wp), PUBLIC ::   rn_bhm_0        !: lateral bilaplacian eddy viscosity          [m4/s] 
     44                                         !! If nn_ahm_ijk_t = 32 a time and space varying Smagorinsky viscosity 
     45                                         !! will be computed. 
     46   REAL(wp), PUBLIC ::   rn_csmc         !: Smagorinsky constant of proportionality  
     47   REAL(wp), PUBLIC ::   rn_minfac       !: Multiplicative factor of theorectical minimum Smagorinsky viscosity 
     48   REAL(wp), PUBLIC ::   rn_maxfac       !: Multiplicative factor of theorectical maximum Smagorinsky viscosity 
    4449 
    4550   LOGICAL , PUBLIC ::   l_ldfdyn_time   !: flag for time variation of the lateral eddy viscosity coef. 
    4651 
    4752   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ahmt, ahmf   !: eddy diffusivity coef. at U- and V-points   [m2/s or m4/s] 
     53   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   dtensq       !: horizontal tension squared         (Smagorinsky only) 
     54   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   dshesq       !: horizontal shearing strain squared (Smagorinsky only) 
     55   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   esqt, esqf   !: Square of the local gridscale (e1e2/(e1+e2))**2            
    4856 
    4957   REAL(wp) ::   r1_12   = 1._wp / 12._wp    ! =1/12 
    5058   REAL(wp) ::   r1_4    = 0.25_wp           ! =1/4 
     59   REAL(wp) ::   r1_8    = 0.125_wp          ! =1/8 
    5160   REAL(wp) ::   r1_288  = 1._wp / 288._wp   ! =1/( 12^2 * 2 ) 
    5261 
    … …  
    7988      !!                  = 31    = F(i,j,k,t) = F(local velocity) (  |u|e  /12   laplacian operator 
    8089      !!                                                           or |u|e^3/12 bilaplacian operator ) 
    81       !!---------------------------------------------------------------------- 
    82       INTEGER  ::   jk                ! dummy loop indices 
     90      !!                  = 32    = F(i,j,k,t) = F(local deformation rate and gridscale) (D and L) (Smagorinsky)   
     91      !!                                                           (   L^2|D|      laplacian operator 
     92      !!                                                           or  L^4|D|/8  bilaplacian operator ) 
     93      !!---------------------------------------------------------------------- 
     94      INTEGER  ::   ji, jj, jk        ! dummy loop indices 
    8395      INTEGER  ::   ierr, inum, ios   ! local integer 
    8496      REAL(wp) ::   zah0              ! local scalar 
    … …  
    8698      NAMELIST/namdyn_ldf/ ln_dynldf_lap, ln_dynldf_blp,                  & 
    8799         &                 ln_dynldf_lev, ln_dynldf_hor, ln_dynldf_iso,   & 
    88          &                 nn_ahm_ijk_t , rn_ahm_0, rn_ahm_b, rn_bhm_0 
     100         &                 nn_ahm_ijk_t , rn_ahm_0, rn_ahm_b, rn_bhm_0,   & 
     101         &                 rn_csmc      , rn_minfac, rn_maxfac 
    89102      !!---------------------------------------------------------------------- 
    90103      ! 
    … …  
    115128         WRITE(numout,*) '      coefficients :' 
    116129         WRITE(numout,*) '         type of time-space variation         nn_ahm_ijk_t  = ', nn_ahm_ijk_t 
    117          WRITE(numout,*) '         lateral laplacian eddy viscosity     rn_ahm_0_lap  = ', rn_ahm_0, ' m2/s' 
     130         WRITE(numout,*) '         lateral laplacian eddy viscosity     rn_ahm_0      = ', rn_ahm_0, ' m2/s' 
    118131         WRITE(numout,*) '         background viscosity (iso case)      rn_ahm_b      = ', rn_ahm_b, ' m2/s' 
    119          WRITE(numout,*) '         lateral bilaplacian eddy viscosity   rn_ahm_0_blp  = ', rn_bhm_0, ' m4/s' 
     132         WRITE(numout,*) '         lateral bilaplacian eddy viscosity   rn_bhm_0      = ', rn_bhm_0, ' m4/s' 
     133         WRITE(numout,*) '      smagorinsky settings (nn_ahm_ijk_t  = 32) :' 
     134         WRITE(numout,*) '         Smagorinsky coefficient              rn_csmc       = ', rn_csmc 
     135         WRITE(numout,*) '         factor multiplier for theorectical lower limit for ' 
     136         WRITE(numout,*) '           Smagorinsky eddy visc (def. 1.0)   rn_minfac    = ', rn_minfac 
     137         WRITE(numout,*) '         factor multiplier for theorectical lower upper for ' 
     138         WRITE(numout,*) '           Smagorinsky eddy visc (def. 1.0)   rn_maxfac    = ', rn_maxfac 
    120139      ENDIF 
    121140 
    … …  
    208227            l_ldfdyn_time = .TRUE.     ! will be calculated by call to ldf_dyn routine in step.F90 
    209228            ! 
     229         CASE(  32  )       !==  time varying 3D field  ==! 
     230            IF(lwp) WRITE(numout,*) '          momentum mixing coef. = F( latitude, longitude, depth , time )' 
     231            IF(lwp) WRITE(numout,*) '             proportional to the local deformation rate and gridscale (Smagorinsky)' 
     232            IF(lwp) WRITE(numout,*) '                                                             : L^2|D| or L^4|D|/8' 
     233            ! 
     234            l_ldfdyn_time = .TRUE.     ! will be calculated by call to ldf_dyn routine in step.F90 
     235            ! 
     236            ! allocate arrays used in ldf_dyn.  
     237            ALLOCATE( dtensq(jpi,jpj) , dshesq(jpi,jpj) , esqt(jpi,jpj) ,  esqf(jpi,jpj) , STAT=ierr ) 
     238            IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'ldf_dyn_init: failed to allocate Smagorinsky arrays') 
     239            ! 
     240            ! Set local gridscale values 
     241            DO jj = 2, jpjm1 
     242               DO ji = fs_2, fs_jpim1 
     243                  esqt(ji,jj) = ( e1e2t(ji,jj) /( e1t(ji,jj) + e2t(ji,jj) ) )**2  
     244                  esqf(ji,jj) = ( e1e2f(ji,jj) /( e1f(ji,jj) + e2f(ji,jj) ) )**2  
     245               END DO 
     246            END DO 
     247            ! 
    210248         CASE DEFAULT 
    211249            CALL ctl_stop('ldf_dyn_init: wrong choice for nn_ahm_ijk_t, the type of space-time variation of ahm') 
    … …  
    232270      !!    nn_ahm_ijk_t = 31    ahmt, ahmf = F(i,j,k,t) = F(local velocity)  
    233271      !!                         ( |u|e /12  or  |u|e^3/12 for laplacian or bilaplacian operator ) 
    234       !!                BLP case : sqrt of the eddy coef, since bilaplacian is en re-entrant laplacian 
    235272      !! 
    236       !! ** action  :    ahmt, ahmf   update at each time step 
     273      !!    nn_ahm_ijk_t = 32    ahmt, ahmf = F(i,j,k,t) = F(local deformation rate and gridscale) (D and L) (Smagorinsky)   
     274      !!                         ( L^2|D|    or  L^4|D|/8  for laplacian or bilaplacian operator ) 
     275      !! 
     276      !! ** note    :    in BLP cases the sqrt of the eddy coef is returned, since bilaplacian is en re-entrant laplacian 
     277      !! ** action  :    ahmt, ahmf   updated at each time step 
    237278      !!---------------------------------------------------------------------- 
    238279      INTEGER, INTENT(in) ::   kt   ! time step index 
    … …  
    240281      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    241282      REAL(wp) ::   zu2pv2_ij_p1, zu2pv2_ij, zu2pv2_ij_m1, zetmax, zefmax   ! local scalar 
     283      REAL(wp) ::   zcmsmag, zstabf_lo, zstabf_up, zdelta, zdb              ! local scalar 
    242284      !!---------------------------------------------------------------------- 
    243285      ! 
    … …  
    248290      CASE(  31  )       !==  time varying 3D field  ==!   = F( local velocity ) 
    249291         ! 
    250          IF( ln_dynldf_lap   ) THEN          !  laplacian operator : |u| e /12 = |u/144| e 
     292         IF( ln_dynldf_lap   ) THEN        ! laplacian operator : |u| e /12 = |u/144| e 
    251293            DO jk = 1, jpkm1 
    252294               DO jj = 2, jpjm1 
    … …  
    280322         CALL lbc_lnk( ahmt, 'T', 1. )   ;   CALL lbc_lnk( ahmf, 'F', 1. ) 
    281323         ! 
     324         ! 
     325      CASE(  32  )       !==  time varying 3D field  ==!   = F( local deformation rate and gridscale ) (Smagorinsky) 
     326         ! 
     327         IF( ln_dynldf_lap .OR. ln_dynldf_blp  ) THEN        ! laplacian operator : (C_smag/pi)^2 L^2 |D| 
     328            ! 
     329            zcmsmag = (rn_csmc/rpi)**2                                              ! (C_smag/pi)^2 
     330            zstabf_lo  = rn_minfac * rn_minfac / ( 2._wp * 4._wp * zcmsmag )        ! lower limit stability factor scaling 
     331            zstabf_up  = rn_maxfac / ( 4._wp * zcmsmag * 2._wp * rdt )              ! upper limit stability factor scaling 
     332            IF( ln_dynldf_blp ) zstabf_lo = ( 16._wp / 9._wp ) * zstabf_lo          ! provide |U|L^3/12 lower limit instead  
     333            !                                                                       ! of |U|L^3/16 in blp case 
     334            DO jk = 1, jpkm1 
     335               ! 
     336               DO jj = 2, jpj 
     337                  DO ji = 2, jpi 
     338                     zdb = ( (  ub(ji,jj,jk) * r1_e2u(ji,jj) -  ub(ji-1,jj,jk) * r1_e2u(ji-1,jj) )  & 
     339                          &                  * r1_e1t(ji,jj) * e2t(ji,jj)                           & 
     340                          & - ( vb(ji,jj,jk) * r1_e1v(ji,jj) -  vb(ji,jj-1,jk) * r1_e1v(ji,jj-1) )  & 
     341                          &                  * r1_e2t(ji,jj) * e1t(ji,jj)    ) * tmask(ji,jj,jk) 
     342                     dtensq(ji,jj) = zdb*zdb 
     343                  END DO 
     344               END DO 
     345               ! 
     346               DO jj = 1, jpjm1 
     347                  DO ji = 1, jpim1 
     348                     zdb = ( (  ub(ji,jj+1,jk) * r1_e1u(ji,jj+1) -  ub(ji,jj,jk) * r1_e1u(ji,jj) )  & 
     349                          &                    * r1_e2f(ji,jj)   * e1f(ji,jj)                       & 
     350                          & + ( vb(ji+1,jj,jk) * r1_e2v(ji+1,jj) -  vb(ji,jj,jk) * r1_e2v(ji,jj) )  & 
     351                          &                    * r1_e1f(ji,jj)   * e2f(ji,jj)  ) * fmask(ji,jj,jk) 
     352                     dshesq(ji,jj) = zdb*zdb 
     353                  END DO 
     354               END DO 
     355               ! 
     356               DO jj = 2, jpjm1 
     357                  DO ji = fs_2, fs_jpim1 
     358                     ! 
     359                     zu2pv2_ij_p1 = ub(ji  ,jj+1,jk) * ub(ji  ,jj+1,jk) + vb(ji+1,jj  ,jk) * vb(ji+1,jj  ,jk) 
     360                     zu2pv2_ij    = ub(ji  ,jj  ,jk) * ub(ji  ,jj  ,jk) + vb(ji  ,jj  ,jk) * vb(ji  ,jj  ,jk) 
     361                     zu2pv2_ij_m1 = ub(ji-1,jj  ,jk) * ub(ji-1,jj  ,jk) + vb(ji  ,jj-1,jk) * vb(ji  ,jj-1,jk) 
     362                                                     ! T-point value 
     363                     zdelta         = zcmsmag * esqt(ji,jj)                                        ! L^2 * (C_smag/pi)^2 
     364                     ahmt(ji,jj,jk) = zdelta * sqrt(          dtensq(ji,jj)   +                        & 
     365                                     &               r1_4 * ( dshesq(ji,jj)   + dshesq(ji,jj-1)   +    & 
     366                                     &                        dshesq(ji-1,jj) + dshesq(ji-1,jj-1) ) ) 
     367                     ahmt(ji,jj,jk) =   MAX( ahmt(ji,jj,jk),   & 
     368                                     &   SQRT( (zu2pv2_ij + zu2pv2_ij_m1) * zdelta * zstabf_lo ) ) ! Impose lower limit == minfac  * |U|L/2 
     369                     ahmt(ji,jj,jk) = MIN( ahmt(ji,jj,jk), zdelta * zstabf_up )                    ! Impose upper limit == maxfac  * L^2/(4*2dt) 
     370                                                     ! F-point value 
     371                     zdelta         = zcmsmag * esqf(ji,jj)                                        ! L^2 * (C_smag/pi)^2 
     372                     ahmf(ji,jj,jk) = zdelta * sqrt(          dshesq(ji,jj)   +                        & 
     373                                     &               r1_4 * ( dtensq(ji,jj)   + dtensq(ji,jj+1)   +    & 
     374                                     &                        dtensq(ji+1,jj) + dtensq(ji+1,jj+1) ) ) 
     375                     ahmf(ji,jj,jk) =   MAX( ahmf(ji,jj,jk),   & 
     376                                     &   SQRT( (zu2pv2_ij + zu2pv2_ij_p1) * zdelta * zstabf_lo ) ) ! Impose lower limit == minfac  * |U|L/2 
     377                     ahmf(ji,jj,jk) = MIN( ahmf(ji,jj,jk), zdelta * zstabf_up )                    ! Impose upper limit == maxfac  * L^2/(4*2dt) 
     378                     ! 
     379                  END DO 
     380               END DO 
     381            END DO 
     382            ! 
     383         ENDIF 
     384         ! 
     385         IF( ln_dynldf_blp ) THEN      ! bilaplacian operator : sqrt( (C_smag/pi)^2 L^4 |D|/8) 
     386                                       !                      = sqrt( A_lap_smag L^2/8 ) 
     387                                       ! stability limits already applied to laplacian values 
     388                                       ! effective default limits are |U|L^3/12 < B_hm < L^4/(32*2dt) 
     389            ! 
     390            DO jk = 1, jpkm1 
     391               ! 
     392               DO jj = 2, jpjm1 
     393                  DO ji = fs_2, fs_jpim1 
     394                     ! 
     395                     ahmt(ji,jj,jk) = sqrt( r1_8 * esqt(ji,jj) * ahmt(ji,jj,jk) ) 
     396                     ahmf(ji,jj,jk) = sqrt( r1_8 * esqf(ji,jj) * ahmf(ji,jj,jk) ) 
     397                     ! 
     398                  END DO 
     399               END DO 
     400            END DO 
     401            ! 
     402         ENDIF 
     403         ! 
     404         CALL lbc_lnk( ahmt, 'T', 1. )   ;   CALL lbc_lnk( ahmf, 'F', 1. ) 
     405         ! 
    282406      END SELECT 
    283407      ! 

  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90

    r6140 r7412  
    77   !!                 !  05-2008  (S. Alderson)  Modified for Interpolation in memory from input grid to model grid 
    88   !!                 !  10-2013  (D. Delrosso, P. Oddo)  suppression of land point prior to interpolation 
     9   !!                 !  12-2015  (J. Harle) Adding BDY on-the-fly interpolation 
    910   !!---------------------------------------------------------------------- 
    1011 
    … …  
    6768      INTEGER                         ::   nreclast     ! last record to be read in the current file 
    6869      CHARACTER(len = 256)            ::   lsmname      ! current name of the NetCDF mask file acting as a key 
     70      INTEGER                         ::   igrd         ! grid type for bdy data 
     71      INTEGER                         ::   ibdy         ! bdy set id number 
    6972   END TYPE FLD 
    7073 
    … …  
    114117CONTAINS 
    115118 
    116    SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset ) 
     119   SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset, jpk_bdy, fvl ) 
    117120      !!--------------------------------------------------------------------- 
    118121      !!                    ***  ROUTINE fld_read  *** 
    … …  
    135138      !                                                     !   kt_offset = +1 => fields at "after"  time level 
    136139      !                                                     !   etc. 
     140      INTEGER  , INTENT(in   ), OPTIONAL     ::   jpk_bdy   ! number of vertical levels in the BDY data 
     141      LOGICAL  , INTENT(in   ), OPTIONAL     ::   fvl   ! number of vertical levels in the BDY data 
     142      !! 
    137143      INTEGER  ::   itmp         ! local variable 
    138144      INTEGER  ::   imf          ! size of the structure sd 
    … …  
    171177         DO jf = 1, imf  
    172178            IF( PRESENT(map) ) imap = map(jf) 
    173             CALL fld_init( kn_fsbc, sd(jf), imap )  ! read each before field (put them in after as they will be swapped) 
     179               IF( PRESENT(jpk_bdy) ) THEN 
     180                  CALL fld_init( kn_fsbc, sd(jf), imap, jpk_bdy, fvl )  ! read each before field (put them in after as they will be swapped) 
     181               ELSE 
     182                  CALL fld_init( kn_fsbc, sd(jf), imap )  ! read each before field (put them in after as they will be swapped) 
     183               ENDIF 
    174184         END DO 
    175185         IF( lwp ) CALL wgt_print()                ! control print 
    … …  
    263273 
    264274               ! read after data 
    265                CALL fld_get( sd(jf), imap ) 
    266  
     275               IF( PRESENT(jpk_bdy) ) THEN 
     276                  CALL fld_get( sd(jf), imap, jpk_bdy, fvl) 
     277               ELSE 
     278                  CALL fld_get( sd(jf), imap ) 
     279               ENDIF 
    267280            ENDIF   ! read new data? 
    268281         END DO                                    ! --- end loop over field --- ! 
    … …  
    302315 
    303316 
    304    SUBROUTINE fld_init( kn_fsbc, sdjf, map ) 
     317   SUBROUTINE fld_init( kn_fsbc, sdjf, map , jpk_bdy, fvl) 
    305318      !!--------------------------------------------------------------------- 
    306319      !!                    ***  ROUTINE fld_init  *** 
    … …  
    309322      !!               - if time interpolation, read before data  
    310323      !!---------------------------------------------------------------------- 
    311       INTEGER  , INTENT(in   ) ::   kn_fsbc   ! sbc computation period (in time step)  
    312       TYPE(FLD), INTENT(inout) ::   sdjf      ! input field related variables 
    313       TYPE(MAP_POINTER),INTENT(in) ::   map   ! global-to-local mapping indices 
     324      INTEGER  , INTENT(in   ) ::   kn_fsbc      ! sbc computation period (in time step)  
     325      TYPE(FLD), INTENT(inout) ::   sdjf         ! input field related variables 
     326      TYPE(MAP_POINTER),INTENT(in) ::   map      ! global-to-local mapping indices 
     327      INTEGER  , INTENT(in), OPTIONAL :: jpk_bdy ! number of vertical levels in the BDY data 
     328      LOGICAL  , INTENT(in), OPTIONAL :: fvl     ! number of vertical levels in the BDY data 
    314329      !! 
    315330      LOGICAL :: llprevyr              ! are we reading previous year  file? 
    … …  
    405420         ENDIF 
    406421         ! 
    407          CALL fld_get( sdjf, map )         ! read before data in after arrays(as we will swap it later) 
     422         ! read before data in after arrays(as we will swap it later) 
     423         IF( PRESENT(jpk_bdy) ) THEN 
     424            CALL fld_get( sdjf, map, jpk_bdy, fvl ) 
     425         ELSE 
     426            CALL fld_get( sdjf, map ) 
     427         ENDIF 
    408428         ! 
    409429         clfmt = "('fld_init : time-interpolation for ', a, ' read previous record = ', i6, ' at time = ', f7.2, ' days')" 
    … …  
    581601 
    582602 
    583    SUBROUTINE fld_get( sdjf, map ) 
     603   SUBROUTINE fld_get( sdjf, map, jpk_bdy, fvl ) 
    584604      !!--------------------------------------------------------------------- 
    585605      !!                    ***  ROUTINE fld_get  *** 
    … …  
    589609      TYPE(FLD)        , INTENT(inout) ::   sdjf   ! input field related variables 
    590610      TYPE(MAP_POINTER), INTENT(in   ) ::   map    ! global-to-local mapping indices 
     611      INTEGER  , INTENT(in), OPTIONAL  ::   jpk_bdy ! number of vertical levels in the bdy data 
     612      LOGICAL  , INTENT(in), OPTIONAL  ::   fvl     ! number of vertical levels in the bdy data 
    591613      ! 
    592614      INTEGER ::   ipk      ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 
    … …  
    601623      ! 
    602624      IF( ASSOCIATED(map%ptr) ) THEN 
    603          IF( sdjf%ln_tint ) THEN   ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map ) 
    604          ELSE                      ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,:  ), sdjf%nrec_a(1), map ) 
    605          ENDIF 
     625         IF( PRESENT(jpk_bdy) ) THEN 
     626            IF( sdjf%ln_tint ) THEN   ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2),                & 
     627                                                        sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl ) 
     628            ELSE                      ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,:  ),                & 
     629                                                        sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl ) 
     630            ENDIF 
     631         ELSE 
     632            IF( sdjf%ln_tint ) THEN   ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2), sdjf%nrec_a(1), map ) 
     633            ELSE                      ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,:  ), sdjf%nrec_a(1), map ) 
     634            ENDIF 
     635         ENDIF         
    606636      ELSE IF( LEN(TRIM(sdjf%wgtname)) > 0 ) THEN 
    607637         CALL wgt_list( sdjf, iw ) 
    … …  
    658688   END SUBROUTINE fld_get 
    659689 
    660  
    661    SUBROUTINE fld_map( num, clvar, dta, nrec, map ) 
     690   SUBROUTINE fld_map( num, clvar, dta, nrec, map, igrd, ibdy, jpk_bdy, fvl ) 
    662691      !!--------------------------------------------------------------------- 
    663692      !!                    ***  ROUTINE fld_map  *** 
    … …  
    666695      !!                using a general mapping (for open boundaries) 
    667696      !!---------------------------------------------------------------------- 
    668 #if defined key_bdy 
    669       USE bdy_oce, ONLY:  dta_global, dta_global2         ! workspace to read in global data arrays 
    670 #endif  
     697 
     698      USE bdy_oce, ONLY: ln_bdy, idx_bdy, dta_global, dta_global_z, dta_global_dz, dta_global2, dta_global2_z, dta_global2_dz                 ! workspace to read in global data arrays 
     699 
    671700      INTEGER                   , INTENT(in ) ::   num     ! stream number 
    672701      CHARACTER(LEN=*)          , INTENT(in ) ::   clvar   ! variable name 
    673       REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   dta   ! output field on model grid (2 dimensional) 
     702      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   dta     ! output field on model grid (2 dimensional) 
    674703      INTEGER                   , INTENT(in ) ::   nrec    ! record number to read (ie time slice) 
    675704      TYPE(MAP_POINTER)         , INTENT(in ) ::   map     ! global-to-local mapping indices 
     705      INTEGER  , INTENT(in), OPTIONAL         ::   igrd, ibdy, jpk_bdy  ! grid type, set number and number of vertical levels in the bdy data 
     706      LOGICAL  , INTENT(in), OPTIONAL         ::   fvl     ! grid type, set number and number of vertical levels in the bdy data 
     707      INTEGER                                 ::   jpkm1_bdy! number of vertical levels in the bdy data minus 1 
    676708      !! 
    677709      INTEGER                                 ::   ipi      ! length of boundary data on local process 
    … …  
    682714      INTEGER                                 ::   ib, ik, ji, jj   ! loop counters 
    683715      INTEGER                                 ::   ierr 
    684       REAL(wp), POINTER, DIMENSION(:,:,:)     ::   dta_read  ! work space for global data 
     716      REAL(wp)                                ::   fv          ! fillvalue  
     717      REAL(wp), POINTER, DIMENSION(:,:,:)     ::   dta_read    ! work space for global data 
     718      REAL(wp), POINTER, DIMENSION(:,:,:)     ::   dta_read_z  ! work space for global data 
     719      REAL(wp), POINTER, DIMENSION(:,:,:)     ::   dta_read_dz ! work space for global data 
    685720      !!--------------------------------------------------------------------- 
    686721      ! 
    … …  
    692727      ilendta = iom_file(num)%dimsz(1,idvar) 
    693728 
    694 #if defined key_bdy 
    695       ipj = iom_file(num)%dimsz(2,idvar) 
    696       IF( map%ll_unstruc) THEN   ! unstructured open boundary data file 
    697          dta_read => dta_global 
    698       ELSE                       ! structured open boundary data file 
    699          dta_read => dta_global2 
     729      IF ( ln_bdy ) THEN 
     730         ipj = iom_file(num)%dimsz(2,idvar) 
     731         IF( map%ll_unstruc) THEN   ! unstructured open boundary data file 
     732            dta_read => dta_global 
     733            IF( PRESENT(jpk_bdy) ) THEN 
     734               IF( jpk_bdy>0 ) THEN 
     735                  dta_read_z => dta_global_z 
     736                  dta_read_dz => dta_global_dz 
     737                  jpkm1_bdy = jpk_bdy-1 
     738               ENDIF 
     739            ENDIF 
     740         ELSE                       ! structured open boundary file 
     741            dta_read => dta_global2 
     742            IF( PRESENT(jpk_bdy) ) THEN 
     743               IF( jpk_bdy>0 ) THEN 
     744                  dta_read_z => dta_global2_z 
     745                  dta_read_dz => dta_global2_dz 
     746                  jpkm1_bdy = jpk_bdy-1 
     747               ENDIF 
     748            ENDIF 
     749         ENDIF 
    700750      ENDIF 
    701 #endif 
    702751 
    703752      IF(lwp) WRITE(numout,*) 'Dim size for ',        TRIM(clvar),' is ', ilendta 
    … …  
    705754      ! 
    706755      SELECT CASE( ipk ) 
    707       CASE(1)        ;   CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1    ), nrec ) 
    708       CASE DEFAULT   ;   CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec ) 
     756      CASE(1)        ;    
     757      CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1    ), nrec ) 
     758         IF ( map%ll_unstruc) THEN ! unstructured open boundary data file 
     759            DO ib = 1, ipi 
     760              DO ik = 1, ipk 
     761                dta(ib,1,ik) =  dta_read(map%ptr(ib),1,ik) 
     762              END DO 
     763            END DO 
     764         ELSE ! we assume that this is a structured open boundary file 
     765            DO ib = 1, ipi 
     766               jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 
     767               ji=map%ptr(ib)-(jj-1)*ilendta 
     768               DO ik = 1, ipk 
     769                  dta(ib,1,ik) =  dta_read(ji,jj,ik) 
     770               END DO 
     771            END DO 
     772         ENDIF 
     773 
     774      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
     775      ! Do we include something here to adjust barotropic velocities ! 
     776      ! in case of a depth difference between bdy files and          ! 
     777      ! bathymetry in the case ln_full_vel = .false. and jpk_bdy>0?  ! 
     778      ! [as the enveloping and parital cells could change H]         ! 
     779      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
     780 
     781      CASE DEFAULT   ; 
     782 
     783      IF( PRESENT(jpk_bdy) .AND. jpk_bdy>0 ) THEN       ! boundary data not on model grid: vertical interpolation 
     784         CALL iom_getatt(num, '_FillValue', fv, cdvar=clvar ) 
     785         dta_read(:,:,:) = -ABS(fv) 
     786         dta_read_z(:,:,:) = 0._wp 
     787         dta_read_dz(:,:,:) = 0._wp 
     788         CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:jpk_bdy), nrec ) 
     789         SELECT CASE( igrd )                   
     790            CASE(1) 
     791               CALL iom_get ( num, jpdom_unknown, 'gdept', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 
     792               CALL iom_get ( num, jpdom_unknown, 'e3t',  dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 
     793            CASE(2)   
     794               CALL iom_get ( num, jpdom_unknown, 'gdepu', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 
     795               CALL iom_get ( num, jpdom_unknown, 'e3u',  dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 
     796            CASE(3) 
     797               CALL iom_get ( num, jpdom_unknown, 'gdepv', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 
     798               CALL iom_get ( num, jpdom_unknown, 'e3v',  dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 
     799         END SELECT 
     800 
     801      IF ( ln_bdy ) &  
     802         CALL fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta) 
     803 
     804      ELSE ! boundary data assumed to be on model grid 
     805          
     806         CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec )                     
     807         IF ( map%ll_unstruc) THEN ! unstructured open boundary data file 
     808            DO ib = 1, ipi 
     809              DO ik = 1, ipk 
     810                dta(ib,1,ik) =  dta_read(map%ptr(ib),1,ik) 
     811              END DO 
     812            END DO 
     813         ELSE ! we assume that this is a structured open boundary file 
     814            DO ib = 1, ipi 
     815               jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 
     816               ji=map%ptr(ib)-(jj-1)*ilendta 
     817               DO ik = 1, ipk 
     818                  dta(ib,1,ik) =  dta_read(ji,jj,ik) 
     819               END DO 
     820            END DO 
     821         ENDIF 
     822      ENDIF ! PRESENT(jpk_bdy) 
    709823      END SELECT 
    710       ! 
    711       IF( map%ll_unstruc ) THEN  ! unstructured open boundary data file 
     824 
     825   END SUBROUTINE fld_map 
     826    
     827   SUBROUTINE fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta) 
     828 
     829      !!--------------------------------------------------------------------- 
     830      !!                    ***  ROUTINE fld_bdy_interp  *** 
     831      !! 
     832      !! ** Purpose :   on the fly vertical interpolation to allow the use of 
     833      !!                boundary data from non-native vertical grid 
     834      !!---------------------------------------------------------------------- 
     835      USE bdy_oce, ONLY:  idx_bdy         ! indexing for map <-> ij transformation 
     836 
     837      REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in )     ::   dta_read      ! work space for global data 
     838      REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in )     ::   dta_read_z    ! work space for global data 
     839      REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in )     ::   dta_read_dz   ! work space for global data 
     840      REAL(wp) , INTENT(in)                                ::   fv            ! fillvalue and alternative -ABS(fv) 
     841      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   dta                        ! output field on model grid (2 dimensional) 
     842      TYPE(MAP_POINTER)         , INTENT(in ) ::   map                        ! global-to-local mapping indices 
     843      LOGICAL  , INTENT(in), OPTIONAL         ::   fvl                        ! grid type, set number and number of vertical levels in the bdy data 
     844      INTEGER  , INTENT(in)                   ::   igrd, ibdy, jpk_bdy        ! number of levels in bdy data 
     845      INTEGER  , INTENT(in)                   ::   ilendta                    ! length of data in file 
     846      !! 
     847      INTEGER                                 ::   ipi                        ! length of boundary data on local process 
     848      INTEGER                                 ::   ipj                        ! length of dummy dimension ( = 1 ) 
     849      INTEGER                                 ::   ipk                        ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk ) 
     850      INTEGER                                 ::   jpkm1_bdy                  ! number of levels in bdy data minus 1 
     851      INTEGER                                 ::   ib, ik, ikk                ! loop counters 
     852      INTEGER                                 ::   ji, jj, zij, zjj           ! temporary indices 
     853      REAL(wp)                                ::   zl, zi, zh                 ! tmp variable for current depth and interpolation factor 
     854      REAL(wp)                                ::   fv_alt, ztrans, ztrans_new ! fillvalue and alternative -ABS(fv) 
     855      CHARACTER (LEN=10)                      ::   ibstr 
     856      !!--------------------------------------------------------------------- 
     857      
     858 
     859      ipi       = SIZE( dta, 1 ) 
     860      ipj       = SIZE( dta_read, 2 ) 
     861      ipk       = SIZE( dta, 3 ) 
     862      jpkm1_bdy = jpk_bdy-1 
     863       
     864      fv_alt = -ABS(fv)  ! set _FillValue < 0 as we make use of MAXVAL and MAXLOC later 
     865      DO ib = 1, ipi