Changeset 12489 for NEMO/trunk/src/OCE/SBC

Timestamp:

2020-02-28T16:55:11+01:00 (4 years ago)

Author:

davestorkey

Message:

Preparation for new timestepping scheme #2390.
Main changes:

1. Initial euler timestep now handled in stp and not in TRA/DYN routines.
2. Renaming of all timestep parameters. In summary, the namelist parameter is now rn_Dt and the current timestep is rDt (and rDt_ice, rDt_trc etc).
3. Renaming of a few miscellaneous parameters, eg. atfp -> rn_atfp (namelist parameter used everywhere) and rau0 -> rho0.

This version gives bit-comparable results to the previous version of the trunk.

Location:

NEMO/trunk/src/OCE/SBC

Files:

• fldread.F90 (modified) (1 diff)
• sbcapr.F90 (modified) (2 diffs)
• sbcblk.F90 (modified) (1 diff)
• sbcblk_skin_coare.F90 (modified) (3 diffs)
• sbcblk_skin_ecmwf.F90 (modified) (1 diff)
• sbccpl.F90 (modified) (4 diffs)
• sbcdcy.F90 (modified) (2 diffs)
• sbcfwb.F90 (modified) (1 diff)
• sbcice_cice.F90 (modified) (4 diffs)
• sbcmod.F90 (modified) (3 diffs)
• sbcrnf.F90 (modified) (4 diffs)

NEMO/trunk/src/OCE/SBC/fldread.F90

@@ -172 +172 @@
       ! Note that all varibles starting by nsec_* are shifted time by +1/2 time step to be centrered
       IF( PRESENT(kit) ) THEN   ! ignore kn_fsbc in this case
-         isecsbc = nsec_year + nsec1jan000 + NINT( (     REAL(      kit,wp) + zt_offset ) * rdt / REAL(nn_baro,wp) )
+         isecsbc = nsec_year + nsec1jan000 + NINT( (     REAL(      kit,wp) + zt_offset ) * rn_Dt / REAL(nn_e,wp) )
       ELSE                      ! middle of sbc time step
          ! note: we use kn_fsbc-1 because nsec_year is defined at the middle of the current time step
-         isecsbc = nsec_year + nsec1jan000 + NINT( ( 0.5*REAL(kn_fsbc-1,wp) + zt_offset ) * rdt )
+         isecsbc = nsec_year + nsec1jan000 + NINT( ( 0.5*REAL(kn_fsbc-1,wp) + zt_offset ) * rn_Dt )
       ENDIF
       imf = SIZE( sd )

NEMO/trunk/src/OCE/SBC/sbcapr.F90

@@ -36 +36 @@
    
    REAL(wp) ::   tarea                ! whole domain mean masked ocean surface
-   REAL(wp) ::   r1_grau              ! = 1.e0 / (grav * rau0)
+   REAL(wp) ::   r1_grau              ! = 1.e0 / (grav * rho0)
    
    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_apr   ! structure of input fields (file informations, fields read)
@@ -98 +98 @@
       ENDIF
       !
-      r1_grau = 1.e0 / (grav * rau0)               !* constant for optimization
+      r1_grau = 1.e0 / (grav * rho0)               !* constant for optimization
       !
       !                                            !* control check

NEMO/trunk/src/OCE/SBC/sbcblk.F90

@@ -278 +278 @@
             ENDIF
             !
-            IF( sf(jfpr)%freqh > 0. .AND. MOD( NINT(3600. * sf(jfpr)%freqh), nn_fsbc * NINT(rdt) ) /= 0 )   &
-               &  CALL ctl_warn( 'sbc_blk_init: sbcmod timestep rdt*nn_fsbc is NOT a submultiple of atmospheric forcing frequency.',   &
-               &                 '               This is not ideal. You should consider changing either rdt or nn_fsbc value...' )
+            IF( sf(jfpr)%freqh > 0. .AND. MOD( NINT(3600. * sf(jfpr)%freqh), nn_fsbc * NINT(rn_Dt) ) /= 0 )   &
+               &  CALL ctl_warn( 'sbc_blk_init: sbcmod timestep rn_Dt*nn_fsbc is NOT a submultiple of atmospheric forcing frequency.',   &
+               &                 '               This is not ideal. You should consider changing either rn_Dt or nn_fsbc value...' )
          ENDIF
       END DO

NEMO/trunk/src/OCE/SBC/sbcblk_skin_coare.F90

@@ -199 +199 @@
          ! Okay test on updated absorbed flux:
          !#LB: remove??? has a strong influence !!!
-         IF( (.NOT. l_exit).AND.(Qnt_ac(ji,jj) + zQabs*rdt <= 0._wp) ) THEN
+         IF( (.NOT. l_exit).AND.(Qnt_ac(ji,jj) + zQabs*rn_Dt <= 0._wp) ) THEN
             l_exit       = .TRUE.
             l_destroy_wl = .TRUE.
@@ -211 +211 @@
             ! 2/ Regardless of WL formed (dT==0 or dT>0), we are in the process to initiate one or warm further it !
 
-            ztac = Tau_ac(ji,jj) + MAX(.002_wp , pTau(ji,jj))*rdt      ! updated momentum integral
+            ztac = Tau_ac(ji,jj) + MAX(.002_wp , pTau(ji,jj))*rn_Dt      ! updated momentum integral
             !PRINT *, '#LBD: updated value for Tac=',  REAL(ztac,4)
 
@@ -218 +218 @@
             DO jl = 1, 5
                zQabs = frac_solar_abs(zHwl)*pQsw(ji,jj) + pQnsol(ji,jj)
-               zqac  = Qnt_ac(ji,jj) + zQabs*rdt ! updated heat absorbed
+               zqac  = Qnt_ac(ji,jj) + zQabs*rn_Dt ! updated heat absorbed
                IF( zqac <= 0._wp ) EXIT
                zHwl = MAX( MIN( Hwl_max , zcd1*ztac/SQRT(zqac)) , 0.1_wp ) ! Warm-layer depth

NEMO/trunk/src/OCE/SBC/sbcblk_skin_ecmwf.F90

@@ -214 +214 @@
          zcst2 = zcst1 / ( 5._wp*zHwl*zusw2 )  !OR: zcst2 = zcst1*rNuwl0 / ( 5._wp*zHwl*zusw2 ) ???
 
-         zcst0 = rdt * (rNuwl0 + 1._wp) / zHwl
+         zcst0 = rn_Dt * (rNuwl0 + 1._wp) / zHwl
 
          ZA = zcst0 * zQabs / ( rNuwl0 * zRhoCp_w )

NEMO/trunk/src/OCE/SBC/sbccpl.F90

@@ -193 +193 @@
 
    REAL(wp) ::   rpref = 101000._wp   ! reference atmospheric pressure[N/m2] 
-   REAL(wp) ::   r1_grau              ! = 1.e0 / (grav * rau0) 
+   REAL(wp) ::   r1_grau              ! = 1.e0 / (grav * rho0) 
 
    INTEGER , ALLOCATABLE, SAVE, DIMENSION(:) ::   nrcvinfo           ! OASIS info argument
@@ -1123 +1123 @@
       !                                                      ! Receive all the atmos. fields (including ice information)
       !                                                      ! ======================================================= !
-      isec = ( kt - nit000 ) * NINT( rdt )                      ! date of exchanges
+      isec = ( kt - nit000 ) * NINT( rn_Dt )                      ! date of exchanges
       DO jn = 1, jprcv                                          ! received fields sent by the atmosphere
          IF( srcv(jn)%laction )   CALL cpl_rcv( jn, isec, frcv(jn)%z3, xcplmask(:,:,1:nn_cplmodel), nrcvinfo(jn) )
@@ -1250 +1250 @@
           IF( kt /= nit000 )   ssh_ibb(:,:) = ssh_ib(:,:)    !* Swap of ssh_ib fields 
 
-          r1_grau = 1.e0 / (grav * rau0)               !* constant for optimization 
+          r1_grau = 1.e0 / (grav * rho0)               !* constant for optimization 
           ssh_ib(:,:) = - ( frcv(jpr_mslp)%z3(:,:,1) - rpref ) * r1_grau    ! equivalent ssh (inverse barometer) 
           apr   (:,:) =     frcv(jpr_mslp)%z3(:,:,1)                         !atmospheric pressure 
@@ -2100 +2100 @@
       !!----------------------------------------------------------------------
       !
-      isec = ( kt - nit000 ) * NINT( rdt )        ! date of exchanges
+      isec = ( kt - nit000 ) * NINT( rn_Dt )        ! date of exchanges
 
       zfr_l(:,:) = 1.- fr_i(:,:)

NEMO/trunk/src/OCE/SBC/sbcdcy.F90

@@ -86 +86 @@
       ! --------------
       ! When are we during the day (from 0 to 1)
-      zlo = ( REAL(nsec_day, wp) - 0.5_wp * rdt ) / rday
-      zup = zlo + ( REAL(nn_fsbc, wp)     * rdt ) / rday
+      zlo = ( REAL(nsec_day, wp) - 0.5_wp * rn_Dt ) / rday
+      zup = zlo + ( REAL(nn_fsbc, wp)     * rn_Dt ) / rday
       !
       IF( nday_qsr == -1 ) THEN       ! first time step only
@@ -251 +251 @@
          END_2D
          !
-         ztmp = rday / ( rdt * REAL(nn_fsbc, wp) )
+         ztmp = rday / ( rn_Dt * REAL(nn_fsbc, wp) )
          rscal(:,:) = rscal(:,:) * ztmp
          !

NEMO/trunk/src/OCE/SBC/sbcfwb.F90

@@ -126 +126 @@
          ENDIF   
          !                                         ! Update fwfold if new year start
-         ikty = 365 * 86400 / rdt                  !!bug  use of 365 days leap year or 360d year !!!!!!!
+         ikty = 365 * 86400 / rn_Dt                  !!bug  use of 365 days leap year or 360d year !!!!!!!
          IF( MOD( kt, ikty ) == 0 ) THEN
             a_fwb_b = a_fwb                           ! mean sea level taking into account the ice+snow
                                                       ! sum over the global domain
-            a_fwb   = glob_sum( 'sbcfwb', e1e2t(:,:) * ( ssh(:,:,Kmm) + snwice_mass(:,:) * r1_rau0 ) )
+            a_fwb   = glob_sum( 'sbcfwb', e1e2t(:,:) * ( ssh(:,:,Kmm) + snwice_mass(:,:) * r1_rho0 ) )
             a_fwb   = a_fwb * 1.e+3 / ( area * rday * 365. )     ! convert in Kg/m3/s = mm/s
 !!gm        !                                                      !!bug 365d year 

NEMO/trunk/src/OCE/SBC/sbcice_cice.F90

@@ -13 +13 @@
    USE dom_oce         ! ocean space and time domain
    USE domvvl
-   USE phycst, only : rcp, rau0, r1_rau0, rhos, rhoi
+   USE phycst, only : rcp, rho0, r1_rho0, rhos, rhoi
    USE in_out_manager  ! I/O manager
    USE iom, ONLY : iom_put,iom_use              ! I/O manager library !!Joakim edit
@@ -228 +228 @@
       IF( .NOT.ln_rstart ) THEN
          IF( ln_ice_embd ) THEN            ! embedded sea-ice: deplete the initial ssh below sea-ice area
-            ssh(:,:,Kmm) = ssh(:,:,Kmm) - snwice_mass(:,:) * r1_rau0
-            ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_rau0
+            ssh(:,:,Kmm) = ssh(:,:,Kmm) - snwice_mass(:,:) * r1_rho0
+            ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_rho0
 
 !!gm This should be put elsewhere....   (same remark for limsbc)
@@ -417 +417 @@
 ! Freezing/melting potential
 ! Calculated over NEMO leapfrog timestep (hence 2*dt)
-      nfrzmlt(:,:) = rau0 * rcp * e3t_m(:,:) * ( Tocnfrz-sst_m(:,:) ) / ( 2.0*dt )
+      nfrzmlt(:,:) = rho0 * rcp * e3t_m(:,:) * ( Tocnfrz-sst_m(:,:) ) / ( 2.0*dt )
 
       ztmp(:,:) = nfrzmlt(:,:)
@@ -450 +450 @@
          zintb = REAL( nn_fsbc + 1 ) / REAL( nn_fsbc ) * 0.5_wp
           !
-         zpice(:,:) = ssh_m(:,:) + (  zintn * snwice_mass(:,:) +  zintb * snwice_mass_b(:,:)  ) * r1_rau0
+         zpice(:,:) = ssh_m(:,:) + (  zintn * snwice_mass(:,:) +  zintb * snwice_mass_b(:,:)  ) * r1_rho0
           !
          !

NEMO/trunk/src/OCE/SBC/sbcmod.F90

@@ -187 +187 @@
       !
       IF( .NOT.ln_usr ) THEN     ! the model calendar needs some specificities (except in user defined case)
-         IF( MOD( rday , rdt ) /= 0. )   CALL ctl_stop( 'the time step must devide the number of second of in a day' )
+         IF( MOD( rday , rn_Dt ) /= 0. )   CALL ctl_stop( 'the time step must devide the number of second of in a day' )
          IF( MOD( rday , 2.  ) /= 0. )   CALL ctl_stop( 'the number of second of in a day must be an even number'    )
-         IF( MOD( rdt  , 2.  ) /= 0. )   CALL ctl_stop( 'the time step (in second) must be an even number'           )
+         IF( MOD( rn_Dt  , 2.  ) /= 0. )   CALL ctl_stop( 'the time step (in second) must be an even number'           )
       ENDIF
       !                       !**  check option consistency
@@ -309 +309 @@
       !     SAS time-step has to be declared in OASIS (mandatory) -> nn_fsbc has to be modified accordingly
       IF( nn_components /= jp_iam_nemo ) THEN
-         IF( nn_components == jp_iam_opa )   nn_fsbc = cpl_freq('O_SFLX') / NINT(rdt)
-         IF( nn_components == jp_iam_sas )   nn_fsbc = cpl_freq('I_SFLX') / NINT(rdt)
+         IF( nn_components == jp_iam_opa )   nn_fsbc = cpl_freq('O_SFLX') / NINT(rn_Dt)
+         IF( nn_components == jp_iam_sas )   nn_fsbc = cpl_freq('I_SFLX') / NINT(rn_Dt)
          !
          IF(lwp)THEN
@@ -331 +331 @@
       ENDIF
       !
-      IF( MOD( rday, REAL(nn_fsbc, wp) * rdt ) /= 0 )   &
+      IF( MOD( rday, REAL(nn_fsbc, wp) * rn_Dt ) /= 0 )   &
          &  CALL ctl_warn( 'sbc_init : nn_fsbc is NOT a multiple of the number of time steps in a day' )
       !
-      IF( ln_dm2dc .AND. NINT(rday) / ( nn_fsbc * NINT(rdt) ) < 8  )   &
+      IF( ln_dm2dc .AND. NINT(rday) / ( nn_fsbc * NINT(rn_Dt) ) < 8  )   &
          &   CALL ctl_warn( 'sbc_init : diurnal cycle for qsr: the sampling of the diurnal cycle is too small...' )
       !

NEMO/trunk/src/OCE/SBC/sbcrnf.F90

@@ -137 +137 @@
          !                                                           ! set temperature & salinity content of runoffs
          IF( ln_rnf_tem ) THEN                                       ! use runoffs temperature data
-            rnf_tsc(:,:,jp_tem) = ( sf_t_rnf(1)%fnow(:,:,1) ) * rnf(:,:) * r1_rau0
+            rnf_tsc(:,:,jp_tem) = ( sf_t_rnf(1)%fnow(:,:,1) ) * rnf(:,:) * r1_rho0
             CALL eos_fzp( sss_m(:,:), ztfrz(:,:) )
             WHERE( sf_t_rnf(1)%fnow(:,:,1) == -999._wp )             ! if missing data value use SST as runoffs temperature
-               rnf_tsc(:,:,jp_tem) = sst_m(:,:) * rnf(:,:) * r1_rau0
+               rnf_tsc(:,:,jp_tem) = sst_m(:,:) * rnf(:,:) * r1_rho0
             END WHERE
          ELSE                                                        ! use SST as runoffs temperature
             !CEOD River is fresh water so must at least be 0 unless we consider ice
-            rnf_tsc(:,:,jp_tem) = MAX( sst_m(:,:), 0.0_wp ) * rnf(:,:) * r1_rau0
+            rnf_tsc(:,:,jp_tem) = MAX( sst_m(:,:), 0.0_wp ) * rnf(:,:) * r1_rho0
          ENDIF
          !                                                           ! use runoffs salinity data
-         IF( ln_rnf_sal )   rnf_tsc(:,:,jp_sal) = ( sf_s_rnf(1)%fnow(:,:,1) ) * rnf(:,:) * r1_rau0
+         IF( ln_rnf_sal )   rnf_tsc(:,:,jp_sal) = ( sf_s_rnf(1)%fnow(:,:,1) ) * rnf(:,:) * r1_rho0
          !                                                           ! else use S=0 for runoffs (done one for all in the init)
                                          CALL iom_put( 'runoffs'     , rnf(:,:)                         )   ! output runoff mass flux
-         IF( iom_use('hflx_rnf_cea') )   CALL iom_put( 'hflx_rnf_cea', rnf_tsc(:,:,jp_tem) * rau0 * rcp )   ! output runoff sensible heat (W/m2)
+         IF( iom_use('hflx_rnf_cea') )   CALL iom_put( 'hflx_rnf_cea', rnf_tsc(:,:,jp_tem) * rho0 * rcp )   ! output runoff sensible heat (W/m2)
       ENDIF
       !
@@ -210 +210 @@
             DO_2D_11_11
                DO jk = 1, nk_rnf(ji,jj)
-                  phdivn(ji,jj,jk) = phdivn(ji,jj,jk) - ( rnf(ji,jj) + rnf_b(ji,jj) ) * zfact * r1_rau0 / h_rnf(ji,jj)
+                  phdivn(ji,jj,jk) = phdivn(ji,jj,jk) - ( rnf(ji,jj) + rnf_b(ji,jj) ) * zfact * r1_rho0 / h_rnf(ji,jj)
                END DO
             END_2D
@@ -221 +221 @@
                !                          ! apply the runoff input flow
                DO jk = 1, nk_rnf(ji,jj)
-                  phdivn(ji,jj,jk) = phdivn(ji,jj,jk) - ( rnf(ji,jj) + rnf_b(ji,jj) ) * zfact * r1_rau0 / h_rnf(ji,jj)
+                  phdivn(ji,jj,jk) = phdivn(ji,jj,jk) - ( rnf(ji,jj) + rnf_b(ji,jj) ) * zfact * r1_rho0 / h_rnf(ji,jj)
                END DO
             END_2D
@@ -227 +227 @@
       ELSE                       !==   runoff put only at the surface   ==!
          h_rnf (:,:)   = e3t (:,:,1,Kmm)        ! update h_rnf to be depth of top box
-         phdivn(:,:,1) = phdivn(:,:,1) - ( rnf(:,:) + rnf_b(:,:) ) * zfact * r1_rau0 / e3t(:,:,1,Kmm)
+         phdivn(:,:,1) = phdivn(:,:,1) - ( rnf(:,:) + rnf_b(:,:) ) * zfact * r1_rho0 / e3t(:,:,1,Kmm)
       ENDIF
       !