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 5868 for branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC – NEMO

Ignore:
Timestamp:
2015-11-06T19:07:02+01:00 (9 years ago)
Author:
jchanut
Message:

Free surface simplification #1620. Step 1: suppress filtered free surface

Location:
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC
Files:
2 deleted
13 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/ASM/asmbkg.F90

    r5836 r5868  
    3434   USE zdfmxl             ! Mixed layer depth 
    3535   USE dom_oce, ONLY :   ndastp 
    36    USE sol_oce, ONLY :   gcx   ! Solver variables defined in memory 
    3736   USE in_out_manager     ! I/O manager 
    3837   USE iom                ! I/O module 
     
    114113            CALL iom_rstput( kt, nitbkg_r, inum, 'en'     , en                ) 
    115114#endif 
    116             CALL iom_rstput( kt, nitbkg_r, inum, 'gcx'    , gcx               ) 
    117115            ! 
    118116            CALL iom_close( inum ) 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90

    r4689 r5868  
    3535   PRIVATE 
    3636 
    37    PUBLIC   bdy_dyn     ! routine called in dynspg_flt (if lk_dynspg_flt) or  
    38                         ! dyn_nxt (if lk_dynspg_ts or lk_dynspg_exp) 
     37   PUBLIC   bdy_dyn    ! routine called in dyn_nxt 
    3938 
    4039#  include "domzgr_substitute.h90" 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90

    r5836 r5868  
    1010   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge 
    1111   !!---------------------------------------------------------------------- 
    12 #if   defined key_bdy   &&   defined key_dynspg_flt 
     12#if   defined key_bdy 
    1313   !!---------------------------------------------------------------------- 
    14    !!   'key_bdy'            AND      unstructured open boundary conditions 
    15    !!   'key_dynspg_flt'                              filtered free surface 
     14   !!   'key_bdy'      unstructured open boundary conditions 
    1615   !!---------------------------------------------------------------------- 
    1716   USE oce             ! ocean dynamics and tracers  
     
    3029   PRIVATE 
    3130 
    32    PUBLIC bdy_vol        ! routine called by dynspg_flt.h90 
     31   PUBLIC bdy_vol       
    3332 
    3433   !! * Substitutions 
     
    4544      !!                      ***  ROUTINE bdyvol  *** 
    4645      !! 
    47       !! ** Purpose :   This routine is called in dynspg_flt to control  
    48       !!      the volume of the system. A correction velocity is calculated 
     46      !! ** Purpose :   This routine controls the volume of the system.  
     47      !!      A correction velocity is calculated 
    4948      !!      to correct the total transport through the unstructured OBC.  
    5049      !!      The total depth used is constant (H0) to be consistent with the  
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r5836 r5868  
    3030   USE zdf_oce         ! ocean vertical physics 
    3131   USE ldftra          ! lateral physics: eddy diffusivity coef. 
    32    USE sol_oce         ! solver variables 
    3332   USE sbc_oce         ! Surface boundary condition: ocean fields 
    3433   USE sbc_ice         ! Surface boundary condition: ice fields 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r5836 r5868  
    3535   USE dtauvd          ! data: U & V current             (dta_uvd routine) 
    3636   USE domvvl          ! varying vertical mesh 
    37    USE dynspg_oce      ! pressure gradient schemes 
    38    USE dynspg_flt      ! filtered free surface 
    39    USE sol_oce         ! ocean solver variables 
    4037   ! 
    4138   USE in_out_manager  ! I/O manager 
     
    133130      ENDIF 
    134131      ! 
    135       IF( lk_agrif ) THEN                  ! read free surface arrays in restart file 
    136          IF( ln_rstart ) THEN 
    137             IF( lk_dynspg_flt )  THEN      ! read or initialize the following fields 
    138                !                           ! gcx, gcxb for agrif_opa_init 
    139                IF( sol_oce_alloc()  > 0 )   CALL ctl_stop('agrif sol_oce_alloc: allocation of arrays failed') 
    140                CALL flt_rst( nit000, 'READ' ) 
    141             ENDIF 
    142          ENDIF                             ! explicit case not coded yet with AGRIF 
    143       ENDIF 
    144       ! 
    145132      !  
    146133      ! Initialize "now" and "before" barotropic velocities: 
     
    445432      !!                 p=integral [ rau*g dz ] 
    446433      !!---------------------------------------------------------------------- 
    447       USE dynspg          ! surface pressure gradient             (dyn_spg routine) 
    448434      USE divhor          ! hor. divergence                       (div_hor routine) 
    449435      USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    450436      ! 
    451437      INTEGER ::   ji, jj, jk        ! dummy loop indices 
    452       INTEGER ::   indic             ! ??? 
    453438      REAL(wp) ::   zmsv, zphv, zmsu, zphu, zalfg     ! temporary scalars 
    454439      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprn 
     
    517502      vb(:,:,:) = vn(:,:,:) 
    518503       
    519       ! WARNING !!!!! 
    520       ! after initializing u and v, we need to calculate the initial streamfunction bsf. 
    521       ! Otherwise, only the trend will be computed and the model will blow up (inconsistency). 
    522       ! to do that, we call dyn_spg with a special trick: 
    523       ! we fill ua and va with the velocities divided by dt, and the streamfunction will be brought to the 
    524       ! right value assuming the velocities have been set up in one time step. 
    525       ! we then set bsfd to zero (first guess for next step is d(psi)/dt = 0.) 
    526       !  sets up s false trend to calculate the barotropic streamfunction. 
    527  
    528       ua(:,:,:) = ub(:,:,:) / rdt 
    529       va(:,:,:) = vb(:,:,:) / rdt 
    530  
    531       ! calls dyn_spg. we assume euler time step, starting from rest. 
    532       indic = 0 
    533       CALL dyn_spg( nit000, indic )       ! surface pressure gradient 
    534       ! 
    535       ! the new velocity is ua*rdt 
    536       ! 
    537       CALL lbc_lnk( ua, 'U', -1. ) 
    538       CALL lbc_lnk( va, 'V', -1. ) 
    539  
    540       ub(:,:,:) = ua(:,:,:) * rdt 
    541       vb(:,:,:) = va(:,:,:) * rdt 
    542       ua(:,:,:) = 0.e0 
    543       va(:,:,:) = 0.e0 
    544       un(:,:,:) = ub(:,:,:) 
    545       vn(:,:,:) = vb(:,:,:) 
    546504      ! 
    547505!!gm  Check  here call to div_hor should not be necessary 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90

    r5643 r5868  
    7878      !!             (lk_vvl=T), the leap-frog is applied on thickness weighted 
    7979      !!             velocity. 
    80       !!             Note also that in filtered free surface (lk_dynspg_flt=T), 
    81       !!             the time stepping has already been done in dynspg module 
    8280      !! 
    8381      !!              * Apply lateral boundary conditions on after velocity  
     
    10199      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    102100      INTEGER  ::   iku, ikv     ! local integers 
    103 #if ! defined key_dynspg_flt 
    104101      REAL(wp) ::   z2dt         ! temporary scalar 
    105 #endif 
    106102      REAL(wp) ::   zue3a, zue3n, zue3b, zuf, zec      ! local scalars 
    107103      REAL(wp) ::   zve3a, zve3n, zve3b, zvf, z1_2dt   !   -      - 
     
    120116         IF(lwp) WRITE(numout,*) '~~~~~~~' 
    121117      ENDIF 
    122  
    123 #if defined key_dynspg_flt 
    124       ! 
    125       ! Next velocity :   Leap-frog time stepping already done in dynspg_flt.F routine 
    126       ! ------------- 
    127  
    128       ! Update after velocity on domain lateral boundaries      (only local domain required) 
    129       ! -------------------------------------------------- 
    130       CALL lbc_lnk( ua, 'U', -1. )         ! local domain boundaries 
    131       CALL lbc_lnk( va, 'V', -1. )  
    132       ! 
    133 #else 
    134118 
    135119# if defined key_dynspg_exp 
     
    201185      CALL Agrif_dyn( kt )             !* AGRIF zoom boundaries 
    202186# endif 
    203 #endif 
    204187 
    205188      IF( l_trddyn ) THEN             ! prepare the atf trend computation + some diagnostics 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90

    r5836 r5868  
    2121   USE dynspg_exp     ! surface pressure gradient     (dyn_spg_exp routine) 
    2222   USE dynspg_ts      ! surface pressure gradient     (dyn_spg_ts  routine) 
    23    USE dynspg_flt     ! surface pressure gradient     (dyn_spg_flt routine) 
    2423   USE dynadv         ! dynamics: vector invariant versus flux form 
    2524   USE dynhpg, ONLY: ln_dynhpg_imp 
     
    3231   USE in_out_manager ! I/O manager 
    3332   USE lib_mpp        ! MPP library 
    34    USE solver         ! solver initialization 
    3533   USE wrk_nemo       ! Memory Allocation 
    3634   USE timing         ! Timing 
     
    5553CONTAINS 
    5654 
    57    SUBROUTINE dyn_spg( kt, kindic ) 
     55   SUBROUTINE dyn_spg( kt ) 
    5856      !!---------------------------------------------------------------------- 
    5957      !!                  ***  ROUTINE dyn_spg  *** 
     
    8078      ! 
    8179      INTEGER, INTENT(in   ) ::   kt       ! ocean time-step index 
    82       INTEGER, INTENT(  out) ::   kindic   ! solver flag 
    8380      ! 
    8481      INTEGER  ::   ji, jj, jk                             ! dummy loop indices 
     
    171168      CASE (  0 )   ;   CALL dyn_spg_exp( kt )              ! explicit 
    172169      CASE (  1 )   ;   CALL dyn_spg_ts ( kt )              ! time-splitting 
    173       CASE (  2 )   ;   CALL dyn_spg_flt( kt, kindic )      ! filtered 
    174170      !                                                     
    175171      END SELECT 
    176172      !                     
    177173      IF( l_trddyn )   THEN                      ! save the surface pressure gradient trends for further diagnostics 
    178          SELECT CASE ( nspg ) 
    179          CASE ( 0, 1 ) 
    180             ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    181             ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    182          CASE( 2 ) 
    183             z2dt = 2. * rdt 
    184             IF( neuler == 0 .AND. kt == nit000 )   z2dt = rdt 
    185             ztrdu(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) / z2dt - ztrdu(:,:,:) 
    186             ztrdv(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) / z2dt - ztrdv(:,:,:) 
    187          END SELECT 
     174         ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
     175         ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    188176         CALL trd_dyn( ztrdu, ztrdv, jpdyn_spg, kt ) 
    189177         ! 
     
    217205         WRITE(numout,*) '     Explicit free surface                  lk_dynspg_exp = ', lk_dynspg_exp 
    218206         WRITE(numout,*) '     Free surface with time splitting       lk_dynspg_ts  = ', lk_dynspg_ts 
    219          WRITE(numout,*) '     Filtered free surface cst volume       lk_dynspg_flt = ', lk_dynspg_flt 
    220207      ENDIF 
    221208 
     
    233220      IF(lk_dynspg_exp)   ioptio = ioptio + 1 
    234221      IF(lk_dynspg_ts )   ioptio = ioptio + 1 
    235       IF(lk_dynspg_flt)   ioptio = ioptio + 1 
    236222      ! 
    237223      IF(  ioptio > 1  .OR. ( ioptio == 0 .AND. .NOT. lk_c1d ) )   & 
     
    242228      IF( lk_dynspg_exp)   nspg =  0 
    243229      IF( lk_dynspg_ts )   nspg =  1 
    244       IF( lk_dynspg_flt)   nspg =  2 
     230 
    245231      ! 
    246232      IF(lwp) THEN 
     
    248234         IF( nspg ==  0 )   WRITE(numout,*) '     explicit free surface' 
    249235         IF( nspg ==  1 )   WRITE(numout,*) '     free surface with time splitting scheme' 
    250          IF( nspg ==  2 )   WRITE(numout,*) '     filtered free surface' 
    251       ENDIF 
    252  
    253 #if defined key_dynspg_flt 
    254       CALL solver_init( nit000 )   ! Elliptic solver initialisation 
    255 #endif 
     236      ENDIF 
    256237      !               ! Control of hydrostatic pressure choice 
    257238      IF( lk_dynspg_ts .AND. ln_dynhpg_imp )   CALL ctl_stop( 'Semi-implicit hpg not compatible with time splitting' ) 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_oce.F90

    r5836 r5868  
    2626#else 
    2727   LOGICAL, PUBLIC, PARAMETER ::   lk_dynspg_ts  = .FALSE.  !: Free surface with time splitting flag 
    28 #endif 
    29 #if   defined key_dynspg_flt 
    30    LOGICAL, PUBLIC, PARAMETER ::   lk_dynspg_flt = .TRUE.   !: Filtered free surface cst volume flag 
    31 #else 
    32    LOGICAL, PUBLIC, PARAMETER ::   lk_dynspg_flt = .FALSE.  !: Filtered free surface cst volume flag 
    3328#endif 
    3429  !                                                                         !!! Time splitting scheme (key_dynspg_ts)  
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/TRD/trd_oce.F90

    r5215 r5868  
    7171   INTEGER, PUBLIC, PARAMETER ::   jpdyn_bfri = 12     !: implicit bottom friction (ln_bfrimp=.TRUE.) 
    7272   INTEGER, PUBLIC, PARAMETER ::   jpdyn_ken  = 13     !: use for calculation of KE 
    73    INTEGER, PUBLIC, PARAMETER ::   jpdyn_spgflt  = 14  !: filter contribution to surface pressure gradient (spg_flt) 
    74    INTEGER, PUBLIC, PARAMETER ::   jpdyn_spgexp  = 15  !: explicit contribution to surface pressure gradient (spg_flt) 
    7573   ! 
    7674   !!---------------------------------------------------------------------- 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/TRD/trddyn.F90

    r5215 r5868  
    113113      CASE( jpdyn_spg )   ;   CALL iom_put( "utrd_spg", putrd )    ! surface pressure gradient 
    114114                              CALL iom_put( "vtrd_spg", pvtrd ) 
    115       CASE( jpdyn_spgexp );   CALL iom_put( "utrd_spgexp", putrd ) ! surface pressure gradient (explicit) 
    116                               CALL iom_put( "vtrd_spgexp", pvtrd ) 
    117       CASE( jpdyn_spgflt );   CALL iom_put( "utrd_spgflt", putrd ) ! surface pressure gradient (filtered) 
    118                               CALL iom_put( "vtrd_spgflt", pvtrd ) 
    119115      CASE( jpdyn_pvo )   ;   CALL iom_put( "utrd_pvo", putrd )    ! planetary vorticity 
    120116                              CALL iom_put( "vtrd_pvo", pvtrd ) 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90

    r5836 r5868  
    120120         CASE( jpdyn_hpg )   ;   CALL iom_put( "ketrd_hpg", zke )    ! hydrostatic pressure gradient 
    121121         CASE( jpdyn_spg )   ;   CALL iom_put( "ketrd_spg", zke )    ! surface pressure gradient 
    122          CASE( jpdyn_spgexp );   CALL iom_put( "ketrd_spgexp", zke ) ! surface pressure gradient (explicit) 
    123          CASE( jpdyn_spgflt );   CALL iom_put( "ketrd_spgflt", zke ) ! surface pressure gradient (filter) 
    124122         CASE( jpdyn_pvo )   ;   CALL iom_put( "ketrd_pvo", zke )    ! planetary vorticity 
    125123         CASE( jpdyn_rvo )   ;   CALL iom_put( "ketrd_rvo", zke )    ! relative  vorticity     (or metric term) 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/step.F90

    r5836 r5868  
    214214#endif 
    215215                                  CALL dyn_hpg( kstp )         ! horizontal gradient of Hydrostatic pressure 
    216                                   CALL dyn_spg( kstp, indic )  ! surface pressure gradient 
     216                                  CALL dyn_spg( kstp )         ! surface pressure gradient 
    217217 
    218218                                  ua_sv(:,:,:) = ua(:,:,:)     ! Save trends (barotropic trend has been fully updated at this stage) 
     
    329329                               CALL dyn_bfr( kstp )         ! bottom friction 
    330330                               CALL dyn_zdf( kstp )         ! vertical diffusion 
    331                                CALL dyn_spg( kstp, indic )  ! surface pressure gradient 
     331                               CALL dyn_spg( kstp )         ! surface pressure gradient 
    332332      ENDIF 
    333333                               CALL dyn_nxt( kstp )         ! lateral velocity at next time step 
  • branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/stpctl.F90

    r3294 r5868  
    1616   USE oce             ! ocean dynamics and tracers variables 
    1717   USE dom_oce         ! ocean space and time domain variables  
    18    USE sol_oce         ! ocean space and time domain variables  
    1918   USE in_out_manager  ! I/O manager 
    2019   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2120   USE lib_mpp         ! distributed memory computing 
    22    USE dynspg_oce      ! pressure gradient schemes  
    2321   USE c1d             ! 1D vertical configuration 
    2422 
     
    4341      !! ** Method  : - Save the time step in numstp 
    4442      !!              - Print it each 50 time steps 
    45       !!              - Print solver statistics in numsol  
    46       !!              - Stop the run IF problem for the solver ( indec < 0 ) 
     43      !!              - Stop the run IF problem ( indic < 0 ) 
    4744      !! 
    4845      !! ** Actions :   'time.step' file containing the last ocean time-step 
     
    5047      !!---------------------------------------------------------------------- 
    5148      INTEGER, INTENT( in ) ::   kt         ! ocean time-step index 
    52       INTEGER, INTENT( inout ) ::   kindic  ! indicator of solver convergence 
     49      INTEGER, INTENT( inout ) ::   kindic  ! error indicator 
    5350      !! 
    5451      INTEGER  ::   ji, jj, jk              ! dummy loop indices 
     
    143140      IF( lk_c1d )  RETURN          ! No log file in case of 1D vertical configuration 
    144141 
    145       ! log file (solver or ssh statistics) 
    146       ! -------- 
    147       IF( lk_dynspg_flt ) THEN      ! elliptic solver statistics (if required) 
    148          ! 
    149          IF(lwp) WRITE(numsol,9200) kt, niter, res, SQRT(epsr)/eps       ! Solver 
    150          ! 
    151          IF( kindic < 0 .AND. zsmin > 0.e0 .AND. zumax <= 20.e0 ) THEN   ! create a abort file if problem found  
    152             IF(lwp) THEN 
    153                WRITE(numout,*) ' stpctl: the elliptic solver DO not converge or explode' 
    154                WRITE(numout,*) ' ====== ' 
    155                WRITE(numout,9200) kt, niter, res, sqrt(epsr)/eps 
    156                WRITE(numout,*) 
    157                WRITE(numout,*) ' stpctl: output of last fields' 
    158                WRITE(numout,*) ' ======  ' 
    159             ENDIF 
    160          ENDIF 
    161          ! 
    162       ELSE                                   !* ssh statistics (and others...) 
    163          IF( kt == nit000 .AND. lwp ) THEN   ! open ssh statistics file (put in solver.stat file) 
    164             CALL ctl_opn( numsol, 'solver.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
    165          ENDIF 
    166          ! 
    167          zssh2 = SUM( sshn(:,:) * sshn(:,:) * tmask_i(:,:) ) 
    168          IF( lk_mpp )   CALL mpp_sum( zssh2 )      ! sum over the global domain 
    169          ! 
    170          IF(lwp) WRITE(numsol,9300) kt, zssh2, zumax, zsmin      ! ssh statistics 
    171          ! 
     142      ! log file (ssh statistics) 
     143      ! --------                                   !* ssh statistics (and others...) 
     144      IF( kt == nit000 .AND. lwp ) THEN   ! open ssh statistics file (put in solver.stat file) 
     145         CALL ctl_opn( numsol, 'solver.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
    172146      ENDIF 
     147      ! 
     148      zssh2 = SUM( sshn(:,:) * sshn(:,:) * tmask_i(:,:) ) 
     149      IF( lk_mpp )   CALL mpp_sum( zssh2 )      ! sum over the global domain 
     150      ! 
     151      IF(lwp) WRITE(numsol,9300) kt, zssh2, zumax, zsmin      ! ssh statistics 
     152      ! 
    173153 
    1741549200  FORMAT('it:', i8, ' iter:', i4, ' r: ',e16.10, ' b: ',e16.10 ) 
Note: See TracChangeset for help on using the changeset viewer.