Changeset 593

Timestamp:

2007-02-09T10:40:34+01:00 (17 years ago)

Author:

opalod

Message:

nemo_v2_update_001 : CT : - add non linear free surface (variable volume) with new cpp key key_vvl

Location:

trunk/NEMO/OPA_SRC

Files:

• eosbn2.F90 (modified) (1 diff)
• istate.F90 (modified) (6 diffs)
• oce.F90 (modified) (1 diff)
• par_EEL_R5.h90 (modified) (4 diffs)
• restart.F90 (modified) (3 diffs)
• step.F90 (modified) (2 diffs)

trunk/NEMO/OPA_SRC/eosbn2.F90

@@ -36 +36 @@
    PUBLIC eos        ! called by step.F90, inidtr.F90, tranpc.F90 and intgrd.F90
    PUBLIC bn2        ! called by step.F90
+   PUBLIC eos_init   ! called by step.F90
 
    !! * Share module variables

trunk/NEMO/OPA_SRC/istate.F90

@@ -39 +39 @@
    USE zpshde          ! partial step: hor. derivative (zps_hde routine)
    USE eosbn2          ! equation of state            (eos bn2 routine)
+   USE domvvl          ! varying vertical mesh
+   USE dynspg_oce      ! pressure gradient schemes
+   USE dynspg_flt      ! pressure gradient schemes
+   USE dynspg_exp      ! pressure gradient schemes
+   USE dynspg_ts       ! pressure gradient schemes
    
    IMPLICIT NONE
@@ -83 +88 @@
          !                                       ! Initialization of ocean to zero
          !     before fields       !       now fields          
-         ;   ub   (:,:,:) = 0.e0   ;   un   (:,:,:) = 0.e0   
-         ;   vb   (:,:,:) = 0.e0   ;   vn   (:,:,:) = 0.e0    
-         ;   rotb (:,:,:) = 0.e0   ;   rotn (:,:,:) = 0.e0  
-         ;   hdivb(:,:,:) = 0.e0   ;   hdivn(:,:,:) = 0.e0  
+         ;   ub   (:,:,:) = 0.e0   ;   un   (:,:,:) = 0.e0   ; sshb(:,:) = 0.e0
+         ;   vb   (:,:,:) = 0.e0   ;   vn   (:,:,:) = 0.e0   ; sshn(:,:) = 0.e0
+         ;   rotb (:,:,:) = 0.e0   ;   rotn (:,:,:) = 0.e0
+         ;   hdivb(:,:,:) = 0.e0   ;   hdivn(:,:,:) = 0.e0
          !
          IF( cp_cfg == 'eel' ) THEN
@@ -126 +131 @@
       ENDIF
 
+      IF( lk_vvl ) THEN
+         ! read free surface arrays in restart file
+         IF( ln_rstart ) THEN
+            IF( lk_dynspg_flt )   CALL flt_rst( nit000, 'READ' )      ! read or initialize the following fields
+            !                                                         ! gcx, gcxb, sshb, sshn
+            IF( lk_dynspg_ts  )   CALL ts_rst ( nit000, 'READ' )      ! read or initialize the following fields
+            !                                                         ! sshb, sshn, sshb_b, sshn_b, un_b, vn_b
+            IF( lk_dynspg_exp )   CALL exp_rst( nit000, 'READ' )      ! read or initialize the following fields
+            !                                                         ! sshb, sshn
+         ENDIF
+         !
+         IF( .NOT. lk_dynspg_flt ) sshbb(:,:) = sshb(:,:)
+         !
+         CALL dom_vvl               ! ssh init and vertical grid update
+
+         CALL eos_init              ! usefull to get the equation state type neos parameter
+
+         CALL bn2( tb, sb, rn2 )    ! before Brunt Vaissala frequency
+
+         IF( .NOT. ln_rstart ) CALL wzv( nit000 ) 
+
+      ENDIF
 
       !                                       ! Vertical velocity
       !                                       ! -----------------
-      CALL wzv( nit000 )                         ! from horizontal divergence
+
+      IF( .NOT. lk_vvl )    CALL wzv( nit000 )                         ! from horizontal divergence
       !
    END SUBROUTINE istate_init
@@ -213 +241 @@
       INTEGER  ::   ijloc
       REAL(wp) ::   zh1, zh2, zslope, zcst, zfcor   ! temporary scalars
-      REAL(wp) ::   zt1  = 12._wp,               &  ! surface temperature value (EEL R5)
-         &          zt2  =  2._wp,               &  ! bottom  temperature value (EEL R5)
-         &          zsal = 35.5_wp,              &  ! constant salinity (EEL R2, R5 and R6)
+      REAL(wp) ::   zt1  = 15._wp,               &  ! surface temperature value (EEL R5)
+         &          zt2  =  5._wp,               &  ! bottom  temperature value (EEL R5)
+         &          zsal = 35.0_wp,              &  ! constant salinity (EEL R2, R5 and R6)
          &          zueel = 0.1_wp                  ! constant uniform zonal velocity (EEL R5)
 # if ! defined key_dynspg_rl
@@ -240 +268 @@
 
             DO jk = 1, jpk
-               tn(:,:,jk) = ( zslope * fsdept(:,:,jk) + zcst  ) * tmask(:,:,jk)
+               tn(:,:,jk) = ( zt2 + zt1 * exp( - fsdept(:,:,jk) / 1000 ) ) * tmask(:,:,jk)
                tb(:,:,jk) = tn(:,:,jk)
             END DO
@@ -291 +319 @@
             sshn(:,:) = sshb(:,:)                   ! set now ssh to the before value
 
-            IF( nn_rstssh /= 1 ) THEN  
-               nn_rstssh = 1                           ! hand-made initilization of ssh 
-               CALL ctl_warn( 'istate_eel: force nn_rstssh = 1' )
+            IF( nn_rstssh /= 0 ) THEN  
+               nn_rstssh = 0                           ! hand-made initilization of ssh 
+               CALL ctl_warn( 'istate_eel: force nn_rstssh = 0' )
             ENDIF
 

trunk/NEMO/OPA_SRC/oce.F90

@@ -61 +61 @@
    !! ------------
    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &   !:
-      sshb, sshn             !: before, now sea surface height (meters)
+      sshb , sshn ,        &  !: before, now sea surface height (meters)
+      sshu , sshv ,        &  !: sea surface height at u- and v- point
+      sshbb, ssha             !: before before sea surface height at t-point
 
 #if defined key_dynspg_rl   ||   defined key_esopa

trunk/NEMO/OPA_SRC/par_EEL_R5.h90

@@ -22 +22 @@
 
       ! data size              !!! * size of all the input files
-      jpidta  = 202   ,      &  !: first horizontal dimension > or = to jpi
-      jpjdta  = 104   ,      &  !: second                     > or = to jpj
-      jpkdta  =  40   ,      &  !: number of levels           > or = to jpk
+      jpidta  =  66   ,      &  !: first horizontal dimension > or = to jpi
+      jpjdta  =  66   ,      &  !: second                     > or = to jpj
+      jpkdta  =  31   ,      &  !: number of levels           > or = to jpk
 
       ! total domain size      !!! * full domain *
@@ -35 +35 @@
 
       ! Domain characteristics
-      jperio  =   0   ,      &  !: lateral cond. type (between 0 and 6)
+      jperio  =   1   ,      &  !: lateral cond. type (between 0 and 6)
       jpisl   =   1   ,      &  !: number of islands
       jpnisl  = jpiglo          !: maximum number of points per island
@@ -67 +67 @@
       ::     &  !:
       ppglam0  =    0.0_wp,   &  !: longitude of first raw and column T-point (jphgr_msh = 1)
-      ppgphi0  =   35.0_wp,   &  !: latitude  of first raw and column T-point (jphgr_msh = 1)
+      ppgphi0  = 43.436430714_wp,   &  !: latitude  of first raw and column T-point (jphgr_msh = 1)
       !                          ! latitude for the Coriolis or Beta parameter (jphgr_msh = 2 or 3)
       ppe1_deg = pp_not_used,   &  !: zonal      grid-spacing (degrees)
       ppe2_deg = pp_not_used,   &  !: meridional grid-spacing (degrees)
       !
-      ppe1_m   = 5000.0_wp,   &  !: zonal      grid-spacing (meters)
-      ppe2_m   = 5000.0_wp       !: meridional grid-spacing (meters)
+      ppe1_m   = 8000.0_wp,   &  !: zonal      grid-spacing (meters)
+      ppe2_m   = 8000.0_wp       !: meridional grid-spacing (meters)
    !!
    !!  Coefficients associated with the vertical coordinate system
    !!
    REAL(wp), PARAMETER  ::       &   !:
-      &     ppsur = pp_to_be_computed       ,  &  !: Computed in domzgr, set ppdzmin, pphmax below
-      &     ppa0  = pp_to_be_computed       ,  &  !:
-      &     ppa1  = pp_to_be_computed       ,  &  !:
+      &     ppsur = -4762.96143546300_wp    ,  &  !: Computed in domzgr, set ppdzmin, pphmax below
+      &     ppa0  =   255.58049070440_wp    ,  &  !:
+      &     ppa1  =   245.58132232490_wp    ,  &  !:
       !
-      &     ppkth = 20._wp                  ,  &  !: (non dimensional): gives the approximate
+      &     ppkth =    21.43336197938_wp    ,  &  !: (non dimensional): gives the approximate
       !                                           !: layer number above which  stretching will
       !                                           !: be maximum. Usually of order jpk/2.
-      &     ppacr = 16.00000000000_wp             !: (non dimensional): stretching factor
+      &     ppacr =     3.00000000000_wp          !: (non dimensional): stretching factor
       !                                           !: for the grid. The highest zacr, the smallest
       !                                           !: the stretching.
@@ -94 +94 @@
    !!
    REAL(wp), PARAMETER ::        &  !:
-      &     ppdzmin = 120._wp               ,  &  !: (meters) vertical thickness of the top layer
-      &     pphmax  = 4000._wp                    !: (meters) Maximum depth of the ocean gdepw(jpk)
+      &     ppdzmin = pp_not_used           ,  &  !: (meters) vertical thickness of the top layer
+      &     pphmax  = pp_not_used                 !: (meters) Maximum depth of the ocean gdepw(jpk)
    !!---------------------------------------------------------------------

trunk/NEMO/OPA_SRC/restart.F90

@@ -127 +127 @@
       CALL iom_rstput( kt, nitrst, numrow, 'un'     , un      )
       CALL iom_rstput( kt, nitrst, numrow, 'vn'     , vn      )
+      IF( lk_vvl ) CALL iom_rstput( kt, nitrst, numrow, 'wn'     , wn      )
       CALL iom_rstput( kt, nitrst, numrow, 'tn'     , tn      )
       CALL iom_rstput( kt, nitrst, numrow, 'sn'     , sn      )
@@ -147 +148 @@
 #endif
 
-      IF( nn_dynhpg_rst == 1 ) THEN
+      IF( nn_dynhpg_rst == 1 .OR. lk_vvl ) THEN
          CALL iom_rstput( kt, nitrst, numrow, 'rhd' , rhd  )
          CALL iom_rstput( kt, nitrst, numrow, 'rhop', rhop )
@@ -271 +272 @@
       CALL iom_get( numror, jpdom_local, 'un'   , un    )        ! now    i-component velocity
       CALL iom_get( numror, jpdom_local, 'vn'   , vn    )        ! now    j-component velocity
+      IF( lk_vvl ) CALL iom_get( numror, jpdom_local, 'wn'   , wn    )        ! now    k-component velocity
       CALL iom_get( numror, jpdom_local, 'tn'   , tn    )        ! now    temperature
       CALL iom_get( numror, jpdom_local, 'sn'   , sn    )        ! now    salinity

trunk/NEMO/OPA_SRC/step.F90

@@ -119 +119 @@
    USE cpl             ! exchanges in coupled mode        (cpl_stp routine)
    USE prtctl          ! Print control                    (prt_ctl routine)
+   USE domvvl          ! variable volume                  (dom_vvl routine)
 
 #if defined key_agrif
@@ -385 +386 @@
                                CALL dyn_nxt( kstp )           ! lateral velocity at next time step
 
+      IF( lk_vvl )             CALL dom_vvl                   ! vertical mesh at next time step
+
 
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>