Changeset 12489 for NEMO/trunk/src
- Timestamp:
- 2020-02-28T16:55:11+01:00 (3 years ago)
- Location:
- NEMO/trunk/src
- Files:
-
- 147 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ABL/abl.F90
r11305 r12489 11 11 USE dom_oce, ONLY: e1t, e1u, e1v, e1f ! scale factors for horizontal grid 12 12 USE dom_oce, ONLY: e2t, e2u, e2v, e2f ! 13 USE dom_oce, ONLY: r dt! oceanic time-step13 USE dom_oce, ONLY: rn_Dt ! oceanic time-step 14 14 USE sbc_oce, ONLY: ght_abl, ghw_abl, e3t_abl, e3w_abl, jpka ! scale factors and altitudes of ABL grid points in the vertical 15 15 -
NEMO/trunk/src/ABL/ablmod.F90
r12353 r12489 152 152 DO jk = 3, jpkam1 153 153 DO ji = 1, jpi ! vector opt. 154 z_elem_a( ji, jk ) = - r dt_abl * Avt_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal155 z_elem_c( ji, jk ) = - r dt_abl * Avt_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal154 z_elem_a( ji, jk ) = - rDt_abl * Avt_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal 155 z_elem_c( ji, jk ) = - rDt_abl * Avt_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal 156 156 z_elem_b( ji, jk ) = e3t_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) ! diagonal 157 157 END DO … … 161 161 ! Neumann at the bottom 162 162 z_elem_a( ji, 2 ) = 0._wp 163 z_elem_c( ji, 2 ) = - r dt_abl * Avt_abl( ji, jj, 2 ) / e3w_abl( 2 )163 z_elem_c( ji, 2 ) = - rDt_abl * Avt_abl( ji, jj, 2 ) / e3w_abl( 2 ) 164 164 ! Homogeneous Neumann at the top 165 z_elem_a( ji, jpka ) = - r dt_abl * Avt_abl( ji, jj, jpka ) / e3w_abl( jpka )165 z_elem_a( ji, jpka ) = - rDt_abl * Avt_abl( ji, jj, jpka ) / e3w_abl( jpka ) 166 166 z_elem_c( ji, jpka ) = 0._wp 167 167 z_elem_b( ji, jpka ) = e3t_abl( jpka ) - z_elem_a( ji, jpka ) … … 184 184 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * psen_ice(ji,jj) * ptm_su(ji,jj) 185 185 #endif 186 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + r dt_abl * zztmp1187 tq_abl ( ji, jj, 2 , nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2 , nt_n, jtra ) + r dt_abl * zztmp2186 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rDt_abl * zztmp1 187 tq_abl ( ji, jj, 2 , nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2 , nt_n, jtra ) + rDt_abl * zztmp2 188 188 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) 189 189 END DO … … 196 196 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pevp_ice(ji, jj) * pssq_ice(ji, jj) 197 197 #endif 198 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + r dt_abl * zztmp1199 tq_abl ( ji, jj, 2 , nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2 , nt_n, jtra ) + r dt_abl * zztmp2198 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rDt_abl * zztmp1 199 tq_abl ( ji, jj, 2 , nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2 , nt_n, jtra ) + rDt_abl * zztmp2 200 200 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) 201 201 END DO … … 242 242 ! Advance u_abl & v_abl to time n+1 243 243 DO_2D_11_11 244 zcff = ( fft_abl(ji,jj) * r dt_abl )*( fft_abl(ji,jj) * rdt_abl ) ! (f dt)**2244 zcff = ( fft_abl(ji,jj) * rDt_abl )*( fft_abl(ji,jj) * rDt_abl ) ! (f dt)**2 245 245 246 246 u_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( & 247 247 & (1._wp-gamma_Cor*(1._wp-gamma_Cor)*zcff)*u_abl( ji, jj, jk, nt_n ) & 248 & + r dt_abl * fft_abl(ji, jj) * v_abl ( ji , jj , jk, nt_n ) ) &248 & + rDt_abl * fft_abl(ji, jj) * v_abl ( ji , jj , jk, nt_n ) ) & 249 249 & / (1._wp + gamma_Cor*gamma_Cor*zcff) 250 250 251 251 v_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( & 252 252 & (1._wp-gamma_Cor*(1._wp-gamma_Cor)*zcff)*v_abl( ji, jj, jk, nt_n ) & 253 & - r dt_abl * fft_abl(ji, jj) * u_abl ( ji , jj, jk, nt_n ) ) &253 & - rDt_abl * fft_abl(ji, jj) * u_abl ( ji , jj, jk, nt_n ) ) & 254 254 & / (1._wp + gamma_Cor*gamma_Cor*zcff) 255 255 END_2D … … 264 264 DO ji = 1, jpi 265 265 u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_a ) & 266 & - r dt_abl * e3t_abl(jk) * fft_abl(ji , jj) * pgv_dta(ji ,jj ,jk)266 & - rDt_abl * e3t_abl(jk) * fft_abl(ji , jj) * pgv_dta(ji ,jj ,jk) 267 267 v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_a ) & 268 & + r dt_abl * e3t_abl(jk) * fft_abl(ji, jj ) * pgu_dta(ji ,jj ,jk)268 & + rDt_abl * e3t_abl(jk) * fft_abl(ji, jj ) * pgu_dta(ji ,jj ,jk) 269 269 END DO 270 270 END DO … … 277 277 DO jk = 1, jpka 278 278 DO ji = 1, jpi 279 u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_a ) - r dt_abl * e3t_abl(jk) * pgu_dta(ji,jj,jk)280 v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_a ) - r dt_abl * e3t_abl(jk) * pgv_dta(ji,jj,jk)279 u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_a ) - rDt_abl * e3t_abl(jk) * pgu_dta(ji,jj,jk) 280 v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_a ) - rDt_abl * e3t_abl(jk) * pgv_dta(ji,jj,jk) 281 281 ENDDO 282 282 ENDDO … … 295 295 DO jk = 3, jpkam1 296 296 DO ji = 1, jpi 297 z_elem_a( ji, jk ) = - r dt_abl * Avm_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal298 z_elem_c( ji, jk ) = - r dt_abl * Avm_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal297 z_elem_a( ji, jk ) = - rDt_abl * Avm_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal 298 z_elem_c( ji, jk ) = - rDt_abl * Avm_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal 299 299 z_elem_b( ji, jk ) = e3t_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) ! diagonal 300 300 END DO … … 304 304 !++ Surface boundary condition 305 305 z_elem_a( ji, 2 ) = 0._wp 306 z_elem_c( ji, 2 ) = - r dt_abl * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 )306 z_elem_c( ji, 2 ) = - rDt_abl * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) 307 307 ! 308 308 zztmp1 = pcd_du(ji, jj) … … 313 313 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pcd_du_ice(ji, jj) * zzice 314 314 #endif 315 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + r dt_abl * zztmp1316 u_abl( ji, jj, 2, nt_a ) = u_abl( ji, jj, 2, nt_a ) + r dt_abl * zztmp2315 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rDt_abl * zztmp1 316 u_abl( ji, jj, 2, nt_a ) = u_abl( ji, jj, 2, nt_a ) + rDt_abl * zztmp2 317 317 318 318 !++ Top Neumann B.C. 319 !z_elem_a( ji, jpka ) = - 0.5_wp * r dt_abl * ( Avm_abl( ji, jj, jpka )+ Avm_abl( ji+1, jj, jpka ) ) / e3w_abl( jpka )319 !z_elem_a( ji, jpka ) = - 0.5_wp * rDt_abl * ( Avm_abl( ji, jj, jpka )+ Avm_abl( ji+1, jj, jpka ) ) / e3w_abl( jpka ) 320 320 !z_elem_c( ji, jpka ) = 0._wp 321 321 !z_elem_b( ji, jpka ) = e3t_abl( jpka ) - z_elem_a( ji, jpka ) … … 362 362 DO jk = 3, jpkam1 363 363 DO ji = 1, jpi 364 z_elem_a( ji, jk ) = -r dt_abl * Avm_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal365 z_elem_c( ji, jk ) = -r dt_abl * Avm_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal364 z_elem_a( ji, jk ) = -rDt_abl * Avm_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal 365 z_elem_c( ji, jk ) = -rDt_abl * Avm_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal 366 366 z_elem_b( ji, jk ) = e3t_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) ! diagonal 367 367 END DO … … 371 371 !++ Surface boundary condition 372 372 z_elem_a( ji, 2 ) = 0._wp 373 z_elem_c( ji, 2 ) = - r dt_abl * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 )373 z_elem_c( ji, 2 ) = - rDt_abl * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) 374 374 ! 375 375 zztmp1 = pcd_du(ji, jj) … … 380 380 zztmp2 = zztmp2 * pfrac_oce(ji,jj) + (1._wp - pfrac_oce(ji,jj)) * pcd_du_ice(ji, jj) * zzice 381 381 #endif 382 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + r dt_abl * zztmp1383 v_abl( ji, jj, 2, nt_a ) = v_abl( ji, jj, 2, nt_a ) + r dt_abl * zztmp2382 z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rDt_abl * zztmp1 383 v_abl( ji, jj, 2, nt_a ) = v_abl( ji, jj, 2, nt_a ) + rDt_abl * zztmp2 384 384 !++ Top Neumann B.C. 385 !z_elem_a( ji, jpka ) = -r dt_abl * Avm_abl( ji, jj, jpka ) / e3w_abl( jpka )385 !z_elem_a( ji, jpka ) = -rDt_abl * Avm_abl( ji, jj, jpka ) / e3w_abl( jpka ) 386 386 !z_elem_c( ji, jpka ) = 0._wp 387 387 !z_elem_b( ji, jpka ) = e3t_abl( jpka ) - z_elem_a( ji, jpka ) … … 436 436 zcff2 = jp_alp3_dyn * zsig**3 + jp_alp2_dyn * zsig**2 & 437 437 & + jp_alp1_dyn * zsig + jp_alp0_dyn 438 zcff = (1._wp-zmsk) + zmsk * zcff2 * r dt ! zcff = 1 for masked points439 ! rdt = rdt_abl / nn_fsbc438 zcff = (1._wp-zmsk) + zmsk * zcff2 * rn_Dt ! zcff = 1 for masked points 439 ! rn_Dt = rDt_abl / nn_fsbc 440 440 zcff = zcff * rest_eq(ji,jj) 441 441 z_cft( ji, jj, jk ) = zcff … … 460 460 zcff2 = jp_alp3_tra * zsig**3 + jp_alp2_tra * zsig**2 & 461 461 & + jp_alp1_tra * zsig + jp_alp0_tra 462 zcff = (1._wp-zmsk) + zmsk * zcff2 * r dt ! zcff = 1 for masked points463 ! rdt = rdt_abl / nn_fsbc462 zcff = (1._wp-zmsk) + zmsk * zcff2 * rn_Dt ! zcff = 1 for masked points 463 ! rn_Dt = rDt_abl / nn_fsbc 464 464 !z_cft( ji, jj, jk ) = zcff 465 465 tq_abl( ji, jj, jk, nt_a, jp_ta ) = (1._wp - zcff ) * tq_abl( ji, jj, jk, nt_a, jp_ta ) & … … 688 688 zbuoy = - Avt_abl( ji, jj, jk ) * zbn2( ji, jj, jk ) 689 689 690 z_elem_a( ji, jk ) = - 0.5_wp * r dt_abl * rn_Sch * ( Avm_abl( ji, jj, jk )+Avm_abl( ji, jj, jk-1 ) ) / e3t_abl( jk ) ! lower-diagonal691 z_elem_c( ji, jk ) = - 0.5_wp * r dt_abl * rn_Sch * ( Avm_abl( ji, jj, jk )+Avm_abl( ji, jj, jk+1 ) ) / e3t_abl( jk+1 ) ! upper-diagonal690 z_elem_a( ji, jk ) = - 0.5_wp * rDt_abl * rn_Sch * ( Avm_abl( ji, jj, jk )+Avm_abl( ji, jj, jk-1 ) ) / e3t_abl( jk ) ! lower-diagonal 691 z_elem_c( ji, jk ) = - 0.5_wp * rDt_abl * rn_Sch * ( Avm_abl( ji, jj, jk )+Avm_abl( ji, jj, jk+1 ) ) / e3t_abl( jk+1 ) ! upper-diagonal 692 692 IF( (zbuoy + zshear) .gt. 0.) THEN ! Patankar trick to avoid negative values of TKE 693 693 z_elem_b( ji, jk ) = e3w_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) & 694 & + e3w_abl(jk) * r dt_abl * rn_Ceps * sqrt(tke_abl( ji, jj, jk, nt_n )) / mxl_abl(ji,jj,jk) ! diagonal695 tke_abl( ji, jj, jk, nt_a ) = e3w_abl(jk) * ( tke_abl( ji, jj, jk, nt_n ) + r dt_abl * ( zbuoy + zshear ) ) ! right-hand-side694 & + e3w_abl(jk) * rDt_abl * rn_Ceps * sqrt(tke_abl( ji, jj, jk, nt_n )) / mxl_abl(ji,jj,jk) ! diagonal 695 tke_abl( ji, jj, jk, nt_a ) = e3w_abl(jk) * ( tke_abl( ji, jj, jk, nt_n ) + rDt_abl * ( zbuoy + zshear ) ) ! right-hand-side 696 696 ELSE 697 697 z_elem_b( ji, jk ) = e3w_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) & 698 & + e3w_abl(jk) * r dt_abl * rn_Ceps * sqrt(tke_abl( ji, jj, jk, nt_n )) / mxl_abl(ji,jj,jk) & ! diagonal699 & - e3w_abl(jk) * r dt_abl * zbuoy700 tke_abl( ji, jj, jk, nt_a ) = e3w_abl(jk) * ( tke_abl( ji, jj, jk, nt_n ) + r dt_abl * zshear ) ! right-hand-side698 & + e3w_abl(jk) * rDt_abl * rn_Ceps * sqrt(tke_abl( ji, jj, jk, nt_n )) / mxl_abl(ji,jj,jk) & ! diagonal 699 & - e3w_abl(jk) * rDt_abl * zbuoy 700 tke_abl( ji, jj, jk, nt_a ) = e3w_abl(jk) * ( tke_abl( ji, jj, jk, nt_n ) + rDt_abl * zshear ) ! right-hand-side 701 701 END IF 702 702 END DO -
NEMO/trunk/src/ABL/par_abl.F90
r11858 r12489 76 76 REAL(wp), PUBLIC, PARAMETER :: gamma_Cor = 0.55_wp 77 77 ! ABL timestep 78 REAL(wp), PUBLIC :: r dt_abl78 REAL(wp), PUBLIC :: rDt_abl 79 79 80 80 !!--------------------------------------------------------------------- -
NEMO/trunk/src/ABL/sbcabl.F90
r12353 r12489 202 202 203 203 ! ABL timestep 204 r dt_abl = nn_fsbc * rdt204 rDt_abl = nn_fsbc * rn_Dt 205 205 206 206 ! Check parameters for dynamics 207 207 zcff = ( jp_alp3_dyn * jp_bmin**3 + jp_alp2_dyn * jp_bmin**2 & 208 & + jp_alp1_dyn * jp_bmin + jp_alp0_dyn ) * r dt_abl208 & + jp_alp1_dyn * jp_bmin + jp_alp0_dyn ) * rDt_abl 209 209 zcff1 = ( jp_alp3_dyn * jp_bmax**3 + jp_alp2_dyn * jp_bmax**2 & 210 & + jp_alp1_dyn * jp_bmax + jp_alp0_dyn ) * r dt_abl210 & + jp_alp1_dyn * jp_bmax + jp_alp0_dyn ) * rDt_abl 211 211 IF(lwp) THEN 212 212 IF(nn_dyn_restore > 0) THEN … … 225 225 ! Check parameters for active tracers 226 226 zcff = ( jp_alp3_tra * jp_bmin**3 + jp_alp2_tra * jp_bmin**2 & 227 & + jp_alp1_tra * jp_bmin + jp_alp0_tra ) * r dt_abl227 & + jp_alp1_tra * jp_bmin + jp_alp0_tra ) * rDt_abl 228 228 zcff1 = ( jp_alp3_tra * jp_bmax**3 + jp_alp2_tra * jp_bmax**2 & 229 & + jp_alp1_tra * jp_bmax + jp_alp0_tra ) * r dt_abl229 & + jp_alp1_tra * jp_bmax + jp_alp0_tra ) * rDt_abl 230 230 IF(lwp) THEN 231 231 WRITE(numout,*) ' ABL Minimum value for tracers restoring = ',zcff -
NEMO/trunk/src/ICE/ice.F90
r11627 r12489 150 150 REAL(wp), PUBLIC :: rn_ecc !: eccentricity of the elliptical yield curve 151 151 INTEGER , PUBLIC :: nn_nevp !: number of iterations for subcycling 152 REAL(wp), PUBLIC :: rn_relast !: ratio => telast/r dt_ice (1/3 or 1/9 depending on nb of subcycling nevp)152 REAL(wp), PUBLIC :: rn_relast !: ratio => telast/rDt_ice (1/3 or 1/9 depending on nb of subcycling nevp) 153 153 ! 154 154 ! !!** ice-advection namelist (namdyn_adv) ** … … 207 207 ! !!** some other parameters 208 208 INTEGER , PUBLIC :: kt_ice !: iteration number 209 REAL(wp), PUBLIC :: r dt_ice !: ice time step210 REAL(wp), PUBLIC :: r1_ rdtice !: = 1. / rdt_ice209 REAL(wp), PUBLIC :: rDt_ice !: ice time step 210 REAL(wp), PUBLIC :: r1_Dt_ice !: = 1. / rDt_ice 211 211 REAL(wp), PUBLIC :: r1_nlay_i !: 1 / nlay_i 212 212 REAL(wp), PUBLIC :: r1_nlay_s !: 1 / nlay_s -
NEMO/trunk/src/ICE/icecor.F90
r12377 r12489 86 86 IF ( nn_icesal == 2 ) THEN ! salinity must stay in bounds [Simin,Simax] ! 87 87 ! !----------------------------------------------------- 88 zzc = rhoi * r1_ rdtice88 zzc = rhoi * r1_Dt_ice 89 89 DO jl = 1, jpl 90 90 DO_2D_11_11 … … 123 123 ! 124 124 IF( ln_icediachk .OR. iom_use('hfxdhc') ) THEN 125 diag_heat(:,:) = - SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_ rdtice & ! W.m-2126 & - SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_ rdtice127 diag_sice(:,:) = SUM( sv_i(:,:,:) - sv_i_b(:,:,:) , dim=3 ) * r1_ rdtice * rhoi128 diag_vice(:,:) = SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_ rdtice * rhoi129 diag_vsnw(:,:) = SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_ rdtice * rhos125 diag_heat(:,:) = - SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_Dt_ice & ! W.m-2 126 & - SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_Dt_ice 127 diag_sice(:,:) = SUM( sv_i(:,:,:) - sv_i_b(:,:,:) , dim=3 ) * r1_Dt_ice * rhoi 128 diag_vice(:,:) = SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_Dt_ice * rhoi 129 diag_vsnw(:,:) = SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_Dt_ice * rhos 130 130 ENDIF 131 131 ! ! concentration tendency (dynamics) 132 132 IF( iom_use('afxdyn') .OR. iom_use('afxthd') .OR. iom_use('afxtot') ) THEN 133 zafx(:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_ rdtice133 zafx(:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_Dt_ice 134 134 CALL iom_put( 'afxdyn' , zafx ) 135 135 ENDIF … … 137 137 CASE( 2 ) !--- thermo trend diagnostics & ice aging 138 138 ! 139 oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * r dt_ice ! ice natural aging incrementation139 oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * rDt_ice ! ice natural aging incrementation 140 140 ! 141 141 IF( ln_icediachk .OR. iom_use('hfxdhc') ) THEN 142 142 diag_heat(:,:) = diag_heat(:,:) & 143 & - SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_ rdtice &144 & - SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_ rdtice143 & - SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_Dt_ice & 144 & - SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_Dt_ice 145 145 diag_sice(:,:) = diag_sice(:,:) & 146 & + SUM( sv_i(:,:,:) - sv_i_b(:,:,:) , dim=3 ) * r1_ rdtice * rhoi146 & + SUM( sv_i(:,:,:) - sv_i_b(:,:,:) , dim=3 ) * r1_Dt_ice * rhoi 147 147 diag_vice(:,:) = diag_vice(:,:) & 148 & + SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_ rdtice * rhoi148 & + SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_Dt_ice * rhoi 149 149 diag_vsnw(:,:) = diag_vsnw(:,:) & 150 & + SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_ rdtice * rhos150 & + SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_Dt_ice * rhos 151 151 CALL iom_put ( 'hfxdhc' , diag_heat ) 152 152 ENDIF 153 153 ! ! concentration tendency (total + thermo) 154 154 IF( iom_use('afxdyn') .OR. iom_use('afxthd') .OR. iom_use('afxtot') ) THEN 155 zafx(:,:) = zafx(:,:) + SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_ rdtice156 CALL iom_put( 'afxthd' , SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_ rdtice )155 zafx(:,:) = zafx(:,:) + SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_Dt_ice 156 CALL iom_put( 'afxthd' , SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_Dt_ice ) 157 157 CALL iom_put( 'afxtot' , zafx ) 158 158 ENDIF -
NEMO/trunk/src/ICE/icectl.F90
r12377 r12489 104 104 105 105 ! -- mass diag -- ! 106 zdiag_mass = ( glob_sum( 'icectl', SUM( v_i * rhoi + v_s * rhos, dim=3 ) * e1e2t ) - pdiag_v ) * r1_ rdtice &106 zdiag_mass = ( glob_sum( 'icectl', SUM( v_i * rhoi + v_s * rhos, dim=3 ) * e1e2t ) - pdiag_v ) * r1_Dt_ice & 107 107 & + glob_sum( 'icectl', ( wfx_bog + wfx_bom + wfx_sum + wfx_sni + wfx_opw + wfx_res + wfx_dyn + & 108 108 & wfx_lam + wfx_pnd + wfx_snw_sni + wfx_snw_sum + wfx_snw_dyn + wfx_snw_sub + & … … 111 111 ! 112 112 ! -- salt diag -- ! 113 zdiag_salt = ( glob_sum( 'icectl', SUM( sv_i * rhoi , dim=3 ) * e1e2t ) - pdiag_s ) * r1_ rdtice &113 zdiag_salt = ( glob_sum( 'icectl', SUM( sv_i * rhoi , dim=3 ) * e1e2t ) - pdiag_s ) * r1_Dt_ice & 114 114 & + glob_sum( 'icectl', ( sfx_bri + sfx_bog + sfx_bom + sfx_sum + sfx_sni + & 115 115 & sfx_opw + sfx_res + sfx_dyn + sfx_sub + sfx_lam ) * e1e2t ) & … … 118 118 ! -- heat diag -- ! 119 119 zdiag_heat = ( glob_sum( 'icectl', ( SUM(SUM(e_i, dim=4), dim=3) + SUM(SUM(e_s, dim=4), dim=3) ) * e1e2t ) - pdiag_t & 120 & ) * r1_ rdtice &120 & ) * r1_Dt_ice & 121 121 & + glob_sum( 'icectl', ( hfx_sum + hfx_bom + hfx_bog + hfx_dif + hfx_opw + hfx_snw & 122 122 & - hfx_thd - hfx_dyn - hfx_res - hfx_sub - hfx_spr ) * e1e2t ) & … … 141 141 ! check conservation issues 142 142 IF( ABS(zdiag_mass) > zchk_m * rn_icechk_glo * zarea ) & 143 & WRITE(numout,*) cd_routine,' : violation mass cons. [kg] = ',zdiag_mass * r dt_ice143 & WRITE(numout,*) cd_routine,' : violation mass cons. [kg] = ',zdiag_mass * rDt_ice 144 144 IF( ABS(zdiag_salt) > zchk_s * rn_icechk_glo * zarea ) & 145 & WRITE(numout,*) cd_routine,' : violation salt cons. [g] = ',zdiag_salt * r dt_ice145 & WRITE(numout,*) cd_routine,' : violation salt cons. [g] = ',zdiag_salt * rDt_ice 146 146 IF( ABS(zdiag_heat) > zchk_t * rn_icechk_glo * zarea ) & 147 & WRITE(numout,*) cd_routine,' : violation heat cons. [J] = ',zdiag_heat * r dt_ice147 & WRITE(numout,*) cd_routine,' : violation heat cons. [J] = ',zdiag_heat * rDt_ice 148 148 ! check negative values 149 149 IF( zdiag_vmin < 0. ) WRITE(numout,*) cd_routine,' : violation v_i < 0 = ',zdiag_vmin … … 160 160 ! it does not mean UM is not conservative (it is checked with above prints) => update (09/2019): same for Prather now 161 161 !IF( ln_adv_Pra .AND. ABS(zvtrp) > zchk_m * rn_icechk_glo * zarea .AND. cd_routine == 'icedyn_adv' ) & 162 ! & WRITE(numout,*) cd_routine,' : violation adv scheme [kg] = ',zvtrp * r dt_ice162 ! & WRITE(numout,*) cd_routine,' : violation adv scheme [kg] = ',zvtrp * rDt_ice 163 163 ENDIF 164 164 ! … … 201 201 IF( lwp ) THEN 202 202 IF( ABS(zdiag_mass) > zchk_m * rn_icechk_glo * zarea ) & 203 & WRITE(numout,*) cd_routine,' : violation mass cons. [kg] = ',zdiag_mass * r dt_ice203 & WRITE(numout,*) cd_routine,' : violation mass cons. [kg] = ',zdiag_mass * rDt_ice 204 204 IF( ABS(zdiag_salt) > zchk_s * rn_icechk_glo * zarea ) & 205 & WRITE(numout,*) cd_routine,' : violation salt cons. [g] = ',zdiag_salt * r dt_ice206 !!IF( ABS(zdiag_heat) > zchk_t * rn_icechk_glo * zarea ) WRITE(numout,*) cd_routine,' : violation heat cons. [J] = ',zdiag_heat * r dt_ice205 & WRITE(numout,*) cd_routine,' : violation salt cons. [g] = ',zdiag_salt * rDt_ice 206 !!IF( ABS(zdiag_heat) > zchk_t * rn_icechk_glo * zarea ) WRITE(numout,*) cd_routine,' : violation heat cons. [J] = ',zdiag_heat * rDt_ice 207 207 ENDIF 208 208 ! … … 250 250 251 251 ! -- mass diag -- ! 252 zdiag_mass = ( SUM( v_i * rhoi + v_s * rhos, dim=3 ) - pdiag_v ) * r1_ rdtice &252 zdiag_mass = ( SUM( v_i * rhoi + v_s * rhos, dim=3 ) - pdiag_v ) * r1_Dt_ice & 253 253 & + ( wfx_bog + wfx_bom + wfx_sum + wfx_sni + wfx_opw + wfx_res + wfx_dyn + wfx_lam + wfx_pnd + & 254 254 & wfx_snw_sni + wfx_snw_sum + wfx_snw_dyn + wfx_snw_sub + wfx_ice_sub + wfx_spr ) & … … 257 257 ! 258 258 ! -- salt diag -- ! 259 zdiag_salt = ( SUM( sv_i * rhoi , dim=3 ) - pdiag_s ) * r1_ rdtice &259 zdiag_salt = ( SUM( sv_i * rhoi , dim=3 ) - pdiag_s ) * r1_Dt_ice & 260 260 & + ( sfx_bri + sfx_bog + sfx_bom + sfx_sum + sfx_sni + sfx_opw + sfx_res + sfx_dyn + sfx_sub + sfx_lam ) & 261 261 & - pdiag_fs … … 263 263 ! 264 264 ! -- heat diag -- ! 265 zdiag_heat = ( SUM( SUM( e_i, dim=4 ), dim=3 ) + SUM( SUM( e_s, dim=4 ), dim=3 ) - pdiag_t ) * r1_ rdtice &265 zdiag_heat = ( SUM( SUM( e_i, dim=4 ), dim=3 ) + SUM( SUM( e_s, dim=4 ), dim=3 ) - pdiag_t ) * r1_Dt_ice & 266 266 & + ( hfx_sum + hfx_bom + hfx_bog + hfx_dif + hfx_opw + hfx_snw & 267 267 & - hfx_thd - hfx_dyn - hfx_res - hfx_sub - hfx_spr ) & … … 455 455 DO jl = 1, jpl 456 456 DO_2D_11_11 457 IF ( ( ( ABS( o_i(ji,jj,jl) ) > r dt_ice ) .OR. &457 IF ( ( ( ABS( o_i(ji,jj,jl) ) > rDt_ice ) .OR. & 458 458 ( ABS( o_i(ji,jj,jl) ) < 0._wp) ) .AND. & 459 459 ( a_i(ji,jj,jl) > 0._wp ) ) THEN … … 651 651 WRITE(numout,*) ' hfx_res : ', hfx_res(ji,jj) 652 652 WRITE(numout,*) ' qsb_ice_bot : ', qsb_ice_bot(ji,jj) 653 WRITE(numout,*) ' qlead : ', qlead(ji,jj) * r1_ rdtice653 WRITE(numout,*) ' qlead : ', qlead(ji,jj) * r1_Dt_ice 654 654 WRITE(numout,*) 655 655 WRITE(numout,*) ' - Salt fluxes at bottom interface ***' -
NEMO/trunk/src/ICE/icedia.F90
r12377 r12489 109 109 ! ---------------------------! 110 110 ! they must be kept outside an IF(iom_use) because of the call to dia_rst below 111 z_frc_volbot = r1_r au0 * glob_sum( 'icedia', -( wfx_ice(:,:) + wfx_snw(:,:) + wfx_err_sub(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater flux ice/snow-ocean112 z_frc_voltop = r1_r au0 * glob_sum( 'icedia', -( wfx_sub(:,:) + wfx_spr(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater flux ice/snow-atm113 z_frc_sal = r1_r au0 * glob_sum( 'icedia', - sfx(:,:) * e1e2t(:,:) ) * 1.e-9 ! salt fluxes ice/snow-ocean111 z_frc_volbot = r1_rho0 * glob_sum( 'icedia', -( wfx_ice(:,:) + wfx_snw(:,:) + wfx_err_sub(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater flux ice/snow-ocean 112 z_frc_voltop = r1_rho0 * glob_sum( 'icedia', -( wfx_sub(:,:) + wfx_spr(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater flux ice/snow-atm 113 z_frc_sal = r1_rho0 * glob_sum( 'icedia', - sfx(:,:) * e1e2t(:,:) ) * 1.e-9 ! salt fluxes ice/snow-ocean 114 114 z_frc_tembot = glob_sum( 'icedia', qt_oce_ai(:,:) * e1e2t(:,:) ) * 1.e-20 ! heat on top of ocean (and below ice) 115 115 z_frc_temtop = glob_sum( 'icedia', qt_atm_oi(:,:) * e1e2t(:,:) ) * 1.e-20 ! heat on top of ice-coean 116 116 ! 117 frc_voltop = frc_voltop + z_frc_voltop * r dt_ice ! km3118 frc_volbot = frc_volbot + z_frc_volbot * r dt_ice ! km3119 frc_sal = frc_sal + z_frc_sal * r dt_ice ! km3*pss120 frc_temtop = frc_temtop + z_frc_temtop * r dt_ice ! 1.e20 J121 frc_tembot = frc_tembot + z_frc_tembot * r dt_ice ! 1.e20 J117 frc_voltop = frc_voltop + z_frc_voltop * rDt_ice ! km3 118 frc_volbot = frc_volbot + z_frc_volbot * rDt_ice ! km3 119 frc_sal = frc_sal + z_frc_sal * rDt_ice ! km3*pss 120 frc_temtop = frc_temtop + z_frc_temtop * rDt_ice ! 1.e20 J 121 frc_tembot = frc_tembot + z_frc_tembot * rDt_ice ! 1.e20 J 122 122 123 123 CALL iom_put( 'ibgfrcvoltop' , frc_voltop ) ! vol forcing ice/snw-atm (km3 equivalent ocean water) … … 128 128 129 129 IF( iom_use('ibgfrchfxtop') .OR. iom_use('ibgfrchfxbot') ) THEN 130 CALL iom_put( 'ibgfrchfxtop' , frc_temtop * z1_e1e2 * 1.e-20 * kt*r dt ) ! heat on top of ice/snw/ocean (W/m2)131 CALL iom_put( 'ibgfrchfxbot' , frc_tembot * z1_e1e2 * 1.e-20 * kt*r dt ) ! heat on top of ocean(below ice) (W/m2)130 CALL iom_put( 'ibgfrchfxtop' , frc_temtop * z1_e1e2 * 1.e-20 * kt*rn_Dt ) ! heat on top of ice/snw/ocean (W/m2) 131 CALL iom_put( 'ibgfrchfxbot' , frc_tembot * z1_e1e2 * 1.e-20 * kt*rn_Dt ) ! heat on top of ocean(below ice) (W/m2) 132 132 ENDIF 133 133 … … 137 137 IF( iom_use('ibgvolume') .OR. iom_use('ibgsaltco') .OR. iom_use('ibgheatco') .OR. iom_use('ibgheatfx') ) THEN 138 138 139 zdiff_vol = r1_r au0 * glob_sum( 'icedia', ( rhoi*vt_i(:,:) + rhos*vt_s(:,:) - vol_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater trend (km3)140 zdiff_sal = r1_r au0 * glob_sum( 'icedia', ( rhoi*st_i(:,:) - sal_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! salt content trend (km3*pss)139 zdiff_vol = r1_rho0 * glob_sum( 'icedia', ( rhoi*vt_i(:,:) + rhos*vt_s(:,:) - vol_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! freshwater trend (km3) 140 zdiff_sal = r1_rho0 * glob_sum( 'icedia', ( rhoi*st_i(:,:) - sal_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-9 ! salt content trend (km3*pss) 141 141 zdiff_tem = glob_sum( 'icedia', ( et_i(:,:) + et_s(:,:) - tem_loc_ini(:,:) ) * e1e2t(:,:) ) * 1.e-20 ! heat content trend (1.e20 J) 142 142 ! + SUM( qevap_ice * a_i_b, dim=3 ) !! clem: I think this term should not be there (but needs a check) -
NEMO/trunk/src/ICE/icedyn_adv.F90
r12377 r12489 93 93 ! diagnostics 94 94 !------------ 95 diag_trp_ei(:,:) = SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_ rdtice96 diag_trp_es(:,:) = SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_ rdtice97 diag_trp_sv(:,:) = SUM( sv_i(:,:,:) - sv_i_b(:,:,:) , dim=3 ) * r1_ rdtice98 diag_trp_vi(:,:) = SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_ rdtice99 diag_trp_vs(:,:) = SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_ rdtice95 diag_trp_ei(:,:) = SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_Dt_ice 96 diag_trp_es(:,:) = SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_Dt_ice 97 diag_trp_sv(:,:) = SUM( sv_i(:,:,:) - sv_i_b(:,:,:) , dim=3 ) * r1_Dt_ice 98 diag_trp_vi(:,:) = SUM( v_i (:,:,:) - v_i_b (:,:,:) , dim=3 ) * r1_Dt_ice 99 diag_trp_vs(:,:) = SUM( v_s (:,:,:) - v_s_b (:,:,:) , dim=3 ) * r1_Dt_ice 100 100 IF( iom_use('icemtrp') ) CALL iom_put( 'icemtrp' , diag_trp_vi * rhoi ) ! ice mass transport 101 101 IF( iom_use('snwmtrp') ) CALL iom_put( 'snwmtrp' , diag_trp_vs * rhos ) ! snw mass transport -
NEMO/trunk/src/ICE/icedyn_adv_pra.F90
r12377 r12489 122 122 ! Note: the advection split is applied at the next time-step in order to avoid blocking global comm. 123 123 ! this should not affect too much the stability 124 zcflnow(1) = MAXVAL( ABS( pu_ice(:,:) ) * r dt_ice * r1_e1u(:,:) )125 zcflnow(1) = MAX( zcflnow(1), MAXVAL( ABS( pv_ice(:,:) ) * r dt_ice * r1_e2v(:,:) ) )124 zcflnow(1) = MAXVAL( ABS( pu_ice(:,:) ) * rDt_ice * r1_e1u(:,:) ) 125 zcflnow(1) = MAX( zcflnow(1), MAXVAL( ABS( pv_ice(:,:) ) * rDt_ice * r1_e2v(:,:) ) ) 126 126 127 127 ! non-blocking global communication send zcflnow and receive zcflprv … … 131 131 ELSE ; icycle = 1 132 132 ENDIF 133 zdt = r dt_ice / REAL(icycle)133 zdt = rDt_ice / REAL(icycle) 134 134 135 135 ! --- transport --- ! … … 687 687 IF ( pv_i(ji,jj,jl) > 0._wp ) THEN 688 688 ! 689 zvs_excess = MAX( 0._wp, pv_s(ji,jj,jl) - pv_i(ji,jj,jl) * (r au0-rhoi) * r1_rhos )689 zvs_excess = MAX( 0._wp, pv_s(ji,jj,jl) - pv_i(ji,jj,jl) * (rho0-rhoi) * r1_rhos ) 690 690 ! 691 691 IF( zvs_excess > 0._wp ) THEN ! snow-ice interface deplets below the ocean surface -
NEMO/trunk/src/ICE/icedyn_adv_umx.F90
r12377 r12489 128 128 ! Note: the advection split is applied at the next time-step in order to avoid blocking global comm. 129 129 ! this should not affect too much the stability 130 zcflnow(1) = MAXVAL( ABS( pu_ice(:,:) ) * r dt_ice * r1_e1u(:,:) )131 zcflnow(1) = MAX( zcflnow(1), MAXVAL( ABS( pv_ice(:,:) ) * r dt_ice * r1_e2v(:,:) ) )130 zcflnow(1) = MAXVAL( ABS( pu_ice(:,:) ) * rDt_ice * r1_e1u(:,:) ) 131 zcflnow(1) = MAX( zcflnow(1), MAXVAL( ABS( pv_ice(:,:) ) * rDt_ice * r1_e2v(:,:) ) ) 132 132 133 133 ! non-blocking global communication send zcflnow and receive zcflprv … … 137 137 ELSE ; icycle = 1 138 138 ENDIF 139 zdt = r dt_ice / REAL(icycle)139 zdt = rDt_ice / REAL(icycle) 140 140 141 141 ! --- transport --- ! … … 1505 1505 IF ( pv_i(ji,jj,jl) > 0._wp ) THEN 1506 1506 ! 1507 zvs_excess = MAX( 0._wp, pv_s(ji,jj,jl) - pv_i(ji,jj,jl) * (r au0-rhoi) * r1_rhos )1507 zvs_excess = MAX( 0._wp, pv_s(ji,jj,jl) - pv_i(ji,jj,jl) * (rho0-rhoi) * r1_rhos ) 1508 1508 ! 1509 1509 IF( zvs_excess > 0._wp ) THEN ! snow-ice interface deplets below the ocean surface -
NEMO/trunk/src/ICE/icedyn_rdgrft.F90
r12377 r12489 250 250 ELSE 251 251 iterate_ridging = 1 252 zdivu (ji) = zfac * r1_ rdtice252 zdivu (ji) = zfac * r1_Dt_ice 253 253 closing_net(ji) = MAX( 0._wp, -zdivu(ji) ) 254 254 opning (ji) = MAX( 0._wp, zdivu(ji) ) … … 455 455 DO jl = 1, jpl 456 456 DO ji = 1, npti 457 zfac = apartf(ji,jl) * closing_gross(ji) * r dt_ice457 zfac = apartf(ji,jl) * closing_gross(ji) * rDt_ice 458 458 IF( zfac > pa_i(ji,jl) .AND. apartf(ji,jl) /= 0._wp ) THEN 459 closing_gross(ji) = pa_i(ji,jl) / apartf(ji,jl) * r1_ rdtice459 closing_gross(ji) = pa_i(ji,jl) / apartf(ji,jl) * r1_Dt_ice 460 460 ENDIF 461 461 END DO … … 467 467 ! Reduce the opening rate in proportion 468 468 DO ji = 1, npti 469 zfac = pato_i(ji) + ( opning(ji) - apartf(ji,0) * closing_gross(ji) ) * r dt_ice469 zfac = pato_i(ji) + ( opning(ji) - apartf(ji,0) * closing_gross(ji) ) * rDt_ice 470 470 IF( zfac < 0._wp ) THEN ! would lead to negative ato_i 471 opning(ji) = apartf(ji,0) * closing_gross(ji) - pato_i(ji) * r1_ rdtice471 opning(ji) = apartf(ji,0) * closing_gross(ji) - pato_i(ji) * r1_Dt_ice 472 472 ELSEIF( zfac > zasum(ji) ) THEN ! would lead to ato_i > asum 473 opning(ji) = apartf(ji,0) * closing_gross(ji) + ( zasum(ji) - pato_i(ji) ) * r1_ rdtice473 opning(ji) = apartf(ji,0) * closing_gross(ji) + ( zasum(ji) - pato_i(ji) ) * r1_Dt_ice 474 474 ENDIF 475 475 END DO … … 515 515 !-------------------------------------------------------- 516 516 DO ji = 1, npti 517 ato_i_1d(ji) = MAX( 0._wp, ato_i_1d(ji) + ( opning(ji) - apartf(ji,0) * closing_gross(ji) ) * r dt_ice )517 ato_i_1d(ji) = MAX( 0._wp, ato_i_1d(ji) + ( opning(ji) - apartf(ji,0) * closing_gross(ji) ) * rDt_ice ) 518 518 END DO 519 519 … … 533 533 534 534 ! area of ridging / rafting ice (airdg1) and of new ridge (airdg2) 535 airdg1 = aridge(ji,jl1) * closing_gross(ji) * r dt_ice536 airft1 = araft (ji,jl1) * closing_gross(ji) * r dt_ice535 airdg1 = aridge(ji,jl1) * closing_gross(ji) * rDt_ice 536 airft1 = araft (ji,jl1) * closing_gross(ji) * rDt_ice 537 537 538 538 airdg2(ji) = airdg1 * hi_hrdg(ji,jl1) … … 575 575 576 576 ! Ice-ocean exchanges associated with ice porosity 577 wfx_dyn_1d(ji) = wfx_dyn_1d(ji) - vsw * rhoi * r1_ rdtice ! increase in ice volume due to seawater frozen in voids578 sfx_dyn_1d(ji) = sfx_dyn_1d(ji) - vsw * sss_1d(ji) * rhoi * r1_ rdtice579 hfx_dyn_1d(ji) = hfx_dyn_1d(ji) + ersw(ji) * r1_ rdtice ! > 0 [W.m-2]577 wfx_dyn_1d(ji) = wfx_dyn_1d(ji) - vsw * rhoi * r1_Dt_ice ! increase in ice volume due to seawater frozen in voids 578 sfx_dyn_1d(ji) = sfx_dyn_1d(ji) - vsw * sss_1d(ji) * rhoi * r1_Dt_ice 579 hfx_dyn_1d(ji) = hfx_dyn_1d(ji) + ersw(ji) * r1_Dt_ice ! > 0 [W.m-2] 580 580 581 581 ! Put the snow lost by ridging into the ocean 582 582 ! Note that esrdg > 0; the ocean must cool to melt snow. If the ocean temp = Tf already, new ice must grow. 583 583 wfx_snw_dyn_1d(ji) = wfx_snw_dyn_1d(ji) + ( rhos * vsrdg(ji) * ( 1._wp - rn_fsnwrdg ) & ! fresh water source for ocean 584 & + rhos * vsrft(ji) * ( 1._wp - rn_fsnwrft ) ) * r1_ rdtice584 & + rhos * vsrft(ji) * ( 1._wp - rn_fsnwrft ) ) * r1_Dt_ice 585 585 586 586 ! virtual salt flux to keep salinity constant 587 587 IF( nn_icesal /= 2 ) THEN 588 588 sirdg2(ji) = sirdg2(ji) - vsw * ( sss_1d(ji) - s_i_1d(ji) ) ! ridge salinity = s_i 589 sfx_bri_1d(ji) = sfx_bri_1d(ji) + sss_1d(ji) * vsw * rhoi * r1_ rdtice & ! put back sss_m into the ocean590 & - s_i_1d(ji) * vsw * rhoi * r1_ rdtice ! and get s_i from the ocean589 sfx_bri_1d(ji) = sfx_bri_1d(ji) + sss_1d(ji) * vsw * rhoi * r1_Dt_ice & ! put back sss_m into the ocean 590 & - s_i_1d(ji) * vsw * rhoi * r1_Dt_ice ! and get s_i from the ocean 591 591 ENDIF 592 592 … … 611 611 IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 612 612 ! Compute ridging /rafting fractions 613 afrdg = aridge(ji,jl1) * closing_gross(ji) * r dt_ice * z1_ai(ji)614 afrft = araft (ji,jl1) * closing_gross(ji) * r dt_ice * z1_ai(ji)613 afrdg = aridge(ji,jl1) * closing_gross(ji) * rDt_ice * z1_ai(ji) 614 afrft = araft (ji,jl1) * closing_gross(ji) * rDt_ice * z1_ai(ji) 615 615 ! Compute ridging /rafting ice and new ridges for es 616 616 esrdg(ji,jk) = ze_s_2d (ji,jk,jl1) * afrdg … … 618 618 ! Put the snow lost by ridging into the ocean 619 619 hfx_dyn_1d(ji) = hfx_dyn_1d(ji) + ( - esrdg(ji,jk) * ( 1._wp - rn_fsnwrdg ) & ! heat sink for ocean (<0, W.m-2) 620 & - esrft(ji,jk) * ( 1._wp - rn_fsnwrft ) ) * r1_ rdtice620 & - esrft(ji,jk) * ( 1._wp - rn_fsnwrft ) ) * r1_Dt_ice 621 621 ! 622 622 ! Remove energy of new ridge to each category jl1 … … 632 632 IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 633 633 ! Compute ridging /rafting fractions 634 afrdg = aridge(ji,jl1) * closing_gross(ji) * r dt_ice * z1_ai(ji)635 afrft = araft (ji,jl1) * closing_gross(ji) * r dt_ice * z1_ai(ji)634 afrdg = aridge(ji,jl1) * closing_gross(ji) * rDt_ice * z1_ai(ji) 635 afrft = araft (ji,jl1) * closing_gross(ji) * rDt_ice * z1_ai(ji) 636 636 ! Compute ridging ice and new ridges for ei 637 637 eirdg(ji,jk) = ze_i_2d (ji,jk,jl1) * afrdg + ersw(ji) * r1_nlay_i -
NEMO/trunk/src/ICE/icedyn_rhg_evp.F90
r12377 r12489 116 116 INTEGER :: jter ! local integers 117 117 ! 118 REAL(wp) :: zrhoco ! r au0 * rn_cio118 REAL(wp) :: zrhoco ! rho0 * rn_cio 119 119 REAL(wp) :: zdtevp, z1_dtevp ! time step for subcycling 120 120 REAL(wp) :: ecc2, z1_ecc2 ! square of yield ellipse eccenticity … … 213 213 ! 1) define some variables and initialize arrays 214 214 !------------------------------------------------------------------------------! 215 zrhoco = r au0 * rn_cio215 zrhoco = rho0 * rn_cio 216 216 217 217 ! ecc2: square of yield ellipse eccenticrity … … 220 220 221 221 ! Time step for subcycling 222 zdtevp = r dt_ice / REAL( nn_nevp )222 zdtevp = rDt_ice / REAL( nn_nevp ) 223 223 z1_dtevp = 1._wp / zdtevp 224 224 225 225 ! alpha parameters (Bouillon 2009) 226 226 IF( .NOT. ln_aEVP ) THEN 227 zalph1 = ( 2._wp * rn_relast * r dt_ice ) * z1_dtevp227 zalph1 = ( 2._wp * rn_relast * rDt_ice ) * z1_dtevp 228 228 zalph2 = zalph1 * z1_ecc2 229 229 -
NEMO/trunk/src/ICE/iceistate.F90
r12399 r12489 374 374 IF( ln_ice_embd ) THEN ! embedded sea-ice: deplete the initial ssh below sea-ice area 375 375 ! 376 ssh(:,:,Kmm) = ssh(:,:,Kmm) - snwice_mass(:,:) * r1_r au0377 ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_r au0376 ssh(:,:,Kmm) = ssh(:,:,Kmm) - snwice_mass(:,:) * r1_rho0 377 ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_rho0 378 378 ! 379 379 IF( .NOT.ln_linssh ) THEN -
NEMO/trunk/src/ICE/icestp.F90
r12377 r12489 338 338 IF( ln_bdy .AND. ln_icediachk ) CALL ctl_warn('par_init: online conservation check does not work with BDY') 339 339 ! 340 r dt_ice = REAL(nn_fsbc) * rdt !--- sea-ice timestep and its inverse341 r1_ rdtice = 1._wp / rdt_ice340 rDt_ice = REAL(nn_fsbc) * rn_Dt !--- sea-ice timestep and its inverse 341 r1_Dt_ice = 1._wp / rDt_ice 342 342 IF(lwp) WRITE(numout,*) 343 IF(lwp) WRITE(numout,*) ' ice timestep r dt_ice = nn_fsbc*rdt = ', rdt_ice343 IF(lwp) WRITE(numout,*) ' ice timestep rDt_ice = nn_fsbc*rn_Dt = ', rDt_ice 344 344 ! 345 345 r1_nlay_i = 1._wp / REAL( nlay_i, wp ) !--- inverse of nlay_i and nlay_s -
NEMO/trunk/src/ICE/icethd.F90
r12377 r12489 116 116 ELSE ! if no ice dynamics => transmit directly the atmospheric stress to the ocean 117 117 DO_2D_00_00 118 zfric(ji,jj) = r1_r au0 * SQRT( 0.5_wp * &118 zfric(ji,jj) = r1_rho0 * SQRT( 0.5_wp * & 119 119 & ( utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj) & 120 120 & + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) ) * tmask(ji,jj,1) … … 136 136 ! 137 137 ! --- Energy received in the lead from atm-oce exchanges, zqld is defined everywhere (J.m-2) --- ! 138 zqld = tmask(ji,jj,1) * r dt_ice * &138 zqld = tmask(ji,jj,1) * rDt_ice * & 139 139 & ( ( 1._wp - at_i_b(ji,jj) ) * qsr_oce(ji,jj) * frq_m(ji,jj) + & 140 140 & ( 1._wp - at_i_b(ji,jj) ) * qns_oce(ji,jj) + qemp_oce(ji,jj) ) 141 141 142 142 ! --- Energy needed to bring ocean surface layer until its freezing (mostly<0 but >0 if supercooling, J.m-2) --- ! 143 zqfr = r au0 * rcp * e3t_m(ji,jj) * ( t_bo(ji,jj) - ( sst_m(ji,jj) + rt0 ) ) * tmask(ji,jj,1) ! both < 0 (t_bo < sst) and > 0 (t_bo > sst)143 zqfr = rho0 * rcp * e3t_m(ji,jj) * ( t_bo(ji,jj) - ( sst_m(ji,jj) + rt0 ) ) * tmask(ji,jj,1) ! both < 0 (t_bo < sst) and > 0 (t_bo > sst) 144 144 zqfr_neg = MIN( zqfr , 0._wp ) ! only < 0 145 145 146 146 ! --- Sensible ocean-to-ice heat flux (mostly>0 but <0 if supercooling, W/m2) 147 147 zfric_u = MAX( SQRT( zfric(ji,jj) ), zfric_umin ) 148 qsb_ice_bot(ji,jj) = rswitch * r au0 * rcp * zch * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ! W.m-2149 150 qsb_ice_bot(ji,jj) = rswitch * MIN( qsb_ice_bot(ji,jj), - zqfr_neg * r1_ rdtice / MAX( at_i(ji,jj), epsi10 ) )148 qsb_ice_bot(ji,jj) = rswitch * rho0 * rcp * zch * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ! W.m-2 149 150 qsb_ice_bot(ji,jj) = rswitch * MIN( qsb_ice_bot(ji,jj), - zqfr_neg * r1_Dt_ice / MAX( at_i(ji,jj), epsi10 ) ) 151 151 ! upper bound for qsb_ice_bot: the heat retrieved from the ocean must be smaller than the heat necessary to reach 152 152 ! the freezing point, so that we do not have SST < T_freeze … … 154 154 155 155 !-- Energy Budget of the leads (J.m-2), source of ice growth in open water. Must be < 0 to form ice 156 qlead(ji,jj) = MIN( 0._wp , zqld - ( qsb_ice_bot(ji,jj) * at_i(ji,jj) * r dt_ice ) - zqfr )156 qlead(ji,jj) = MIN( 0._wp , zqld - ( qsb_ice_bot(ji,jj) * at_i(ji,jj) * rDt_ice ) - zqfr ) 157 157 158 158 ! If there is ice and leads are warming => transfer energy from the lead budget and use it for bottom melting 159 159 ! If the grid cell is fully covered by ice (no leads) => transfer energy from the lead budget to the ice bottom budget 160 160 IF( ( zqld >= 0._wp .AND. at_i(ji,jj) > 0._wp ) .OR. at_i(ji,jj) >= (1._wp - epsi10) ) THEN 161 fhld (ji,jj) = rswitch * zqld * r1_ rdtice / MAX( at_i(ji,jj), epsi10 ) ! divided by at_i since this is (re)multiplied by a_i in icethd_dh.F90161 fhld (ji,jj) = rswitch * zqld * r1_Dt_ice / MAX( at_i(ji,jj), epsi10 ) ! divided by at_i since this is (re)multiplied by a_i in icethd_dh.F90 162 162 qlead(ji,jj) = 0._wp 163 163 ELSE … … 185 185 ! Third step in iceupdate.F90 : heat from ice-ocean mass exchange (zf_mass) + solar 186 186 qt_oce_ai(:,:) = ( 1._wp - at_i_b(:,:) ) * qns_oce(:,:) + qemp_oce(:,:) & ! Non solar heat flux received by the ocean 187 & - qlead(:,:) * r1_ rdtice & ! heat flux taken from the ocean where there is open water ice formation187 & - qlead(:,:) * r1_Dt_ice & ! heat flux taken from the ocean where there is open water ice formation 188 188 & - at_i (:,:) * qsb_ice_bot(:,:) & ! heat flux taken by sensible flux 189 189 & - at_i (:,:) * fhld (:,:) ! heat flux taken during bottom growth/melt -
NEMO/trunk/src/ICE/icethd_da.F90
r12377 r12489 128 128 zwlat = zm1 * ( MAX( 0._wp, sst_1d(ji) - ( t_bo_1d(ji) - rt0 ) ) )**zm2 ! Melt speed rate [m/s] 129 129 ! 130 zda_tot(ji) = MIN( zwlat * zperi * r dt_ice, at_i_1d(ji) ) ! sea ice concentration decrease (>0)130 zda_tot(ji) = MIN( zwlat * zperi * rDt_ice, at_i_1d(ji) ) ! sea ice concentration decrease (>0) 131 131 132 132 ! --- Distribute reduction among ice categories and calculate associated ice-ocean fluxes --- ! … … 137 137 138 138 ! Contribution to salt flux 139 sfx_lam_1d(ji) = sfx_lam_1d(ji) + rhoi * h_i_1d(ji) * zda * s_i_1d(ji) * r1_ rdtice139 sfx_lam_1d(ji) = sfx_lam_1d(ji) + rhoi * h_i_1d(ji) * zda * s_i_1d(ji) * r1_Dt_ice 140 140 141 141 ! Contribution to heat flux into the ocean [W.m-2], (<0) 142 hfx_thd_1d(ji) = hfx_thd_1d(ji) - zda * r1_ rdtice * ( h_i_1d(ji) * r1_nlay_i * SUM( e_i_1d(ji,1:nlay_i) ) &142 hfx_thd_1d(ji) = hfx_thd_1d(ji) - zda * r1_Dt_ice * ( h_i_1d(ji) * r1_nlay_i * SUM( e_i_1d(ji,1:nlay_i) ) & 143 143 + h_s_1d(ji) * r1_nlay_s * SUM( e_s_1d(ji,1:nlay_s) ) ) 144 144 145 145 ! Contribution to mass flux 146 wfx_lam_1d(ji) = wfx_lam_1d(ji) + zda * r1_ rdtice * ( rhoi * h_i_1d(ji) + rhos * h_s_1d(ji) )146 wfx_lam_1d(ji) = wfx_lam_1d(ji) + zda * r1_Dt_ice * ( rhoi * h_i_1d(ji) + rhos * h_s_1d(ji) ) 147 147 148 148 ! new concentration -
NEMO/trunk/src/ICE/icethd_dh.F90
r10786 r12489 76 76 REAL(wp) :: zgrr ! bottom growth rate 77 77 REAL(wp) :: zt_i_new ! bottom formation temperature 78 REAL(wp) :: z1_rho ! 1/(rhos+r au0-rhoi)78 REAL(wp) :: z1_rho ! 1/(rhos+rho0-rhoi) 79 79 80 80 REAL(wp) :: zQm ! enthalpy exchanged with the ocean (J/m2), >0 towards the ocean … … 130 130 ! 131 131 DO ji = 1, npti 132 zq_top(ji) = MAX( 0._wp, qml_ice_1d(ji) * r dt_ice )132 zq_top(ji) = MAX( 0._wp, qml_ice_1d(ji) * rDt_ice ) 133 133 END DO 134 134 ! … … 138 138 zdum = qns_ice_1d(ji) + qsr_ice_1d(ji) - qtr_ice_top_1d(ji) - qcn_ice_top_1d(ji) 139 139 qml_ice_1d(ji) = zdum * MAX( 0._wp , SIGN( 1._wp, t_su_1d(ji) - rt0 ) ) 140 zq_top(ji) = MAX( 0._wp, qml_ice_1d(ji) * r dt_ice )140 zq_top(ji) = MAX( 0._wp, qml_ice_1d(ji) * rDt_ice ) 141 141 END DO 142 142 ! … … 145 145 DO ji = 1, npti 146 146 zf_tt(ji) = qcn_ice_bot_1d(ji) + qsb_ice_bot_1d(ji) + fhld_1d(ji) 147 zq_bot(ji) = MAX( 0._wp, zf_tt(ji) * r dt_ice )147 zq_bot(ji) = MAX( 0._wp, zf_tt(ji) * rDt_ice ) 148 148 END DO 149 149 … … 172 172 DO ji = 1, npti 173 173 IF( t_s_1d(ji,jk) > rt0 ) THEN 174 hfx_res_1d (ji) = hfx_res_1d (ji) + e_s_1d(ji,jk) * zh_s(ji,jk) * a_i_1d(ji) * r1_ rdtice ! heat flux to the ocean [W.m-2], < 0175 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhos * zh_s(ji,jk) * a_i_1d(ji) * r1_ rdtice ! mass flux174 hfx_res_1d (ji) = hfx_res_1d (ji) + e_s_1d(ji,jk) * zh_s(ji,jk) * a_i_1d(ji) * r1_Dt_ice ! heat flux to the ocean [W.m-2], < 0 175 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhos * zh_s(ji,jk) * a_i_1d(ji) * r1_Dt_ice ! mass flux 176 176 ! updates 177 177 dh_s_mlt(ji) = dh_s_mlt(ji) - zh_s(ji,jk) … … 193 193 ! 194 194 ! --- precipitation --- 195 zdh_s_pre (ji) = zsnw(ji) * sprecip_1d(ji) * r dt_ice * r1_rhos / at_i_1d(ji) ! thickness change195 zdh_s_pre (ji) = zsnw(ji) * sprecip_1d(ji) * rDt_ice * r1_rhos / at_i_1d(ji) ! thickness change 196 196 zqprec (ji) = - qprec_ice_1d(ji) ! enthalpy of the precip (>0, J.m-3) 197 197 ! 198 hfx_spr_1d(ji) = hfx_spr_1d(ji) + zdh_s_pre(ji) * a_i_1d(ji) * zqprec(ji) * r1_ rdtice ! heat flux from snow precip (>0, W.m-2)199 wfx_spr_1d(ji) = wfx_spr_1d(ji) - rhos * a_i_1d(ji) * zdh_s_pre(ji) * r1_ rdtice ! mass flux, <0198 hfx_spr_1d(ji) = hfx_spr_1d(ji) + zdh_s_pre(ji) * a_i_1d(ji) * zqprec(ji) * r1_Dt_ice ! heat flux from snow precip (>0, W.m-2) 199 wfx_spr_1d(ji) = wfx_spr_1d(ji) - rhos * a_i_1d(ji) * zdh_s_pre(ji) * r1_Dt_ice ! mass flux, <0 200 200 201 201 ! --- melt of falling snow --- … … 203 203 zdeltah (ji,1) = - rswitch * zq_top(ji) / MAX( zqprec(ji) , epsi20 ) ! thickness change 204 204 zdeltah (ji,1) = MAX( - zdh_s_pre(ji), zdeltah(ji,1) ) ! bound melting 205 hfx_snw_1d (ji) = hfx_snw_1d (ji) - zdeltah(ji,1) * a_i_1d(ji) * zqprec(ji) * r1_ rdtice ! heat used to melt snow (W.m-2, >0)206 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,1) * r1_ rdtice ! snow melting only = water into the ocean (then without snow precip), >0205 hfx_snw_1d (ji) = hfx_snw_1d (ji) - zdeltah(ji,1) * a_i_1d(ji) * zqprec(ji) * r1_Dt_ice ! heat used to melt snow (W.m-2, >0) 206 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,1) * r1_Dt_ice ! snow melting only = water into the ocean (then without snow precip), >0 207 207 208 208 ! updates available heat + precipitations after melting … … 243 243 zdh_s_mel(ji) = zdh_s_mel(ji) + zdeltah(ji,jk) 244 244 245 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,jk) * a_i_1d(ji) * e_s_1d (ji,jk) * r1_ rdtice ! heat used to melt snow(W.m-2, >0)246 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,jk) * r1_ rdtice ! snow melting only = water into the ocean (then without snow precip)245 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,jk) * a_i_1d(ji) * e_s_1d (ji,jk) * r1_Dt_ice ! heat used to melt snow(W.m-2, >0) 246 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,jk) * r1_Dt_ice ! snow melting only = water into the ocean (then without snow precip) 247 247 248 248 ! updates available heat + thickness … … 264 264 IF( evap_ice_1d(ji) > 0._wp ) THEN 265 265 ! 266 zdh_s_sub (ji) = MAX( - h_s_1d(ji) , - evap_ice_1d(ji) * r1_rhos * r dt_ice )267 zevap_rema(ji) = evap_ice_1d(ji) * r dt_ice + zdh_s_sub(ji) * rhos ! remaining evap in kg.m-2 (used for ice melting later on)266 zdh_s_sub (ji) = MAX( - h_s_1d(ji) , - evap_ice_1d(ji) * r1_rhos * rDt_ice ) 267 zevap_rema(ji) = evap_ice_1d(ji) * rDt_ice + zdh_s_sub(ji) * rhos ! remaining evap in kg.m-2 (used for ice melting later on) 268 268 zdeltah (ji,1) = MAX( zdh_s_sub(ji), - zdh_s_pre(ji) ) 269 269 270 270 hfx_sub_1d (ji) = hfx_sub_1d(ji) + & ! Heat flux by sublimation [W.m-2], < 0 (sublimate snow that had fallen, then pre-existing snow) 271 271 & ( zdeltah(ji,1) * zqprec(ji) + ( zdh_s_sub(ji) - zdeltah(ji,1) ) * e_s_1d(ji,1) ) & 272 & * a_i_1d(ji) * r1_ rdtice273 wfx_snw_sub_1d(ji) = wfx_snw_sub_1d(ji) - rhos * a_i_1d(ji) * zdh_s_sub(ji) * r1_ rdtice ! Mass flux by sublimation272 & * a_i_1d(ji) * r1_Dt_ice 273 wfx_snw_sub_1d(ji) = wfx_snw_sub_1d(ji) - rhos * a_i_1d(ji) * zdh_s_sub(ji) * r1_Dt_ice ! Mass flux by sublimation 274 274 275 275 ! new snow thickness … … 328 328 zfmdt = - rhoi * zdeltah(ji,jk) ! Recompute mass flux [kg/m2, >0] 329 329 330 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_ rdtice ! Heat flux to the ocean [W.m-2], <0330 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_Dt_ice ! Heat flux to the ocean [W.m-2], <0 331 331 ! ice enthalpy zEi is "sent" to the ocean 332 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_ rdtice ! Salt flux332 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_Dt_ice ! Salt flux 333 333 ! using s_i_1d and not sz_i_1d(jk) is ok 334 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_ rdtice ! Mass flux334 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_Dt_ice ! Mass flux 335 335 336 336 ELSE !-- Surface melting … … 354 354 zQm = zfmdt * zEw ! Energy of the melt water sent to the ocean [J/m2, <0] 355 355 356 sfx_sum_1d(ji) = sfx_sum_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_ rdtice ! Salt flux >0356 sfx_sum_1d(ji) = sfx_sum_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_Dt_ice ! Salt flux >0 357 357 ! using s_i_1d and not sz_i_1d(jk) is ok) 358 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_ rdtice ! Heat flux [W.m-2], < 0359 hfx_sum_1d(ji) = hfx_sum_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_ rdtice ! Heat flux used in this process [W.m-2], > 0358 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_Dt_ice ! Heat flux [W.m-2], < 0 359 hfx_sum_1d(ji) = hfx_sum_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_Dt_ice ! Heat flux used in this process [W.m-2], > 0 360 360 ! 361 wfx_sum_1d(ji) = wfx_sum_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_ rdtice ! Mass flux361 wfx_sum_1d(ji) = wfx_sum_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_Dt_ice ! Mass flux 362 362 363 363 END IF … … 369 369 dh_i_sub(ji) = dh_i_sub(ji) + zdum 370 370 371 sfx_sub_1d(ji) = sfx_sub_1d(ji) - rhoi * a_i_1d(ji) * zdum * s_i_1d(ji) * r1_ rdtice ! Salt flux >0371 sfx_sub_1d(ji) = sfx_sub_1d(ji) - rhoi * a_i_1d(ji) * zdum * s_i_1d(ji) * r1_Dt_ice ! Salt flux >0 372 372 ! clem: flux is sent to the ocean for simplicity 373 373 ! but salt should remain in the ice except 374 374 ! if all ice is melted. => must be corrected 375 hfx_sub_1d(ji) = hfx_sub_1d(ji) + zdum * e_i_1d(ji,jk) * a_i_1d(ji) * r1_ rdtice ! Heat flux [W.m-2], < 0376 377 wfx_ice_sub_1d(ji) = wfx_ice_sub_1d(ji) - rhoi * a_i_1d(ji) * zdum * r1_ rdtice ! Mass flux > 0375 hfx_sub_1d(ji) = hfx_sub_1d(ji) + zdum * e_i_1d(ji,jk) * a_i_1d(ji) * r1_Dt_ice ! Heat flux [W.m-2], < 0 376 377 wfx_ice_sub_1d(ji) = wfx_ice_sub_1d(ji) - rhoi * a_i_1d(ji) * zdum * r1_Dt_ice ! Mass flux > 0 378 378 379 379 ! update remaining mass flux … … 400 400 ! remaining "potential" evap is sent to ocean 401 401 DO ji = 1, npti 402 wfx_err_sub_1d(ji) = wfx_err_sub_1d(ji) - zevap_rema(ji) * a_i_1d(ji) * r1_ rdtice ! <=0 (net evap for the ocean in kg.m-2.s-1)402 wfx_err_sub_1d(ji) = wfx_err_sub_1d(ji) - zevap_rema(ji) * a_i_1d(ji) * r1_Dt_ice ! <=0 (net evap for the ocean in kg.m-2.s-1) 403 403 END DO 404 404 … … 428 428 !--- zswi12 if 2.0e-8 < dh/dt < 3.6e-7 429 429 !--- zswi2 if dh/dt > 3.6e-7 430 zgrr = MIN( 1.0e-3, MAX ( dh_i_bog(ji) * r1_ rdtice , epsi10 ) )430 zgrr = MIN( 1.0e-3, MAX ( dh_i_bog(ji) * r1_Dt_ice , epsi10 ) ) 431 431 zswi2 = MAX( 0._wp , SIGN( 1._wp , zgrr - 3.6e-7 ) ) 432 432 zswi12 = MAX( 0._wp , SIGN( 1._wp , zgrr - 2.0e-8 ) ) * ( 1.0 - zswi2 ) … … 448 448 zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0) 449 449 450 dh_i_bog(ji) = r dt_ice * MAX( 0._wp , zf_tt(ji) / ( zdE * rhoi ) )450 dh_i_bog(ji) = rDt_ice * MAX( 0._wp , zf_tt(ji) / ( zdE * rhoi ) ) 451 451 452 452 END DO … … 454 454 zfmdt = - rhoi * dh_i_bog(ji) ! Mass flux x time step (kg/m2, < 0) 455 455 456 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_ rdtice ! Heat flux to the ocean [W.m-2], >0457 hfx_bog_1d(ji) = hfx_bog_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_ rdtice ! Heat flux used in this process [W.m-2], <0458 459 sfx_bog_1d(ji) = sfx_bog_1d(ji) - rhoi * a_i_1d(ji) * dh_i_bog(ji) * s_i_new(ji) * r1_ rdtice ! Salt flux, <0460 461 wfx_bog_1d(ji) = wfx_bog_1d(ji) - rhoi * a_i_1d(ji) * dh_i_bog(ji) * r1_ rdtice ! Mass flux, <0456 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_Dt_ice ! Heat flux to the ocean [W.m-2], >0 457 hfx_bog_1d(ji) = hfx_bog_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_Dt_ice ! Heat flux used in this process [W.m-2], <0 458 459 sfx_bog_1d(ji) = sfx_bog_1d(ji) - rhoi * a_i_1d(ji) * dh_i_bog(ji) * s_i_new(ji) * r1_Dt_ice ! Salt flux, <0 460 461 wfx_bog_1d(ji) = wfx_bog_1d(ji) - rhoi * a_i_1d(ji) * dh_i_bog(ji) * r1_Dt_ice ! Mass flux, <0 462 462 463 463 ! update heat content (J.m-2) and layer thickness … … 490 490 zfmdt = - zdeltah(ji,jk) * rhoi ! Mass flux x time step > 0 491 491 492 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_ rdtice ! Heat flux to the ocean [W.m-2], <0492 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_1d(ji) * zEi * r1_Dt_ice ! Heat flux to the ocean [W.m-2], <0 493 493 ! ice enthalpy zEi is "sent" to the ocean 494 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_ rdtice ! Salt flux494 sfx_res_1d(ji) = sfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_Dt_ice ! Salt flux 495 495 ! using s_i_1d and not sz_i_1d(jk) is ok 496 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_ rdtice ! Mass flux496 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_Dt_ice ! Mass flux 497 497 498 498 ! update heat content (J.m-2) and layer thickness … … 520 520 zQm = zfmdt * zEw ! Heat exchanged with ocean 521 521 522 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_ rdtice ! Heat flux to the ocean [W.m-2], <0523 hfx_bom_1d(ji) = hfx_bom_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_ rdtice ! Heat used in this process [W.m-2], >0524 525 sfx_bom_1d(ji) = sfx_bom_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_ rdtice ! Salt flux522 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_Dt_ice ! Heat flux to the ocean [W.m-2], <0 523 hfx_bom_1d(ji) = hfx_bom_1d(ji) - zfmdt * a_i_1d(ji) * zdE * r1_Dt_ice ! Heat used in this process [W.m-2], >0 524 525 sfx_bom_1d(ji) = sfx_bom_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * s_i_1d(ji) * r1_Dt_ice ! Salt flux 526 526 ! using s_i_1d and not sz_i_1d(jk) is ok 527 527 528 wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_ rdtice ! Mass flux528 wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoi * a_i_1d(ji) * zdeltah(ji,jk) * r1_Dt_ice ! Mass flux 529 529 530 530 ! update heat content (J.m-2) and layer thickness … … 556 556 557 557 zq_rema(ji) = zq_rema(ji) + zdeltah(ji,1) * e_s_1d(ji,1) ! update available heat (J.m-2) 558 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,1) * a_i_1d(ji) * e_s_1d(ji,1) * r1_ rdtice ! Heat used to melt snow, W.m-2 (>0)559 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,1) * r1_ rdtice ! Mass flux558 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,1) * a_i_1d(ji) * e_s_1d(ji,1) * r1_Dt_ice ! Heat used to melt snow, W.m-2 (>0) 559 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) - rhos * a_i_1d(ji) * zdeltah(ji,1) * r1_Dt_ice ! Mass flux 560 560 dh_s_mlt(ji) = dh_s_mlt(ji) + zdeltah(ji,1) 561 561 ! 562 562 ! Remaining heat flux (W.m-2) is sent to the ocean heat budget 563 qt_oce_ai_1d(ji) = qt_oce_ai_1d(ji) + ( zq_rema(ji) * a_i_1d(ji) ) * r1_ rdtice563 qt_oce_ai_1d(ji) = qt_oce_ai_1d(ji) + ( zq_rema(ji) * a_i_1d(ji) ) * r1_Dt_ice 564 564 565 565 IF( ln_icectl .AND. zq_rema(ji) < 0. .AND. lwp ) WRITE(numout,*) 'ALERTE zq_rema <0 = ', zq_rema(ji) … … 571 571 ! When snow load excesses Archimede's limit, snow-ice interface goes down under sea-level, 572 572 ! flooding of seawater transforms snow into ice dh_snowice is positive for the ice 573 z1_rho = 1._wp / ( rhos+r au0-rhoi )573 z1_rho = 1._wp / ( rhos+rho0-rhoi ) 574 574 DO ji = 1, npti 575 575 ! 576 dh_snowice(ji) = MAX( 0._wp , ( rhos * h_s_1d(ji) + (rhoi-r au0) * h_i_1d(ji) ) * z1_rho )576 dh_snowice(ji) = MAX( 0._wp , ( rhos * h_s_1d(ji) + (rhoi-rho0) * h_i_1d(ji) ) * z1_rho ) 577 577 578 578 h_i_1d(ji) = h_i_1d(ji) + dh_snowice(ji) … … 584 584 zQm = zfmdt * zEw 585 585 586 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_ rdtice ! Heat flux587 588 sfx_sni_1d(ji) = sfx_sni_1d(ji) + sss_1d(ji) * a_i_1d(ji) * zfmdt * r1_ rdtice ! Salt flux586 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_1d(ji) * zEw * r1_Dt_ice ! Heat flux 587 588 sfx_sni_1d(ji) = sfx_sni_1d(ji) + sss_1d(ji) * a_i_1d(ji) * zfmdt * r1_Dt_ice ! Salt flux 589 589 590 590 ! Case constant salinity in time: virtual salt flux to keep salinity constant 591 591 IF( nn_icesal /= 2 ) THEN 592 sfx_bri_1d(ji) = sfx_bri_1d(ji) - sss_1d (ji) * a_i_1d(ji) * zfmdt * r1_ rdtice & ! put back sss_m into the ocean593 & - s_i_1d(ji) * a_i_1d(ji) * dh_snowice(ji) * rhoi * r1_ rdtice ! and get rn_icesal from the ocean592 sfx_bri_1d(ji) = sfx_bri_1d(ji) - sss_1d (ji) * a_i_1d(ji) * zfmdt * r1_Dt_ice & ! put back sss_m into the ocean 593 & - s_i_1d(ji) * a_i_1d(ji) * dh_snowice(ji) * rhoi * r1_Dt_ice ! and get rn_icesal from the ocean 594 594 ENDIF 595 595 596 596 ! Mass flux: All snow is thrown in the ocean, and seawater is taken to replace the volume 597 wfx_sni_1d(ji) = wfx_sni_1d(ji) - a_i_1d(ji) * dh_snowice(ji) * rhoi * r1_ rdtice598 wfx_snw_sni_1d(ji) = wfx_snw_sni_1d(ji) + a_i_1d(ji) * dh_snowice(ji) * rhos * r1_ rdtice597 wfx_sni_1d(ji) = wfx_sni_1d(ji) - a_i_1d(ji) * dh_snowice(ji) * rhoi * r1_Dt_ice 598 wfx_snw_sni_1d(ji) = wfx_snw_sni_1d(ji) + a_i_1d(ji) * dh_snowice(ji) * rhos * r1_Dt_ice 599 599 600 600 ! update heat content (J.m-2) and layer thickness … … 618 618 ! mass & energy loss to the ocean 619 619 hfx_res_1d(ji) = hfx_res_1d(ji) + ( 1._wp - rswitch ) * & 620 & ( e_s_1d(ji,jk) * h_s_1d(ji) * r1_nlay_s * a_i_1d(ji) * r1_ rdtice ) ! heat flux to the ocean [W.m-2], < 0620 & ( e_s_1d(ji,jk) * h_s_1d(ji) * r1_nlay_s * a_i_1d(ji) * r1_Dt_ice ) ! heat flux to the ocean [W.m-2], < 0 621 621 wfx_res_1d(ji) = wfx_res_1d(ji) + ( 1._wp - rswitch ) * & 622 & ( rhos * h_s_1d(ji) * r1_nlay_s * a_i_1d(ji) * r1_ rdtice ) ! mass flux622 & ( rhos * h_s_1d(ji) * r1_nlay_s * a_i_1d(ji) * r1_Dt_ice ) ! mass flux 623 623 ! update energy (mass is updated in the next loop) 624 624 e_s_1d(ji,jk) = rswitch * e_s_1d(ji,jk) -
NEMO/trunk/src/ICE/icethd_do.F90
r12377 r12489 141 141 ! Physical constants 142 142 zhicrit = 0.04 ! frazil ice thickness 143 ztwogp = 2. * r au0 / ( grav * 0.3 * ( rau0 - rhoi ) ) ! reduced grav143 ztwogp = 2. * rho0 / ( grav * 0.3 * ( rho0 - rhoi ) ) ! reduced grav 144 144 zsqcd = 1.0 / SQRT( 1.3 * zcai ) ! 1/SQRT(airdensity*drag) 145 145 zgamafr = 0.03 … … 289 289 290 290 ! Contribution to heat flux to the ocean [W.m-2], >0 291 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * zEw * r1_ rdtice291 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * zEw * r1_Dt_ice 292 292 ! Total heat flux used in this process [W.m-2] 293 hfx_opw_1d(ji) = hfx_opw_1d(ji) - zfmdt * zdE * r1_ rdtice293 hfx_opw_1d(ji) = hfx_opw_1d(ji) - zfmdt * zdE * r1_Dt_ice 294 294 ! mass flux 295 wfx_opw_1d(ji) = wfx_opw_1d(ji) - zv_newice(ji) * rhoi * r1_ rdtice295 wfx_opw_1d(ji) = wfx_opw_1d(ji) - zv_newice(ji) * rhoi * r1_Dt_ice 296 296 ! salt flux 297 sfx_opw_1d(ji) = sfx_opw_1d(ji) - zv_newice(ji) * rhoi * zs_newice(ji) * r1_ rdtice297 sfx_opw_1d(ji) = sfx_opw_1d(ji) - zv_newice(ji) * rhoi * zs_newice(ji) * r1_Dt_ice 298 298 END DO 299 299 -
NEMO/trunk/src/ICE/icethd_ent.F90
r10069 r12489 129 129 ! then we should not (* a_i) again but not important since this is just to check that remap error is ~0 130 130 DO ji = 1, npti 131 hfx_err_rem_1d(ji) = hfx_err_rem_1d(ji) + a_i_1d(ji) * r1_ rdtice * &131 hfx_err_rem_1d(ji) = hfx_err_rem_1d(ji) + a_i_1d(ji) * r1_Dt_ice * & 132 132 & ( SUM( qnew(ji,1:nlay_i) ) * zhnew(ji) - SUM( eh_i_old(ji,0:nlay_i+1) ) ) 133 133 END DO -
NEMO/trunk/src/ICE/icethd_pnd.F90
r12377 r12489 165 165 ! melt pond mass flux (<0) 166 166 IF( zdv_mlt > 0._wp ) THEN 167 zfac = zfr_mlt * zdv_mlt * rhow * r1_ rdtice167 zfac = zfr_mlt * zdv_mlt * rhow * r1_Dt_ice 168 168 wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac 169 169 ! -
NEMO/trunk/src/ICE/icethd_sal.F90
r12377 r12489 68 68 CASE( 2 ) ! time varying salinity with linear profile ! 69 69 ! !---------------------------------------------! 70 z1_time_gd = 1._wp / rn_time_gd * r dt_ice71 z1_time_fl = 1._wp / rn_time_fl * r dt_ice70 z1_time_gd = 1._wp / rn_time_gd * rDt_ice 71 z1_time_fl = 1._wp / rn_time_fl * rDt_ice 72 72 ! 73 73 DO ji = 1, npti … … 98 98 99 99 ! Salt flux 100 sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoi * a_i_1d(ji) * h_i_1d(ji) * ( zs_i_fl + zs_i_gd ) * r1_ rdtice100 sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoi * a_i_1d(ji) * h_i_1d(ji) * ( zs_i_fl + zs_i_gd ) * r1_Dt_ice 101 101 ENDIF 102 102 END DO -
NEMO/trunk/src/ICE/icethd_zdf_bl99.F90
r12396 r12489 320 320 DO ji = 1, npti 321 321 zcpi = rcpi + zgamma * sz_i_1d(ji,jk) / MAX( ( t_i_1d(ji,jk) - rt0 ) * ( ztiold(ji,jk) - rt0 ), epsi10 ) 322 zeta_i(ji,jk) = r dt_ice * r1_rhoi * z1_h_i(ji) / MAX( epsi10, zcpi )322 zeta_i(ji,jk) = rDt_ice * r1_rhoi * z1_h_i(ji) / MAX( epsi10, zcpi ) 323 323 END DO 324 324 END DO … … 326 326 DO jk = 1, nlay_s 327 327 DO ji = 1, npti 328 zeta_s(ji,jk) = r dt_ice * r1_rhos * r1_rcpi * z1_h_s(ji)328 zeta_s(ji,jk) = rDt_ice * r1_rhos * r1_rcpi * z1_h_s(ji) 329 329 END DO 330 330 END DO … … 826 826 IF( t_su_1d(ji) < rt0 ) THEN ! case T_su < 0degC 827 827 zhfx_err = ( qns_ice_1d(ji) + qsr_ice_1d(ji) - zradtr_i(ji,nlay_i) - qcn_ice_bot_1d(ji) & 828 & + zdq * r1_ rdtice ) * a_i_1d(ji)828 & + zdq * r1_Dt_ice ) * a_i_1d(ji) 829 829 ELSE ! case T_su = 0degC 830 830 zhfx_err = ( qcn_ice_top_1d(ji) + qtr_ice_top_1d(ji) - zradtr_i(ji,nlay_i) - qcn_ice_bot_1d(ji) & 831 & + zdq * r1_ rdtice ) * a_i_1d(ji)831 & + zdq * r1_Dt_ice ) * a_i_1d(ji) 832 832 ENDIF 833 833 … … 835 835 836 836 zhfx_err = ( qcn_ice_top_1d(ji) + qtr_ice_top_1d(ji) - zradtr_i(ji,nlay_i) - qcn_ice_bot_1d(ji) & 837 & + zdq * r1_ rdtice ) * a_i_1d(ji)837 & + zdq * r1_Dt_ice ) * a_i_1d(ji) 838 838 839 839 ENDIF … … 843 843 ! 844 844 ! hfx_dif = Heat flux diagnostic of sensible heat used to warm/cool ice in W.m-2 845 hfx_dif_1d(ji) = hfx_dif_1d(ji) - zdq * r1_ rdtice * a_i_1d(ji)845 hfx_dif_1d(ji) = hfx_dif_1d(ji) - zdq * r1_Dt_ice * a_i_1d(ji) 846 846 ! 847 847 END DO -
NEMO/trunk/src/ICE/iceupdate.F90
r12377 r12489 171 171 snwice_mass (ji,jj) = tmask(ji,jj,1) * ( rhos * vt_s(ji,jj) + rhoi * vt_i(ji,jj) ) 172 172 ! ! time evolution of snow+ice mass 173 snwice_fmass (ji,jj) = ( snwice_mass(ji,jj) - snwice_mass_b(ji,jj) ) * r1_ rdtice173 snwice_fmass (ji,jj) = ( snwice_mass(ji,jj) - snwice_mass_b(ji,jj) ) * r1_Dt_ice 174 174 175 175 END_2D … … 329 329 ENDIF 330 330 331 zrhoco = r au0 * rn_cio331 zrhoco = rho0 * rn_cio 332 332 ! 333 333 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step) -
NEMO/trunk/src/ICE/icevar.F90
r12377 r12489 488 488 DO_3D_11_11( 1, nlay_i ) 489 489 ! update exchanges with ocean 490 hfx_res(ji,jj) = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_i(ji,jj,jk,jl) * r1_ rdtice ! W.m-2 <0490 hfx_res(ji,jj) = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_i(ji,jj,jk,jl) * r1_Dt_ice ! W.m-2 <0 491 491 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * zswitch(ji,jj) 492 492 t_i(ji,jj,jk,jl) = t_i(ji,jj,jk,jl) * zswitch(ji,jj) + rt0 * ( 1._wp - zswitch(ji,jj) ) … … 495 495 DO_3D_11_11( 1, nlay_s ) 496 496 ! update exchanges with ocean 497 hfx_res(ji,jj) = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_s(ji,jj,jk,jl) * r1_ rdtice ! W.m-2 <0497 hfx_res(ji,jj) = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_s(ji,jj,jk,jl) * r1_Dt_ice ! W.m-2 <0 498 498 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * zswitch(ji,jj) 499 499 t_s(ji,jj,jk,jl) = t_s(ji,jj,jk,jl) * zswitch(ji,jj) + rt0 * ( 1._wp - zswitch(ji,jj) ) … … 505 505 DO_2D_11_11 506 506 ! update exchanges with ocean 507 sfx_res(ji,jj) = sfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * sv_i(ji,jj,jl) * rhoi * r1_ rdtice508 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_i (ji,jj,jl) * rhoi * r1_ rdtice509 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl) * rhos * r1_ rdtice507 sfx_res(ji,jj) = sfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * sv_i(ji,jj,jl) * rhoi * r1_Dt_ice 508 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_i (ji,jj,jl) * rhoi * r1_Dt_ice 509 wfx_res(ji,jj) = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl) * rhos * r1_Dt_ice 510 510 ! 511 511 a_i (ji,jj,jl) = a_i (ji,jj,jl) * zswitch(ji,jj) … … 717 717 !! ** Purpose : compute the equivalent ssh in lead when sea ice is embedded 718 718 !! 719 !! ** Method : ssh_lead = ssh + (Mice + Msnow) / r au0719 !! ** Method : ssh_lead = ssh + (Mice + Msnow) / rho0 720 720 !! 721 721 !! ** Reference : Jean-Michel Campin, John Marshall, David Ferreira, … … 747 747 zintb = REAL( nn_fsbc + 1 ) / REAL( nn_fsbc ) * 0.5_wp 748 748 ! 749 zsnwiceload(:,:) = ( zintn * psnwice_mass(:,:) + zintb * psnwice_mass_b(:,:) ) * r1_r au0749 zsnwiceload(:,:) = ( zintn * psnwice_mass(:,:) + zintb * psnwice_mass_b(:,:) ) * r1_rho0 750 750 ! 751 751 ELSE … … 937 937 ! In case snow load is in excess that would lead to transformation from snow to ice 938 938 ! Then, transfer the snow excess into the ice (different from icethd_dh) 939 zdh = MAX( 0._wp, ( rhos * ph_s(ji,jl) + ( rhoi - r au0 ) * ph_i(ji,jl) ) * r1_rau0 )939 zdh = MAX( 0._wp, ( rhos * ph_s(ji,jl) + ( rhoi - rho0 ) * ph_i(ji,jl) ) * r1_rho0 ) 940 940 ! recompute h_i, h_s avoiding out of bounds values 941 941 ph_i(ji,jl) = MIN( hi_max(jl), ph_i(ji,jl) + zdh ) -
NEMO/trunk/src/ICE/icewri.F90
r12377 r12489 87 87 ! Standard outputs 88 88 !----------------- 89 zrho1 = ( r au0 - rhoi ) * r1_rau0 ; zrho2 = rhos * r1_rau089 zrho1 = ( rho0 - rhoi ) * r1_rho0 ; zrho2 = rhos * r1_rho0 90 90 ! masks 91 91 CALL iom_put( 'icemask' , zmsk00 ) ! ice mask 0% -
NEMO/trunk/src/NST/agrif_oce_sponge.F90
r12377 r12489 439 439 440 440 !* set relaxation time scale 441 IF( neuler == 0 .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_tra / ( rdt )442 ELSE ; ztrelax = rn_trelax_tra / (2._wp * r dt )441 IF( l_1st_euler .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_tra / ( rn_Dt ) 442 ELSE ; ztrelax = rn_trelax_tra / (2._wp * rn_Dt ) 443 443 ENDIF 444 444 … … 596 596 #endif 597 597 !* set relaxation time scale 598 IF( neuler == 0 .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_dyn / ( rdt )599 ELSE ; ztrelax = rn_trelax_dyn / (2._wp * r dt )598 IF( l_1st_euler .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_dyn / ( rn_Dt ) 599 ELSE ; ztrelax = rn_trelax_dyn / (2._wp * rn_Dt ) 600 600 ENDIF 601 601 ! … … 772 772 # endif 773 773 !* set relaxation time scale 774 IF( neuler == 0 .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_dyn / ( rdt )775 ELSE ; ztrelax = rn_trelax_dyn / (2._wp * r dt )774 IF( l_1st_euler .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_dyn / ( rn_Dt ) 775 ELSE ; ztrelax = rn_trelax_dyn / (2._wp * rn_Dt ) 776 776 ENDIF 777 777 ! -
NEMO/trunk/src/NST/agrif_oce_update.F90
r12377 r12489 256 256 ! 2) BEFORE fields: 257 257 !------------------ 258 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0) )) THEN258 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler) )) THEN 259 259 ! 260 260 ! Vertical scale factor interpolations … … 351 351 ENDDO 352 352 353 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN353 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 354 354 ! Add asselin part 355 355 DO jn = 1,jpts … … 361 361 ztnu = tabres_child(ji,jj,jk,jn) * e3t(ji,jj,jk,Kmm_a) 362 362 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 363 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) &363 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 364 364 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 365 365 ENDIF … … 381 381 END DO 382 382 ! 383 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN383 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 384 384 ts(i1:i2,j1:j2,1:jpkm1,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,1:jpkm1,1:jpts,Kmm_a) 385 385 ENDIF … … 422 422 ENDDO 423 423 !< jc tmp 424 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN424 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 425 425 ! Add asselin part 426 426 DO jn = 1,jpts … … 432 432 ztnu = tabres(ji,jj,jk,jn) 433 433 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 434 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) &434 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 435 435 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 436 436 ENDIF … … 452 452 END DO 453 453 ! 454 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN454 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 455 455 ts(i1:i2,j1:j2,k1:k2,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,k1:k2,1:jpts,Kmm_a) 456 456 ENDIF … … 551 551 DO jj=j1,j2 552 552 DO ji=i1,i2 553 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN ! Add asselin part553 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part 554 554 zub = uu(ji,jj,jk,Kbb_a) * e3u(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 555 555 zuno = uu(ji,jj,jk,Kmm_a) * e3u(ji,jj,jk,Krhs_a) 556 556 zunu = tabres_child(ji,jj,jk) * e3u(ji,jj,jk,Kmm_a) 557 uu(ji,jj,jk,Kbb_a) = ( zub + atfp * ( zunu - zuno) ) &557 uu(ji,jj,jk,Kbb_a) = ( zub + rn_atfp * ( zunu - zuno) ) & 558 558 & * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a) 559 559 ENDIF … … 564 564 END DO 565 565 ! 566 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN566 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 567 567 uu(i1:i2,j1:j2,1:jpkm1,Kbb_a) = uu(i1:i2,j1:j2,1:jpkm1,Kmm_a) 568 568 ENDIF … … 597 597 tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e2u(ji,jj) 598 598 ! 599 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN ! Add asselin part599 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part 600 600 zub = uu(ji,jj,jk,Kbb_a) * e3u(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 601 601 zuno = uu(ji,jj,jk,Kmm_a) * e3u(ji,jj,jk,Krhs_a) 602 602 zunu = tabres(ji,jj,jk,1) 603 uu(ji,jj,jk,Kbb_a) = ( zub + atfp * ( zunu - zuno) ) &603 uu(ji,jj,jk,Kbb_a) = ( zub + rn_atfp * ( zunu - zuno) ) & 604 604 & * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a) 605 605 ENDIF … … 610 610 END DO 611 611 ! 612 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN612 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 613 613 uu(i1:i2,j1:j2,k1:k2,Kbb_a) = uu(i1:i2,j1:j2,k1:k2,Kmm_a) 614 614 ENDIF … … 751 751 DO jj=j1,j2 752 752 DO ji=i1,i2 753 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN ! Add asselin part753 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part 754 754 zvb = vv(ji,jj,jk,Kbb_a) * e3v(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 755 755 zvno = vv(ji,jj,jk,Kmm_a) * e3v(ji,jj,jk,Krhs_a) 756 756 zvnu = tabres_child(ji,jj,jk) * e3v(ji,jj,jk,Kmm_a) 757 vv(ji,jj,jk,Kbb_a) = ( zvb + atfp * ( zvnu - zvno) ) &757 vv(ji,jj,jk,Kbb_a) = ( zvb + rn_atfp * ( zvnu - zvno) ) & 758 758 & * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a) 759 759 ENDIF … … 764 764 END DO 765 765 ! 766 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN766 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 767 767 vv(i1:i2,j1:j2,1:jpkm1,Kbb_a) = vv(i1:i2,j1:j2,1:jpkm1,Kmm_a) 768 768 ENDIF … … 801 801 tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e1v(ji,jj) 802 802 ! 803 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN ! Add asselin part803 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part 804 804 zvb = vv(ji,jj,jk,Kbb_a) * e3v(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 805 805 zvno = vv(ji,jj,jk,Kmm_a) * e3v(ji,jj,jk,Krhs_a) 806 806 zvnu = tabres(ji,jj,jk,1) 807 vv(ji,jj,jk,Kbb_a) = ( zvb + atfp * ( zvnu - zvno) ) &807 vv(ji,jj,jk,Kbb_a) = ( zvb + rn_atfp * ( zvnu - zvno) ) & 808 808 & * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a) 809 809 ENDIF … … 814 814 END DO 815 815 ! 816 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN816 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 817 817 vv(i1:i2,j1:j2,k1:k2,Kbb_a) = vv(i1:i2,j1:j2,k1:k2,Kmm_a) 818 818 ENDIF … … 907 907 ! Update barotropic velocities: 908 908 IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 909 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN ! Add asselin part909 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part 910 910 zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a) 911 uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + atfp * zcorr * umask(ji,jj,1)911 uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + rn_atfp * zcorr * umask(ji,jj,1) 912 912 END IF 913 913 ENDIF … … 928 928 END DO 929 929 ! 930 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN930 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 931 931 uu_b(i1:i2,j1:j2,Kbb_a) = uu_b(i1:i2,j1:j2,Kmm_a) 932 932 ENDIF … … 973 973 ! Update barotropic velocities: 974 974 IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 975 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN ! Add asselin part975 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part 976 976 zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a) 977 vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + atfp * zcorr * vmask(ji,jj,1)977 vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + rn_atfp * zcorr * vmask(ji,jj,1) 978 978 END IF 979 979 ENDIF … … 994 994 END DO 995 995 ! 996 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN996 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 997 997 vv_b(i1:i2,j1:j2,Kbb_a) = vv_b(i1:i2,j1:j2,Kmm_a) 998 998 ENDIF … … 1021 1021 END DO 1022 1022 ELSE 1023 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN1023 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 1024 1024 DO jj=j1,j2 1025 1025 DO ji=i1,i2 1026 1026 ssh(ji,jj,Kbb_a) = ssh(ji,jj,Kbb_a) & 1027 & + atfp * ( tabres(ji,jj) - ssh(ji,jj,Kmm_a) ) * tmask(ji,jj,1)1027 & + rn_atfp * ( tabres(ji,jj) - ssh(ji,jj,Kmm_a) ) * tmask(ji,jj,1) 1028 1028 END DO 1029 1029 END DO … … 1036 1036 END DO 1037 1037 ! 1038 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN1038 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 1039 1039 ssh(i1:i2,j1:j2,Kbb_a) = ssh(i1:i2,j1:j2,Kmm_a) 1040 1040 ENDIF … … 1117 1117 IF (western_side) THEN 1118 1118 DO jj=j1,j2 1119 zcor = r dt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj))1119 zcor = rn_Dt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj)) 1120 1120 ssh(i1 ,jj,Kmm_a) = ssh(i1 ,jj,Kmm_a) + zcor 1121 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) ssh(i1 ,jj,Kbb_a) = ssh(i1 ,jj,Kbb_a) +atfp * zcor1121 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(i1 ,jj,Kbb_a) = ssh(i1 ,jj,Kbb_a) + rn_atfp * zcor 1122 1122 END DO 1123 1123 ENDIF 1124 1124 IF (eastern_side) THEN 1125 1125 DO jj=j1,j2 1126 zcor = - r dt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj))1126 zcor = - rn_Dt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj)) 1127 1127 ssh(i2+1,jj,Kmm_a) = ssh(i2+1,jj,Kmm_a) + zcor 1128 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) ssh(i2+1,jj,Kbb_a) = ssh(i2+1,jj,Kbb_a) +atfp * zcor1128 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(i2+1,jj,Kbb_a) = ssh(i2+1,jj,Kbb_a) + rn_atfp * zcor 1129 1129 END DO 1130 1130 ENDIF … … 1205 1205 IF (southern_side) THEN 1206 1206 DO ji=i1,i2 1207 zcor = r dt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1))1207 zcor = rn_Dt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1)) 1208 1208 ssh(ji,j1 ,Kmm_a) = ssh(ji,j1 ,Kmm_a) + zcor 1209 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) ssh(ji,j1 ,Kbb_a) = ssh(ji,j1,Kbb_a) +atfp * zcor1209 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(ji,j1 ,Kbb_a) = ssh(ji,j1,Kbb_a) + rn_atfp * zcor 1210 1210 END DO 1211 1211 ENDIF 1212 1212 IF (northern_side) THEN 1213 1213 DO ji=i1,i2 1214 zcor = - r dt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2))1214 zcor = - rn_Dt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2)) 1215 1215 ssh(ji,j2+1,Kmm_a) = ssh(ji,j2+1,Kmm_a) + zcor 1216 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) ssh(ji,j2+1,Kbb_a) = ssh(ji,j2+1,Kbb_a) +atfp * zcor1216 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(ji,j2+1,Kbb_a) = ssh(ji,j2+1,Kbb_a) + rn_atfp * zcor 1217 1217 END DO 1218 1218 ENDIF … … 1359 1359 ! hdiv(i1:i2,j1:j2,1:jpkm1) = e3t(i1:i2,j1:j2,1:jpkm1,Kbb_a) 1360 1360 1361 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0) )) THEN1361 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler) )) THEN 1362 1362 DO jk = 1, jpkm1 1363 1363 DO jj=j1,j2 1364 1364 DO ji=i1,i2 1365 1365 e3t(ji,jj,jk,Kbb_a) = e3t(ji,jj,jk,Kbb_a) & 1366 & + atfp * ( ptab(ji,jj,jk) - e3t(ji,jj,jk,Kmm_a) )1366 & + rn_atfp * ( ptab(ji,jj,jk) - e3t(ji,jj,jk,Kmm_a) ) 1367 1367 END DO 1368 1368 END DO … … 1422 1422 END DO 1423 1423 ! 1424 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN1424 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 1425 1425 e3t (i1:i2,j1:j2,1:jpk,Kbb_a) = e3t (i1:i2,j1:j2,1:jpk,Kmm_a) 1426 1426 e3w (i1:i2,j1:j2,1:jpk,Kbb_a) = e3w (i1:i2,j1:j2,1:jpk,Kmm_a) -
NEMO/trunk/src/NST/agrif_top_sponge.F90
r12377 r12489 137 137 138 138 !* set relaxation time scale 139 IF( neuler == 0 .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_tra / ( rdt )140 ELSE ; ztrelax = rn_trelax_tra / (2._wp * r dt )139 IF( l_1st_euler .AND. lk_agrif_fstep ) THEN ; ztrelax = rn_trelax_tra / ( rn_Dt ) 140 ELSE ; ztrelax = rn_trelax_tra / (2._wp * rn_Dt ) 141 141 ENDIF 142 142 -
NEMO/trunk/src/NST/agrif_top_update.F90
r12377 r12489 125 125 ENDDO 126 126 ! 127 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN127 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 128 128 ! Add asselin part 129 129 DO jn = 1,jptra … … 135 135 ztnu = tabres_child(ji,jj,jk,jn) * e3t(ji,jj,jk,Kmm_a) 136 136 ztno = tr(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 137 tr(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) &137 tr(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 138 138 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 139 139 ENDIF … … 155 155 END DO 156 156 ! 157 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN157 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 158 158 tr(i1:i2,j1:j2,1:jpkm1,1:jptra,Kbb_a) = tr(i1:i2,j1:j2,1:jpkm1,1:jptra,Kmm_a) 159 159 ENDIF … … 199 199 ENDDO 200 200 !< jc tmp 201 IF (.NOT.(lk_agrif_fstep.AND.( neuler==0))) THEN201 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 202 202 ! Add asselin part 203 203 DO jn = n1,n2 … … 209 209 ztnu = tabres(ji,jj,jk,jn) 210 210 ztno = tr(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 211 tr(ji,jj,jk,jn,Kbb_a) = ( ztb + atfp * ( ztnu - ztno) ) &211 tr(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 212 212 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 213 213 ENDIF … … 229 229 END DO 230 230 ! 231 IF (( neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN231 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 232 232 tr(i1:i2,j1:j2,k1:k2,n1:n2,Kbb_a) = tr(i1:i2,j1:j2,k1:k2,n1:n2,Kmm_a) 233 233 ENDIF -
NEMO/trunk/src/NST/agrif_user.F90
r12377 r12489 202 202 203 203 ! Check time steps 204 IF( NINT(Agrif_Rhot()) * NINT(r dt) .NE. Agrif_Parent(rdt) ) THEN205 WRITE(cl_check1,*) NINT(Agrif_Parent(r dt))206 WRITE(cl_check2,*) NINT(r dt)207 WRITE(cl_check3,*) NINT(Agrif_Parent(r dt)/Agrif_Rhot())204 IF( NINT(Agrif_Rhot()) * NINT(rn_Dt) .NE. Agrif_Parent(rn_Dt) ) THEN 205 WRITE(cl_check1,*) NINT(Agrif_Parent(rn_Dt)) 206 WRITE(cl_check2,*) NINT(rn_Dt) 207 WRITE(cl_check3,*) NINT(Agrif_Parent(rn_Dt)/Agrif_Rhot()) 208 208 CALL ctl_stop( 'Incompatible time step between ocean grids', & 209 209 & 'parent grid value : '//cl_check1 , & … … 613 613 IF( check_namelist ) THEN 614 614 ! Check time steps 615 IF( NINT(Agrif_Rhot()) * NINT(r dt) .NE. Agrif_Parent(rdt) ) THEN616 WRITE(cl_check1,*) Agrif_Parent(r dt)617 WRITE(cl_check2,*) r dt618 WRITE(cl_check3,*) r dt*Agrif_Rhot()615 IF( NINT(Agrif_Rhot()) * NINT(rn_Dt) .NE. Agrif_Parent(rn_Dt) ) THEN 616 WRITE(cl_check1,*) Agrif_Parent(rn_Dt) 617 WRITE(cl_check2,*) rn_Dt 618 WRITE(cl_check3,*) rn_Dt*Agrif_Rhot() 619 619 CALL ctl_stop( 'incompatible time step between grids', & 620 620 & 'parent grid value : '//cl_check1 , & -
NEMO/trunk/src/OCE/ASM/asminc.F90
r12377 r12489 487 487 ENDIF 488 488 ! 489 IF(lwp) WRITE(numout,*) ' ==>>> Euler time step switch is ', neuler489 IF(lwp) WRITE(numout,*) ' ==>>> Euler time step switch is ', l_1st_euler 490 490 ! 491 491 IF( lk_asminc ) THEN !== data assimilation ==! … … 534 534 ! 535 535 it = kt - nit000 + 1 536 zincwgt = wgtiau(it) / r dt ! IAU weight for the current time step536 zincwgt = wgtiau(it) / rn_Dt ! IAU weight for the current time step 537 537 ! 538 538 IF(lwp) THEN … … 577 577 IF ( kt == nitdin_r ) THEN 578 578 ! 579 neuler = 0! Force Euler forward step579 l_1st_euler = .TRUE. ! Force Euler forward step 580 580 ! 581 581 ! Initialize the now fields with the background + increment … … 651 651 ! 652 652 it = kt - nit000 + 1 653 zincwgt = wgtiau(it) / r dt ! IAU weight for the current time step653 zincwgt = wgtiau(it) / rn_Dt ! IAU weight for the current time step 654 654 ! 655 655 IF(lwp) THEN … … 677 677 IF ( kt == nitdin_r ) THEN 678 678 ! 679 neuler = 0! Force Euler forward step679 l_1st_euler = .TRUE. ! Force Euler forward step 680 680 ! 681 681 ! Initialize the now fields with the background + increment … … 722 722 ! 723 723 it = kt - nit000 + 1 724 zincwgt = wgtiau(it) / r dt ! IAU weight for the current time step724 zincwgt = wgtiau(it) / rn_Dt ! IAU weight for the current time step 725 725 ! 726 726 IF(lwp) THEN … … 753 753 IF ( kt == nitdin_r ) THEN 754 754 ! 755 neuler = 0! Force Euler forward step755 l_1st_euler = .TRUE. ! Force Euler forward step 756 756 ! 757 757 ssh(:,:,Kmm) = ssh_bkg(:,:) + ssh_bkginc(:,:) ! Initialize the now fields the background + increment … … 841 841 it = kt - nit000 + 1 842 842 zincwgt = wgtiau(it) ! IAU weight for the current time step 843 ! note this is not a tendency so should not be divided by r dt (as with the tracer and other increments)843 ! note this is not a tendency so should not be divided by rn_Dt (as with the tracer and other increments) 844 844 ! 845 845 IF(lwp) THEN … … 876 876 #if defined key_cice && defined key_asminc 877 877 ! Sea-ice : CICE case. Pass ice increment tendency into CICE 878 ndaice_da(:,:) = seaice_bkginc(:,:) * zincwgt / r dt878 ndaice_da(:,:) = seaice_bkginc(:,:) * zincwgt / rn_Dt 879 879 #endif 880 880 ! … … 896 896 IF ( kt == nitdin_r ) THEN 897 897 ! 898 neuler = 0! Force Euler forward step898 l_1st_euler = 0 ! Force Euler forward step 899 899 ! 900 900 ! Sea-ice : SI3 case … … 926 926 #if defined key_cice && defined key_asminc 927 927 ! Sea-ice : CICE case. Pass ice increment tendency into CICE 928 ndaice_da(:,:) = seaice_bkginc(:,:) / r dt928 ndaice_da(:,:) = seaice_bkginc(:,:) / rn_Dt 929 929 #endif 930 930 IF ( .NOT. PRESENT(kindic) ) THEN … … 959 959 ! ! fwf : ice formation and melting 960 960 ! 961 ! zfons = ( -nfresh_da(ji,jj)*soce + nfsalt_da(ji,jj) )*r dt961 ! zfons = ( -nfresh_da(ji,jj)*soce + nfsalt_da(ji,jj) )*rn_Dt 962 962 ! 963 963 ! ! change salinity down to mixed layer depth … … 1000 1000 ! 1001 1001 ! ! ! salt exchanges at the ice/ocean interface 1002 ! ! zpmess = zfons / r dt_ice ! rdt_ice is ice timestep1002 ! ! zpmess = zfons / rDt_ice ! rDt_ice is ice timestep 1003 1003 ! ! 1004 1004 ! !! Adjust fsalt. A +ve fsalt means adding salt to ocean -
NEMO/trunk/src/OCE/BDY/bdyice.F90
r11536 r12489 179 179 180 180 ! Then, a) transfer the snow excess into the ice (different from icethd_dh) 181 zdh = MAX( 0._wp, ( rhos * h_s(ji,jj,jl) + ( rhoi - r au0 ) * h_i(ji,jj,jl) ) * r1_rau0 )181 zdh = MAX( 0._wp, ( rhos * h_s(ji,jj,jl) + ( rhoi - rho0 ) * h_i(ji,jj,jl) ) * r1_rho0 ) 182 182 ! Or, b) transfer all the snow into ice (if incoming ice is likely to melt as it comes into a warmer environment) 183 183 !zdh = MAX( 0._wp, h_s(ji,jj,jl) * rhos / rhoi ) -
NEMO/trunk/src/OCE/BDY/bdylib.F90
r12377 r12489 240 240 ! Centred derivative is calculated as average of "left" and "right" derivatives for 241 241 ! this reason. 242 ! Note no r dt factor in expression for zdt because it cancels in the expressions for242 ! Note no rn_Dt factor in expression for zdt because it cancels in the expressions for 243 243 ! zrx and zry. 244 244 zdt = phia(iibm1 ,ijbm1 ) - phib(iibm1 ,ijbm1 ) … … 259 259 zout = sign( 1., zrx ) 260 260 zout = 0.5*( zout + abs(zout) ) 261 zwgt = 2.*r dt*( (1.-zout) * idx%nbd(jb,igrd) + zout * idx%nbdout(jb,igrd) )261 zwgt = 2.*rn_Dt*( (1.-zout) * idx%nbd(jb,igrd) + zout * idx%nbdout(jb,igrd) ) 262 262 ! only apply radiation on outflow points 263 263 if( ll_npo ) then !! NPO version !! … … 425 425 zout = sign( 1., zrx ) 426 426 zout = 0.5*( zout + abs(zout) ) 427 zwgt = 2.*r dt*( (1.-zout) * idx%nbd(jb,igrd) + zout * idx%nbdout(jb,igrd) )427 zwgt = 2.*rn_Dt*( (1.-zout) * idx%nbd(jb,igrd) + zout * idx%nbdout(jb,igrd) ) 428 428 ! only apply radiation on outflow points 429 429 if( ll_npo ) then !! NPO version !! -
NEMO/trunk/src/OCE/BDY/bdytides.F90
r12377 r12489 297 297 ! Absolute time from model initialization: 298 298 IF( PRESENT(kit) ) THEN 299 z_arg = ( REAL(kt, wp) + ( REAL(kit, wp) + zt_offset - 1. ) / REAL(nn_ baro, wp) ) * rdt299 z_arg = ( REAL(kt, wp) + ( REAL(kit, wp) + zt_offset - 1. ) / REAL(nn_e, wp) ) * rn_Dt 300 300 ELSE 301 z_arg = ( REAL(kt, wp) + zt_offset ) * r dt301 z_arg = ( REAL(kt, wp) + zt_offset ) * rn_Dt 302 302 ENDIF 303 303 304 304 ! Linear ramp on tidal component at open boundaries 305 305 zramp = 1. 306 IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - REAL(nit000,wp)*r dt)/(rn_tide_ramp_dt*rday),0.),1.)306 IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - REAL(nit000,wp)*rn_Dt)/(rn_tide_ramp_dt*rday),0.),1.) 307 307 308 308 DO ib_bdy = 1,nb_bdy … … 319 319 ! We refresh nodal factors every day below 320 320 ! This should be done somewhere else 321 IF ( ( nsec_day == NINT(0.5_wp * r dt) .OR. kt==nit000 ) .AND. lk_first_btstp ) THEN322 ! 323 kt_tide = kt - NINT((REAL(nsec_day,wp) - 0.5_wp * r dt)/rdt)321 IF ( ( nsec_day == NINT(0.5_wp * rn_Dt) .OR. kt==nit000 ) .AND. lk_first_btstp ) THEN 322 ! 323 kt_tide = kt - NINT((REAL(nsec_day,wp) - 0.5_wp * rn_Dt)/rn_Dt) 324 324 ! 325 325 IF(lwp) THEN … … 333 333 ! 334 334 ENDIF 335 zoff = REAL(-kt_tide,wp) * r dt ! time offset relative to nodal factor computation time335 zoff = REAL(-kt_tide,wp) * rn_Dt ! time offset relative to nodal factor computation time 336 336 ! 337 337 ! If time splitting, initialize arrays from slow varying open boundary data: -
NEMO/trunk/src/OCE/BDY/bdyvol.F90
r12377 r12489 77 77 ! Calculate the cumulate surface Flux z_cflxemp (m3/s) over all the domain 78 78 ! ----------------------------------------------------------------------- 79 IF ( kc == 1 ) z_cflxemp = glob_sum( 'bdyvol', ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) * bdytmask(:,:) * e1e2t(:,:) ) / r au079 IF ( kc == 1 ) z_cflxemp = glob_sum( 'bdyvol', ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) * bdytmask(:,:) * e1e2t(:,:) ) / rho0 80 80 81 81 ! Compute bdy surface each cycle if non linear free surface -
NEMO/trunk/src/OCE/DIA/dia25h.F90
r12377 r12489 140 140 ! ----------------- 141 141 ! Define frequency of summing to create 25 h mean 142 IF( MOD( 3600,NINT(r dt) ) == 0 ) THEN143 i_steps = 3600/NINT(r dt)142 IF( MOD( 3600,NINT(rn_Dt) ) == 0 ) THEN 143 i_steps = 3600/NINT(rn_Dt) 144 144 ELSE 145 CALL ctl_stop('STOP', 'dia_wri_tide: timestep must give MOD(3600,r dt) = 0 otherwise no hourly values are possible')145 CALL ctl_stop('STOP', 'dia_wri_tide: timestep must give MOD(3600,rn_Dt) = 0 otherwise no hourly values are possible') 146 146 ENDIF 147 147 -
NEMO/trunk/src/OCE/DIA/diaar5.F90
r12377 r12489 103 103 END DO 104 104 CALL iom_put( 'volcello' , zrhd(:,:,:) ) ! WARNING not consistent with CMIP DR where volcello is at ca. 2000 105 CALL iom_put( 'masscello' , r au0 * e3t(:,:,:,Kmm) * tmask(:,:,:) ) ! ocean mass105 CALL iom_put( 'masscello' , rho0 * e3t(:,:,:,Kmm) * tmask(:,:,:) ) ! ocean mass 106 106 ENDIF 107 107 ! … … 181 181 CALL iom_put( 'sshsteric', zssh_steric ) 182 182 ! ! ocean bottom pressure 183 zztmp = r au0 * grav * 1.e-4_wp ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa183 zztmp = rho0 * grav * 1.e-4_wp ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa 184 184 zbotpres(:,:) = zztmp * ( zbotpres(:,:) + ssh(:,:,Kmm) + thick0(:,:) ) 185 185 CALL iom_put( 'botpres', zbotpres ) … … 213 213 ztemp = glob_sum( 'diaar5', ztsn(:,:,1,jp_tem) ) 214 214 zsal = glob_sum( 'diaar5', ztsn(:,:,1,jp_sal) ) 215 zmass = r au0 * ( zarho + zvol )215 zmass = rho0 * ( zarho + zvol ) 216 216 ! 217 217 CALL iom_put( 'masstot', zmass ) … … 251 251 z2d(:,:) = 0._wp 252 252 DO_3D_11_11( 1, jpkm1 ) 253 z2d(ji,jj) = z2d(ji,jj) + r au0 * e3t(ji,jj,jk,Kmm) * ztpot(ji,jj,jk)253 z2d(ji,jj) = z2d(ji,jj) + rho0 * e3t(ji,jj,jk,Kmm) * ztpot(ji,jj,jk) 254 254 END_3D 255 255 CALL iom_put( 'tosmint_pot', z2d ) … … 285 285 ELSE 286 286 DO_3D_11_11( 1, jpk ) 287 zpe(ji,jj) = zpe(ji,jj) + avt(ji,jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * r au0 * e3w(ji,jj,jk,Kmm)287 zpe(ji,jj) = zpe(ji,jj) + avt(ji,jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * rho0 * e3w(ji,jj,jk,Kmm) 288 288 END_3D 289 289 ENDIF … … 325 325 CALL lbc_lnk( 'diaar5', z2d, 'U', -1. ) 326 326 IF( cptr == 'adv' ) THEN 327 IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , r au0_rcp * z2d ) ! advective heat transport in i-direction328 IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , r au0 * z2d ) ! advective salt transport in i-direction327 IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rho0_rcp * z2d ) ! advective heat transport in i-direction 328 IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rho0 * z2d ) ! advective salt transport in i-direction 329 329 ENDIF 330 330 IF( cptr == 'ldf' ) THEN 331 IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , r au0_rcp * z2d ) ! diffusive heat transport in i-direction332 IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , r au0 * z2d ) ! diffusive salt transport in i-direction331 IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rho0_rcp * z2d ) ! diffusive heat transport in i-direction 332 IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rho0 * z2d ) ! diffusive salt transport in i-direction 333 333 ENDIF 334 334 ! … … 339 339 CALL lbc_lnk( 'diaar5', z2d, 'V', -1. ) 340 340 IF( cptr == 'adv' ) THEN 341 IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , r au0_rcp * z2d ) ! advective heat transport in j-direction342 IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , r au0 * z2d ) ! advective salt transport in j-direction341 IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rho0_rcp * z2d ) ! advective heat transport in j-direction 342 IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rho0 * z2d ) ! advective salt transport in j-direction 343 343 ENDIF 344 344 IF( cptr == 'ldf' ) THEN 345 IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , r au0_rcp * z2d ) ! diffusive heat transport in j-direction346 IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , r au0 * z2d ) ! diffusive salt transport in j-direction345 IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rho0_rcp * z2d ) ! diffusive heat transport in j-direction 346 IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rho0 * z2d ) ! diffusive salt transport in j-direction 347 347 ENDIF 348 348 -
NEMO/trunk/src/OCE/DIA/diacfl.F90
r12377 r12489 52 52 ! 53 53 INTEGER :: ji, jj, jk ! dummy loop indices 54 REAL(wp) :: z 2dt, zCu_max, zCv_max, zCw_max! local scalars54 REAL(wp) :: zCu_max, zCv_max, zCw_max ! local scalars 55 55 INTEGER , DIMENSION(3) :: iloc_u , iloc_v , iloc_w , iloc ! workspace 56 56 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zCu_cfl, zCv_cfl, zCw_cfl ! workspace … … 59 59 IF( ln_timing ) CALL timing_start('dia_cfl') 60 60 ! 61 ! ! setup timestep multiplier to account for initial Eulerian timestep62 IF( neuler == 0 .AND. kt == nit000 ) THEN ; z2dt = rdt63 ELSE ; z2dt = rdt * 2._wp64 ENDIF65 !66 !67 61 DO_3D_11_11( 1, jpk ) 68 zCu_cfl(ji,jj,jk) = ABS( uu(ji,jj,jk,Kmm) ) * z2dt / e1u (ji,jj) ! for i-direction69 zCv_cfl(ji,jj,jk) = ABS( vv(ji,jj,jk,Kmm) ) * z2dt / e2v (ji,jj) ! for j-direction70 zCw_cfl(ji,jj,jk) = ABS( ww(ji,jj,jk) ) * z2dt / e3w(ji,jj,jk,Kmm) ! for k-direction62 zCu_cfl(ji,jj,jk) = ABS( uu(ji,jj,jk,Kmm) ) * rDt / e1u (ji,jj) ! for i-direction 63 zCv_cfl(ji,jj,jk) = ABS( vv(ji,jj,jk,Kmm) ) * rDt / e2v (ji,jj) ! for j-direction 64 zCw_cfl(ji,jj,jk) = ABS( ww(ji,jj,jk) ) * rDt / e3w(ji,jj,jk,Kmm) ! for k-direction 71 65 END_3D 72 66 ! … … 118 112 WRITE(numcfl,*) '******************************************' 119 113 WRITE(numcfl,FMT='(3x,a12,6x,f7.4,1x,i4,1x,i4,1x,i4)') 'Run Max Cu', rCu_max, nCu_loc(1), nCu_loc(2), nCu_loc(3) 120 WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', z2dt/rCu_max114 WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', rDt/rCu_max 121 115 WRITE(numcfl,*) '******************************************' 122 116 WRITE(numcfl,FMT='(3x,a12,6x,f7.4,1x,i4,1x,i4,1x,i4)') 'Run Max Cv', rCv_max, nCv_loc(1), nCv_loc(2), nCv_loc(3) 123 WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', z2dt/rCv_max117 WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', rDt/rCv_max 124 118 WRITE(numcfl,*) '******************************************' 125 119 WRITE(numcfl,FMT='(3x,a12,6x,f7.4,1x,i4,1x,i4,1x,i4)') 'Run Max Cw', rCw_max, nCw_loc(1), nCw_loc(2), nCw_loc(3) 126 WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', z2dt/rCw_max120 WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', rDt/rCw_max 127 121 CLOSE( numcfl ) 128 122 ! … … 131 125 WRITE(numout,*) 'dia_cfl : Maximum Courant number information for the run ' 132 126 WRITE(numout,*) '~~~~~~~' 133 WRITE(numout,*) ' Max Cu = ', rCu_max, ' at (i,j,k) = (',nCu_loc(1),nCu_loc(2),nCu_loc(3),') => dt/C = ', z2dt/rCu_max134 WRITE(numout,*) ' Max Cv = ', rCv_max, ' at (i,j,k) = (',nCv_loc(1),nCv_loc(2),nCv_loc(3),') => dt/C = ', z2dt/rCv_max135 WRITE(numout,*) ' Max Cw = ', rCw_max, ' at (i,j,k) = (',nCw_loc(1),nCw_loc(2),nCw_loc(3),') => dt/C = ', z2dt/rCw_max127 WRITE(numout,*) ' Max Cu = ', rCu_max, ' at (i,j,k) = (',nCu_loc(1),nCu_loc(2),nCu_loc(3),') => dt/C = ', rDt/rCu_max 128 WRITE(numout,*) ' Max Cv = ', rCv_max, ' at (i,j,k) = (',nCv_loc(1),nCv_loc(2),nCv_loc(3),') => dt/C = ', rDt/rCv_max 129 WRITE(numout,*) ' Max Cw = ', rCw_max, ' at (i,j,k) = (',nCw_loc(1),nCw_loc(2),nCw_loc(3),') => dt/C = ', rDt/rCw_max 136 130 ENDIF 137 131 ! -
NEMO/trunk/src/OCE/DIA/diadct.F90
r12377 r12489 676 676 zsn = interp(Kmm,k%I,k%J,jk,'V',ts(:,:,:,jp_sal,Kmm) ) 677 677 zrhop = interp(Kmm,k%I,k%J,jk,'V',rhop) 678 zrhoi = interp(Kmm,k%I,k%J,jk,'V',rhd*r au0+rau0)678 zrhoi = interp(Kmm,k%I,k%J,jk,'V',rhd*rho0+rho0) 679 679 zsshn = 0.5*( ssh(k%I,k%J,Kmm) + ssh(k%I,k%J+1,Kmm) ) * vmask(k%I,k%J,1) 680 680 CASE(2,3) … … 682 682 zsn = interp(Kmm,k%I,k%J,jk,'U',ts(:,:,:,jp_sal,Kmm) ) 683 683 zrhop = interp(Kmm,k%I,k%J,jk,'U',rhop) 684 zrhoi = interp(Kmm,k%I,k%J,jk,'U',rhd*r au0+rau0)684 zrhoi = interp(Kmm,k%I,k%J,jk,'U',rhd*rho0+rho0) 685 685 zsshn = 0.5*( ssh(k%I,k%J,Kmm) + ssh(k%I+1,k%J,Kmm) ) * umask(k%I,k%J,1) 686 686 END SELECT … … 849 849 zsn = interp(Kmm,k%I,k%J,jk,'V',ts(:,:,:,jp_sal,Kmm) ) 850 850 zrhop = interp(Kmm,k%I,k%J,jk,'V',rhop) 851 zrhoi = interp(Kmm,k%I,k%J,jk,'V',rhd*r au0+rau0)851 zrhoi = interp(Kmm,k%I,k%J,jk,'V',rhd*rho0+rho0) 852 852 853 853 CASE(2,3) … … 855 855 zsn = interp(Kmm,k%I,k%J,jk,'U',ts(:,:,:,jp_sal,Kmm) ) 856 856 zrhop = interp(Kmm,k%I,k%J,jk,'U',rhop) 857 zrhoi = interp(Kmm,k%I,k%J,jk,'U',rhd*r au0+rau0)857 zrhoi = interp(Kmm,k%I,k%J,jk,'U',rhd*rho0+rho0) 858 858 zsshn = 0.5*( ssh(k%I,k%J,Kmm) + ssh(k%I+1,k%J,Kmm) ) * umask(k%I,k%J,1) 859 859 END SELECT -
NEMO/trunk/src/OCE/DIA/diadetide.F90
r12377 r12489 9 9 USE in_out_manager , ONLY : lwp, numout 10 10 USE iom , ONLY : iom_put 11 USE dom_oce , ONLY : r dt, nsec_day11 USE dom_oce , ONLY : rn_Dt, nsec_day 12 12 USE phycst , ONLY : rpi 13 13 USE tide_mod … … 100 100 zwght = 0.0_wp 101 101 DO jn = 1, ndiadetide 102 ztmp = ( tdiadetide(jn) - REAL( nsec_day, KIND=wp ) ) / r dt102 ztmp = ( tdiadetide(jn) - REAL( nsec_day, KIND=wp ) ) / rn_Dt 103 103 IF ( ( ztmp < 0.5_wp ).AND.( ztmp >= -0.5_wp ) ) THEN 104 104 zwght = zwght + 1.0_wp / REAL( ndiadetide, KIND=wp ) -
NEMO/trunk/src/OCE/DIA/diahsb.F90
r12377 r12489 91 91 ! 1 - Trends due to forcing ! 92 92 ! ------------------------- ! 93 z_frc_trd_v = r1_r au0 * glob_sum( 'diahsb', - ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) * surf(:,:) ) ! volume fluxes93 z_frc_trd_v = r1_rho0 * glob_sum( 'diahsb', - ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) * surf(:,:) ) ! volume fluxes 94 94 z_frc_trd_t = glob_sum( 'diahsb', sbc_tsc(:,:,jp_tem) * surf(:,:) ) ! heat fluxes 95 95 z_frc_trd_s = glob_sum( 'diahsb', sbc_tsc(:,:,jp_sal) * surf(:,:) ) ! salt fluxes … … 101 101 & + glob_sum( 'diahsb', ( risf_cav_tsc(:,:,jp_tem) + risf_par_tsc(:,:,jp_tem) ) * surf(:,:) ) 102 102 ! ! Add penetrative solar radiation 103 IF( ln_traqsr ) z_frc_trd_t = z_frc_trd_t + r1_r au0_rcp * glob_sum( 'diahsb', qsr (:,:) * surf(:,:) )103 IF( ln_traqsr ) z_frc_trd_t = z_frc_trd_t + r1_rho0_rcp * glob_sum( 'diahsb', qsr (:,:) * surf(:,:) ) 104 104 ! ! Add geothermal heat flux 105 105 IF( ln_trabbc ) z_frc_trd_t = z_frc_trd_t + glob_sum( 'diahsb', qgh_trd0(:,:) * surf(:,:) ) … … 121 121 ENDIF 122 122 123 frc_v = frc_v + z_frc_trd_v * r dt124 frc_t = frc_t + z_frc_trd_t * r dt125 frc_s = frc_s + z_frc_trd_s * r dt123 frc_v = frc_v + z_frc_trd_v * rn_Dt 124 frc_t = frc_t + z_frc_trd_t * rn_Dt 125 frc_s = frc_s + z_frc_trd_s * rn_Dt 126 126 ! ! Advection flux through fixed surface (z=0) 127 127 IF( ln_linssh ) THEN 128 frc_wn_t = frc_wn_t + z_wn_trd_t * r dt129 frc_wn_s = frc_wn_s + z_wn_trd_s * r dt128 frc_wn_t = frc_wn_t + z_wn_trd_t * rn_Dt 129 frc_wn_s = frc_wn_s + z_wn_trd_s * rn_Dt 130 130 ENDIF 131 131 … … 197 197 198 198 CALL iom_put( 'bgfrcvol' , frc_v * 1.e-9 ) ! vol - surface forcing (km3) 199 CALL iom_put( 'bgfrctem' , frc_t * r au0 * rcp * 1.e-20 ) ! hc - surface forcing (1.e20 J)200 CALL iom_put( 'bgfrchfx' , frc_t * r au0 * rcp / & ! hc - surface forcing (W/m2)201 & ( surf_tot * kt * r dt ) )199 CALL iom_put( 'bgfrctem' , frc_t * rho0 * rcp * 1.e-20 ) ! hc - surface forcing (1.e20 J) 200 CALL iom_put( 'bgfrchfx' , frc_t * rho0 * rcp / & ! hc - surface forcing (W/m2) 201 & ( surf_tot * kt * rn_Dt ) ) 202 202 CALL iom_put( 'bgfrcsal' , frc_s * 1.e-9 ) ! sc - surface forcing (psu*km3) 203 203 … … 205 205 CALL iom_put( 'bgtemper' , zdiff_hc / zvol_tot ) ! Temperature drift (C) 206 206 CALL iom_put( 'bgsaline' , zdiff_sc / zvol_tot ) ! Salinity drift (PSU) 207 CALL iom_put( 'bgheatco' , zdiff_hc * 1.e-20 * r au0 * rcp ) ! Heat content drift (1.e20 J)208 CALL iom_put( 'bgheatfx' , zdiff_hc * r au0 * rcp / & ! Heat flux drift (W/m2)209 & ( surf_tot * kt * r dt ) )207 CALL iom_put( 'bgheatco' , zdiff_hc * 1.e-20 * rho0 * rcp ) ! Heat content drift (1.e20 J) 208 CALL iom_put( 'bgheatfx' , zdiff_hc * rho0 * rcp / & ! Heat flux drift (W/m2) 209 & ( surf_tot * kt * rn_Dt ) ) 210 210 CALL iom_put( 'bgsaltco' , zdiff_sc * 1.e-9 ) ! Salt content drift (psu*km3) 211 211 CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9 ) ! volume ssh drift (km3) … … 225 225 CALL iom_put( 'bgtemper' , zdiff_hc1 / zvol_tot) ! Heat content drift (C) 226 226 CALL iom_put( 'bgsaline' , zdiff_sc1 / zvol_tot) ! Salt content drift (PSU) 227 CALL iom_put( 'bgheatco' , zdiff_hc1 * 1.e-20 * r au0 * rcp ) ! Heat content drift (1.e20 J)228 CALL iom_put( 'bgheatfx' , zdiff_hc1 * r au0 * rcp / & ! Heat flux drift (W/m2)229 & ( surf_tot * kt * r dt ) )227 CALL iom_put( 'bgheatco' , zdiff_hc1 * 1.e-20 * rho0 * rcp ) ! Heat content drift (1.e20 J) 228 CALL iom_put( 'bgheatfx' , zdiff_hc1 * rho0 * rcp / & ! Heat flux drift (W/m2) 229 & ( surf_tot * kt * rn_Dt ) ) 230 230 CALL iom_put( 'bgsaltco' , zdiff_sc1 * 1.e-9 ) ! Salt content drift (psu*km3) 231 231 CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9 ) ! volume ssh drift (km3) -
NEMO/trunk/src/OCE/DIA/diahth.F90
r12377 r12489 261 261 zzdep = 300. 262 262 CALL dia_hth_htc( Kmm, zzdep, ts(:,:,:,jp_tem,Kmm), htc3 ) 263 CALL iom_put( 'hc300', r au0_rcp * htc3 ) ! vertically integrated heat content (J/m2)263 CALL iom_put( 'hc300', rho0_rcp * htc3 ) ! vertically integrated heat content (J/m2) 264 264 ENDIF 265 265 ! … … 270 270 zzdep = 700. 271 271 CALL dia_hth_htc( Kmm, zzdep, ts(:,:,:,jp_tem,Kmm), htc7 ) 272 CALL iom_put( 'hc700', r au0_rcp * htc7 ) ! vertically integrated heat content (J/m2)272 CALL iom_put( 'hc700', rho0_rcp * htc7 ) ! vertically integrated heat content (J/m2) 273 273 274 274 ENDIF … … 280 280 zzdep = 2000. 281 281 CALL dia_hth_htc( Kmm, zzdep, ts(:,:,:,jp_tem,Kmm), htc20 ) 282 CALL iom_put( 'hc2000', r au0_rcp * htc20 ) ! vertically integrated heat content (J/m2)282 CALL iom_put( 'hc2000', rho0_rcp * htc20 ) ! vertically integrated heat content (J/m2) 283 283 ENDIF 284 284 ! -
NEMO/trunk/src/OCE/DIA/dianam.F90
r10068 r12489 72 72 73 73 IF( llfsec .OR. kfreq < 0 ) THEN ; inbsec = kfreq ! output frequency already in seconds 74 ELSE ; inbsec = kfreq * NINT( r dt ) ! from time-step to seconds74 ELSE ; inbsec = kfreq * NINT( rn_Dt ) ! from time-step to seconds 75 75 ENDIF 76 76 iddss = NINT( rday ) ! number of seconds in 1 day … … 116 116 ! date of the beginning and the end of the run 117 117 118 zdrun = r dt / rday * REAL( nitend - nit000, wp ) ! length of the run in days119 zjul = fjulday - r dt / rday118 zdrun = rn_Dt / rday * REAL( nitend - nit000, wp ) ! length of the run in days 119 zjul = fjulday - rn_Dt / rday 120 120 CALL ju2ymds( zjul , iyear1, imonth1, iday1, zsec1 ) ! year/month/day of the beginning of run 121 121 CALL ju2ymds( zjul + zdrun, iyear2, imonth2, iday2, zsec2 ) ! year/month/day of the end of run -
NEMO/trunk/src/OCE/DIA/diaptr.F90
r12377 r12489 50 50 51 51 REAL(wp) :: rc_sv = 1.e-6_wp ! conversion from m3/s to Sverdrup 52 REAL(wp) :: rc_pwatt = 1.e-15_wp ! conversion from W to PW (further x r au0 x Cp)53 REAL(wp) :: rc_ggram = 1.e-9_wp ! conversion from g to Gg (further x r au0)52 REAL(wp) :: rc_pwatt = 1.e-15_wp ! conversion from W to PW (further x rho0 x Cp) 53 REAL(wp) :: rc_ggram = 1.e-9_wp ! conversion from g to Gg (further x rho0) 54 54 55 55 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk ! T-point basin interior masks … … 346 346 IF( dia_ptr_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_ptr_init : unable to allocate arrays' ) 347 347 348 rc_pwatt = rc_pwatt * r au0_rcp ! conversion from K.s-1 to PetaWatt349 rc_ggram = rc_ggram * r au0 ! conversion from m3/s to Gg/s348 rc_pwatt = rc_pwatt * rho0_rcp ! conversion from K.s-1 to PetaWatt 349 rc_ggram = rc_ggram * rho0 ! conversion from m3/s to Gg/s 350 350 351 351 IF( lk_mpp ) CALL mpp_ini_znl( numout ) ! Define MPI communicator for zonal sum -
NEMO/trunk/src/OCE/DIA/diawri.F90
r12377 r12489 173 173 174 174 IF ( iom_use("taubot") ) THEN ! bottom stress 175 zztmp = r au0 * 0.25175 zztmp = rho0 * 0.25 176 176 z2d(:,:) = 0._wp 177 177 DO_2D_00_00 … … 212 212 IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN ! vertical mass transport & its square value 213 213 ! Caution: in the VVL case, it only correponds to the baroclinic mass transport. 214 z2d(:,:) = r au0 * e1e2t(:,:)214 z2d(:,:) = rho0 * e1e2t(:,:) 215 215 DO jk = 1, jpk 216 216 z3d(:,:,jk) = ww(:,:,jk) * z2d(:,:) … … 249 249 z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm) * tmask(ji,jj,jk) 250 250 END_3D 251 CALL iom_put( "heatc", r au0_rcp * z2d ) ! vertically integrated heat content (J/m2)251 CALL iom_put( "heatc", rho0_rcp * z2d ) ! vertically integrated heat content (J/m2) 252 252 ENDIF 253 253 … … 257 257 z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) * tmask(ji,jj,jk) 258 258 END_3D 259 CALL iom_put( "saltc", r au0 * z2d ) ! vertically integrated salt content (PSU*kg/m2)259 CALL iom_put( "saltc", rho0 * z2d ) ! vertically integrated salt content (PSU*kg/m2) 260 260 ENDIF 261 261 ! … … 279 279 z2d(:,:) = 0.e0 280 280 DO jk = 1, jpkm1 281 z3d(:,:,jk) = r au0 * uu(:,:,jk,Kmm) * e2u(:,:) * e3u(:,:,jk,Kmm) * umask(:,:,jk)281 z3d(:,:,jk) = rho0 * uu(:,:,jk,Kmm) * e2u(:,:) * e3u(:,:,jk,Kmm) * umask(:,:,jk) 282 282 z2d(:,:) = z2d(:,:) + z3d(:,:,jk) 283 283 END DO … … 308 308 z3d(:,:,jpk) = 0.e0 309 309 DO jk = 1, jpkm1 310 z3d(:,:,jk) = r au0 * vv(:,:,jk,Kmm) * e1v(:,:) * e3v(:,:,jk,Kmm) * vmask(:,:,jk)310 z3d(:,:,jk) = rho0 * vv(:,:,jk,Kmm) * e1v(:,:) * e3v(:,:,jk,Kmm) * vmask(:,:,jk) 311 311 END DO 312 312 CALL iom_put( "v_masstr", z3d ) ! mass transport in j-direction … … 337 337 END_3D 338 338 CALL lbc_lnk( 'diawri', z2d, 'T', -1. ) 339 CALL iom_put( "tosmint", r au0 * z2d ) ! Vertical integral of temperature339 CALL iom_put( "tosmint", rho0 * z2d ) ! Vertical integral of temperature 340 340 ENDIF 341 341 IF( iom_use("somint") ) THEN … … 345 345 END_3D 346 346 CALL lbc_lnk( 'diawri', z2d, 'T', -1. ) 347 CALL iom_put( "somint", r au0 * z2d ) ! Vertical integral of salinity347 CALL iom_put( "somint", rho0 * z2d ) ! Vertical integral of salinity 348 348 ENDIF 349 349 … … 432 432 clop = "x" ! no use of the mask value (require less cpu time and otherwise the model crashes) 433 433 #if defined key_diainstant 434 zsto = nn_write * r dt434 zsto = nn_write * rn_Dt 435 435 clop = "inst("//TRIM(clop)//")" 436 436 #else 437 zsto=r dt437 zsto=rn_Dt 438 438 clop = "ave("//TRIM(clop)//")" 439 439 #endif 440 zout = nn_write * r dt441 zmax = ( nitend - nit000 + 1 ) * r dt440 zout = nn_write * rn_Dt 441 zmax = ( nitend - nit000 + 1 ) * rn_Dt 442 442 443 443 ! Define indices of the horizontal output zoom and vertical limit storage … … 460 460 461 461 ! Compute julian date from starting date of the run 462 CALL ymds2ju( nyear, nmonth, nday, r dt, zjulian )462 CALL ymds2ju( nyear, nmonth, nday, rn_Dt, zjulian ) 463 463 zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment 464 464 IF(lwp)WRITE(numout,*) … … 482 482 CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit 483 483 & iimi, iima-iimi+1, ijmi, ijma-ijmi+1, & 484 & nit000-1, zjulian, r dt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set )484 & nit000-1, zjulian, rn_Dt, nh_T, nid_T, domain_id=nidom, snc4chunks=snc4set ) 485 485 CALL histvert( nid_T, "deptht", "Vertical T levels", & ! Vertical grid: gdept 486 486 & "m", ipk, gdept_1d, nz_T, "down" ) … … 518 518 CALL histbeg( clhstnam, jpi, glamu, jpj, gphiu, & ! Horizontal grid: glamu and gphiu 519 519 & iimi, iima-iimi+1, ijmi, ijma-ijmi+1, & 520 & nit000-1, zjulian, r dt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set )520 & nit000-1, zjulian, rn_Dt, nh_U, nid_U, domain_id=nidom, snc4chunks=snc4set ) 521 521 CALL histvert( nid_U, "depthu", "Vertical U levels", & ! Vertical grid: gdept 522 522 & "m", ipk, gdept_1d, nz_U, "down" ) … … 531 531 CALL histbeg( clhstnam, jpi, glamv, jpj, gphiv, & ! Horizontal grid: glamv and gphiv 532 532 & iimi, iima-iimi+1, ijmi, ijma-ijmi+1, & 533 & nit000-1, zjulian, r dt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set )533 & nit000-1, zjulian, rn_Dt, nh_V, nid_V, domain_id=nidom, snc4chunks=snc4set ) 534 534 CALL histvert( nid_V, "depthv", "Vertical V levels", & ! Vertical grid : gdept 535 535 & "m", ipk, gdept_1d, nz_V, "down" ) … … 544 544 CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit 545 545 & iimi, iima-iimi+1, ijmi, ijma-ijmi+1, & 546 & nit000-1, zjulian, r dt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set )546 & nit000-1, zjulian, rn_Dt, nh_W, nid_W, domain_id=nidom, snc4chunks=snc4set ) 547 547 CALL histvert( nid_W, "depthw", "Vertical W levels", & ! Vertical grid: gdepw 548 548 & "m", ipk, gdepw_1d, nz_W, "down" ) … … 554 554 CALL histbeg( clhstnam, jpi, glamt, jpj, gphit, & ! Horizontal grid: glamt and gphit 555 555 & iimi, iima-iimi+1, ijmi, ijma-ijmi+1, & 556 & nit000-1, zjulian, r dt, nh_A, nid_A, domain_id=nidom, snc4chunks=snc4set )556 & nit000-1, zjulian, rn_Dt, nh_A, nid_A, domain_id=nidom, snc4chunks=snc4set ) 557 557 CALL histvert( nid_A, "ght_abl", "Vertical T levels", & ! Vertical grid: gdept 558 558 & "m", ipka, ght_abl(2:jpka), nz_A, "up" ) -
NEMO/trunk/src/OCE/DIU/diu_coolskin.F90
r12377 r12489 67 67 68 68 69 SUBROUTINE diurnal_sst_coolskin_step(psqflux, pstauflux, psrho, rdt)69 SUBROUTINE diurnal_sst_coolskin_step(psqflux, pstauflux, psrho, pDt) 70 70 !!---------------------------------------------------------------------- 71 71 !! *** ROUTINE diurnal_sst_takaya_step *** … … 81 81 REAL(wp), INTENT(IN), DIMENSION(jpi,jpj) :: pstauflux ! Wind stress (kg/ m s^2) 82 82 REAL(wp), INTENT(IN), DIMENSION(jpi,jpj) :: psrho ! Water density (kg/m^3) 83 REAL(wp), INTENT(IN) :: rdt ! Time-step83 REAL(wp), INTENT(IN) :: pDt ! Time-step 84 84 85 85 ! Local variables -
NEMO/trunk/src/OCE/DIU/diu_layers.F90
r12377 r12489 39 39 ! Cool skin 40 40 41 CALL diurnal_sst_coolskin_step( qns, taum, rhop(:,:,1), r dt)41 CALL diurnal_sst_coolskin_step( qns, taum, rhop(:,:,1), rn_Dt) 42 42 43 43 CALL iom_put( "sst_wl" , x_dsst ) ! warm layer (write out before update below). … … 45 45 46 46 ! Diurnal warm layer model 47 CALL diurnal_sst_takaya_step( kstp, qsr, qns, taum, rhop(:,:,1), r dt)47 CALL diurnal_sst_takaya_step( kstp, qsr, qns, taum, rhop(:,:,1), rn_Dt) 48 48 49 49 END SUBROUTINE diurnal_layers -
NEMO/trunk/src/OCE/DOM/daymod.F90
r12377 r12489 20 20 !! ------------------------------- 21 21 !! sbcmod assume that the time step is dividing the number of second of 22 !! in a day, i.e. ===> MOD( rday, r dt ) == 022 !! in a day, i.e. ===> MOD( rday, rn_Dt ) == 0 23 23 !! except when user defined forcing is used (see sbcmod.F90) 24 24 !!---------------------------------------------------------------------- … … 73 73 ! 74 74 ! max number of seconds between each restart 75 IF( REAL( nitend - nit000 + 1 ) * r dt > REAL( HUGE( nsec1jan000 ) ) ) THEN75 IF( REAL( nitend - nit000 + 1 ) * rn_Dt > REAL( HUGE( nsec1jan000 ) ) ) THEN 76 76 CALL ctl_stop( 'The number of seconds between each restart exceeds the integer 4 max value: 2^31-1. ', & 77 77 & 'You must do a restart at higher frequency (or remove this stop and recompile the code in I8)' ) … … 79 79 nsecd = NINT( rday ) 80 80 nsecd05 = NINT( 0.5 * rday ) 81 ndt = NINT( r dt )82 ndt05 = NINT( 0.5 * r dt )81 ndt = NINT( rn_Dt ) 82 ndt05 = NINT( 0.5 * rn_Dt ) 83 83 84 84 IF( .NOT. l_offline ) CALL day_rst( nit000, 'READ' ) … … 239 239 nsec_monday = nsec_monday + ndt 240 240 nsec_day = nsec_day + ndt 241 adatrj = adatrj + r dt / rday242 fjulday = fjulday + r dt / rday241 adatrj = adatrj + rn_Dt / rday 242 fjulday = fjulday + rn_Dt / rday 243 243 IF( ABS(fjulday - REAL(NINT(fjulday),wp)) < zprec ) fjulday = REAL(NINT(fjulday),wp) ! avoid truncation error 244 244 IF( ABS(adatrj - REAL(NINT(adatrj ),wp)) < zprec ) adatrj = REAL(NINT(adatrj ),wp) ! avoid truncation error … … 309 309 !! In both those options, the exact duration of the experiment 310 310 !! since the beginning (cumulated duration of all previous restart runs) 311 !! is not stored in the restart and is assumed to be (nit000-1)*r dt.311 !! is not stored in the restart and is assumed to be (nit000-1)*rn_Dt. 312 312 !! This is valid is the time step has remained constant. 313 313 !! … … 379 379 isecond = ( nhour * NINT(rhhmm) + nminute ) * NINT(rmmss) 380 380 IF( isecond - ndt05 .lt. 0 ) ndastp = ndastp - 1 ! Start hour is specified in the namelist (default 0) 381 adatrj = ( REAL( nit000-1, wp ) * r dt ) / rday381 adatrj = ( REAL( nit000-1, wp ) * rn_Dt ) / rday 382 382 ! note this is wrong if time step has changed during run 383 383 ENDIF … … 389 389 isecond = ( nhour * NINT(rhhmm) + nminute ) * NINT(rmmss) 390 390 IF( isecond - ndt05 .LT. 0 ) ndastp = ndastp - 1 ! Start hour is specified in the namelist (default 0) 391 adatrj = ( REAL( nit000-1, wp ) * r dt ) / rday391 adatrj = ( REAL( nit000-1, wp ) * rn_Dt ) / rday 392 392 ENDIF 393 393 IF( ABS(adatrj - REAL(NINT(adatrj),wp)) < 0.1 / rday ) adatrj = REAL(NINT(adatrj),wp) ! avoid truncation error -
NEMO/trunk/src/OCE/DOM/dom_oce.F90
r12377 r12489 33 33 LOGICAL , PUBLIC :: ln_linssh !: =T linear free surface ==>> model level are fixed in time 34 34 LOGICAL , PUBLIC :: ln_meshmask !: =T create a mesh-mask file (mesh_mask.nc) 35 REAL(wp), PUBLIC :: rn_ rdt!: time step for the dynamics and tracer35 REAL(wp), PUBLIC :: rn_Dt !: time step for the dynamics and tracer 36 36 REAL(wp), PUBLIC :: rn_atfp !: asselin time filter parameter 37 INTEGER , PUBLIC :: nn_euler !: =0 start with forward time step or not (=1)37 LOGICAL , PUBLIC :: ln_1st_euler !: =T start with forward time step or not (=F) 38 38 LOGICAL , PUBLIC :: ln_crs !: Apply grid coarsening to dynamical model output or online passive tracers 39 39 … … 49 49 LOGICAL, PUBLIC :: ln_bt_auto !: Set number of barotropic iterations automatically 50 50 INTEGER, PUBLIC :: nn_bt_flt !: Filter choice 51 INTEGER, PUBLIC :: nn_ baro !: Number of barotropic iterations during one baroclinic step (rdt)51 INTEGER, PUBLIC :: nn_e !: Number of barotropic iterations during one baroclinic step (rn_Dt) 52 52 REAL(wp), PUBLIC :: rn_bt_cmax !: Maximum allowed courant number (used if ln_bt_auto=T) 53 53 REAL(wp), PUBLIC :: rn_bt_alpha !: Time stepping diffusion parameter 54 54 55 55 56 ! !! old non-DOCTOR names still used in the model57 REAL(wp), PUBLIC :: atfp !: asselin time filter parameter58 REAL(wp), PUBLIC :: rdt !: time step for the dynamics and tracer59 60 56 ! !!! associated variables 61 INTEGER , PUBLIC :: neuler !: restart euler forward option (0=Euler) 62 REAL(wp), PUBLIC :: r2dt !: = 2*rdt except at nit000 (=rdt) if neuler=0 57 LOGICAL , PUBLIC :: l_1st_euler !: Euler 1st time-step flag (=T if ln_restart=F or ln_1st_euler=T) 58 REAL(wp), PUBLIC :: rDt, r1_Dt !: Current model timestep and reciprocal 59 !: rDt = 2 * rn_Dt if leapfrog and l_1st_euler = F 60 !: = rn_Dt if leapfrog and l_1st_euler = T 61 !: = rn_Dt if RK3 63 62 64 63 !!---------------------------------------------------------------------- -
NEMO/trunk/src/OCE/DOM/domain.F90
r12377 r12489 287 287 & nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl , & 288 288 & nn_it000, nn_itend , nn_date0 , nn_time0 , nn_leapy , nn_istate , & 289 & nn_stock, nn_write , ln_mskland , ln_clobber , nn_chunksz, nn_euler ,&289 & nn_stock, nn_write , ln_mskland , ln_clobber , nn_chunksz, ln_1st_euler , & 290 290 & ln_cfmeta, ln_xios_read, nn_wxios 291 NAMELIST/namdom/ ln_linssh, rn_ rdt, rn_atfp, ln_crs, ln_meshmask291 NAMELIST/namdom/ ln_linssh, rn_Dt, rn_atfp, ln_crs, ln_meshmask 292 292 #if defined key_netcdf4 293 293 NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip … … 317 317 WRITE(numout,*) ' restart output directory cn_ocerst_outdir= ', TRIM( cn_ocerst_outdir ) 318 318 WRITE(numout,*) ' restart logical ln_rstart = ', ln_rstart 319 WRITE(numout,*) ' start with forward time step nn_euler = ', nn_euler319 WRITE(numout,*) ' start with forward time step ln_1st_euler = ', ln_1st_euler 320 320 WRITE(numout,*) ' control of time step nn_rstctl = ', nn_rstctl 321 321 WRITE(numout,*) ' number of the first time step nn_it000 = ', nn_it000 … … 353 353 nleapy = nn_leapy 354 354 ninist = nn_istate 355 neuler = nn_euler356 IF( neuler == 1.AND. .NOT. ln_rstart ) THEN355 l_1st_euler = ln_1st_euler 356 IF( .NOT. l_1st_euler .AND. .NOT. ln_rstart ) THEN 357 357 IF(lwp) WRITE(numout,*) 358 358 IF(lwp) WRITE(numout,*)' ==>>> Start from rest (ln_rstart=F)' 359 IF(lwp) WRITE(numout,*)' an Euler initial time step is used : nn_euler is forced to 0'360 neuler = 0359 IF(lwp) WRITE(numout,*)' an Euler initial time step is used : l_1st_euler is forced to .true. ' 360 l_1st_euler = .true. 361 361 ENDIF 362 362 ! ! control of output frequency … … 408 408 WRITE(numout,*) ' linear free surface (=T) ln_linssh = ', ln_linssh 409 409 WRITE(numout,*) ' create mesh/mask file ln_meshmask = ', ln_meshmask 410 WRITE(numout,*) ' ocean time step rn_ rdt = ', rn_rdt410 WRITE(numout,*) ' ocean time step rn_Dt = ', rn_Dt 411 411 WRITE(numout,*) ' asselin time filter parameter rn_atfp = ', rn_atfp 412 412 WRITE(numout,*) ' online coarsening of dynamical fields ln_crs = ', ln_crs 413 413 ENDIF 414 414 ! 415 ! ! conversion DOCTOR names into model names (this should disappear soon)416 atfp = rn_atfp417 r dt = rn_rdt415 !! Initialise current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3 416 rDt = 2._wp * rn_Dt 417 r1_Dt = 1._wp / rDt 418 418 419 419 IF( TRIM(Agrif_CFixed()) == '0' ) THEN -
NEMO/trunk/src/OCE/DOM/domvvl.F90
r12377 r12489 235 235 IF( ln_vvl_ztilde_as_zstar ) THEN ! z-star emulation using z-tile 236 236 frq_rst_e3t(:,:) = 0._wp !Ignore namelist settings 237 frq_rst_hdv(:,:) = 1._wp / r dt237 frq_rst_hdv(:,:) = 1._wp / rn_Dt 238 238 ENDIF 239 239 IF ( ln_vvl_zstar_at_eqtor ) THEN ! use z-star in vicinity of the Equator … … 247 247 ! values inside the equatorial band (ztilde as zstar) 248 248 frq_rst_e3t(ji,jj) = 0.0_wp 249 frq_rst_hdv(ji,jj) = 1.0_wp / r dt249 frq_rst_hdv(ji,jj) = 1.0_wp / rn_Dt 250 250 ELSE ! transition band (2.5 to 6 degrees N/S) 251 251 ! ! (linearly transition from z-tilde to z-star) … … 253 253 & * ( 1.0_wp - COS( rad*(ABS(gphit(ji,jj))-2.5_wp) & 254 254 & * 180._wp / 3.5_wp ) ) 255 frq_rst_hdv(ji,jj) = (1.0_wp / r dt) &256 & + ( frq_rst_hdv(ji,jj)-(1.e0_wp / r dt) )*0.5_wp &255 frq_rst_hdv(ji,jj) = (1.0_wp / rn_Dt) & 256 & + ( frq_rst_hdv(ji,jj)-(1.e0_wp / rn_Dt) )*0.5_wp & 257 257 & * ( 1._wp - COS( rad*(ABS(gphit(ji,jj))-2.5_wp) & 258 258 & * 180._wp / 3.5_wp ) ) … … 264 264 ij0 = 128 ; ij1 = 135 ; 265 265 frq_rst_e3t( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.0_wp 266 frq_rst_hdv( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1.e0_wp / r dt266 frq_rst_hdv( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1.e0_wp / rn_Dt 267 267 ENDIF 268 268 ENDIF … … 319 319 INTEGER :: ji, jj, jk ! dummy loop indices 320 320 INTEGER , DIMENSION(3) :: ijk_max, ijk_min ! temporary integers 321 REAL(wp) :: z 2dt, z_tmin, z_tmax! local scalars321 REAL(wp) :: z_tmin, z_tmax ! local scalars 322 322 LOGICAL :: ll_do_bclinic ! local logical 323 323 REAL(wp), DIMENSION(jpi,jpj) :: zht, z_scale, zwu, zwv, zhdiv … … 373 373 IF( kt > nit000 ) THEN 374 374 DO jk = 1, jpkm1 375 hdiv_lf(:,:,jk) = hdiv_lf(:,:,jk) - r dt * frq_rst_hdv(:,:) &375 hdiv_lf(:,:,jk) = hdiv_lf(:,:,jk) - rn_Dt * frq_rst_hdv(:,:) & 376 376 & * ( hdiv_lf(:,:,jk) - e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) ) 377 377 END DO … … 430 430 ! 4 - Time stepping of baroclinic scale factors 431 431 ! --------------------------------------------- 432 ! Leapfrog time stepping433 ! ~~~~~~~~~~~~~~~~~~~~~~434 IF( neuler == 0 .AND. kt == nit000 ) THEN435 z2dt = rdt436 ELSE437 z2dt = 2.0_wp * rdt438 ENDIF439 432 CALL lbc_lnk( 'domvvl', tilde_e3t_a(:,:,:), 'T', 1._wp ) 440 tilde_e3t_a(:,:,:) = tilde_e3t_b(:,:,:) + z2dt * tmask(:,:,:) * tilde_e3t_a(:,:,:)433 tilde_e3t_a(:,:,:) = tilde_e3t_b(:,:,:) + rDt * tmask(:,:,:) * tilde_e3t_a(:,:,:) 441 434 442 435 ! Maximum deformation control … … 624 617 ! - ML - e3(t/u/v)_b are allready computed in dynnxt. 625 618 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 626 IF( neuler == 0 .AND. kt == nit000) THEN619 IF( l_1st_euler ) THEN 627 620 tilde_e3t_b(:,:,:) = tilde_e3t_n(:,:,:) 628 621 ELSE 629 622 tilde_e3t_b(:,:,:) = tilde_e3t_n(:,:,:) & 630 & + atfp * ( tilde_e3t_b(:,:,:) - 2.0_wp * tilde_e3t_n(:,:,:) + tilde_e3t_a(:,:,:) )623 & + rn_atfp * ( tilde_e3t_b(:,:,:) - 2.0_wp * tilde_e3t_n(:,:,:) + tilde_e3t_a(:,:,:) ) 631 624 ENDIF 632 625 tilde_e3t_n(:,:,:) = tilde_e3t_a(:,:,:) … … 821 814 e3t(:,:,:,Kbb) = e3t_0(:,:,:) 822 815 END WHERE 823 IF( neuler == 0) THEN816 IF( l_1st_euler ) THEN 824 817 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 825 818 ENDIF … … 827 820 IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart files' 828 821 IF(lwp) write(numout,*) 'e3t_n set equal to e3t_b.' 829 IF(lwp) write(numout,*) ' neuler is forced to 0'822 IF(lwp) write(numout,*) 'l_1st_euler is forced to true' 830 823 CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 831 824 e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 832 neuler = 0825 l_1st_euler = .true. 833 826 ELSE IF( id2 > 0 ) THEN 834 827 IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kbb) not found in restart files' 835 828 IF(lwp) write(numout,*) 'e3t_b set equal to e3t_n.' 836 IF(lwp) write(numout,*) ' neuler is forced to 0'829 IF(lwp) write(numout,*) 'l_1st_euler is forced to true' 837 830 CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 838 831 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 839 neuler = 0832 l_1st_euler = .true. 840 833 ELSE 841 834 IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart file' 842 835 IF(lwp) write(numout,*) 'Compute scale factor from sshn' 843 IF(lwp) write(numout,*) ' neuler is forced to 0'836 IF(lwp) write(numout,*) 'l_1st_euler is forced to true' 844 837 DO jk = 1, jpk 845 838 e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) & … … 848 841 END DO 849 842 e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 850 neuler = 0843 l_1st_euler = .true. 851 844 ENDIF 852 845 ! ! ----------- ! … … 1008 1001 WRITE(numout,*) ' rn_rst_e3t = 0.e0' 1009 1002 WRITE(numout,*) ' hard-wired : z-tilde cutoff frequency of low-pass filter (days)' 1010 WRITE(numout,*) ' rn_lf_cutoff = 1.0/r dt'1003 WRITE(numout,*) ' rn_lf_cutoff = 1.0/rn_Dt' 1011 1004 ELSE 1012 1005 WRITE(numout,*) ' z-tilde to zstar restoration timescale (days) rn_rst_e3t = ', rn_rst_e3t -
NEMO/trunk/src/OCE/DOM/istate.F90
r12377 r12489 92 92 ! ! --------------- 93 93 numror = 0 ! define numror = 0 -> no restart file to read 94 neuler = 0! Set time-step indicator at nit000 (euler forward)94 l_1st_euler = .true. ! Set time-step indicator at nit000 (euler forward) 95 95 CALL day_init ! model calendar (using both namelist and restart infos) 96 96 ! ! Initialization of ocean to zero -
NEMO/trunk/src/OCE/DOM/phycst.F90
r10068 r12489 39 39 REAL(wp), PUBLIC :: rt0 = 273.15_wp !: freezing point of fresh water [Kelvin] 40 40 41 REAL(wp), PUBLIC :: r au0 !: volumic mass of reference [kg/m3]42 REAL(wp), PUBLIC :: r1_r au0 !: = 1. / rau0 [m3/kg]41 REAL(wp), PUBLIC :: rho0 !: volumic mass of reference [kg/m3] 42 REAL(wp), PUBLIC :: r1_rho0 !: = 1. / rho0 [m3/kg] 43 43 REAL(wp), PUBLIC :: rcp !: ocean specific heat [J/Kelvin] 44 44 REAL(wp), PUBLIC :: r1_rcp !: = 1. / rcp [Kelvin/J] 45 REAL(wp), PUBLIC :: r au0_rcp !: = rau0 * rcp46 REAL(wp), PUBLIC :: r1_r au0_rcp !: = 1. / ( rau0 * rcp )45 REAL(wp), PUBLIC :: rho0_rcp !: = rho0 * rcp 46 REAL(wp), PUBLIC :: r1_rho0_rcp !: = 1. / ( rho0 * rcp ) 47 47 48 48 REAL(wp), PUBLIC :: emic = 0.97_wp !: emissivity of snow or ice (not used?) -
NEMO/trunk/src/OCE/DYN/dynatf.F90
r12377 r12489 87 87 !! arrays to start the next time step: 88 88 !! (puu(Kmm),pvv(Kmm)) = (puu(Kmm),pvv(Kmm)) 89 !! + atfp [ (puu(Kbb),pvv(Kbb)) + (puu(Kaa),pvv(Kaa)) - 2 (puu(Kmm),pvv(Kmm)) ]89 !! + rn_atfp [ (puu(Kbb),pvv(Kbb)) + (puu(Kaa),pvv(Kaa)) - 2 (puu(Kmm),pvv(Kmm)) ] 90 90 !! Note that with flux form advection and non linear free surface, 91 91 !! the time filter is applied on thickness weighted velocity. … … 157 157 ! 158 158 IF( l_trddyn ) THEN ! prepare the atf trend computation + some diagnostics 159 z1_2dt = 1._wp / (2. * rdt) ! Euler or leap-frog time step160 IF( neuler == 0 .AND. kt == nit000 ) z1_2dt = 1._wp / rdt161 159 ! 162 160 ! ! Kinetic energy and Conversion … … 164 162 ! 165 163 IF( ln_dyn_trd ) THEN ! 3D output: total momentum trends 166 zua(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) * z1_2dt167 zva(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) * z1_2dt164 zua(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) * r1_Dt 165 zva(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) * r1_Dt 168 166 CALL iom_put( "utrd_tot", zua ) ! total momentum trends, except the asselin time filter 169 167 CALL iom_put( "vtrd_tot", zva ) … … 178 176 ! ------------------------------------------ 179 177 180 IF( .NOT. ( neuler == 0 .AND. kt == nit000 )) THEN !* Leap-Frog : Asselin time filter178 IF( .NOT. l_1st_euler ) THEN !* Leap-Frog : Asselin time filter 181 179 ! ! =============! 182 180 IF( ln_linssh ) THEN ! Fixed volume ! 183 181 ! ! =============! 184 182 DO_3D_11_11( 1, jpkm1 ) 185 puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) )186 pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) )183 puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + rn_atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) ) 184 pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + rn_atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) ) 187 185 END_3D 188 186 ! ! ================! … … 193 191 ALLOCATE( ze3t_f(jpi,jpj,jpk), zwfld(jpi,jpj) ) 194 192 DO jk = 1, jpkm1 195 ze3t_f(:,:,jk) = pe3t(:,:,jk,Kmm) + atfp * ( pe3t(:,:,jk,Kbb) - 2._wp * pe3t(:,:,jk,Kmm) + pe3t(:,:,jk,Kaa) )193 ze3t_f(:,:,jk) = pe3t(:,:,jk,Kmm) + rn_atfp * ( pe3t(:,:,jk,Kbb) - 2._wp * pe3t(:,:,jk,Kmm) + pe3t(:,:,jk,Kaa) ) 196 194 END DO 197 195 ! Add volume filter correction: compatibility with tracer advection scheme 198 196 ! => time filter + conservation correction 199 zcoef = atfp * rdt * r1_rau0197 zcoef = rn_atfp * rn_Dt * r1_rho0 200 198 zwfld(:,:) = emp_b(:,:) - emp(:,:) 201 199 IF ( ln_rnf ) zwfld(:,:) = zwfld(:,:) - ( rnf_b(:,:) - rnf(:,:) ) … … 209 207 ! to manage rnf, isf and possibly in the futur icb, tide water glacier (...) 210 208 ! ...(kt, coef, ktop, kbot, hz, fwf_b, fwf) 211 IF ( ln_isf ) CALL isf_dynatf( kt, Kmm, ze3t_f, atfp * rdt )209 IF ( ln_isf ) CALL isf_dynatf( kt, Kmm, ze3t_f, rn_atfp * rn_Dt ) 212 210 ! 213 211 pe3t(:,:,1:jpkm1,Kmm) = ze3t_f(:,:,1:jpkm1) ! filtered scale factor at T-points … … 218 216 CALL dom_vvl_interpol( pe3t(:,:,:,Kmm), pe3v(:,:,:,Kmm), 'V' ) 219 217 DO_3D_11_11( 1, jpkm1 ) 220 puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) )221 pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) )218 puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + rn_atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) ) 219 pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + rn_atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) ) 222 220 END_3D 223 221 ! … … 236 234 zve3b = pe3v(ji,jj,jk,Kbb) * pvv(ji,jj,jk,Kbb) 237 235 ! 238 puu(ji,jj,jk,Kmm) = ( zue3n + atfp * ( zue3b - 2._wp * zue3n + zue3a ) ) / ze3u_f(ji,jj,jk)239 pvv(ji,jj,jk,Kmm) = ( zve3n + atfp * ( zve3b - 2._wp * zve3n + zve3a ) ) / ze3v_f(ji,jj,jk)236 puu(ji,jj,jk,Kmm) = ( zue3n + rn_atfp * ( zue3b - 2._wp * zue3n + zue3a ) ) / ze3u_f(ji,jj,jk) 237 pvv(ji,jj,jk,Kmm) = ( zve3n + rn_atfp * ( zve3b - 2._wp * zve3n + zve3a ) ) / ze3v_f(ji,jj,jk) 240 238 END_3D 241 239 pe3u(:,:,1:jpkm1,Kmm) = ze3u_f(:,:,1:jpkm1) … … 263 261 ENDIF 264 262 ! 265 ENDIF ! neuler /= 0263 ENDIF ! .NOT. l_1st_euler 266 264 ! 267 265 ! Set "now" and "before" barotropic velocities for next time step: -
NEMO/trunk/src/OCE/DYN/dynspg.F90
r12377 r12489 67 67 !! ln_apr_dyn=T : the atmospheric pressure forcing is applied 68 68 !! as the gradient of the inverse barometer ssh: 69 !! apgu = - 1/r au0 di[apr] = 0.5*grav di[ssh_ib+ssh_ibb]70 !! apgv = - 1/r au0 dj[apr] = 0.5*grav dj[ssh_ib+ssh_ibb]71 !! Note that as all external forcing a time averaging over a two r dt69 !! apgu = - 1/rho0 di[apr] = 0.5*grav di[ssh_ib+ssh_ibb] 70 !! apgv = - 1/rho0 dj[apr] = 0.5*grav dj[ssh_ib+ssh_ibb] 71 !! Note that as all external forcing a time averaging over a two rn_Dt 72 72 !! period is used to prevent the divergence of odd and even time step. 73 73 !!---------------------------------------------------------------------- … … 78 78 ! 79 79 INTEGER :: ji, jj, jk ! dummy loop indices 80 REAL(wp) :: z2dt, zg_2, zintp, zgr au0r, zld ! local scalars80 REAL(wp) :: z2dt, zg_2, zintp, zgrho0r, zld ! local scalars 81 81 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpice 82 82 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv … … 114 114 ! 115 115 ! Update tide potential at the beginning of current time step 116 zt0step = REAL(nsec_day, wp)-0.5_wp*r dt116 zt0step = REAL(nsec_day, wp)-0.5_wp*rn_Dt 117 117 CALL upd_tide(zt0step, Kmm) 118 118 ! … … 134 134 ALLOCATE( zpice(jpi,jpj) ) 135 135 zintp = REAL( MOD( kt-1, nn_fsbc ) ) / REAL( nn_fsbc ) 136 zgr au0r = - grav * r1_rau0137 zpice(:,:) = ( zintp * snwice_mass(:,:) + ( 1.- zintp ) * snwice_mass_b(:,:) ) * zgr au0r136 zgrho0r = - grav * r1_rho0 137 zpice(:,:) = ( zintp * snwice_mass(:,:) + ( 1.- zintp ) * snwice_mass_b(:,:) ) * zgrho0r 138 138 DO_2D_00_00 139 139 spgu(ji,jj) = spgu(ji,jj) + ( zpice(ji+1,jj) - zpice(ji,jj) ) * r1_e1u(ji,jj) … … 183 183 NAMELIST/namdyn_spg/ ln_dynspg_exp , ln_dynspg_ts, & 184 184 & ln_bt_fw, ln_bt_av , ln_bt_auto , & 185 & nn_ baro, rn_bt_cmax, nn_bt_flt, rn_bt_alpha185 & nn_e , rn_bt_cmax, nn_bt_flt, rn_bt_alpha 186 186 !!---------------------------------------------------------------------- 187 187 ! … … 222 222 ! 223 223 IF( nspg == np_TS ) THEN ! split-explicit scheme initialisation 224 CALL dyn_spg_ts_init ! do it first: set nn_ baroused to allocate some arrays later on224 CALL dyn_spg_ts_init ! do it first: set nn_e used to allocate some arrays later on 225 225 ENDIF 226 226 ! -
NEMO/trunk/src/OCE/DYN/dynspg_exp.F90
r12377 r12489 49 49 !! momentum trend the surface pressure gradient : 50 50 !! (uu(rhs),vv(rhs)) = (uu(rhs),vv(rhs)) + (spgu,spgv) 51 !! where spgu = -1/r au0 d/dx(ps) = -g/e1u di( ssh(now) )52 !! spgv = -1/r au0 d/dy(ps) = -g/e2v dj( ssh(now) )51 !! where spgu = -1/rho0 d/dx(ps) = -g/e1u di( ssh(now) ) 52 !! spgv = -1/rho0 d/dy(ps) = -g/e2v dj( ssh(now) ) 53 53 !! 54 54 !! ** Action : (puu(:,:,:,Krhs),pvv(:,:,:,Krhs)) trend of horizontal velocity increased by -
NEMO/trunk/src/OCE/DYN/dynspg_ts.F90
r12377 r12489 72 72 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_adv , vn_adv !: Advection vel. at "now" barocl. step 73 73 ! 74 INTEGER, SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_ baro <= 2.5 nn_baro75 REAL(wp),SAVE :: r dtbt! Barotropic time step74 INTEGER, SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_e <= 2.5 nn_e 75 REAL(wp),SAVE :: rDt_e ! Barotropic time step 76 76 ! 77 77 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: wgtbtp1, wgtbtp2 ! 1st & 2nd weights used in time filtering of barotropic fields … … 102 102 ierr(:) = 0 103 103 ! 104 ALLOCATE( wgtbtp1(3*nn_ baro), wgtbtp2(3*nn_baro), zwz(jpi,jpj), STAT=ierr(1) )104 ALLOCATE( wgtbtp1(3*nn_e), wgtbtp2(3*nn_e), zwz(jpi,jpj), STAT=ierr(1) ) 105 105 IF( ln_dynvor_een .OR. ln_dynvor_eeT ) & 106 106 & ALLOCATE( ftnw(jpi,jpj) , ftne(jpi,jpj) , ftsw(jpi,jpj) , ftse(jpi,jpj), STAT=ierr(2) ) … … 150 150 LOGICAL :: ll_init ! =T : special startup of 2d equations 151 151 INTEGER :: noffset ! local integers : time offset for bdy update 152 REAL(wp) :: r1_ 2dt_b, z1_hu, z1_hv ! local scalars152 REAL(wp) :: r1_Dt_b, z1_hu, z1_hv ! local scalars 153 153 REAL(wp) :: za0, za1, za2, za3 ! - - 154 154 REAL(wp) :: zztmp, zldg ! - - … … 180 180 ! zwdramp = 1._wp / (rn_wdmin2 - rn_wdmin1) ! more general ramp 181 181 ! ! inverse of baroclinic time step 182 IF( kt == nit000 .AND. neuler == 0 ) THEN ; r1_2dt_b = 1._wp / ( rdt ) 183 ELSE ; r1_2dt_b = 1._wp / ( 2._wp * rdt ) 184 ENDIF 182 r1_Dt_b = 1._wp / rDt 185 183 ! 186 184 ll_init = ln_bt_av ! if no time averaging, then no specific restart 187 185 ll_fw_start = .FALSE. 188 186 ! ! time offset in steps for bdy data update 189 IF( .NOT.ln_bt_fw ) THEN ; noffset = - nn_ baro187 IF( .NOT.ln_bt_fw ) THEN ; noffset = - nn_e 190 188 ELSE ; noffset = 0 191 189 ENDIF … … 198 196 IF(lwp) WRITE(numout,*) 199 197 ! 200 IF( neuler == 0) ll_init=.TRUE.201 ! 202 IF( ln_bt_fw .OR. neuler == 0) THEN198 IF( l_1st_euler ) ll_init=.TRUE. 199 ! 200 IF( ln_bt_fw .OR. l_1st_euler ) THEN 203 201 ll_fw_start =.TRUE. 204 202 noffset = 0 … … 209 207 CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 210 208 ! 211 ENDIF 212 ! 213 ! If forward start at previous time step, and centered integration, 214 ! then update averaging weights: 215 IF (.NOT.ln_bt_fw .AND.( neuler==0 .AND. kt==nit000+1 ) ) THEN 216 ll_fw_start=.FALSE. 217 CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 218 ENDIF 219 ! 220 209 ELSEIF( kt == nit000 + 1 ) THEN !* initialisation 2nd time-step 210 ! 211 IF( .NOT.ln_bt_fw ) THEN 212 ! If we did an Euler timestep on the first timestep we need to reset ll_fw_start 213 ! and the averaging weights. We don't have an easy way of telling whether we did 214 ! an Euler timestep on the first timestep (because l_1st_euler is reset to .false. 215 ! at the end of the first timestep) so just do this in all cases. 216 ll_fw_start = .FALSE. 217 CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 218 ENDIF 219 ! 220 ENDIF 221 ! 221 222 ! ----------------------------------------------------------------------------- 222 223 ! Phase 1 : Coupling between general trend and barotropic estimates (1st step) … … 302 303 IF( ln_bt_fw ) THEN ! Add wind forcing 303 304 DO_2D_00_00 304 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_r au0 * utau(ji,jj) * r1_hu(ji,jj,Kmm)305 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_r au0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm)305 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) 306 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm) 306 307 END_2D 307 308 ELSE 308 zztmp = r1_r au0 * r1_2309 zztmp = r1_rho0 * r1_2 309 310 DO_2D_00_00 310 311 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) … … 319 320 ! ! --------------------------------------------------- ! 320 321 IF (ln_bt_fw) THEN ! FORWARD integration: use kt+1/2 fluxes (NOW+1/2) 321 zssh_frc(:,:) = r1_r au0 * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) )322 zssh_frc(:,:) = r1_rho0 * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) 322 323 ELSE ! CENTRED integration: use kt-1/2 + kt+1/2 fluxes (NOW) 323 zztmp = r1_r au0 * r1_2324 zztmp = r1_rho0 * r1_2 324 325 zssh_frc(:,:) = zztmp * ( emp(:,:) + emp_b(:,:) & 325 326 & - rnf(:,:) - rnf_b(:,:) & … … 428 429 ! Update tide potential at the beginning of current time substep 429 430 IF( ln_tide_pot .AND. ln_tide ) THEN 430 zt0substep = REAL(nsec_day, wp) - 0.5_wp*r dt + (jn + noffset - 1) * rdt / REAL(nn_baro, wp)431 zt0substep = REAL(nsec_day, wp) - 0.5_wp*rn_Dt + (jn + noffset - 1) * rn_Dt / REAL(nn_e, wp) 431 432 CALL upd_tide(zt0substep, Kmm) 432 433 END IF … … 494 495 IF( .NOT.Agrif_Root() .AND. ln_bt_fw ) CALL agrif_dyn_ts_flux( jn, zhU, zhV ) 495 496 #endif 496 IF( ln_wd_il ) CALL wad_lmt_bt(zhU, zhV, sshn_e, zssh_frc, r dtbt) !!gm wad_lmt_bt use of lbc_lnk on zhU, zhV497 IF( ln_wd_il ) CALL wad_lmt_bt(zhU, zhV, sshn_e, zssh_frc, rDt_e) !!gm wad_lmt_bt use of lbc_lnk on zhU, zhV 497 498 498 499 IF( ln_wd_dl ) THEN ! un_e and vn_e are set to zero at faces where … … 509 510 DO_2D_00_00 510 511 zhdiv = ( zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1) ) * r1_e1e2t(ji,jj) 511 ssha_e(ji,jj) = ( sshn_e(ji,jj) - r dtbt* ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj)512 ssha_e(ji,jj) = ( sshn_e(ji,jj) - rDt_e * ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj) 512 513 END_2D 513 514 ! … … 599 600 DO_2D_00_00 600 601 ua_e(ji,jj) = ( un_e(ji,jj) & 601 & + r dtbt* ( zu_spg(ji,jj) &602 & + rDt_e * ( zu_spg(ji,jj) & 602 603 & + zu_trd(ji,jj) & 603 604 & + zu_frc(ji,jj) ) & … … 605 606 606 607 va_e(ji,jj) = ( vn_e(ji,jj) & 607 & + r dtbt* ( zv_spg(ji,jj) &608 & + rDt_e * ( zv_spg(ji,jj) & 608 609 & + zv_trd(ji,jj) & 609 610 & + zv_frc(ji,jj) ) & … … 624 625 ! 625 626 ua_e(ji,jj) = ( hu_e (ji,jj) * un_e (ji,jj) & 626 & + r dtbt* ( zhu_bck * zu_spg (ji,jj) & !627 & + rDt_e * ( zhu_bck * zu_spg (ji,jj) & ! 627 628 & + zhup2_e(ji,jj) * zu_trd (ji,jj) & ! 628 629 & + hu(ji,jj,Kmm) * zu_frc (ji,jj) ) ) * z1_hu 629 630 ! 630 631 va_e(ji,jj) = ( hv_e (ji,jj) * vn_e (ji,jj) & 631 & + r dtbt* ( zhv_bck * zv_spg (ji,jj) & !632 & + rDt_e * ( zhv_bck * zv_spg (ji,jj) & ! 632 633 & + zhvp2_e(ji,jj) * zv_trd (ji,jj) & ! 633 634 & + hv(ji,jj,Kmm) * zv_frc (ji,jj) ) ) * z1_hv … … 637 638 IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 638 639 DO_2D_00_00 639 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - r dtbt* zCdU_u(ji,jj) * hur_e(ji,jj))640 va_e(ji,jj) = va_e(ji,jj) /(1.0 - r dtbt* zCdU_v(ji,jj) * hvr_e(ji,jj))640 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rDt_e * zCdU_u(ji,jj) * hur_e(ji,jj)) 641 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rDt_e * zCdU_v(ji,jj) * hvr_e(ji,jj)) 641 642 END_2D 642 643 ENDIF … … 701 702 ! Set advection velocity correction: 702 703 IF (ln_bt_fw) THEN 703 IF( .NOT.( kt == nit000 .AND. neuler==0) ) THEN704 IF( .NOT.( kt == nit000 .AND. l_1st_euler ) ) THEN 704 705 DO_2D_11_11 705 706 zun_save = un_adv(ji,jj) 706 707 zvn_save = vn_adv(ji,jj) 707 708 ! ! apply the previously computed correction 708 un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - atfp * un_bf(ji,jj) )709 vn_adv(ji,jj) = r1_2 * ( vb2_b(ji,jj) + zvn_save - atfp * vn_bf(ji,jj) )709 un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - rn_atfp * un_bf(ji,jj) ) 710 vn_adv(ji,jj) = r1_2 * ( vb2_b(ji,jj) + zvn_save - rn_atfp * vn_bf(ji,jj) ) 710 711 ! ! Update corrective fluxes for next time step 711 un_bf(ji,jj) = atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) )712 vn_bf(ji,jj) = atfp * vn_bf(ji,jj) + ( zvn_save - vb2_b(ji,jj) )712 un_bf(ji,jj) = rn_atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) ) 713 vn_bf(ji,jj) = rn_atfp * vn_bf(ji,jj) + ( zvn_save - vb2_b(ji,jj) ) 713 714 ! ! Save integrated transport for next computation 714 715 ub2_b(ji,jj) = zun_save … … 728 729 IF( ln_dynadv_vec .OR. ln_linssh ) THEN 729 730 DO jk=1,jpkm1 730 puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) ) * r1_ 2dt_b731 pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) ) * r1_ 2dt_b731 puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) ) * r1_Dt_b 732 pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) ) * r1_Dt_b 732 733 END DO 733 734 ELSE … … 744 745 ! 745 746 DO jk=1,jpkm1 746 puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_ 2dt_b747 pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_ 2dt_b747 puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_Dt_b 748 pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_Dt_b 748 749 END DO 749 750 ! Save barotropic velocities not transport: … … 808 809 LOGICAL, INTENT(in) :: ll_fw ! forward time splitting =.true. 809 810 INTEGER, INTENT(inout) :: jpit ! cycle length 810 REAL(wp), DIMENSION(3*nn_ baro), INTENT(inout) :: zwgt1, & ! Primary weights811 REAL(wp), DIMENSION(3*nn_e), INTENT(inout) :: zwgt1, & ! Primary weights 811 812 zwgt2 ! Secondary weights 812 813 … … 820 821 ! Set time index when averaged value is requested 821 822 IF (ll_fw) THEN 822 jic = nn_ baro823 jic = nn_e 823 824 ELSE 824 jic = 2 * nn_ baro825 jic = 2 * nn_e 825 826 ENDIF 826 827 … … 828 829 IF (ll_av) THEN 829 830 ! Define simple boxcar window for primary weights 830 ! (width = nn_ baro, centered around jic)831 ! (width = nn_e, centered around jic) 831 832 SELECT CASE ( nn_bt_flt ) 832 833 CASE( 0 ) ! No averaging … … 834 835 jpit = jic 835 836 836 CASE( 1 ) ! Boxcar, width = nn_ baro837 DO jn = 1, 3*nn_ baro838 za1 = ABS(float(jn-jic))/float(nn_ baro)837 CASE( 1 ) ! Boxcar, width = nn_e 838 DO jn = 1, 3*nn_e 839 za1 = ABS(float(jn-jic))/float(nn_e) 839 840 IF (za1 < 0.5_wp) THEN 840 841 zwgt1(jn) = 1._wp … … 843 844 ENDDO 844 845 845 CASE( 2 ) ! Boxcar, width = 2 * nn_ baro846 DO jn = 1, 3*nn_ baro847 za1 = ABS(float(jn-jic))/float(nn_ baro)846 CASE( 2 ) ! Boxcar, width = 2 * nn_e 847 DO jn = 1, 3*nn_e 848 za1 = ABS(float(jn-jic))/float(nn_e) 848 849 IF (za1 < 1._wp) THEN 849 850 zwgt1(jn) = 1._wp … … 889 890 IF( TRIM(cdrw) == 'READ' ) THEN ! Read/initialise 890 891 ! ! --------------- 891 IF( ln_rstart .AND. ln_bt_fw .AND. ( neuler/=0) ) THEN !* Read the restart file892 IF( ln_rstart .AND. ln_bt_fw .AND. (.NOT.l_1st_euler) ) THEN !* Read the restart file 892 893 CALL iom_get( numror, jpdom_autoglo, 'ub2_b' , ub2_b (:,:), ldxios = lrxios ) 893 894 CALL iom_get( numror, jpdom_autoglo, 'vb2_b' , vb2_b (:,:), ldxios = lrxios ) … … 975 976 976 977 ! Estimate number of iterations to satisfy a max courant number= rn_bt_cmax 977 IF( ln_bt_auto ) nn_ baro = CEILING( rdt / rn_bt_cmax * zcmax)978 IF( ln_bt_auto ) nn_e = CEILING( rn_Dt / rn_bt_cmax * zcmax) 978 979 979 r dtbt = rdt / REAL( nn_baro, wp )980 zcmax = zcmax * r dtbt980 rDt_e = rn_Dt / REAL( nn_e , wp ) 981 zcmax = zcmax * rDt_e 981 982 ! Print results 982 983 IF(lwp) WRITE(numout,*) … … 984 985 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~' 985 986 IF( ln_bt_auto ) THEN 986 IF(lwp) WRITE(numout,*) ' ln_ts_auto =.true. Automatically set nn_ baro'987 IF(lwp) WRITE(numout,*) ' ln_ts_auto =.true. Automatically set nn_e ' 987 988 IF(lwp) WRITE(numout,*) ' Max. courant number allowed: ', rn_bt_cmax 988 989 ELSE 989 IF(lwp) WRITE(numout,*) ' ln_ts_auto=.false.: Use nn_ baro in namelist nn_baro = ', nn_baro990 IF(lwp) WRITE(numout,*) ' ln_ts_auto=.false.: Use nn_e in namelist nn_e = ', nn_e 990 991 ENDIF 991 992 992 993 IF(ln_bt_av) THEN 993 IF(lwp) WRITE(numout,*) ' ln_bt_av =.true. ==> Time averaging over nn_ barotime steps is on '994 IF(lwp) WRITE(numout,*) ' ln_bt_av =.true. ==> Time averaging over nn_e time steps is on ' 994 995 ELSE 995 996 IF(lwp) WRITE(numout,*) ' ln_bt_av =.false. => No time averaging of barotropic variables ' … … 1011 1012 SELECT CASE ( nn_bt_flt ) 1012 1013 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' Dirac' 1013 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = nn_ baro'1014 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = 2*nn_ baro'1014 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = nn_e' 1015 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = 2*nn_e' 1015 1016 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_bt_flt: should 0,1, or 2' ) 1016 1017 END SELECT 1017 1018 ! 1018 1019 IF(lwp) WRITE(numout,*) ' ' 1019 IF(lwp) WRITE(numout,*) ' nn_ baro = ', nn_baro1020 IF(lwp) WRITE(numout,*) ' Barotropic time step [s] is :', r dtbt1020 IF(lwp) WRITE(numout,*) ' nn_e = ', nn_e 1021 IF(lwp) WRITE(numout,*) ' Barotropic time step [s] is :', rDt_e 1021 1022 IF(lwp) WRITE(numout,*) ' Maximum Courant number is :', zcmax 1022 1023 ! … … 1030 1031 ENDIF 1031 1032 IF( zcmax>0.9_wp ) THEN 1032 CALL ctl_stop( 'dynspg_ts ERROR: Maximum Courant number is greater than 0.9: Inc. nn_ baro!' )1033 CALL ctl_stop( 'dynspg_ts ERROR: Maximum Courant number is greater than 0.9: Inc. nn_e !' ) 1033 1034 ENDIF 1034 1035 ! … … 1429 1430 ! 1430 1431 IF( ln_wd_il ) THEN ! W/D : use the "clipped" bottom friction !!gm explain WHY, please ! 1431 zztmp = -1._wp / r dtbt1432 zztmp = -1._wp / rDt_e 1432 1433 DO_2D_00_00 1433 1434 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) * wdrampu(ji,jj) * MAX( & -
NEMO/trunk/src/OCE/DYN/dynzdf.F90
r12377 r12489 92 92 ENDIF 93 93 ENDIF 94 ! !* set time step95 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt = rdt ! = rdt (restart with Euler time stepping)96 ELSEIF( kt <= nit000 + 1 ) THEN ; r2dt = 2. * rdt ! = 2 rdt (leapfrog)97 ENDIF98 !99 94 ! !* explicit top/bottom drag case 100 95 IF( .NOT.ln_drgimp ) CALL zdf_drg_exp( kt, Kmm, puu(:,:,:,Kbb), pvv(:,:,:,Kbb), puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add top/bottom friction trend to (puu(Kaa),pvv(Kaa)) … … 112 107 IF( ln_dynadv_vec .OR. ln_linssh ) THEN ! applied on velocity 113 108 DO jk = 1, jpkm1 114 puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kbb) + r 2dt * puu(:,:,jk,Krhs) ) * umask(:,:,jk)115 pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kbb) + r 2dt * pvv(:,:,jk,Krhs) ) * vmask(:,:,jk)109 puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kbb) + rDt * puu(:,:,jk,Krhs) ) * umask(:,:,jk) 110 pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kbb) + rDt * pvv(:,:,jk,Krhs) ) * vmask(:,:,jk) 116 111 END DO 117 112 ELSE ! applied on thickness weighted velocity 118 113 DO jk = 1, jpkm1 119 114 puu(:,:,jk,Kaa) = ( e3u(:,:,jk,Kbb) * puu(:,:,jk,Kbb) & 120 & + r 2dt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs) ) / e3u(:,:,jk,Kaa) * umask(:,:,jk)115 & + rDt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs) ) / e3u(:,:,jk,Kaa) * umask(:,:,jk) 121 116 pvv(:,:,jk,Kaa) = ( e3v(:,:,jk,Kbb) * pvv(:,:,jk,Kbb) & 122 & + r 2dt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs) ) / e3v(:,:,jk,Kaa) * vmask(:,:,jk)117 & + rDt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs) ) / e3v(:,:,jk,Kaa) * vmask(:,:,jk) 123 118 END DO 124 119 ENDIF … … 138 133 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 139 134 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 140 puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + r 2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua141 pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + r 2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va135 puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 136 pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 142 137 END_2D 143 138 IF( ln_isfcav ) THEN ! Ocean cavities (ISF) … … 147 142 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 148 143 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 149 puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + r 2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua150 pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + r 2dt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va144 puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 145 pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 151 146 END_2D 152 147 END IF … … 156 151 ! 157 152 ! !* Matrix construction 158 zdt = r 2dt * 0.5153 zdt = rDt * 0.5 159 154 IF( ln_zad_Aimp ) THEN !! 160 155 SELECT CASE( nldf_dyn ) … … 232 227 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 233 228 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) ! after scale factor at T-point 234 zwd(ji,jj,iku) = zwd(ji,jj,iku) - r 2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua229 zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 235 230 END_2D 236 231 IF ( ln_isfcav ) THEN ! top friction (always implicit) … … 239 234 iku = miku(ji,jj) ! ocean top level at u- and v-points 240 235 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) ! after scale factor at T-point 241 zwd(ji,jj,iku) = zwd(ji,jj,iku) - r 2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua236 zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 242 237 END_2D 243 238 END IF … … 265 260 DO_2D_00_00 266 261 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa) 267 puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + r 2dt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) ) &268 & / ( ze3ua * r au0 ) * umask(ji,jj,1)262 puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) ) & 263 & / ( ze3ua * rho0 ) * umask(ji,jj,1) 269 264 END_2D 270 265 DO_3D_00_00( 2, jpkm1 ) … … 282 277 ! 283 278 ! !* Matrix construction 284 zdt = r 2dt * 0.5279 zdt = rDt * 0.5 285 280 IF( ln_zad_Aimp ) THEN !! 286 281 SELECT CASE( nldf_dyn ) … … 357 352 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 358 353 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) ! after scale factor at T-point 359 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - r 2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va354 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va 360 355 END_2D 361 356 IF ( ln_isfcav ) THEN … … 363 358 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 364 359 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) ! after scale factor at T-point 365 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - r 2dt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va360 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va 366 361 END_2D 367 362 ENDIF … … 389 384 DO_2D_00_00 390 385 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa) 391 pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + r 2dt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) ) &392 & / ( ze3va * r au0 ) * vmask(ji,jj,1)386 pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) ) & 387 & / ( ze3va * rho0 ) * vmask(ji,jj,1) 393 388 END_2D 394 389 DO_3D_00_00( 2, jpkm1 ) … … 404 399 ! 405 400 IF( l_trddyn ) THEN ! save the vertical diffusive trends for further diagnostics 406 ztrdu(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) / r 2dt - ztrdu(:,:,:)407 ztrdv(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) / r 2dt - ztrdv(:,:,:)401 ztrdu(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) / rDt - ztrdu(:,:,:) 402 ztrdv(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) / rDt - ztrdv(:,:,:) 408 403 CALL trd_dyn( ztrdu, ztrdv, jpdyn_zdf, kt, Kmm ) 409 404 DEALLOCATE( ztrdu, ztrdv ) -
NEMO/trunk/src/OCE/DYN/sshwzv.F90
r12377 r12489 75 75 REAL(wp), DIMENSION(jpi,jpj,jpt), INTENT(inout) :: pssh ! sea-surface height 76 76 ! 77 INTEGER :: jk ! dummy loop indice78 REAL(wp) :: z 2dt, zcoef ! local scalars77 INTEGER :: jk ! dummy loop index 78 REAL(wp) :: zcoef ! local scalar 79 79 REAL(wp), DIMENSION(jpi,jpj) :: zhdiv ! 2D workspace 80 80 !!---------------------------------------------------------------------- … … 88 88 ENDIF 89 89 ! 90 z2dt = 2._wp * rdt ! set time step size (Euler/Leapfrog) 91 IF( neuler == 0 .AND. kt == nit000 ) z2dt = rdt 92 zcoef = 0.5_wp * r1_rau0 90 zcoef = 0.5_wp * r1_rho0 93 91 94 92 ! !------------------------------! … … 96 94 ! !------------------------------! 97 95 IF(ln_wd_il) THEN 98 CALL wad_lmt(pssh(:,:,Kbb), zcoef * (emp_b(:,:) + emp(:,:)), z2dt, Kmm, uu, vv )96 CALL wad_lmt(pssh(:,:,Kbb), zcoef * (emp_b(:,:) + emp(:,:)), rDt, Kmm, uu, vv ) 99 97 ENDIF 100 98 … … 109 107 ! compute the vertical velocity which can be used to compute the non-linear terms of the momentum equations. 110 108 ! 111 pssh(:,:,Kaa) = ( pssh(:,:,Kbb) - z2dt * ( zcoef * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) ) ) * ssmask(:,:)109 pssh(:,:,Kaa) = ( pssh(:,:,Kbb) - rDt * ( zcoef * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) ) ) * ssmask(:,:) 112 110 ! 113 111 #if defined key_agrif … … 152 150 ! 153 151 INTEGER :: ji, jj, jk ! dummy loop indices 154 REAL(wp) :: z1_2dt ! local scalars155 152 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zhdiv 156 153 !!---------------------------------------------------------------------- … … 168 165 ! ! Now Vertical Velocity ! 169 166 ! !------------------------------! 170 z1_2dt = 1. / ( 2. * rdt ) ! set time step size (Euler/Leapfrog)171 IF( neuler == 0 .AND. kt == nit000 ) z1_2dt = 1. / rdt172 167 ! 173 168 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! z_tilde and layer cases … … 187 182 ! computation of w 188 183 pww(:,:,jk) = pww(:,:,jk+1) - ( e3t(:,:,jk,Kmm) * hdiv(:,:,jk) + zhdiv(:,:,jk) & 189 & + z1_2dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) ) ) * tmask(:,:,jk)184 & + r1_Dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) ) ) * tmask(:,:,jk) 190 185 END DO 191 186 ! IF( ln_vvl_layer ) pww(:,:,:) = 0.e0 … … 195 190 ! computation of w 196 191 pww(:,:,jk) = pww(:,:,jk+1) - ( e3t(:,:,jk,Kmm) * hdiv(:,:,jk) & 197 & + z1_2dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) ) ) * tmask(:,:,jk)192 & + r1_Dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) ) ) * tmask(:,:,jk) 198 193 END DO 199 194 ENDIF … … 227 222 !! ** Method : - apply Asselin time fiter to now ssh (excluding the forcing 228 223 !! from the filter, see Leclair and Madec 2010) and swap : 229 !! pssh(:,:,Kmm) = pssh(:,:,Kaa) + atfp * ( pssh(:,:,Kbb) -2 pssh(:,:,Kmm) + pssh(:,:,Kaa) )230 !! - atfp * rdt * ( emp_b - emp ) / rau0224 !! pssh(:,:,Kmm) = pssh(:,:,Kaa) + rn_atfp * ( pssh(:,:,Kbb) -2 pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 225 !! - rn_atfp * rn_Dt * ( emp_b - emp ) / rho0 231 226 !! 232 227 !! ** action : - pssh(:,:,Kmm) time filtered … … 249 244 ENDIF 250 245 ! !== Euler time-stepping: no filter, just swap ==! 251 IF ( .NOT.( neuler == 0 .AND. kt == nit000) ) THEN ! Only do time filtering for leapfrog timesteps246 IF ( .NOT.( l_1st_euler ) ) THEN ! Only do time filtering for leapfrog timesteps 252 247 ! ! filtered "now" field 253 pssh(:,:,Kmm) = pssh(:,:,Kmm) + atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) )248 pssh(:,:,Kmm) = pssh(:,:,Kmm) + rn_atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 254 249 IF( .NOT.ln_linssh ) THEN ! "now" <-- with forcing removed 255 zcoef = atfp * rdt * r1_rau0250 zcoef = rn_atfp * rn_Dt * r1_rho0 256 251 pssh(:,:,Kmm) = pssh(:,:,Kmm) - zcoef * ( emp_b(:,:) - emp (:,:) & 257 252 & - rnf_b(:,:) + rnf (:,:) & … … 260 255 261 256 ! ice sheet coupling 262 IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) pssh(:,:,Kbb) = pssh(:,:,Kbb) - atfp * rdt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:)257 IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) pssh(:,:,Kbb) = pssh(:,:,Kbb) - rn_atfp * rn_Dt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 263 258 264 259 ENDIF … … 311 306 DO_3D_00_00( 1, jpkm1 ) 312 307 z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 313 ! 2*r dt and not r2dt (for restartability)314 Cu_adv(ji,jj,jk) = 2._wp * r dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) ) &308 ! 2*rn_Dt and not rDt (for restartability) 309 Cu_adv(ji,jj,jk) = 2._wp * rn_Dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) ) & 315 310 & + ( MAX( e2u(ji ,jj)*e3u(ji ,jj,jk,Kmm)*uu(ji ,jj,jk,Kmm) + un_td(ji ,jj,jk), 0._wp ) - & 316 311 & MIN( e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kmm)*uu(ji-1,jj,jk,Kmm) + un_td(ji-1,jj,jk), 0._wp ) ) & … … 324 319 DO_3D_00_00( 1, jpkm1 ) 325 320 z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 326 ! 2*r dt and not r2dt (for restartability)327 Cu_adv(ji,jj,jk) = 2._wp * r dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) ) &321 ! 2*rn_Dt and not rDt (for restartability) 322 Cu_adv(ji,jj,jk) = 2._wp * rn_Dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) ) & 328 323 & + ( MAX( e2u(ji ,jj)*e3u(ji ,jj,jk,Kmm)*uu(ji ,jj,jk,Kmm), 0._wp ) - & 329 324 & MIN( e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kmm)*uu(ji-1,jj,jk,Kmm), 0._wp ) ) & -
NEMO/trunk/src/OCE/DYN/wet_dry.F90
r12377 r12489 270 270 271 271 272 SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, r dtbt)272 SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rDt_e ) 273 273 !!---------------------------------------------------------------------- 274 274 !! *** ROUTINE wad_lmt *** … … 280 280 !! ** Action : - calculate flux limiter and W/D flag 281 281 !!---------------------------------------------------------------------- 282 REAL(wp) , INTENT(in ) :: r dtbt! ocean time-step index282 REAL(wp) , INTENT(in ) :: rDt_e ! ocean time-step index 283 283 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zflxu, zflxv, sshn_e, zssh_frc 284 284 ! … … 299 299 zdepwd = 50._wp ! maximum depth that ocean cells can have W/D processes 300 300 ! 301 z2dt = r dtbt301 z2dt = rDt_e 302 302 ! 303 303 zflxp(:,:) = 0._wp -
NEMO/trunk/src/OCE/FLO/flo4rk.F90
r12377 r12489 130 130 ! computation of Runge-Kutta factor 131 131 DO jfl = 1, jpnfl 132 zrkxfl(jfl,jind) = r dt*zufl(jfl)133 zrkyfl(jfl,jind) = r dt*zvfl(jfl)134 zrkzfl(jfl,jind) = r dt*zwfl(jfl)132 zrkxfl(jfl,jind) = rn_Dt*zufl(jfl) 133 zrkyfl(jfl,jind) = rn_Dt*zvfl(jfl) 134 zrkzfl(jfl,jind) = rn_Dt*zwfl(jfl) 135 135 END DO 136 136 IF( jind /= 4 ) THEN -
NEMO/trunk/src/OCE/FLO/floblk.F90
r12377 r12489 233 233 ! test to know if the "age" of the float is not bigger than the 234 234 ! time step 235 IF( zagenewfl(jfl) > r dt ) THEN236 zttfl(jfl) = (r dt-zagefl(jfl)) / zvol237 zagenewfl(jfl) = r dt235 IF( zagenewfl(jfl) > rn_Dt ) THEN 236 zttfl(jfl) = (rn_Dt-zagefl(jfl)) / zvol 237 zagenewfl(jfl) = rn_Dt 238 238 ENDIF 239 239 … … 340 340 ifin = 1 341 341 DO jfl = 1, jpnfl 342 IF( zagefl(jfl) < r dt ) ifin = 0342 IF( zagefl(jfl) < rn_Dt ) ifin = 0 343 343 tpifl(jfl) = zgifl(jfl) + 0.5 344 344 tpjfl(jfl) = zgjfl(jfl) + 0.5 … … 347 347 ifin = 1 348 348 DO jfl = 1, jpnfl 349 IF( zagefl(jfl) < r dt ) ifin = 0349 IF( zagefl(jfl) < rn_Dt ) ifin = 0 350 350 tpifl(jfl) = zgifl(jfl) + 0.5 351 351 tpjfl(jfl) = zgjfl(jfl) + 0.5 -
NEMO/trunk/src/OCE/FLO/flowri.F90
r12377 r12489 122 122 ztem(jfl) = ts(iafloc,ibfloc,icfl,jp_tem,Kmm) 123 123 zsal (jfl) = ts(iafloc,ibfloc,icfl,jp_sal,Kmm) 124 zrho (jfl) = (rhd(iafloc,ibfloc,icfl)+1)*r au0124 zrho (jfl) = (rhd(iafloc,ibfloc,icfl)+1)*rho0 125 125 126 126 ENDIF … … 142 142 ztem(jfl) = ts(iafloc,ibfloc,icfl,jp_tem,Kmm) 143 143 zsal(jfl) = ts(iafloc,ibfloc,icfl,jp_sal,Kmm) 144 zrho(jfl) = (rhd(iafloc,ibfloc,icfl)+1)*r au0144 zrho(jfl) = (rhd(iafloc,ibfloc,icfl)+1)*rho0 145 145 146 146 ENDIF … … 245 245 !------------------------------- 246 246 irec = INT( (kt-nn_it000+1)/nn_writefl ) +1 247 ztime = ( kt-nn_it000 + 1 ) * r dt247 ztime = ( kt-nn_it000 + 1 ) * rn_Dt 248 248 249 249 CALL flioputv( numflo , 'time_counter', ztime , start=(/irec/) ) -
NEMO/trunk/src/OCE/ICB/icbini.F90
r12472 r12489 60 60 !! - setup either test icebergs or calving file 61 61 !!---------------------------------------------------------------------- 62 REAL(wp), INTENT(in) :: pdt ! iceberg time-step (r dt*nn_fsbc)62 REAL(wp), INTENT(in) :: pdt ! iceberg time-step (rn_Dt*nn_fsbc) 63 63 INTEGER , INTENT(in) :: kt ! time step number 64 64 ! -
NEMO/trunk/src/OCE/ICB/icbtrj.F90
r10068 r12489 74 74 75 75 ! compute end time step date 76 zfjulday = fjulday + r dt / rday * REAL( nitend - nit000 + 1 , wp)76 zfjulday = fjulday + rn_Dt / rday * REAL( nitend - nit000 + 1 , wp) 77 77 IF( ABS(zfjulday - REAL(NINT(zfjulday),wp)) < 0.1 / rday ) zfjulday = REAL(NINT(zfjulday),wp) ! avoid truncation error 78 78 CALL ju2ymds( zfjulday, iyear, imonth, iday, zsec ) -
NEMO/trunk/src/OCE/IOM/iom.F90
r12377 r12489 274 274 ! 275 275 ! set time step length 276 dtime%second = r dt276 dtime%second = rn_Dt 277 277 CALL xios_set_timestep( dtime ) 278 278 ! … … 410 410 IF(cdmdl == "OPA") THEN 411 411 !from restart.F90 412 CALL iom_set_rstw_var_active("r dt")412 CALL iom_set_rstw_var_active("rn_Dt") 413 413 IF ( .NOT. ln_diurnal_only ) THEN 414 414 CALL iom_set_rstw_var_active('ub' ) … … 448 448 449 449 i = 0 450 i = i + 1; fields(i)%vname="r dt"; fields(i)%grid="grid_scalar"450 i = i + 1; fields(i)%vname="rn_Dt"; fields(i)%grid="grid_scalar" 451 451 i = i + 1; fields(i)%vname="un"; fields(i)%grid="grid_N_3D" 452 452 i = i + 1; fields(i)%vname="ub"; fields(i)%grid="grid_N_3D" … … 2358 2358 idx = INDEX(clname,'@startdate@') + INDEX(clname,'@STARTDATE@') 2359 2359 DO WHILE ( idx /= 0 ) 2360 cldate = iom_sdate( fjulday - r dt / rday )2360 cldate = iom_sdate( fjulday - rn_Dt / rday ) 2361 2361 clname = clname(1:idx-1)//TRIM(cldate)//clname(idx+11:LEN_TRIM(clname)) 2362 2362 idx = INDEX(clname,'@startdate@') + INDEX(clname,'@STARTDATE@') … … 2365 2365 idx = INDEX(clname,'@startdatefull@') + INDEX(clname,'@STARTDATEFULL@') 2366 2366 DO WHILE ( idx /= 0 ) 2367 cldate = iom_sdate( fjulday - r dt / rday, ldfull = .TRUE. )2367 cldate = iom_sdate( fjulday - rn_Dt / rday, ldfull = .TRUE. ) 2368 2368 clname = clname(1:idx-1)//TRIM(cldate)//clname(idx+15:LEN_TRIM(clname)) 2369 2369 idx = INDEX(clname,'@startdatefull@') + INDEX(clname,'@STARTDATEFULL@') … … 2372 2372 idx = INDEX(clname,'@enddate@') + INDEX(clname,'@ENDDATE@') 2373 2373 DO WHILE ( idx /= 0 ) 2374 cldate = iom_sdate( fjulday + r dt / rday * REAL( nitend - nit000, wp ), ld24 = .TRUE. )2374 cldate = iom_sdate( fjulday + rn_Dt / rday * REAL( nitend - nit000, wp ), ld24 = .TRUE. ) 2375 2375 clname = clname(1:idx-1)//TRIM(cldate)//clname(idx+9:LEN_TRIM(clname)) 2376 2376 idx = INDEX(clname,'@enddate@') + INDEX(clname,'@ENDDATE@') … … 2379 2379 idx = INDEX(clname,'@enddatefull@') + INDEX(clname,'@ENDDATEFULL@') 2380 2380 DO WHILE ( idx /= 0 ) 2381 cldate = iom_sdate( fjulday + r dt / rday * REAL( nitend - nit000, wp ), ld24 = .TRUE., ldfull = .TRUE. )2381 cldate = iom_sdate( fjulday + rn_Dt / rday * REAL( nitend - nit000, wp ), ld24 = .TRUE., ldfull = .TRUE. ) 2382 2382 clname = clname(1:idx-1)//TRIM(cldate)//clname(idx+13:LEN_TRIM(clname)) 2383 2383 idx = INDEX(clname,'@enddatefull@') + INDEX(clname,'@ENDDATEFULL@') -
NEMO/trunk/src/OCE/IOM/restart.F90
r12377 r12489 144 144 !!---------------------------------------------------------------------- 145 145 IF(lwxios) CALL iom_swap( cwxios_context ) 146 CALL iom_rstput( kt, nitrst, numrow, 'rdt' , r dt , ldxios = lwxios) ! dynamics time step146 CALL iom_rstput( kt, nitrst, numrow, 'rdt' , rn_Dt , ldxios = lwxios) ! dynamics time step 147 147 CALL iom_delay_rst( 'WRITE', 'OCE', numrow ) ! save only ocean delayed global communication variables 148 148 … … 247 247 IF( iom_varid( numror, 'rdt', ldstop = .FALSE. ) > 0 ) THEN 248 248 CALL iom_get( numror, 'rdt', zrdt, ldxios = lrxios ) 249 IF( zrdt /= rdt ) neuler = 0 249 IF( zrdt /= rn_Dt ) THEN 250 IF(lwp) WRITE( numout,*) 251 IF(lwp) WRITE( numout,*) 'rst_read: rdt not equal to the read one' 252 IF(lwp) WRITE( numout,*) 253 IF(lwp) WRITE( numout,*) ' ==>>> forced euler first time-step' 254 l_1st_euler = .TRUE. 255 ENDIF 250 256 ENDIF 251 257 … … 256 262 IF ( ln_diurnal_only ) THEN 257 263 IF(lwp) WRITE( numout, * ) & 258 & "rst_read:- ln_diurnal_only set, setting rhop=r au0"259 rhop = r au0264 & "rst_read:- ln_diurnal_only set, setting rhop=rho0" 265 rhop = rho0 260 266 CALL iom_get( numror, jpdom_autoglo, 'tn' , w3d, ldxios = lrxios ) 261 267 ts(:,:,1,jp_tem,Kmm) = w3d(:,:,1) … … 270 276 CALL iom_get( numror, jpdom_autoglo, 'sshb' ,ssh(:,: ,Kbb), ldxios = lrxios ) 271 277 ELSE 272 neuler = 0278 l_1st_euler = .TRUE. ! before field not found, forced euler 1st time-step 273 279 ENDIF 274 280 ! … … 284 290 ENDIF 285 291 ! 286 IF( neuler == 0 ) THEN ! Euler restart (neuler=0)292 IF( l_1st_euler ) THEN ! Euler restart 287 293 ts (:,:,:,:,Kbb) = ts (:,:,:,:,Kmm) ! all before fields set to now values 288 294 uu (:,:,: ,Kbb) = uu (:,:,: ,Kmm) -
NEMO/trunk/src/OCE/ISF/isfcav.F90
r12343 r12489 24 24 USE oce , ONLY: ts ! ocean tracers 25 25 USE par_oce , ONLY: jpi,jpj ! ocean space and time domain 26 USE phycst , ONLY: grav,r au0,rau0_rcp,r1_rau0_rcp ! physical constants26 USE phycst , ONLY: grav,rho0,rho0_rcp,r1_rho0_rcp ! physical constants 27 27 USE eosbn2 , ONLY: ln_teos10 ! use ln_teos10 or not 28 28 ! … … 85 85 ! 86 86 ! initialisation 87 IF (TRIM(cn_gammablk) == 'vel_stab' ) zqoce_b (:,:) = ptsc(:,:,jp_tem) * r au0_rcp ! last time step total heat fluxes (to speed up convergence)87 IF (TRIM(cn_gammablk) == 'vel_stab' ) zqoce_b (:,:) = ptsc(:,:,jp_tem) * rho0_rcp ! last time step total heat fluxes (to speed up convergence) 88 88 ! 89 89 ! compute ice shelf melting … … 142 142 ! 143 143 ! set temperature content 144 ptsc(:,:,jp_tem) = - zqh(:,:) * r1_r au0_rcp144 ptsc(:,:,jp_tem) = - zqh(:,:) * r1_rho0_rcp 145 145 ! 146 146 ! write restart variables (qoceisf, qhcisf, fwfisf for now and before) … … 215 215 risf_lamb1 =-0.0564_wp 216 216 risf_lamb2 = 0.0773_wp 217 risf_lamb3 =-7.8633e-8 * grav * r au0217 risf_lamb3 =-7.8633e-8 * grav * rho0 218 218 ELSE ! linearisation from table 4 (Asay-Davis et al., 2015) 219 219 risf_lamb1 =-0.0573_wp 220 220 risf_lamb2 = 0.0832_wp 221 risf_lamb3 =-7.5300e-8 * grav * r au0221 risf_lamb3 =-7.5300e-8 * grav * rho0 222 222 ENDIF 223 223 -
NEMO/trunk/src/OCE/ISF/isfcavmlt.F90
r12340 r12489 17 17 18 18 USE dom_oce ! ocean space and time domain 19 USE phycst , ONLY: rcp, r au0, rau0_rcp ! physical constants19 USE phycst , ONLY: rcp, rho0, rho0_rcp ! physical constants 20 20 USE eosbn2 , ONLY: eos_fzp ! equation of state 21 21 … … 161 161 ! 162 162 ! compute ocean-ice heat flux and then derive fwf assuming that ocean heat flux equal latent heat 163 pqfwf(:,:) = - pgt(:,:) * r au0_rcp * zthd(:,:) / rLfusisf ! fresh water flux ( > 0 out )163 pqfwf(:,:) = - pgt(:,:) * rho0_rcp * zthd(:,:) / rLfusisf ! fresh water flux ( > 0 out ) 164 164 pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocea-ice flux ( > 0 out ) 165 165 pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux ( > 0 out ) … … 213 213 ! 214 214 ! compute coeficient to solve the 2nd order equation 215 zeps1 = r au0_rcp * pgt(ji,jj)216 zeps2 = rLfusisf * r au0 * pgs(ji,jj)215 zeps1 = rho0_rcp * pgt(ji,jj) 216 zeps2 = rLfusisf * rho0 * pgs(ji,jj) 217 217 zeps3 = rhoisf * rcpisf * rkappa / MAX(risfdep(ji,jj),zeps) 218 218 zeps4 = risf_lamb2 + risf_lamb3 * risfdep(ji,jj) … … 238 238 ! 239 239 ! compute the upward water and heat flux (eq. 24 and eq. 26) 240 pqfwf(ji,jj) = r au0 * pgs(ji,jj) * ( zsfrz - pstbl(ji,jj) ) / MAX(zsfrz,zeps) ! fresh water flux (> 0 out)241 pqoce(ji,jj) = r au0_rcp * pgt(ji,jj) * zthd (ji,jj) ! ocean-ice heat flux (> 0 out)240 pqfwf(ji,jj) = rho0 * pgs(ji,jj) * ( zsfrz - pstbl(ji,jj) ) / MAX(zsfrz,zeps) ! fresh water flux (> 0 out) 241 pqoce(ji,jj) = rho0_rcp * pgt(ji,jj) * zthd (ji,jj) ! ocean-ice heat flux (> 0 out) 242 242 pqhc (ji,jj) = rcp * pqfwf(ji,jj) * ztfrz(ji,jj) ! heat content flux (> 0 out) 243 243 ! -
NEMO/trunk/src/OCE/ISF/isfcpl.F90
r12353 r12489 68 68 ! 69 69 ! start on an euler time step 70 neuler = 070 l_1st_euler = .TRUE. 71 71 ! 72 72 ! allocation and initialisation to 0 … … 502 502 ! compute run length 503 503 nstp_iscpl = nitend - nit000 + 1 504 rdt_iscpl = nstp_iscpl * rn_ rdt504 rdt_iscpl = nstp_iscpl * rn_Dt 505 505 z1_rdtiscpl = 1._wp / rdt_iscpl 506 506 -
NEMO/trunk/src/OCE/ISF/isfdynatf.F90
r12372 r12489 13 13 USE isf_oce 14 14 15 USE phycst , ONLY: r1_r au0 ! physical constant15 USE phycst , ONLY: r1_rho0 ! physical constant 16 16 USE dom_oce, ONLY: tmask, ssmask, ht, e3t, r1_e1e2t ! time and space domain 17 17 … … 39 39 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pe3t_f ! time filtered scale factor to be corrected 40 40 ! 41 REAL(wp) , INTENT(in ) :: pcoef ! atfp * rdt * r1_rau041 REAL(wp) , INTENT(in ) :: pcoef ! rn_atfp * rn_Dt * r1_rho0 42 42 !!-------------------------------------------------------------------- 43 43 INTEGER :: jk ! loop index … … 70 70 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: pfrac, phtbl ! fraction of bottom cell included in tbl, tbl thickness 71 71 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: pfwf , pfwf_b ! now/before fwf 72 REAL(wp), INTENT(in ) :: pcoef ! atfp * rdt * r1_rau072 REAL(wp), INTENT(in ) :: pcoef ! rn_atfp * rn_Dt * r1_rho0 73 73 !!---------------------------------------------------------------------- 74 74 INTEGER :: ji,jj,jk … … 77 77 ! 78 78 ! compute fwf conservation correction 79 zfwfinc(:,:) = pcoef * ( pfwf_b(:,:) - pfwf(:,:) ) / ( ht(:,:) + 1._wp - ssmask(:,:) ) * r1_r au079 zfwfinc(:,:) = pcoef * ( pfwf_b(:,:) - pfwf(:,:) ) / ( ht(:,:) + 1._wp - ssmask(:,:) ) * r1_rho0 80 80 ! 81 81 ! add the increment -
NEMO/trunk/src/OCE/ISF/isfhdiv.F90
r12340 r12489 16 16 17 17 USE dom_oce ! time and space domain 18 USE phycst , ONLY: r1_r au0 ! physical constant18 USE phycst , ONLY: r1_rho0 ! physical constant 19 19 USE in_out_manager ! 20 20 … … 96 96 ! 97 97 ! compute integrated divergence correction 98 zhdiv(:,:) = 0.5_wp * ( pfwf(:,:) + pfwf_b(:,:) ) * r1_r au0 / phtbl(:,:)98 zhdiv(:,:) = 0.5_wp * ( pfwf(:,:) + pfwf_b(:,:) ) * r1_rho0 / phtbl(:,:) 99 99 ! 100 100 ! update divergence at each level affected by ice shelf top boundary layer -
NEMO/trunk/src/OCE/ISF/isfpar.F90
r12077 r12489 24 24 USE dom_oce , ONLY: bathy ! ocean space and time domain 25 25 USE par_oce , ONLY: jpi,jpj ! ocean space and time domain 26 USE phycst , ONLY: r1_r au0_rcp ! physical constants26 USE phycst , ONLY: r1_rho0_rcp ! physical constants 27 27 ! 28 28 USE in_out_manager ! I/O manager … … 88 88 ! 89 89 ! set temperature content 90 ptsc(:,:,jp_tem) = zqh(:,:) * r1_r au0_rcp90 ptsc(:,:,jp_tem) = zqh(:,:) * r1_rho0_rcp 91 91 ! 92 92 ! write restart variables (qoceisf, qhcisf, fwfisf for now and before) -
NEMO/trunk/src/OCE/ISF/isfparmlt.F90
r12077 r12489 13 13 USE dom_oce ! ocean space and time domain 14 14 USE oce , ONLY: ts ! ocean dynamics and tracers 15 USE phycst , ONLY: rcp, r au0 ! physical constants15 USE phycst , ONLY: rcp, rho0 ! physical constants 16 16 USE eosbn2 , ONLY: eos_fzp ! equation of state 17 17 … … 148 148 ! 149 149 ! 2. ------------Net heat flux and fresh water flux due to the ice shelf 150 pqoce(:,:) = r au0 * rcp * rn_gammat0 * risfLeff(:,:) * e1t(:,:) * ( ztavg(:,:) - ztfrz(:,:) ) * r1_e1e2t(:,:)150 pqoce(:,:) = rho0 * rcp * rn_gammat0 * risfLeff(:,:) * e1t(:,:) * ( ztavg(:,:) - ztfrz(:,:) ) * r1_e1e2t(:,:) 151 151 pqfwf(:,:) = - pqoce(:,:) / rLfusisf ! derived from the latent heat flux 152 152 pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux -
NEMO/trunk/src/OCE/LDF/ldfdyn.F90
r12377 r12489 407 407 zcmsmag = (rn_csmc/rpi)**2 ! (C_smag/pi)^2 408 408 zstabf_lo = rn_minfac * rn_minfac / ( 2._wp * 12._wp * 12._wp * zcmsmag ) ! lower limit stability factor scaling 409 zstabf_up = rn_maxfac / ( 4._wp * zcmsmag * 2._wp * r dt ) ! upper limit stability factor scaling409 zstabf_up = rn_maxfac / ( 4._wp * zcmsmag * 2._wp * rn_Dt ) ! upper limit stability factor scaling 410 410 IF( ln_dynldf_blp ) zstabf_lo = ( 16._wp / 9._wp ) * zstabf_lo ! provide |U|L^3/12 lower limit instead 411 411 ! ! of |U|L^3/16 in blp case -
NEMO/trunk/src/OCE/LDF/ldftra.F90
r12377 r12489 820 820 ! 821 821 IF( iom_use('weiv_masstr') ) THEN ! vertical mass transport & its square value 822 zw2d(:,:) = r au0 * e1e2t(:,:)822 zw2d(:,:) = rho0 * e1e2t(:,:) 823 823 DO jk = 1, jpk 824 824 zw3d(:,:,jk) = zw3d(:,:,jk) * zw2d(:,:) … … 830 830 zw3d(:,:,:) = 0.e0 831 831 DO jk = 1, jpkm1 832 zw3d(:,:,jk) = r au0 * ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) )832 zw3d(:,:,jk) = rho0 * ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) ) 833 833 END DO 834 834 CALL iom_put( "ueiv_masstr", zw3d ) ! mass transport in i-direction 835 835 ENDIF 836 836 ! 837 zztmp = 0.5_wp * r au0 * rcp837 zztmp = 0.5_wp * rho0 * rcp 838 838 IF( iom_use('ueiv_heattr') .OR. iom_use('ueiv_heattr3d') ) THEN 839 839 zw2d(:,:) = 0._wp … … 853 853 zw3d(:,:,:) = 0.e0 854 854 DO jk = 1, jpkm1 855 zw3d(:,:,jk) = r au0 * ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) )855 zw3d(:,:,jk) = rho0 * ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) ) 856 856 END DO 857 857 CALL iom_put( "veiv_masstr", zw3d ) ! mass transport in i-direction -
NEMO/trunk/src/OCE/OBS/diaobs.F90
r12377 r12489 539 539 ENDIF 540 540 541 idaystp = NINT( rday / r dt )541 idaystp = NINT( rday / rn_Dt ) 542 542 543 543 !----------------------------------------------------------------------- … … 774 774 & rday 775 775 USE dom_oce, ONLY : & ! Ocean space and time domain variables 776 & r dt776 & rn_Dt 777 777 778 778 IMPLICIT NONE … … 805 805 !! Compute number of days + number of hours + min since initial time 806 806 !!---------------------------------------------------------------------- 807 zdayfrc = kstp * r dt / rday807 zdayfrc = kstp * rn_Dt / rday 808 808 zdayfrc = zdayfrc - aint(zdayfrc) 809 809 imin = imin + int( zdayfrc * 24 * 60 ) … … 816 816 iday=iday+1 817 817 END DO 818 iday = iday + kstp * r dt / rday818 iday = iday + kstp * rn_Dt / rday 819 819 820 820 !----------------------------------------------------------------------- -
NEMO/trunk/src/OCE/OBS/obs_prep.F90
r12377 r12489 613 613 !! * Modules used 614 614 USE dom_oce, ONLY : & ! Geographical information 615 & r dt615 & rn_Dt 616 616 USE phycst, ONLY : & ! Physical constants 617 617 & rday, & … … 662 662 663 663 ! Intialize the number of time steps per day 664 idaystp = NINT( rday / r dt )664 idaystp = NINT( rday / rn_Dt ) 665 665 666 666 !--------------------------------------------------------------------- … … 732 732 733 733 ! Add in the number of time steps to the observation minute 734 zminstp = rmmss / r dt734 zminstp = rmmss / rn_Dt 735 735 zhoustp = rhhmm * zminstp 736 736 -
NEMO/trunk/src/OCE/SBC/fldread.F90
r12377 r12489 172 172 ! Note that all varibles starting by nsec_* are shifted time by +1/2 time step to be centrered 173 173 IF( PRESENT(kit) ) THEN ! ignore kn_fsbc in this case 174 isecsbc = nsec_year + nsec1jan000 + NINT( ( REAL( kit,wp) + zt_offset ) * r dt / REAL(nn_baro,wp) )174 isecsbc = nsec_year + nsec1jan000 + NINT( ( REAL( kit,wp) + zt_offset ) * rn_Dt / REAL(nn_e,wp) ) 175 175 ELSE ! middle of sbc time step 176 176 ! note: we use kn_fsbc-1 because nsec_year is defined at the middle of the current time step 177 isecsbc = nsec_year + nsec1jan000 + NINT( ( 0.5*REAL(kn_fsbc-1,wp) + zt_offset ) * r dt )177 isecsbc = nsec_year + nsec1jan000 + NINT( ( 0.5*REAL(kn_fsbc-1,wp) + zt_offset ) * rn_Dt ) 178 178 EN