Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • nemo/nemo
  • sparonuz/nemo
  • hatfield/nemo
  • extdevs/nemo
4 results
Show changes
Showing
with 794 additions and 667 deletions
......@@ -84,7 +84,7 @@ CONTAINS
INTEGER :: ji, jj, jl ! loop indices
!!-------------------------------------------------------------------
ALLOCATE( diag_dvpn_mlt(jpi,jpj), diag_dvpn_lid(jpi,jpj), diag_dvpn_drn(jpi,jpj), diag_dvpn_rnf(jpi,jpj) )
ALLOCATE( diag_dvpn_mlt(A2D(0)) , diag_dvpn_lid(A2D(0)) , diag_dvpn_drn(A2D(0)) , diag_dvpn_rnf(A2D(0)) )
ALLOCATE( diag_dvpn_mlt_1d(jpij), diag_dvpn_lid_1d(jpij), diag_dvpn_drn_1d(jpij), diag_dvpn_rnf_1d(jpij) )
!
diag_dvpn_mlt (:,:) = 0._wp ; diag_dvpn_drn (:,:) = 0._wp
......@@ -98,7 +98,7 @@ CONTAINS
at_i(:,:) = SUM( a_i, dim=3 )
!
DO jl = 1, jpl
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF( v_i(ji,jj,jl) < epsi10 .OR. at_i(ji,jj) < epsi10 ) THEN
wfx_pnd (ji,jj) = wfx_pnd(ji,jj) + ( v_ip(ji,jj,jl) + v_il(ji,jj,jl) ) * rhow * r1_Dt_ice
a_ip (ji,jj,jl) = 0._wp
......@@ -115,7 +115,7 @@ CONTAINS
! Identify grid cells with ice
!------------------------------
npti = 0 ; nptidx(:) = 0
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF( at_i(ji,jj) >= epsi10 ) THEN
npti = npti + 1
nptidx( npti ) = (jj - 1) * jpi + ji
......@@ -137,6 +137,8 @@ CONTAINS
END SELECT
ENDIF
! the following fields need to be updated in the halos (done in icethd): a_ip, v_ip, v_il, h_ip, h_il
!------------------------------------
! Diagnostics
!------------------------------------
......@@ -529,7 +531,7 @@ CONTAINS
zv_pnd , & ! volume of meltwater contributing to ponds
zv_mlt ! total amount of meltwater produced
REAL(wp), DIMENSION(jpi,jpj) :: zvolp_ini , & !! total melt pond water available before redistribution and drainage
REAL(wp), DIMENSION(A2D(0)) :: zvolp_ini , & !! total melt pond water available before redistribution and drainage
zvolp , & !! total melt pond water volume
zvolp_res !! remaining melt pond water available after drainage
......@@ -589,7 +591,7 @@ CONTAINS
zvolp(:,:) = 0._wp
DO jl = 1, jpl
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF ( a_i(ji,jj,jl) > epsi10 ) THEN
......@@ -637,7 +639,7 @@ CONTAINS
IF( ln_pnd_lids ) THEN
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF ( at_i(ji,jj) > 0.01 .AND. hm_i(ji,jj) > rn_himin .AND. zvolp_ini(ji,jj) > zvp_min * at_i(ji,jj) ) THEN
......@@ -764,7 +766,7 @@ CONTAINS
DO jl = 1, jpl
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
! ! zap lids on small ponds
! IF ( a_i(ji,jj,jl) > epsi10 .AND. v_ip(ji,jj,jl) < epsi10 &
......@@ -826,7 +828,7 @@ CONTAINS
!!
!!------------------------------------------------------------------
REAL (wp), DIMENSION(jpi,jpj), INTENT(INOUT) :: &
REAL (wp), DIMENSION(A2D(0)), INTENT(INOUT) :: &
zvolp, & ! total available pond water
zdvolp ! remaining meltwater after redistribution
......@@ -865,10 +867,10 @@ CONTAINS
INTEGER :: ji, jj, jk, jl ! loop indices
a_ip(:,:,:) = 0._wp
h_ip(:,:,:) = 0._wp
a_ip(A2D(0),:) = 0._wp
h_ip(A2D(0),:) = 0._wp
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
IF ( at_i(ji,jj) > 0.01 .AND. hm_i(ji,jj) > rn_himin .AND. zvolp(ji,jj) > zvp_min * at_i(ji,jj) ) THEN
......
......@@ -55,7 +55,7 @@ CONTAINS
!!-------------------------------------------------------------------
!! *** ROUTINE ice_update_alloc ***
!!-------------------------------------------------------------------
ALLOCATE( utau_oce(jpi,jpj), vtau_oce(jpi,jpj), tmod_io(jpi,jpj), STAT=ice_update_alloc )
ALLOCATE( utau_oce(A2D(0)), vtau_oce(A2D(0)), tmod_io(A2D(1)), STAT=ice_update_alloc )
!
CALL mpp_sum( 'iceupdate', ice_update_alloc )
IF( ice_update_alloc /= 0 ) CALL ctl_stop( 'STOP', 'ice_update_alloc: failed to allocate arrays' )
......@@ -104,22 +104,29 @@ CONTAINS
! Net heat flux on top of the ice-ocean (W.m-2)
!----------------------------------------------
IF( ln_cndflx ) THEN ! ice-atm interface = conduction (and melting) fluxes
qt_atm_oi(:,:) = ( 1._wp - at_i_b(:,:) ) * ( qns_oce(:,:) + qsr_oce(:,:) ) + qemp_oce(:,:) + &
& SUM( a_i_b * ( qcn_ice + qml_ice + qtr_ice_top ), dim=3 ) + qemp_ice(:,:)
DO_2D( 0, 0, 0, 0 )
qt_atm_oi(ji,jj) = ( 1._wp - at_i_b(ji,jj) ) * ( qns_oce(ji,jj) + qsr_oce(ji,jj) ) + qemp_oce(ji,jj) &
& + SUM( a_i_b(ji,jj,:) * ( qcn_ice(ji,jj,:) + qml_ice(ji,jj,:) + qtr_ice_top(ji,jj,:) ) ) &
& + qemp_ice(ji,jj)
END_2D
ELSE ! ice-atm interface = solar and non-solar fluxes
qt_atm_oi(:,:) = qns_tot(:,:) + qsr_tot(:,:)
DO_2D( 0, 0, 0, 0 )
qt_atm_oi(ji,jj) = qns_tot(ji,jj) + qsr_tot(ji,jj)
END_2D
ENDIF
! --- case we bypass ice thermodynamics --- !
IF( .NOT. ln_icethd ) THEN ! we suppose ice is impermeable => ocean is isolated from atmosphere
qt_atm_oi (:,:) = ( 1._wp - at_i_b(:,:) ) * ( qns_oce(:,:) + qsr_oce(:,:) ) + qemp_oce(:,:)
qt_oce_ai (:,:) = ( 1._wp - at_i_b(:,:) ) * qns_oce(:,:) + qemp_oce(:,:)
emp_ice (:,:) = 0._wp
qemp_ice (:,:) = 0._wp
qevap_ice (:,:,:) = 0._wp
DO_2D( 0, 0, 0, 0 )
qt_atm_oi (ji,jj) = ( 1._wp - at_i_b(ji,jj) ) * ( qns_oce(ji,jj) + qsr_oce(ji,jj) ) + qemp_oce(ji,jj)
qt_oce_ai (ji,jj) = ( 1._wp - at_i_b(ji,jj) ) * qns_oce(ji,jj) + qemp_oce(ji,jj)
emp_ice (ji,jj) = 0._wp
qemp_ice (ji,jj) = 0._wp
qevap_ice (ji,jj,:) = 0._wp
END_2D
ENDIF
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
! Solar heat flux reaching the ocean (max) = zqsr (W.m-2)
!---------------------------------------------------
......@@ -183,20 +190,22 @@ CONTAINS
snwice_mass_b(ji,jj) = snwice_mass(ji,jj) ! save mass from the previous ice time step
! ! new mass per unit area
snwice_mass (ji,jj) = tmask(ji,jj,1) * ( rhos * vt_s(ji,jj) + rhoi * vt_i(ji,jj) + rhow * (vt_ip(ji,jj) + vt_il(ji,jj)) )
! ! time evolution of snow+ice mass
snwice_fmass (ji,jj) = ( snwice_mass(ji,jj) - snwice_mass_b(ji,jj) ) * r1_Dt_ice
END_2D
CALL lbc_lnk( 'iceupdate', snwice_mass, 'T', 1.0_wp, snwice_mass_b, 'T', 1.0_wp ) ! needed for sshwzv and dynspg_ts (lbc on emp is done in sbcmod)
! time evolution of snow+ice mass
snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) * r1_Dt_ice
! Storing the transmitted variables
!----------------------------------
fr_i (:,:) = at_i(:,:) ! Sea-ice fraction
tn_ice(:,:,:) = t_su(:,:,:) ! Ice surface temperature
tn_ice(:,:,:) = t_su(A2D(0),:) ! Ice surface temperature
! Snow/ice albedo (only if sent to coupler, useless in forced mode)
!------------------------------------------------------------------
CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_eff, h_ip, cloud_fra, alb_ice ) ! ice albedo
CALL ice_alb( ln_pnd_alb, t_su(A2D(0),:), h_i(A2D(0),:), h_s(A2D(0),:), a_ip_eff(:,:,:), h_ip(A2D(0),:), cloud_fra(:,:), & ! <<== in
& alb_ice(:,:,:) ) ! ==>> out
!
IF( lrst_ice ) THEN !* write snwice_mass fields in the restart file
CALL update_rst( 'WRITE', kt )
......@@ -216,8 +225,8 @@ CONTAINS
IF( iom_use('sfxopw' ) ) CALL iom_put( 'sfxopw', sfx_opw * 1.e-03 ) ! salt flux from open water formation
IF( iom_use('sfxdyn' ) ) CALL iom_put( 'sfxdyn', sfx_dyn * 1.e-03 ) ! salt flux from ridging rafting
IF( iom_use('sfxbri' ) ) CALL iom_put( 'sfxbri', sfx_bri * 1.e-03 ) ! salt flux from brines
IF( iom_use('sfxres' ) ) CALL iom_put( 'sfxres', sfx_res * 1.e-03 ) ! salt flux from undiagnosed processes
IF( iom_use('sfxsub' ) ) CALL iom_put( 'sfxsub', sfx_sub * 1.e-03 ) ! salt flux from sublimation
IF( iom_use('sfxres' ) ) CALL iom_put( 'sfxres', sfx_res(A2D(0)) * 1.e-03 ) ! salt flux from undiagnosed processes
! --- mass fluxes [kg/m2/s] --- !
CALL iom_put( 'emp_oce', emp_oce ) ! emp over ocean (taking into account the snow blown away from the ice)
......@@ -232,13 +241,13 @@ CONTAINS
CALL iom_put( 'vfxsni' , wfx_sni ) ! mass flux from snow-ice formation
CALL iom_put( 'vfxopw' , wfx_opw ) ! mass flux from growth in open water
CALL iom_put( 'vfxdyn' , wfx_dyn ) ! mass flux from dynamics (ridging)
CALL iom_put( 'vfxres' , wfx_res ) ! mass flux from undiagnosed processes
CALL iom_put( 'vfxpnd' , wfx_pnd ) ! mass flux from melt ponds
CALL iom_put( 'vfxsub' , wfx_ice_sub ) ! mass flux from ice sublimation (ice-atm.)
CALL iom_put( 'vfxsub_err', wfx_err_sub ) ! "excess" of sublimation sent to ocean
CALL iom_put( 'vfxres' , wfx_res(A2D(0)) ) ! mass flux from undiagnosed processes
IF ( iom_use( 'vfxthin' ) ) THEN ! mass flux from ice growth in open water + thin ice (<20cm) => comparable to observations
ALLOCATE( z2d(jpi,jpj) )
ALLOCATE( z2d(A2D(0)) )
WHERE( hm_i(:,:) < 0.2 .AND. hm_i(:,:) > 0. ) ; z2d = wfx_bog
ELSEWHERE ; z2d = 0._wp
END WHERE
......@@ -264,8 +273,8 @@ CONTAINS
IF( iom_use('qtr_ice_top') ) CALL iom_put( 'qtr_ice_top', SUM( qtr_ice_top * a_i_b, dim=3 ) ) ! solar flux transmitted thru ice surface
IF( iom_use('qt_oce' ) ) CALL iom_put( 'qt_oce' , ( qsr_oce + qns_oce ) * ( 1._wp - at_i_b ) + qemp_oce )
IF( iom_use('qt_ice' ) ) CALL iom_put( 'qt_ice' , SUM( ( qns_ice + qsr_ice ) * a_i_b, dim=3 ) + qemp_ice )
IF( iom_use('qt_oce_ai' ) ) CALL iom_put( 'qt_oce_ai' , qt_oce_ai * tmask(:,:,1) ) ! total heat flux at the ocean surface: interface oce-(ice+atm)
IF( iom_use('qt_atm_oi' ) ) CALL iom_put( 'qt_atm_oi' , qt_atm_oi * tmask(:,:,1) ) ! total heat flux at the oce-ice surface: interface atm-(ice+oce)
IF( iom_use('qt_oce_ai' ) ) CALL iom_put( 'qt_oce_ai' , qt_oce_ai * smask0 ) ! total heat flux at the ocean surface: interface oce-(ice+atm)
IF( iom_use('qt_atm_oi' ) ) CALL iom_put( 'qt_atm_oi' , qt_atm_oi * smask0 ) ! total heat flux at the oce-ice surface: interface atm-(ice+oce)
IF( iom_use('qemp_oce' ) ) CALL iom_put( 'qemp_oce' , qemp_oce ) ! Downward Heat Flux from E-P over ocean
IF( iom_use('qemp_ice' ) ) CALL iom_put( 'qemp_ice' , qemp_ice ) ! Downward Heat Flux from E-P over ice
......@@ -282,9 +291,9 @@ CONTAINS
! heat fluxes associated with mass exchange (freeze/melt/precip...)
CALL iom_put ('hfxthd' , hfx_thd ) !
CALL iom_put ('hfxdyn' , hfx_dyn ) !
CALL iom_put ('hfxres' , hfx_res ) !
CALL iom_put ('hfxsub' , hfx_sub ) !
CALL iom_put ('hfxspr' , hfx_spr ) ! Heat flux from snow precip heat content
CALL iom_put ('hfxres' , hfx_res(A2D(0)) ) !
! other heat fluxes
IF( iom_use('hfxsensib' ) ) CALL iom_put( 'hfxsensib' , qsb_ice_bot * at_i_b ) ! Sensible oceanic heat flux
......@@ -314,7 +323,7 @@ CONTAINS
!!
!! ** Action : * at each ice time step (every nn_fsbc time step):
!! - compute the modulus of ice-ocean relative velocity
!! (*rho*Cd) at T-point (C-grid) or I-point (B-grid)
!! (*rho*Cd) at T-point (C-grid)
!! tmod_io = rhoco * | U_ice-U_oce |
!! - update the modulus of stress at ocean surface
!! taum = (1-a) * taum + a * tmod_io * | U_ice-U_oce |
......@@ -325,19 +334,19 @@ CONTAINS
!!
!! NB: - ice-ocean rotation angle no more allowed
!! - here we make an approximation: taum is only computed every ice time step
!! This avoids mutiple average to pass from T -> U,V grids and next from U,V grids
!! to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton...
!! This avoids mutiple average to pass from U,V grids to T grids
!! taum is used in TKE and GLS, which should not be too sensitive to this approximaton...
!!
!! ** Outputs : - utau, vtau : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes
!! ** Outputs : - utau, vtau : surface ocean i- and j-stress (T-pts) updated with ice-ocean fluxes
!! - taum : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes
!!---------------------------------------------------------------------
INTEGER , INTENT(in) :: kt ! ocean time-step index
REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents
!
INTEGER :: ji, jj ! dummy loop indices
REAL(wp) :: zat_u, zutau_ice, zu_t, zmodt ! local scalar
REAL(wp) :: zat_v, zvtau_ice, zv_t, zrhoco ! - -
REAL(wp) :: zflagi ! - -
REAL(wp) :: zutau_ice, zu_t, zmodt ! local scalar
REAL(wp) :: zvtau_ice, zv_t, zrhoco ! - -
REAL(wp) :: zflagi ! - -
!!---------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('iceupdate')
......@@ -350,46 +359,51 @@ CONTAINS
zrhoco = rho0 * rn_cio
!
IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step)
DO_2D( 0, 0, 0, 0 ) !* update the modulus of stress at ocean surface (T-point)
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) !* rhoco * |U_ice-U_oce| at T-point
! ! 2*(U_ice-U_oce) at T-point
zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj) ! u_oce = ssu_m
zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1) ! v_oce = ssv_m
! ! |U_ice-U_oce|^2
zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t )
!
tmod_io(ji,jj) = zrhoco * SQRT( zmodt )
END_2D
IF( nn_hls == 1 ) CALL lbc_lnk( 'iceupdate', tmod_io, 'T', 1._wp )
!
DO_2D( 0, 0, 0, 0 ) !* save the air-ocean stresses at ice time-step
! ! 2*(U_ice-U_oce) at T-point
zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj)
zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1)
! ! |U_ice-U_oce|^2
zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj) ! u_oce = ssu_m
zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1) ! v_oce = ssv_m
! ! |U_ice-U_oce|^2
zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t )
! ! update the ocean stress modulus
taum(ji,jj) = ( 1._wp - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * zrhoco * zmodt
tmod_io(ji,jj) = zrhoco * SQRT( zmodt ) ! rhoco * |U_ice-U_oce| at T-point
!
utau_oce(ji,jj) = utau(ji,jj)
vtau_oce(ji,jj) = vtau(ji,jj)
END_2D
CALL lbc_lnk( 'iceupdate', taum, 'T', 1.0_wp, tmod_io, 'T', 1.0_wp )
!
utau_oce(:,:) = utau(:,:) !* save the air-ocean stresses at ice time-step
vtau_oce(:,:) = vtau(:,:)
!
ENDIF
!
! !== every ocean time-step ==!
IF ( ln_drgice_imp ) THEN
! Save drag with right sign to update top drag in the ocean implicit friction
rCdU_ice(:,:) = -r1_rho0 * tmod_io(:,:) * at_i(:,:) * tmask(:,:,1)
DO_2D( 1, 1, 1, 1 )
rCdU_ice(ji,jj) = -r1_rho0 * tmod_io(ji,jj) * at_i(ji,jj) * tmask(ji,jj,1)
END_2D
zflagi = 0._wp
ELSE
zflagi = 1._wp
ENDIF
!
DO_2D( 0, 0, 0, 0 ) !* update the stress WITHOUT an ice-ocean rotation angle
! ice area at u and v-points
zat_u = ( at_i(ji,jj) * tmask(ji,jj,1) + at_i (ji+1,jj ) * tmask(ji+1,jj ,1) ) &
& / MAX( 1.0_wp , tmask(ji,jj,1) + tmask(ji+1,jj ,1) )
zat_v = ( at_i(ji,jj) * tmask(ji,jj,1) + at_i (ji ,jj+1 ) * tmask(ji ,jj+1,1) ) &
& / MAX( 1.0_wp , tmask(ji,jj,1) + tmask(ji ,jj+1,1) )
! ! linearized quadratic drag formulation
zutau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) )
zvtau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) )
zutau_ice = 0.5_wp * tmod_io(ji,jj) * ( u_ice(ji,jj) + u_ice(ji-1,jj) - pu_oce(ji,jj) - pu_oce(ji-1,jj) )
zvtau_ice = 0.5_wp * tmod_io(ji,jj) * ( v_ice(ji,jj) + v_ice(ji,jj-1) - pv_oce(ji,jj) - pv_oce(ji,jj-1) )
! ! stresses at the ocean surface
utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_ice
vtau(ji,jj) = ( 1._wp - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice
utau(ji,jj) = ( 1._wp - at_i(ji,jj) ) * utau_oce(ji,jj) + at_i(ji,jj) * zutau_ice
vtau(ji,jj) = ( 1._wp - at_i(ji,jj) ) * vtau_oce(ji,jj) + at_i(ji,jj) * zvtau_ice
END_2D
CALL lbc_lnk( 'iceupdate', utau, 'U', -1.0_wp, vtau, 'V', -1.0_wp ) ! lateral boundary condition
!
IF( ln_timing ) CALL timing_stop('iceupdate')
!
......@@ -446,12 +460,12 @@ CONTAINS
CALL iom_get( numrir, jpdom_auto, 'snwice_mass_b', snwice_mass_b )
ELSE ! start from rest
IF(lwp) WRITE(numout,*) ' ==>> previous run without snow-ice mass output then set it'
snwice_mass (:,:) = tmask(:,:,1) * ( rhos * vt_s(:,:) + rhoi * vt_i(:,:) )
snwice_mass (:,:) = tmask(:,:,1) * ( rhos * vt_s(:,:) + rhoi * vt_i(:,:) + rhow * (vt_ip(:,:) + vt_il(:,:)) )
snwice_mass_b(:,:) = snwice_mass(:,:)
ENDIF
ELSE !* Start from rest
IF(lwp) WRITE(numout,*) ' ==>> start from rest: set the snow-ice mass'
snwice_mass (:,:) = tmask(:,:,1) * ( rhos * vt_s(:,:) + rhoi * vt_i(:,:) )
snwice_mass (:,:) = tmask(:,:,1) * ( rhos * vt_s(:,:) + rhoi * vt_i(:,:) + rhow * (vt_ip(:,:) + vt_il(:,:)) )
snwice_mass_b(:,:) = snwice_mass(:,:)
ENDIF
!
......
......@@ -117,70 +117,88 @@ CONTAINS
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z1_at_i, z1_vt_i, z1_vt_s
!!-------------------------------------------------------------------
!
! ! integrated values
vt_i(:,:) = SUM( v_i (:,:,:) , dim=3 )
vt_s(:,:) = SUM( v_s (:,:,:) , dim=3 )
st_i(:,:) = SUM( sv_i(:,:,:) , dim=3 )
at_i(:,:) = SUM( a_i (:,:,:) , dim=3 )
et_s(:,:) = SUM( SUM( e_s (:,:,:,:), dim=4 ), dim=3 )
et_i(:,:) = SUM( SUM( e_i (:,:,:,:), dim=4 ), dim=3 )
! full arrays: vt_i, vt_s, at_i, vt_ip, vt_il, at_ip
! reduced arrays: the rest
!
at_ip(:,:) = SUM( a_ip(:,:,:), dim=3 ) ! melt ponds
vt_ip(:,:) = SUM( v_ip(:,:,:), dim=3 )
vt_il(:,:) = SUM( v_il(:,:,:), dim=3 )
! ! integrated values
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
vt_i(ji,jj) = SUM( v_i (ji,jj,:) )
vt_s(ji,jj) = SUM( v_s (ji,jj,:) )
at_i(ji,jj) = SUM( a_i (ji,jj,:) )
!
at_ip(ji,jj) = SUM( a_ip(ji,jj,:) ) ! melt ponds
vt_ip(ji,jj) = SUM( v_ip(ji,jj,:) )
vt_il(ji,jj) = SUM( v_il(ji,jj,:) )
END_2D
DO_2D( 0, 0, 0, 0 )
st_i(ji,jj) = SUM( sv_i(ji,jj,:) )
et_s(ji,jj) = SUM( SUM( e_s (ji,jj,:,:), dim=2 ) )
et_i(ji,jj) = SUM( SUM( e_i (ji,jj,:,:), dim=2 ) )
END_2D
!
ato_i(:,:) = 1._wp - at_i(:,:) ! open water fraction
!
!!GS: tm_su always needed by ABL over sea-ice
ALLOCATE( z1_at_i(jpi,jpj) )
WHERE( at_i(:,:) > epsi20 ) ; z1_at_i(:,:) = 1._wp / at_i(:,:)
ELSEWHERE ; z1_at_i(:,:) = 0._wp
ALLOCATE( z1_at_i(A2D(0)) )
WHERE( at_i(A2D(0)) > epsi20 ) ; z1_at_i(:,:) = 1._wp / at_i(A2D(0))
ELSEWHERE ; z1_at_i(:,:) = 0._wp
END WHERE
tm_su(:,:) = SUM( t_su(:,:,:) * a_i(:,:,:) , dim=3 ) * z1_at_i(:,:)
WHERE( at_i(:,:)<=epsi20 ) tm_su(:,:) = rt0
DO_2D( 0, 0, 0, 0 )
IF( at_i(ji,jj)<=epsi20 ) THEN
tm_su(ji,jj) = rt0
ELSE
tm_su(ji,jj) = SUM( t_su(ji,jj,:) * a_i(ji,jj,:) ) * z1_at_i(ji,jj)
ENDIF
END_2D
!
! The following fields are calculated for diagnostics and outputs only
! ==> Do not use them for other purposes
IF( kn > 1 ) THEN
!
ALLOCATE( z1_vt_i(jpi,jpj) , z1_vt_s(jpi,jpj) )
WHERE( vt_i(:,:) > epsi20 ) ; z1_vt_i(:,:) = 1._wp / vt_i(:,:)
ELSEWHERE ; z1_vt_i(:,:) = 0._wp
ALLOCATE( z1_vt_i(A2D(0)) , z1_vt_s(A2D(0)) )
WHERE( vt_i(A2D(0)) > epsi20 ) ; z1_vt_i(:,:) = 1._wp / vt_i(A2D(0))
ELSEWHERE ; z1_vt_i(:,:) = 0._wp
END WHERE
WHERE( vt_s(:,:) > epsi20 ) ; z1_vt_s(:,:) = 1._wp / vt_s(:,:)
ELSEWHERE ; z1_vt_s(:,:) = 0._wp
WHERE( vt_s(A2D(0)) > epsi20 ) ; z1_vt_s(:,:) = 1._wp / vt_s(A2D(0))
ELSEWHERE ; z1_vt_s(:,:) = 0._wp
END WHERE
!
! ! mean ice/snow thickness
hm_i(:,:) = vt_i(:,:) * z1_at_i(:,:)
hm_s(:,:) = vt_s(:,:) * z1_at_i(:,:)
!
! ! mean temperature (K), salinity and age
tm_si(:,:) = SUM( t_si(:,:,:) * a_i(:,:,:) , dim=3 ) * z1_at_i(:,:)
om_i (:,:) = SUM( oa_i(:,:,:) , dim=3 ) * z1_at_i(:,:)
sm_i (:,:) = st_i(:,:) * z1_vt_i(:,:)
!
tm_i(:,:) = 0._wp
tm_s(:,:) = 0._wp
DO jl = 1, jpl
DO jk = 1, nlay_i
tm_i(:,:) = tm_i(:,:) + r1_nlay_i * t_i (:,:,jk,jl) * v_i(:,:,jl) * z1_vt_i(:,:)
END DO
DO jk = 1, nlay_s
tm_s(:,:) = tm_s(:,:) + r1_nlay_s * t_s (:,:,jk,jl) * v_s(:,:,jl) * z1_vt_s(:,:)
DO_2D( 0, 0, 0, 0 )
hm_i(ji,jj) = vt_i(ji,jj) * z1_at_i(ji,jj)
hm_s(ji,jj) = vt_s(ji,jj) * z1_at_i(ji,jj)
!
! ! mean temperature (K), salinity and age
tm_si(ji,jj) = SUM( t_si(ji,jj,:) * a_i(ji,jj,:) ) * z1_at_i(ji,jj)
om_i (ji,jj) = SUM( oa_i(ji,jj,:) ) * z1_at_i(ji,jj)
sm_i (ji,jj) = st_i(ji,jj) * z1_vt_i(ji,jj)
!
tm_i(ji,jj) = 0._wp
tm_s(ji,jj) = 0._wp
DO jl = 1, jpl
DO jk = 1, nlay_i
tm_i(ji,jj) = tm_i(ji,jj) + r1_nlay_i * t_i (ji,jj,jk,jl) * v_i(ji,jj,jl) * z1_vt_i(ji,jj)
END DO
DO jk = 1, nlay_s
tm_s(ji,jj) = tm_s(ji,jj) + r1_nlay_s * t_s (ji,jj,jk,jl) * v_s(ji,jj,jl) * z1_vt_s(ji,jj)
END DO
END DO
END DO
!
!
END_2D
! ! put rt0 where there is no ice
WHERE( at_i(:,:)<=epsi20 )
WHERE( at_i(A2D(0)) <= epsi20 )
tm_si(:,:) = rt0
tm_i (:,:) = rt0
tm_s (:,:) = rt0
END WHERE
!
! ! mean melt pond depth
WHERE( at_ip(:,:) > epsi20 ) ; hm_ip(:,:) = vt_ip(:,:) / at_ip(:,:) ; hm_il(:,:) = vt_il(:,:) / at_ip(:,:)
ELSEWHERE ; hm_ip(:,:) = 0._wp ; hm_il(:,:) = 0._wp
WHERE( at_ip(A2D(0)) > epsi20 )
hm_ip(:,:) = vt_ip(A2D(0)) / at_ip(A2D(0))
hm_il(:,:) = vt_il(A2D(0)) / at_ip(A2D(0))
ELSEWHERE
hm_ip(:,:) = 0._wp
hm_il(:,:) = 0._wp
END WHERE
!
DEALLOCATE( z1_vt_i , z1_vt_s )
......@@ -206,7 +224,8 @@ CONTAINS
REAL(wp) :: zlay_i, zlay_s ! - -
REAL(wp), PARAMETER :: zhl_max = 0.015_wp ! pond lid thickness above which the ponds disappear from the albedo calculation
REAL(wp), PARAMETER :: zhl_min = 0.005_wp ! pond lid thickness below which the full pond area is used in the albedo calculation
REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_a_i, z1_v_i, z1_a_ip, za_s_fra
REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_a_i, z1_v_i, z1_a_ip
REAL(wp), DIMENSION(A2D(0),jpl) :: za_s_fra
!!-------------------------------------------------------------------
!!gm Question 2: It is possible to define existence of sea-ice in a common way between
......@@ -247,15 +266,15 @@ CONTAINS
h_ip(:,:,:) = v_ip(:,:,:) * z1_a_ip(:,:,:)
h_il(:,:,:) = v_il(:,:,:) * z1_a_ip(:,:,:)
! !--- melt pond effective area (used for albedo)
a_ip_frac(:,:,:) = a_ip(:,:,:) * z1_a_i(:,:,:)
WHERE ( h_il(:,:,:) <= zhl_min ) ; a_ip_eff(:,:,:) = a_ip_frac(:,:,:) ! lid is very thin. Expose all the pond
ELSEWHERE( h_il(:,:,:) >= zhl_max ) ; a_ip_eff(:,:,:) = 0._wp ! lid is very thick. Cover all the pond up with ice and snow
ELSEWHERE ; a_ip_eff(:,:,:) = a_ip_frac(:,:,:) * & ! lid is in between. Expose part of the pond
& ( zhl_max - h_il(:,:,:) ) / ( zhl_max - zhl_min )
a_ip_frac(:,:,:) = a_ip(A2D(0),:) * z1_a_i(A2D(0),:)
WHERE ( h_il(A2D(0),:) <= zhl_min ) ; a_ip_eff(:,:,:) = a_ip_frac(:,:,:) ! lid is very thin. Expose all the pond
ELSEWHERE( h_il(A2D(0),:) >= zhl_max ) ; a_ip_eff(:,:,:) = 0._wp ! lid is very thick. Cover all the pond up with ice and snow
ELSEWHERE ; a_ip_eff(:,:,:) = a_ip_frac(:,:,:) * & ! lid is in between. Expose part of the pond
& ( zhl_max - h_il(A2D(0),:) ) / ( zhl_max - zhl_min )
END WHERE
!
CALL ice_var_snwfra( h_s, za_s_fra ) ! calculate ice fraction covered by snow
a_ip_eff = MIN( a_ip_eff, 1._wp - za_s_fra ) ! make sure (a_ip_eff + a_s_fra) <= 1
CALL ice_var_snwfra( h_s(A2D(0),:), za_s_fra(:,:,:) ) ! calculate ice fraction covered by snow
a_ip_eff(:,:,:) = MIN( a_ip_eff(:,:,:), 1._wp - za_s_fra(:,:,:) ) ! make sure (a_ip_eff + a_s_fra) <= 1
!
! !--- salinity (with a minimum value imposed everywhere)
IF( nn_icesal == 2 ) THEN
......@@ -300,9 +319,9 @@ CONTAINS
END DO
!
! integrated values
vt_i (:,:) = SUM( v_i , dim=3 )
vt_s (:,:) = SUM( v_s , dim=3 )
at_i (:,:) = SUM( a_i , dim=3 )
vt_i (:,:) = SUM( v_i, dim=3 )
vt_s (:,:) = SUM( v_s, dim=3 )
at_i (:,:) = SUM( a_i, dim=3 )
!
END SUBROUTINE ice_var_glo2eqv
......@@ -538,7 +557,7 @@ CONTAINS
sfx_res(ji,jj) = sfx_res(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * sv_i(ji,jj,jl) * rhoi * r1_Dt_ice
wfx_res(ji,jj) = wfx_res(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * v_i (ji,jj,jl) * rhoi * r1_Dt_ice
wfx_res(ji,jj) = wfx_res(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl) * rhos * r1_Dt_ice
wfx_pnd(ji,jj) = wfx_pnd(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * ( v_ip(ji,jj,jl)+v_il(ji,jj,jl) ) * rhow * r1_Dt_ice
wfx_res(ji,jj) = wfx_res(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * ( v_ip(ji,jj,jl)+v_il(ji,jj,jl) ) * rhow * r1_Dt_ice
!
a_i (ji,jj,jl) = a_i (ji,jj,jl) * zswitch(ji,jj)
v_i (ji,jj,jl) = v_i (ji,jj,jl) * zswitch(ji,jj)
......@@ -640,11 +659,11 @@ CONTAINS
psv_i (ji,jj,jl) = 0._wp
ENDIF
IF( pv_ip(ji,jj,jl) < 0._wp .OR. pv_il(ji,jj,jl) < 0._wp .OR. pa_ip(ji,jj,jl) <= 0._wp ) THEN
wfx_pnd(ji,jj) = wfx_pnd(ji,jj) + pv_il(ji,jj,jl) * rhow * z1_dt
wfx_res(ji,jj) = wfx_res(ji,jj) + pv_il(ji,jj,jl) * rhow * z1_dt
pv_il (ji,jj,jl) = 0._wp
ENDIF
IF( pv_ip(ji,jj,jl) < 0._wp .OR. pa_ip(ji,jj,jl) <= 0._wp ) THEN
wfx_pnd(ji,jj) = wfx_pnd(ji,jj) + pv_ip(ji,jj,jl) * rhow * z1_dt
wfx_res(ji,jj) = wfx_res(ji,jj) + pv_ip(ji,jj,jl) * rhow * z1_dt
pv_ip (ji,jj,jl) = 0._wp
ENDIF
END_2D
......@@ -713,15 +732,19 @@ CONTAINS
!! instead of setting everything to zero as just below
bv_i (:,:,:) = 0._wp
DO jl = 1, jpl
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, nlay_i )
DO_3D( 0, 0, 0, 0, 1, nlay_i )
IF( t_i(ji,jj,jk,jl) < rt0 - epsi10 ) THEN
bv_i(ji,jj,jl) = bv_i(ji,jj,jl) - rTmlt * sz_i(ji,jj,jk,jl) * r1_nlay_i / ( t_i(ji,jj,jk,jl) - rt0 )
ENDIF
END_3D
END DO
WHERE( vt_i(:,:) > epsi20 ) ; bvm_i(:,:) = SUM( bv_i(:,:,:) * v_i(:,:,:) , dim=3 ) / vt_i(:,:)
ELSEWHERE ; bvm_i(:,:) = 0._wp
END WHERE
DO_2D( 0, 0, 0, 0 )
IF( vt_i(ji,jj) > epsi20 ) THEN
bvm_i(ji,jj) = SUM( bv_i(ji,jj,:) * v_i(ji,jj,:) ) / vt_i(ji,jj)
ELSE
bvm_i(ji,jj) = 0._wp
ENDIF
END_2D
!
END SUBROUTINE ice_var_bv
......@@ -1286,8 +1309,8 @@ CONTAINS
!!
!!-------------------------------------------------------------------
SUBROUTINE ice_var_snwfra_3d( ph_s, pa_s_fra )
REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: ph_s ! snow thickness
REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pa_s_fra ! ice fraction covered by snow
REAL(wp), DIMENSION(A2D(0),jpl), INTENT(in ) :: ph_s ! snow thickness
REAL(wp), DIMENSION(A2D(0),jpl), INTENT( out) :: pa_s_fra ! ice fraction covered by snow
IF ( nn_snwfra == 0 ) THEN ! basic 0 or 1 snow cover
WHERE( ph_s > 0._wp ) ; pa_s_fra = 1._wp
ELSEWHERE ; pa_s_fra = 0._wp
......@@ -1344,8 +1367,8 @@ CONTAINS
!!--------------------------------------------------------------------------
!!gm I think it can be usefull to set this as a FUNCTION, not a SUBROUTINE....
SUBROUTINE ice_var_snwblow_2d( pin, pout )
REAL(wp), DIMENSION(:,:), INTENT(in ) :: pin ! previous fraction lead ( 1. - a_i_b )
REAL(wp), DIMENSION(:,:), INTENT(inout) :: pout
REAL(wp), DIMENSION(A2D(0)), INTENT(in ) :: pin ! previous fraction lead ( 1. - a_i_b )
REAL(wp), DIMENSION(A2D(0)), INTENT(inout) :: pout
pout = ( 1._wp - ( pin )**rn_snwblow )
END SUBROUTINE ice_var_snwblow_2d
......
......@@ -26,7 +26,6 @@ MODULE icewri
USE iom ! I/O manager library
USE lib_mpp ! MPP library
USE lib_fortran ! fortran utilities (glob_sum + no signed zero)
USE lbclnk ! lateral boundary conditions (or mpp links)
USE timing ! Timing
IMPLICIT NONE
......@@ -53,10 +52,10 @@ CONTAINS
INTEGER :: ji, jj, jk, jl ! dummy loop indices
REAL(wp) :: z2da, z2db, zrho1, zrho2
REAL(wp) :: zmiss_val ! missing value retrieved from xios
REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace
REAL(wp), DIMENSION(jpi,jpj) :: zmsk00, zmsk05, zmsk15, zmsksn ! O%, 5% and 15% concentration mask and snow mask
REAL(wp), DIMENSION(jpi,jpj,jpl) :: zmsk00l, zmsksnl ! cat masks
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zfast, zalb, zmskalb ! 2D workspace
REAL(wp), DIMENSION(A2D(0)) :: z2d ! 2D workspace
REAL(wp), DIMENSION(A2D(0)) :: zmsk00, zmsk05, zmsk15, zmsksn ! O%, 5% and 15% concentration mask and snow mask
REAL(wp), DIMENSION(A2D(0),jpl) :: zmsk00l, zmsksnl ! cat masks
REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zfast, zalb, zmskalb ! 2D workspace
!
! Global ice diagnostics (SIMIP)
REAL(wp) :: zdiag_area_nh, zdiag_extt_nh, zdiag_volu_nh ! area, extent, volume
......@@ -72,14 +71,14 @@ CONTAINS
CALL ice_var_bv
! tresholds for outputs
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
zmsk00(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice , 0 if no ice
zmsk05(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.05_wp ) ) ! 1 if 5% ice , 0 if less
zmsk15(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.15_wp ) ) ! 1 if 15% ice, 0 if less
zmsksn(ji,jj) = MAX( 0._wp , SIGN( 1._wp , vt_s(ji,jj) - epsi06 ) ) ! 1 if snow , 0 if no snow
END_2D
DO jl = 1, jpl
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
zmsk00l(ji,jj,jl) = MAX( 0._wp , SIGN( 1._wp , a_i(ji,jj,jl) - epsi06 ) )
zmsksnl(ji,jj,jl) = MAX( 0._wp , SIGN( 1._wp , v_s(ji,jj,jl) - epsi06 ) )
END_2D
......@@ -96,102 +95,111 @@ CONTAINS
CALL iom_put( 'icepres' , zmsk00 ) ! Ice presence (1 or 0)
!
! general fields
IF( iom_use('icemass' ) ) CALL iom_put( 'icemass', vt_i * rhoi * zmsk00 ) ! Ice mass per cell area
IF( iom_use('snwmass' ) ) CALL iom_put( 'snwmass', vt_s * rhos * zmsksn ) ! Snow mass per cell area
IF( iom_use('iceconc' ) ) CALL iom_put( 'iceconc', at_i * zmsk00 ) ! ice concentration
IF( iom_use('icevolu' ) ) CALL iom_put( 'icevolu', vt_i * zmsk00 ) ! ice volume = mean ice thickness over the cell
IF( iom_use('icethic' ) ) CALL iom_put( 'icethic', hm_i * zmsk00 ) ! ice thickness
IF( iom_use('snwthic' ) ) CALL iom_put( 'snwthic', hm_s * zmsk00 ) ! snw thickness
IF( iom_use('icebrv' ) ) CALL iom_put( 'icebrv' , bvm_i* 100. * zmsk00 ) ! brine volume
IF( iom_use('iceage' ) ) CALL iom_put( 'iceage' , om_i / rday * zmsk15 + zmiss_val * ( 1._wp - zmsk15 ) ) ! ice age
IF( iom_use('icehnew' ) ) CALL iom_put( 'icehnew', ht_i_new ) ! new ice thickness formed in the leads
IF( iom_use('snwvolu' ) ) CALL iom_put( 'snwvolu', vt_s * zmsksn ) ! snow volume
IF( iom_use('icefrb' ) ) THEN ! Ice freeboard
z2d(:,:) = ( zrho1 * hm_i(:,:) - zrho2 * hm_s(:,:) )
IF( iom_use('icemass' ) ) CALL iom_put( 'icemass', vt_i(A2D(0)) * rhoi * zmsk00 ) ! Ice mass per cell area
IF( iom_use('snwmass' ) ) CALL iom_put( 'snwmass', vt_s(A2D(0)) * rhos * zmsksn ) ! Snow mass per cell area
IF( iom_use('iceconc' ) ) CALL iom_put( 'iceconc', at_i(A2D(0)) * zmsk00 ) ! ice concentration
IF( iom_use('icevolu' ) ) CALL iom_put( 'icevolu', vt_i(A2D(0)) * zmsk00 ) ! ice volume = mean ice thickness over the cell
IF( iom_use('icethic' ) ) CALL iom_put( 'icethic', hm_i(:,:) * zmsk00 ) ! ice thickness
IF( iom_use('snwthic' ) ) CALL iom_put( 'snwthic', hm_s(:,:) * zmsk00 ) ! snw thickness
IF( iom_use('icebrv' ) ) CALL iom_put( 'icebrv' , bvm_i(:,:)* 100. * zmsk00 ) ! brine volume
IF( iom_use('iceage' ) ) CALL iom_put( 'iceage' , om_i(:,:) / rday * zmsk15 + zmiss_val * ( 1._wp - zmsk15 ) ) ! ice age
IF( iom_use('icehnew' ) ) CALL iom_put( 'icehnew', ht_i_new(:,:) ) ! new ice thickness formed in the leads
IF( iom_use('snwvolu' ) ) CALL iom_put( 'snwvolu', vt_s(A2D(0)) * zmsksn ) ! snow volume
IF( iom_use('icefrb' ) ) THEN ! Ice freeboard
z2d(:,:) = zrho1 * hm_i(:,:) - zrho2 * hm_s(:,:)
WHERE( z2d < 0._wp ) z2d = 0._wp
CALL iom_put( 'icefrb' , z2d * zmsk00 )
CALL iom_put( 'icefrb' , z2d * zmsk00 )
ENDIF
! melt ponds
IF( iom_use('iceapnd' ) ) CALL iom_put( 'iceapnd', at_ip * zmsk00 ) ! melt pond total fraction
IF( iom_use('icehpnd' ) ) CALL iom_put( 'icehpnd', hm_ip * zmsk00 ) ! melt pond depth
IF( iom_use('icevpnd' ) ) CALL iom_put( 'icevpnd', vt_ip * zmsk00 ) ! melt pond total volume per unit area
IF( iom_use('icehlid' ) ) CALL iom_put( 'icehlid', hm_il * zmsk00 ) ! melt pond lid depth
IF( iom_use('icevlid' ) ) CALL iom_put( 'icevlid', vt_il * zmsk00 ) ! melt pond lid total volume per unit area
IF( iom_use('iceapnd' ) ) CALL iom_put( 'iceapnd', at_ip(A2D(0)) * zmsk00 ) ! melt pond total fraction
IF( iom_use('icehpnd' ) ) CALL iom_put( 'icehpnd', hm_ip(:,:) * zmsk00 ) ! melt pond depth
IF( iom_use('icevpnd' ) ) CALL iom_put( 'icevpnd', vt_ip(A2D(0)) * zmsk00 ) ! melt pond total volume per unit area
IF( iom_use('icehlid' ) ) CALL iom_put( 'icehlid', hm_il(:,:) * zmsk00 ) ! melt pond lid depth
IF( iom_use('icevlid' ) ) CALL iom_put( 'icevlid', vt_il(A2D(0)) * zmsk00 ) ! melt pond lid total volume per unit area
! salt
IF( iom_use('icesalt' ) ) CALL iom_put( 'icesalt', sm_i * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! mean ice salinity
IF( iom_use('icesalm' ) ) CALL iom_put( 'icesalm', st_i * rhoi * 1.0e-3 * zmsk00 ) ! Mass of salt in sea ice per cell area
IF( iom_use('icesalt' ) ) CALL iom_put( 'icesalt', sm_i(:,:) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! mean ice salinity
IF( iom_use('icesalm' ) ) CALL iom_put( 'icesalm', st_i(:,:) * rhoi * 1.0e-3 * zmsk00 ) ! Mass of salt in sea ice per cell area
! heat
IF( iom_use('icetemp' ) ) CALL iom_put( 'icetemp', ( tm_i - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! ice mean temperature
IF( iom_use('snwtemp' ) ) CALL iom_put( 'snwtemp', ( tm_s - rt0 ) * zmsksn + zmiss_val * ( 1._wp - zmsksn ) ) ! snw mean temperature
IF( iom_use('icettop' ) ) CALL iom_put( 'icettop', ( tm_su - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! temperature at the ice surface
IF( iom_use('icetbot' ) ) CALL iom_put( 'icetbot', ( t_bo - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! temperature at the ice bottom
IF( iom_use('icetsni' ) ) CALL iom_put( 'icetsni', ( tm_si - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! temperature at the snow-ice interface
IF( iom_use('icehc' ) ) CALL iom_put( 'icehc' , -et_i * zmsk00 ) ! ice heat content
IF( iom_use('snwhc' ) ) CALL iom_put( 'snwhc' , -et_s * zmsksn ) ! snow heat content
IF( iom_use('icetemp' ) ) CALL iom_put( 'icetemp', ( tm_i(:,:) - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! ice mean temperature
IF( iom_use('snwtemp' ) ) CALL iom_put( 'snwtemp', ( tm_s(:,:) - rt0 ) * zmsksn + zmiss_val * ( 1._wp - zmsksn ) ) ! snw mean temperature
IF( iom_use('icettop' ) ) CALL iom_put( 'icettop', ( tm_su(:,:) - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! temperature at the ice surface
IF( iom_use('icetbot' ) ) CALL iom_put( 'icetbot', ( t_bo(A2D(0)) - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! temperature at the ice bottom
IF( iom_use('icetsni' ) ) CALL iom_put( 'icetsni', ( tm_si(:,:) - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! temperature at the snow-ice interface
IF( iom_use('icehc' ) ) CALL iom_put( 'icehc' , -et_i(:,:) * zmsk00 ) ! ice heat content
IF( iom_use('snwhc' ) ) CALL iom_put( 'snwhc' , -et_s(:,:) * zmsksn ) ! snow heat content
! momentum
IF( iom_use('uice' ) ) CALL iom_put( 'uice' , u_ice ) ! ice velocity u
IF( iom_use('vice' ) ) CALL iom_put( 'vice' , v_ice ) ! ice velocity v
IF( iom_use('uice' ) ) CALL iom_put( 'uice' , u_ice(A2D(0)) ) ! ice velocity u
IF( iom_use('vice' ) ) CALL iom_put( 'vice' , v_ice(A2D(0)) ) ! ice velocity v
!
IF( iom_use('icevel') .OR. iom_use('fasticepres') ) THEN ! module of ice velocity & fast ice
ALLOCATE( zfast(jpi,jpj) )
IF( iom_use('icevel') .OR. iom_use('fasticepres') ) THEN ! module of ice velocity & fast ice
ALLOCATE( zfast(A2D(0)) )
DO_2D( 0, 0, 0, 0 )
z2da = u_ice(ji,jj) + u_ice(ji-1,jj)
z2db = v_ice(ji,jj) + v_ice(ji,jj-1)
z2d(ji,jj) = 0.5_wp * SQRT( z2da * z2da + z2db * z2db )
END_2D
CALL lbc_lnk( 'icewri', z2d, 'T', 1.0_wp )
CALL iom_put( 'icevel', z2d )
WHERE( z2d(:,:) < 5.e-04_wp .AND. zmsk15(:,:) == 1._wp ) ; zfast(:,:) = 1._wp ! record presence of fast ice
WHERE( z2d(:,:) < 5.e-04_wp .AND. zmsk15(:,:) == 1._wp ) ; zfast(:,:) = 1._wp ! record presence of fast ice
ELSEWHERE ; zfast(:,:) = 0._wp
END WHERE
CALL iom_put( 'fasticepres', zfast )
DEALLOCATE( zfast )
ENDIF
!
IF( iom_use('icealb') .OR. iom_use('albedo') ) THEN ! ice albedo and surface albedo
ALLOCATE( zalb(jpi,jpj), zmskalb(jpi,jpj) )
IF( iom_use('icealb') .OR. iom_use('albedo') ) THEN ! ice albedo and surface albedo
ALLOCATE( zalb(A2D(0)), zmskalb(A2D(0)) )
! ice albedo
WHERE( at_i_b < 1.e-03 )
WHERE( at_i_b(:,:) < 1.e-03 )
zmskalb(:,:) = 0._wp
zalb (:,:) = rn_alb_oce
ELSEWHERE
zmskalb(:,:) = 1._wp
zalb (:,:) = SUM( alb_ice * a_i_b, dim=3 ) / at_i_b
zalb (:,:) = SUM( alb_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) / at_i_b(:,:)
END WHERE
CALL iom_put( 'icealb' , zalb * zmskalb + zmiss_val * ( 1._wp - zmskalb ) )
! ice+ocean albedo
zalb(:,:) = SUM( alb_ice * a_i_b, dim=3 ) + rn_alb_oce * ( 1._wp - at_i_b )
zalb(:,:) = SUM( alb_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) + rn_alb_oce * ( 1._wp - at_i_b(:,:) )
CALL iom_put( 'albedo' , zalb )
DEALLOCATE( zalb, zmskalb )
ENDIF
!
! --- category-dependent fields --- !
IF( iom_use('icemask_cat' ) ) CALL iom_put( 'icemask_cat' , zmsk00l ) ! ice mask 0%
IF( iom_use('iceconc_cat' ) ) CALL iom_put( 'iceconc_cat' , a_i * zmsk00l ) ! area for categories
IF( iom_use('icethic_cat' ) ) CALL iom_put( 'icethic_cat' , h_i * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! thickness for categories
IF( iom_use('snwthic_cat' ) ) CALL iom_put( 'snwthic_cat' , h_s * zmsksnl + zmiss_val * ( 1._wp - zmsksnl ) ) ! snow depth for categories
IF( iom_use('icesalt_cat' ) ) CALL iom_put( 'icesalt_cat' , s_i * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! salinity for categories
IF( iom_use('iceage_cat' ) ) CALL iom_put( 'iceage_cat' , o_i / rday * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice age
IF( iom_use('icetemp_cat' ) ) CALL iom_put( 'icetemp_cat' , ( SUM( t_i, dim=3 ) * r1_nlay_i - rt0 ) &
& * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice temperature
IF( iom_use('snwtemp_cat' ) ) CALL iom_put( 'snwtemp_cat' , ( SUM( t_s, dim=3 ) * r1_nlay_s - rt0 ) &
& * zmsksnl + zmiss_val * ( 1._wp - zmsksnl ) ) ! snow temperature
IF( iom_use('icettop_cat' ) ) CALL iom_put( 'icettop_cat' , ( t_su - rt0 ) * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! surface temperature
IF( iom_use('icebrv_cat' ) ) CALL iom_put( 'icebrv_cat' , bv_i * 100. * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! brine volume
IF( iom_use('iceapnd_cat' ) ) CALL iom_put( 'iceapnd_cat' , a_ip * zmsk00l ) ! melt pond frac for categories
IF( iom_use('icevpnd_cat' ) ) CALL iom_put( 'icevpnd_cat' , v_ip * zmsk00l ) ! melt pond volume for categories
IF( iom_use('icehpnd_cat' ) ) CALL iom_put( 'icehpnd_cat' , h_ip * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! melt pond thickness for categories
IF( iom_use('icehlid_cat' ) ) CALL iom_put( 'icehlid_cat' , h_il * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! melt pond lid thickness for categories
IF( iom_use('iceafpnd_cat') ) CALL iom_put( 'iceafpnd_cat', a_ip_frac * zmsk00l ) ! melt pond frac per ice area for categories
IF( iom_use('iceaepnd_cat') ) CALL iom_put( 'iceaepnd_cat', a_ip_eff * zmsk00l ) ! melt pond effective frac for categories
IF( iom_use('icealb_cat' ) ) CALL iom_put( 'icealb_cat' , alb_ice * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice albedo for categories
IF( iom_use('icemask_cat' ) ) CALL iom_put( 'icemask_cat' , zmsk00l ) ! ice mask 0%
IF( iom_use('iceconc_cat' ) ) CALL iom_put( 'iceconc_cat' , a_i(A2D(0),:) * zmsk00l ) ! area for categories
IF( iom_use('icethic_cat' ) ) CALL iom_put( 'icethic_cat' , h_i(A2D(0),:) * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! thickness for categories
IF( iom_use('snwthic_cat' ) ) CALL iom_put( 'snwthic_cat' , h_s(A2D(0),:) * zmsksnl + zmiss_val &
& * ( 1._wp - zmsksnl ) ) ! snow depth for categories
IF( iom_use('icesalt_cat' ) ) CALL iom_put( 'icesalt_cat' , s_i(A2D(0),:) * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! salinity for categories
IF( iom_use('iceage_cat' ) ) CALL iom_put( 'iceage_cat' , o_i(A2D(0),:) / rday * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! ice age
IF( iom_use('icetemp_cat' ) ) CALL iom_put( 'icetemp_cat' , ( SUM( t_i(A2D(0),:,:), dim=3 ) * r1_nlay_i - rt0 ) &
& * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice temperature
IF( iom_use('snwtemp_cat' ) ) CALL iom_put( 'snwtemp_cat' , ( SUM( t_s(A2D(0),:,:), dim=3 ) * r1_nlay_s - rt0 ) &
& * zmsksnl + zmiss_val * ( 1._wp - zmsksnl ) ) ! snow temperature
IF( iom_use('icettop_cat' ) ) CALL iom_put( 'icettop_cat' , ( t_su(A2D(0),:) - rt0 ) * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! surface temperature
IF( iom_use('icebrv_cat' ) ) CALL iom_put( 'icebrv_cat' , bv_i(:,:,:) * 100. * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! brine volume
IF( iom_use('iceapnd_cat' ) ) CALL iom_put( 'iceapnd_cat' , a_ip(A2D(0),:) * zmsk00l ) ! melt pond frac for categories
IF( iom_use('icevpnd_cat' ) ) CALL iom_put( 'icevpnd_cat' , v_ip(A2D(0),:) * zmsk00l ) ! melt pond volume for categories
IF( iom_use('icehpnd_cat' ) ) CALL iom_put( 'icehpnd_cat' , h_ip(A2D(0),:) * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! melt pond thickness for categories
IF( iom_use('icehlid_cat' ) ) CALL iom_put( 'icehlid_cat' , h_il(A2D(0),:) * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! melt pond lid thickness for categories
IF( iom_use('iceafpnd_cat') ) CALL iom_put( 'iceafpnd_cat', a_ip_frac(:,:,:) * zmsk00l ) ! melt pond frac per ice area for categories
IF( iom_use('iceaepnd_cat') ) CALL iom_put( 'iceaepnd_cat', a_ip_eff(:,:,:) * zmsk00l ) ! melt pond effective frac for categories
IF( iom_use('icealb_cat' ) ) CALL iom_put( 'icealb_cat' , alb_ice(:,:,:) * zmsk00l + zmiss_val &
& * ( 1._wp - zmsk00l ) ) ! ice albedo for categories
!------------------
! Add-ons for SIMIP
!------------------
! trends
IF( iom_use('dmithd') ) CALL iom_put( 'dmithd', - wfx_bog - wfx_bom - wfx_sum - wfx_sni - wfx_opw - wfx_lam - wfx_res ) ! Sea-ice mass change from thermodynamics
IF( iom_use('dmithd') ) CALL iom_put( 'dmithd', - wfx_bog(:,:) - wfx_bom(:,:) - wfx_sum(:,:) - wfx_sni(:,:) &
& - wfx_opw(:,:) - wfx_lam(:,:) - wfx_res(A2D(0)) ) ! Sea-ice mass change from thermodynamics
IF( iom_use('dmidyn') ) CALL iom_put( 'dmidyn', - wfx_dyn + rhoi * diag_trp_vi ) ! Sea-ice mass change from dynamics(kg/m2/s)
IF( iom_use('dmiopw') ) CALL iom_put( 'dmiopw', - wfx_opw ) ! Sea-ice mass change through growth in open water
IF( iom_use('dmibog') ) CALL iom_put( 'dmibog', - wfx_bog ) ! Sea-ice mass change through basal growth
......@@ -211,17 +219,23 @@ CONTAINS
IF( iom_use('NH_icearea') .OR. iom_use('NH_icevolu') .OR. iom_use('NH_iceextt') .OR. &
& iom_use('SH_icearea') .OR. iom_use('SH_icevolu') .OR. iom_use('SH_iceextt') ) THEN
!
WHERE( ff_t(:,:) > 0._wp ) ; z2d(:,:) = 1._wp
ELSEWHERE ; z2d(:,:) = 0.
WHERE( ff_t(A2D(0)) > 0._wp ) ; z2d(:,:) = 1._wp
ELSEWHERE ; z2d(:,:) = 0.
END WHERE
!
IF( iom_use('NH_icearea') ) zdiag_area_nh = glob_sum( 'icewri', at_i * z2d * e1e2t * 1.e-12 )
IF( iom_use('NH_icevolu') ) zdiag_volu_nh = glob_sum( 'icewri', vt_i * z2d * e1e2t * 1.e-12 )
IF( iom_use('NH_iceextt') ) zdiag_extt_nh = glob_sum( 'icewri', z2d * e1e2t * 1.e-12 * zmsk15 )
IF( iom_use('NH_icearea') ) zdiag_area_nh = glob_sum( 'icewri', at_i(A2D(0)) * z2d * e1e2t(A2D(0)) &
& * 1.e-12 )
IF( iom_use('NH_icevolu') ) zdiag_volu_nh = glob_sum( 'icewri', vt_i(A2D(0)) * z2d * e1e2t(A2D(0)) &
& * 1.e-12 )
IF( iom_use('NH_iceextt') ) zdiag_extt_nh = glob_sum( 'icewri', z2d * e1e2t(A2D(0)) &
& * 1.e-12 * zmsk15 )
!
IF( iom_use('SH_icearea') ) zdiag_area_sh = glob_sum( 'icewri', at_i * ( 1._wp - z2d ) * e1e2t * 1.e-12 )
IF( iom_use('SH_icevolu') ) zdiag_volu_sh = glob_sum( 'icewri', vt_i * ( 1._wp - z2d ) * e1e2t * 1.e-12 )
IF( iom_use('SH_iceextt') ) zdiag_extt_sh = glob_sum( 'icewri', ( 1._wp - z2d ) * e1e2t * 1.e-12 * zmsk15 )
IF( iom_use('SH_icearea') ) zdiag_area_sh = glob_sum( 'icewri', at_i(A2D(0)) * ( 1._wp - z2d ) * e1e2t(A2D(0)) &
& * 1.e-12 )
IF( iom_use('SH_icevolu') ) zdiag_volu_sh = glob_sum( 'icewri', vt_i(A2D(0)) * ( 1._wp - z2d ) * e1e2t(A2D(0)) &
& * 1.e-12 )
IF( iom_use('SH_iceextt') ) zdiag_extt_sh = glob_sum( 'icewri', ( 1._wp - z2d ) * e1e2t(A2D(0)) &
& * 1.e-12 * zmsk15 )
!
CALL iom_put( 'NH_icearea' , zdiag_area_nh )
CALL iom_put( 'NH_icevolu' , zdiag_volu_nh )
......@@ -258,25 +272,25 @@ CONTAINS
!
!! The file is open in dia_wri_state (ocean routine)
CALL iom_rstput( 0, 0, kid, 'sithic', hm_i ) ! Ice thickness
CALL iom_rstput( 0, 0, kid, 'siconc', at_i ) ! Ice concentration
CALL iom_rstput( 0, 0, kid, 'sitemp', tm_i - rt0 ) ! Ice temperature
CALL iom_rstput( 0, 0, kid, 'sivelu', u_ice ) ! i-Ice speed
CALL iom_rstput( 0, 0, kid, 'sivelv', v_ice ) ! j-Ice speed
CALL iom_rstput( 0, 0, kid, 'sistru', utau_ice ) ! i-Wind stress over ice
CALL iom_rstput( 0, 0, kid, 'sistrv', vtau_ice ) ! i-Wind stress over ice
CALL iom_rstput( 0, 0, kid, 'sisflx', qsr ) ! Solar flx over ocean
CALL iom_rstput( 0, 0, kid, 'sinflx', qns ) ! NonSolar flx over ocean
CALL iom_rstput( 0, 0, kid, 'snwpre', sprecip ) ! Snow precipitation
CALL iom_rstput( 0, 0, kid, 'sisali', sm_i ) ! Ice salinity
CALL iom_rstput( 0, 0, kid, 'sivolu', vt_i ) ! Ice volume
CALL iom_rstput( 0, 0, kid, 'sidive', divu_i*1.0e8 ) ! Ice divergence
CALL iom_rstput( 0, 0, kid, 'si_amp', at_ip ) ! Melt pond fraction
CALL iom_rstput( 0, 0, kid, 'si_vmp', vt_ip ) ! Melt pond volume
CALL iom_rstput( 0, 0, kid, 'sithicat', h_i ) ! Ice thickness
CALL iom_rstput( 0, 0, kid, 'siconcat', a_i ) ! Ice concentration
CALL iom_rstput( 0, 0, kid, 'sisalcat', s_i ) ! Ice salinity
CALL iom_rstput( 0, 0, kid, 'snthicat', h_s ) ! Snw thickness
CALL iom_rstput( 0, 0, kid, 'sithic', hm_i(A2D(0)) ) ! Ice thickness
CALL iom_rstput( 0, 0, kid, 'siconc', at_i(A2D(0)) ) ! Ice concentration
CALL iom_rstput( 0, 0, kid, 'sitemp', tm_i(A2D(0)) - rt0 ) ! Ice temperature
CALL iom_rstput( 0, 0, kid, 'sivelu', u_ice(A2D(0)) ) ! i-Ice speed
CALL iom_rstput( 0, 0, kid, 'sivelv', v_ice(A2D(0)) ) ! j-Ice speed
CALL iom_rstput( 0, 0, kid, 'sistru', utau_ice(A2D(0)) ) ! i-Wind stress over ice
CALL iom_rstput( 0, 0, kid, 'sistrv', vtau_ice(A2D(0)) ) ! i-Wind stress over ice
CALL iom_rstput( 0, 0, kid, 'sisflx', qsr(A2D(0)) ) ! Solar flx over ocean
CALL iom_rstput( 0, 0, kid, 'sinflx', qns(A2D(0)) ) ! NonSolar flx over ocean
CALL iom_rstput( 0, 0, kid, 'snwpre', sprecip(A2D(0)) ) ! Snow precipitation
CALL iom_rstput( 0, 0, kid, 'sisali', sm_i(A2D(0)) ) ! Ice salinity
CALL iom_rstput( 0, 0, kid, 'sivolu', vt_i(A2D(0)) ) ! Ice volume
CALL iom_rstput( 0, 0, kid, 'sidive', divu_i(A2D(0))*1.0e8 ) ! Ice divergence
CALL iom_rstput( 0, 0, kid, 'si_amp', at_ip(A2D(0)) ) ! Melt pond fraction
CALL iom_rstput( 0, 0, kid, 'si_vmp', vt_ip(A2D(0)) ) ! Melt pond volume
CALL iom_rstput( 0, 0, kid, 'sithicat', h_i(A2D(0),:) ) ! Ice thickness
CALL iom_rstput( 0, 0, kid, 'siconcat', a_i(A2D(0),:) ) ! Ice concentration
CALL iom_rstput( 0, 0, kid, 'sisalcat', s_i(A2D(0),:) ) ! Ice salinity
CALL iom_rstput( 0, 0, kid, 'snthicat', h_s(A2D(0),:) ) ! Snw thickness
END SUBROUTINE ice_wri_state
......
......@@ -112,12 +112,10 @@ CONTAINS
CALL dom_qco_zgr( Kbb_a, Kmm_a )
#endif
#if defined key_si3
CALL lbc_lnk( 'finalize_lbc_for_agrif', a_i, 'T',1._wp, v_i,'T',1._wp, &
& v_s, 'T',1._wp, sv_i,'T',1._wp, oa_i,'T',1._wp, &
& a_ip,'T',1._wp, v_ip,'T',1._wp, v_il,'T',1._wp )
CALL lbc_lnk( 'finalize_lbc_for_agrif', t_su,'T',1._wp )
CALL lbc_lnk( 'finalize_lbc_for_agrif', e_s,'T',1._wp )
CALL lbc_lnk( 'finalize_lbc_for_agrif', e_i,'T',1._wp )
CALL lbc_lnk( 'finalize_lbc_for_agrif', a_i, 'T',1._wp, v_i,'T',1._wp, &
& v_s, 'T',1._wp, sv_i,'T',1._wp, oa_i,'T',1._wp, &
& a_ip,'T',1._wp, v_ip,'T',1._wp, v_il,'T',1._wp, t_su,'T',1._wp )
CALL lbc_lnk( 'finalize_lbc_for_agrif', e_i,'T',1._wp, e_s,'T',1._wp )
CALL lbc_lnk( 'finalize_lbc_for_agrif', u_ice, 'U', -1._wp, v_ice, 'V', -1._wp )
#endif
#if defined key_top
......
......@@ -59,12 +59,10 @@ CONTAINS
Agrif_UseSpecialValue = .TRUE.
CALL Agrif_init_variable(tra_iceini_id,procname=interp_tra_ice)
!
CALL lbc_lnk( 'agrif_istate_ice', a_i,'T',1._wp, v_i,'T',1._wp, &
& v_s,'T',1._wp, sv_i,'T',1._wp, oa_i,'T',1._wp, &
& a_ip,'T',1._wp, v_ip,'T',1._wp, v_il,'T',1._wp )
CALL lbc_lnk( 'agrif_istate_ice', t_su,'T',1._wp )
CALL lbc_lnk( 'agrif_istate_ice', e_s,'T',1._wp )
CALL lbc_lnk( 'agrif_istate_ice', e_i,'T',1._wp )
CALL lbc_lnk( 'agrif_istate_ice', a_i,'T',1._wp, v_i,'T',1._wp, &
& v_s,'T',1._wp, sv_i,'T',1._wp, oa_i,'T',1._wp, &
& a_ip,'T',1._wp, v_ip,'T',1._wp, v_il,'T',1._wp, t_su,'T',1._wp )
CALL lbc_lnk( 'agrif_istate_ice', e_i,'T',1._wp, e_s,'T',1._wp )
!
! Set u_ice, v_ice:
use_sign_north = .TRUE.
......
......@@ -270,7 +270,7 @@ CONTAINS
ibdy2 = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhox() ! halo + land + nbghostcells
!
IF( .NOT.ln_dynspg_ts ) THEN ! Store transport
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
DO jj = 1, jpj
uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
......@@ -278,7 +278,7 @@ CONTAINS
END DO
ENDIF
!
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
zub(ji,:) = 0._wp
zhub(ji,:) = 0._wp
DO jk = 1, jpkm1
......@@ -300,7 +300,7 @@ CONTAINS
END DO
END DO
!
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
zvb(ji,:) = 0._wp
zhvb(ji,:) = 0._wp
DO jk = 1, jpkm1
......@@ -330,14 +330,14 @@ CONTAINS
ibdy2 = jpiglo - ( nn_hls + 2 )
!
IF( .NOT.ln_dynspg_ts ) THEN
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
DO jj = 1, jpj
uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
END DO
END DO
ENDIF
!
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
zub(ji,:) = 0._wp
zhub(ji,:) = 0._wp
DO jk = 1, jpkm1
......@@ -363,14 +363,14 @@ CONTAINS
ibdy2 = jpiglo - ( nn_hls + 1 )
!
IF( .NOT.ln_dynspg_ts ) THEN
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
DO jj = 1, jpj
vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
END DO
END DO
ENDIF
!
DO ji = mi0(ibdy1), mi1(ibdy2)
DO ji = mi0(ibdy1,nn_hls), mi1(ibdy2,nn_hls)
zvb(ji,:) = 0._wp
zhvb(ji,:) = 0._wp
DO jk = 1, jpkm1
......@@ -400,7 +400,7 @@ CONTAINS
jbdy2 = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhoy()
!
IF( .NOT.ln_dynspg_ts ) THEN
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
DO ji = 1, jpi
uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
......@@ -408,7 +408,7 @@ CONTAINS
END DO
ENDIF
!
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
zvb(:,jj) = 0._wp
zhvb(:,jj) = 0._wp
DO jk=1,jpkm1
......@@ -430,7 +430,7 @@ CONTAINS
END DO
END DO
!
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
zub(:,jj) = 0._wp
zhub(:,jj) = 0._wp
DO jk = 1, jpkm1
......@@ -460,14 +460,14 @@ CONTAINS
jbdy2 = jpjglo - ( nn_hls + 2 )
!
IF( .NOT.ln_dynspg_ts ) THEN
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
DO ji = 1, jpi
vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
END DO
END DO
ENDIF
!
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
zvb(:,jj) = 0._wp
zhvb(:,jj) = 0._wp
DO jk=1,jpkm1
......@@ -493,14 +493,14 @@ CONTAINS
jbdy2 = jpjglo - ( nn_hls + 1 )
!
IF( .NOT.ln_dynspg_ts ) THEN
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
DO ji = 1, jpi
uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
END DO
END DO
ENDIF
!
DO jj = mj0(jbdy1), mj1(jbdy2)
DO jj = mj0(jbdy1,nn_hls), mj1(jbdy2,nn_hls)
zub(:,jj) = 0._wp
zhub(:,jj) = 0._wp
DO jk = 1, jpkm1
......@@ -553,7 +553,7 @@ CONTAINS
IF( lk_west ) THEN
istart = nn_hls + 2 ! halo + land + 1
iend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhox() ! halo + land + nbghostcells
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj=1,jpj
va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
......@@ -565,7 +565,7 @@ CONTAINS
IF( lk_east ) THEN
istart = jpiglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhox()
iend = jpiglo - ( nn_hls + 1 )
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj=1,jpj
va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
......@@ -573,7 +573,7 @@ CONTAINS
END DO
istart = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox()
iend = jpiglo - ( nn_hls + 2 )
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj=1,jpj
ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
END DO
......@@ -584,7 +584,7 @@ CONTAINS
IF( lk_south ) THEN
jstart = nn_hls + 2
jend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhoy()
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji=1,jpi
ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
......@@ -597,14 +597,14 @@ CONTAINS
IF( lk_north ) THEN
jstart = jpjglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhoy()
jend = jpjglo - ( nn_hls + 1 )
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji=1,jpi
ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
END DO
END DO
jstart = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy()
jend = jpjglo - ( nn_hls + 2 )
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji=1,jpi
va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
END DO
......@@ -642,7 +642,7 @@ CONTAINS
IF( lk_west ) THEN
istart = nn_hls + 2
iend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhox()
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj=1,jpj
zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
......@@ -654,14 +654,14 @@ CONTAINS
IF( lk_east ) THEN
istart = jpiglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhox()
iend = jpiglo - ( nn_hls + 1 )
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj=1,jpj
zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
END DO
END DO
istart = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox()
iend = jpiglo - ( nn_hls + 2 )
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj=1,jpj
zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
END DO
......@@ -672,7 +672,7 @@ CONTAINS
IF( lk_south ) THEN
jstart = nn_hls + 2
jend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhoy()
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji=1,jpi
zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
......@@ -684,14 +684,14 @@ CONTAINS
IF( lk_north ) THEN
jstart = jpjglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhoy()
jend = jpjglo - ( nn_hls + 1 )
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji=1,jpi
zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
END DO
END DO
jstart = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy()
jend = jpjglo - ( nn_hls + 2 )
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji=1,jpi
zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
END DO
......@@ -801,7 +801,7 @@ CONTAINS
istart = nn_hls + 2 ! halo + land + 1
iend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhox() ! halo + land + nbghostcells
IF (lk_div_cons) iend = istart
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj = 1, jpj
ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
END DO
......@@ -813,7 +813,7 @@ CONTAINS
istart = jpiglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhox() ! halo + land + nbghostcells - 1
iend = jpiglo - ( nn_hls + 1 ) ! halo + land + 1 - 1
IF (lk_div_cons) istart = iend
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj = 1, jpj
ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
END DO
......@@ -825,7 +825,7 @@ CONTAINS
jstart = nn_hls + 2 ! halo + land + 1
jend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhoy() ! halo + land + nbghostcells
IF (lk_div_cons) jend = jstart
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji = 1, jpi
ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
END DO
......@@ -837,7 +837,7 @@ CONTAINS
jstart = jpjglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhoy() ! halo + land + nbghostcells - 1
jend = jpjglo - ( nn_hls + 1 ) ! halo + land + 1 - 1
IF (lk_div_cons) jstart = jend
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji = 1, jpi
ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
END DO
......@@ -870,7 +870,7 @@ CONTAINS
istart = nn_hls + 2 ! halo + land + 1
iend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhox() ! halo + land + nbghostcells
IF (lk_div_cons) iend = istart
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj = 1, jpj
ssha_e(ji,jj) = hbdy(ji,jj)
END DO
......@@ -882,7 +882,7 @@ CONTAINS
istart = jpiglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhox() ! halo + land + nbghostcells - 1
iend = jpiglo - ( nn_hls + 1 ) ! halo + land + 1 - 1
IF (lk_div_cons) istart = iend
DO ji = mi0(istart), mi1(iend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj = 1, jpj
ssha_e(ji,jj) = hbdy(ji,jj)
END DO
......@@ -894,7 +894,7 @@ CONTAINS
jstart = nn_hls + 2 ! halo + land + 1
jend = nn_hls + nbghostcells + nn_shift_bar*Agrif_Rhoy() ! halo + land + nbghostcells
IF (lk_div_cons) jend = jstart
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji = 1, jpi
ssha_e(ji,jj) = hbdy(ji,jj)
END DO
......@@ -906,7 +906,7 @@ CONTAINS
jstart = jpjglo - ( nn_hls + nbghostcells -1 ) - nn_shift_bar*Agrif_Rhoy() ! halo + land + nbghostcells - 1
jend = jpjglo - ( nn_hls + 1 ) ! halo + land + 1 - 1
IF (lk_div_cons) jstart = jend
DO jj = mj0(jstart), mj1(jend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji = 1, jpi
ssha_e(ji,jj) = hbdy(ji,jj)
END DO
......@@ -1551,7 +1551,7 @@ CONTAINS
DO ji=i1,i2
DO jj=j1,j2
IF (utint_stage(ji,jj)==0) THEN
zx = 2._wp*MOD(ABS(mig0(ji)-nbghostcells_x_w), INT(Agrif_Rhox()))/zrhox - 1._wp
zx = 2._wp*MOD(ABS(mig(ji,0)-nbghostcells_x_w), INT(Agrif_Rhox()))/zrhox - 1._wp
ubdy(ji,jj) = ubdy(ji,jj) + 0.25_wp*(1._wp-zx*zx) * ptab(ji,jj) &
& / zrhoy *r1_e2u(ji,jj) * umask(ji,jj,1)
utint_stage(ji,jj) = 1
......@@ -1671,7 +1671,7 @@ CONTAINS
DO ji=i1,i2
DO jj=j1,j2
IF (vtint_stage(ji,jj)==0) THEN
zy = 2._wp*MOD(ABS(mjg0(jj)-nbghostcells_y_s), INT(Agrif_Rhoy()))/zrhoy - 1._wp
zy = 2._wp*MOD(ABS(mjg(jj,0)-nbghostcells_y_s), INT(Agrif_Rhoy()))/zrhoy - 1._wp
vbdy(ji,jj) = vbdy(ji,jj) + 0.25_wp*(1._wp-zy*zy) * ptab(ji,jj) &
& / zrhox * r1_e1v(ji,jj) * vmask(ji,jj,1)
vtint_stage(ji,jj) = 1
......@@ -1755,7 +1755,7 @@ CONTAINS
DO jj = j1, j2
DO ji = i1, i2
IF( ABS( ptab(ji,jj) - glamt(ji,jj) ) > ztst ) THEN
WRITE(numout,*) ' Agrif error for glamt: parent, child, i, j ', ptab(ji,jj), glamt(ji,jj), mig0(ji), mig0(jj)
WRITE(numout,*) ' Agrif error for glamt: parent, child, i, j ', ptab(ji,jj), glamt(ji,jj), mig(ji,0), mjg(jj,0)
! kindic_agr = kindic_agr + 1
ENDIF
END DO
......@@ -1784,7 +1784,7 @@ CONTAINS
DO jj = j1, j2
DO ji = i1, i2
IF( ABS( ptab(ji,jj) - gphit(ji,jj) ) > ztst ) THEN
WRITE(numout,*) ' Agrif error for gphit: parent, child, i, j ', ptab(ji,jj), gphit(ji,jj), mig0(ji), mig0(jj)
WRITE(numout,*) ' Agrif error for gphit: parent, child, i, j ', ptab(ji,jj), gphit(ji,jj), mig(ji,0), mjg(jj,0)
! kindic_agr = kindic_agr + 1
ENDIF
END DO
......@@ -1999,8 +1999,8 @@ CONTAINS
iend = nn_hls + nbghostcells + ispon ! halo + land + nbghostcells + sponge
jstart = nn_hls + 2
jend = jpjglo - nn_hls - 1
DO ji = mi0(istart), mi1(iend)
DO jj = mj0(jstart), mj1(jend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2008,7 +2008,7 @@ CONTAINS
END DO
ENDIF
END DO
DO jj = mj0(jstart), mj1(jend-1)
DO jj = mj0(jstart,nn_hls), mj1(jend-1,nn_hls)
IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2017,8 +2017,8 @@ CONTAINS
ENDIF
END DO
END DO
DO ji = mi0(istart), mi1(iend-1)
DO jj = mj0(jstart), mj1(jend)
DO ji = mi0(istart,nn_hls), mi1(iend-1,nn_hls)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2036,8 +2036,8 @@ CONTAINS
iend = jpiglo - nn_hls - 1 ! halo + land + 1 - 1
jstart = nn_hls + 2
jend = jpjglo - nn_hls - 1
DO ji = mi0(istart), mi1(iend)
DO jj = mj0(jstart), mj1(jend)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2045,7 +2045,7 @@ CONTAINS
END DO
ENDIF
END DO
DO jj = mj0(jstart), mj1(jend-1)
DO jj = mj0(jstart,nn_hls), mj1(jend-1,nn_hls)
IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2054,8 +2054,8 @@ CONTAINS
ENDIF
END DO
END DO
DO ji = mi0(istart), mi1(iend-1)
DO jj = mj0(jstart), mj1(jend)
DO ji = mi0(istart,nn_hls), mi1(iend-1,nn_hls)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2073,8 +2073,8 @@ CONTAINS
jend = nn_hls + nbghostcells + ispon ! halo + land + nbghostcells + sponge
istart = nn_hls + 2
iend = jpiglo - nn_hls - 1
DO jj = mj0(jstart), mj1(jend)
DO ji = mi0(istart), mi1(iend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2082,7 +2082,7 @@ CONTAINS
END DO
ENDIF
END DO
DO ji = mi0(istart), mi1(iend-1)
DO ji = mi0(istart,nn_hls), mi1(iend-1,nn_hls)
IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2091,8 +2091,8 @@ CONTAINS
ENDIF
END DO
END DO
DO jj = mj0(jstart), mj1(jend-1)
DO ji = mi0(istart), mi1(iend)
DO jj = mj0(jstart,nn_hls), mj1(jend-1,nn_hls)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2110,8 +2110,8 @@ CONTAINS
jend = jpjglo - nn_hls - 1 ! halo + land + 1 - 1
istart = nn_hls + 2
iend = jpiglo - nn_hls - 1
DO jj = mj0(jstart), mj1(jend)
DO ji = mi0(istart), mi1(iend)
DO jj = mj0(jstart,nn_hls), mj1(jend,nn_hls)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2119,7 +2119,7 @@ CONTAINS
END DO
ENDIF
END DO
DO ji = mi0(istart), mi1(iend-1)
DO ji = mi0(istart,nn_hls), mi1(iend-1,nn_hls)
IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......@@ -2128,8 +2128,8 @@ CONTAINS
ENDIF
END DO
END DO
DO jj = mj0(jstart), mj1(jend-1)
DO ji = mi0(istart), mi1(iend)
DO jj = mj0(jstart,nn_hls), mj1(jend-1,nn_hls)
DO ji = mi0(istart,nn_hls), mi1(iend,nn_hls)
IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
IF ( .NOT.ln_vert_remap) THEN
DO jk = 1, jpkm1
......
......@@ -161,15 +161,15 @@ CONTAINS
IF( lk_west ) THEN ! --- West --- !
ind1 = nn_hls + nbghostcells ! halo + nbghostcells
ind2 = nn_hls + nbghostcells + ispongearea
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = REAL(ind2 - mig(ji), wp) * z1_ispongearea
ztabramp(ji,jj) = REAL(ind2 - mig(ji,nn_hls), wp) * z1_ispongearea
END DO
END DO
! ghost cells:
ind1 = 1
ind2 = nn_hls + nbghostcells ! halo + nbghostcells
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = 1._wp
END DO
......@@ -178,15 +178,15 @@ CONTAINS
IF( lk_east ) THEN ! --- East --- !
ind1 = jpiglo - ( nn_hls + nbghostcells -1 ) - ispongearea - 1
ind2 = jpiglo - ( nn_hls + nbghostcells -1 ) - 1 ! halo + land + nbghostcells - 1
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mig(ji) - ind1, wp) * z1_ispongearea )
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mig(ji,nn_hls) - ind1, wp) * z1_ispongearea )
END DO
END DO
! ghost cells:
ind1 = jpiglo - ( nn_hls + nbghostcells -1 ) - 1 ! halo + land + nbghostcells - 1
ind2 = jpiglo - 1
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = 1._wp
END DO
......@@ -195,15 +195,15 @@ CONTAINS
IF( lk_south ) THEN ! --- South --- !
ind1 = nn_hls + nbghostcells ! halo + nbghostcells
ind2 = nn_hls + nbghostcells + jspongearea
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(ind2 - mjg(jj), wp) * z1_jspongearea )
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(ind2 - mjg(jj,nn_hls), wp) * z1_jspongearea )
END DO
END DO
! ghost cells:
ind1 = 1
ind2 = nn_hls + nbghostcells ! halo + nbghostcells
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = 1._wp
END DO
......@@ -212,15 +212,15 @@ CONTAINS
IF( lk_north ) THEN ! --- North --- !
ind1 = jpjglo - ( nn_hls + nbghostcells -1 ) - jspongearea - 1
ind2 = jpjglo - ( nn_hls + nbghostcells -1 ) - 1 ! halo + nbghostcells - 1
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mjg(jj) - ind1, wp) * z1_jspongearea )
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mjg(jj,nn_hls) - ind1, wp) * z1_jspongearea )
END DO
END DO
! ghost cells:
ind1 = jpjglo - ( nn_hls + nbghostcells -1 ) ! halo + land + nbghostcells - 1
ind2 = jpjglo
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = 1._wp
END DO
......@@ -294,15 +294,15 @@ CONTAINS
IF( lk_west ) THEN ! --- West --- !
ind1 = nn_hls + nbghostcells + ishift
ind2 = nn_hls + nbghostcells + ishift + ispongearea
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = REAL(ind2 - mig(ji), wp) * z1_ispongearea
ztabramp(ji,jj) = REAL(ind2 - mig(ji,nn_hls), wp) * z1_ispongearea
END DO
END DO
! ghost cells:
ind1 = 1
ind2 = nn_hls + nbghostcells + ishift ! halo + nbghostcells
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = 1._wp
END DO
......@@ -311,15 +311,15 @@ CONTAINS
IF( lk_east ) THEN ! --- East --- !
ind1 = jpiglo - ( nn_hls + nbghostcells -1 + ishift) - ispongearea - 1
ind2 = jpiglo - ( nn_hls + nbghostcells -1 + ishift) - 1 ! halo + nbghostcells - 1
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mig(ji) - ind1, wp) * z1_ispongearea )
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mig(ji,nn_hls) - ind1, wp) * z1_ispongearea )
END DO
END DO
! ghost cells:
ind1 = jpiglo - ( nn_hls + nbghostcells -1 + ishift) - 1 ! halo + nbghostcells - 1
ind2 = jpiglo - 1
DO ji = mi0(ind1), mi1(ind2)
DO ji = mi0(ind1,nn_hls), mi1(ind2,nn_hls)
DO jj = 1, jpj
ztabramp(ji,jj) = 1._wp
END DO
......@@ -328,15 +328,15 @@ CONTAINS
IF( lk_south ) THEN ! --- South --- !
ind1 = nn_hls + nbghostcells + jshift ! halo + nbghostcells
ind2 = nn_hls + nbghostcells + jshift + jspongearea
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(ind2 - mjg(jj), wp) * z1_jspongearea )
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(ind2 - mjg(jj,nn_hls), wp) * z1_jspongearea )
END DO
END DO
! ghost cells:
ind1 = 1
ind2 = nn_hls + nbghostcells + jshift ! halo + land + nbghostcells
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = 1._wp
END DO
......@@ -345,15 +345,15 @@ CONTAINS
IF( lk_north ) THEN ! --- North --- !
ind1 = jpjglo - ( nn_hls + nbghostcells -1 + jshift) - jspongearea - 1
ind2 = jpjglo - ( nn_hls + nbghostcells -1 + jshift) - 1 ! halo + land + nbghostcells - 1
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mjg(jj) - ind1, wp) * z1_jspongearea )
ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL(mjg(jj,nn_hls) - ind1, wp) * z1_jspongearea )
END DO
END DO
! ghost cells:
ind1 = jpjglo - ( nn_hls + nbghostcells -1 + jshift) ! halo + land + nbghostcells - 1
ind2 = jpjglo
DO jj = mj0(ind1), mj1(ind2)
DO jj = mj0(ind1,nn_hls), mj1(ind2,nn_hls)
DO ji = 1, jpi
ztabramp(ji,jj) = 1._wp
END DO
......@@ -730,7 +730,7 @@ CONTAINS
jmax = j2-1
ind1 = jpjglo - ( nn_hls + nbghostcells + 1 ) ! North
DO jj = mj0(ind1), mj1(ind1)
DO jj = mj0(ind1,nn_hls), mj1(ind1,nn_hls)
jmax = MIN(jmax,jj)
END DO
......@@ -858,7 +858,7 @@ CONTAINS
imax = i2 - 1
ind1 = jpiglo - ( nn_hls + nbghostcells + 1 ) ! East
DO ji = mi0(ind1), mi1(ind1)
DO ji = mi0(ind1,nn_hls), mi1(ind1,nn_hls)
imax = MIN(imax,ji)
END DO
......@@ -958,7 +958,7 @@ CONTAINS
jmax = j2-1
ind1 = jpjglo - ( nn_hls + nbghostcells + 1 ) ! North
DO jj = mj0(ind1), mj1(ind1)
DO jj = mj0(ind1,nn_hls), mj1(ind1,nn_hls)
jmax = MIN(jmax,jj)
END DO
......@@ -1025,7 +1025,7 @@ CONTAINS
imax = i2 - 1
ind1 = jpiglo - ( nn_hls + nbghostcells + 1 ) ! East
DO ji = mi0(ind1), mi1(ind1)
DO ji = mi0(ind1,nn_hls), mi1(ind1,nn_hls)
imax = MIN(imax,ji)
END DO
......
......@@ -1893,7 +1893,7 @@ CONTAINS
DO jk=k1,k2-1
IF (ABS((ptab(ji,jj,jk)-e3u_0(ji,jj,jk))*umask(ji,jj,jk)).GE.1.e-6) THEN
kindic_agr = kindic_agr + 1
print *, 'erro u-pt', mig0(ji), mjg0(jj), jk, mbku(ji,jj), ikbot, ptab(ji,jj,jk), e3u_0(ji,jj,jk)
PRINT *, 'erro u-pt', mig(ji,0), mjg(jj,0), jk, mbku(ji,jj), ikbot, ptab(ji,jj,jk), e3u_0(ji,jj,jk)
ENDIF
END DO
ENDIF
......@@ -1933,7 +1933,7 @@ CONTAINS
DO jk=k1,k2-1
IF (ABS((ptab(ji,jj,jk)-e3v_0(ji,jj,jk))*vmask(ji,jj,jk)).GE.1.e-6) THEN
kindic_agr = kindic_agr + 1
print *, 'erro v-pt', mig0(ji), mjg0(jj), mbkv(ji,jj), ptab(ji,jj,jk), e3v_0(ji,jj,jk)
PRINT *, 'erro v-pt', mig(ji,0), mjg(jj,0), mbkv(ji,jj), ptab(ji,jj,jk), e3v_0(ji,jj,jk)
ENDIF
END DO
ENDIF
......
......@@ -1095,8 +1095,8 @@
!!----------------------------------------------------------------------
!
SELECT CASE( i )
CASE(1) ; indglob = mig(indloc)
CASE(2) ; indglob = mjg(indloc)
CASE(1) ; indglob = mig(indloc,nn_hls)
CASE(2) ; indglob = mjg(indloc,nn_hls)
CASE DEFAULT ; indglob = indloc
END SELECT
!
......@@ -1115,10 +1115,10 @@
INTEGER, INTENT(out) :: jmin, jmax
!!----------------------------------------------------------------------
!
imin = mig( 1 )
jmin = mjg( 1 )
imax = mig(jpi)
jmax = mjg(jpj)
imin = mig( 1 ,nn_hls)
jmin = mjg( 1 ,nn_hls)
imax = mig(jpi,nn_hls)
jmax = mjg(jpj,nn_hls)
!
END SUBROUTINE Agrif_get_proc_info
......
......@@ -491,10 +491,10 @@ CONTAINS
! Find lenght of boundaries and rim on local mpi domain
!------------------------------------------------------
!
iwe = mig(1)
ies = mig(jpi)
iso = mjg(1)
ino = mjg(jpj)
iwe = mig( 1,nn_hls)
ies = mig(jpi,nn_hls)
iso = mjg( 1,nn_hls)
ino = mjg(jpj,nn_hls)
!
DO ib_bdy = 1, nb_bdy
DO igrd = 1, jpbgrd
......@@ -554,8 +554,8 @@ CONTAINS
& nbrdta(ib,igrd,ib_bdy) == ir ) THEN
!
icount = icount + 1
idx_bdy(ib_bdy)%nbi(icount,igrd) = nbidta(ib,igrd,ib_bdy) - mig(1) + 1 ! global to local indexes
idx_bdy(ib_bdy)%nbj(icount,igrd) = nbjdta(ib,igrd,ib_bdy) - mjg(1) + 1 ! global to local indexes
idx_bdy(ib_bdy)%nbi(icount,igrd) = nbidta(ib,igrd,ib_bdy) - mig(1,nn_hls) + 1 ! global to local indexes
idx_bdy(ib_bdy)%nbj(icount,igrd) = nbjdta(ib,igrd,ib_bdy) - mjg(1,nn_hls) + 1 ! global to local indexes
idx_bdy(ib_bdy)%nbr(icount,igrd) = nbrdta(ib,igrd,ib_bdy)
idx_bdy(ib_bdy)%nbmap(icount,igrd) = ib
ENDIF
......@@ -1014,7 +1014,7 @@ CONTAINS
DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
IF( mig0(ii) > 2 .AND. mig0(ii) < Ni0glo-2 .AND. mjg0(ij) > 2 .AND. mjg0(ij) < Nj0glo-2 ) THEN
IF( mig(ii,0) > 2 .AND. mig(ii,0) < Ni0glo-2 .AND. mjg(ij,0) > 2 .AND. mjg(ij,0) < Nj0glo-2 ) THEN
WRITE(ctmp1,*) ' Orlanski is not safe when the open boundaries are on the interior of the computational domain'
CALL ctl_stop( ctmp1 )
END IF
......@@ -1090,7 +1090,7 @@ CONTAINS
! This error check only works if you are using the bdyXmask arrays (which are set to 0 on rims)
IF( i_offset == 1 .and. zefl + zwfl == 2._wp ) THEN
icount = icount + 1
IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii),mjg(ij)
IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii,nn_hls),mjg(ij,nn_hls)
ELSE
ztmp(ii,ij) = -zwfl + zefl
ENDIF
......@@ -1130,7 +1130,7 @@ CONTAINS
znfl = zmask(ii,ij+j_offset )
! This error check only works if you are using the bdyXmask arrays (which are set to 0 on rims)
IF( j_offset == 1 .and. znfl + zsfl == 2._wp ) THEN
IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii),mjg(ij)
IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii,nn_hls),mjg(ij,nn_hls)
icount = icount + 1
ELSE
ztmp(ii,ij) = -zsfl + znfl
......@@ -1594,8 +1594,8 @@ CONTAINS
ztestmask(1:2)=0.
DO ji = 1, jpi
DO jj = 1, jpj
IF( mig0(ji) == jpiwob(ib) .AND. mjg0(jj) == jpjwdt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mig0(ji) == jpiwob(ib) .AND. mjg0(jj) == jpjwft(ib) ) ztestmask(2) = tmask(ji,jj,1)
IF( mig(ji,0) == jpiwob(ib) .AND. mjg(jj,0) == jpjwdt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mig(ji,0) == jpiwob(ib) .AND. mjg(jj,0) == jpjwft(ib) ) ztestmask(2) = tmask(ji,jj,1)
END DO
END DO
CALL mpp_sum( 'bdyini', ztestmask, 2 ) ! sum over the global domain
......@@ -1630,8 +1630,8 @@ CONTAINS
ztestmask(1:2)=0.
DO ji = 1, jpi
DO jj = 1, jpj
IF( mig0(ji) == jpieob(ib)+1 .AND. mjg0(jj) == jpjedt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mig0(ji) == jpieob(ib)+1 .AND. mjg0(jj) == jpjeft(ib) ) ztestmask(2) = tmask(ji,jj,1)
IF( mig(ji,0) == jpieob(ib)+1 .AND. mjg(jj,0) == jpjedt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mig(ji,0) == jpieob(ib)+1 .AND. mjg(jj,0) == jpjeft(ib) ) ztestmask(2) = tmask(ji,jj,1)
END DO
END DO
CALL mpp_sum( 'bdyini', ztestmask, 2 ) ! sum over the global domain
......@@ -1666,8 +1666,8 @@ CONTAINS
ztestmask(1:2)=0.
DO ji = 1, jpi
DO jj = 1, jpj
IF( mjg0(jj) == jpjsob(ib) .AND. mig0(ji) == jpisdt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mjg0(jj) == jpjsob(ib) .AND. mig0(ji) == jpisft(ib) ) ztestmask(2) = tmask(ji,jj,1)
IF( mjg(jj,0) == jpjsob(ib) .AND. mig(ji,0) == jpisdt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mjg(jj,0) == jpjsob(ib) .AND. mig(ji,0) == jpisft(ib) ) ztestmask(2) = tmask(ji,jj,1)
END DO
END DO
CALL mpp_sum( 'bdyini', ztestmask, 2 ) ! sum over the global domain
......@@ -1688,8 +1688,8 @@ CONTAINS
ztestmask(1:2)=0.
DO ji = 1, jpi
DO jj = 1, jpj
IF( mjg0(jj) == jpjnob(ib)+1 .AND. mig0(ji) == jpindt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mjg0(jj) == jpjnob(ib)+1 .AND. mig0(ji) == jpinft(ib) ) ztestmask(2) = tmask(ji,jj,1)
IF( mjg(jj,0) == jpjnob(ib)+1 .AND. mig(ji,0) == jpindt(ib) ) ztestmask(1) = tmask(ji,jj,1)
IF( mjg(jj,0) == jpjnob(ib)+1 .AND. mig(ji,0) == jpinft(ib) ) ztestmask(2) = tmask(ji,jj,1)
END DO
END DO
CALL mpp_sum( 'bdyini', ztestmask, 2 ) ! sum over the global domain
......
......@@ -211,8 +211,8 @@ CONTAINS
! sbc fields
CALL crs_dom_ope( ssh(:,:,Kmm) , 'VOL', 'T', tmask, sshn_crs , p_e12=e1e2t, p_e3=ze3t , psgn=1.0_wp )
CALL crs_dom_ope( utau , 'SUM', 'U', umask, utau_crs , p_e12=e2u , p_surf_crs=e2u_crs , psgn=1.0_wp )
CALL crs_dom_ope( vtau , 'SUM', 'V', vmask, vtau_crs , p_e12=e1v , p_surf_crs=e1v_crs , psgn=1.0_wp )
CALL crs_dom_ope( utau , 'SUM', 'T', tmask, utau_crs , p_e12=e2t , p_surf_crs=e2t_crs , psgn=1.0_wp ) !clem tau: check psgn ??
CALL crs_dom_ope( vtau , 'SUM', 'T', tmask, vtau_crs , p_e12=e1t , p_surf_crs=e1t_crs , psgn=1.0_wp ) !
CALL crs_dom_ope( wndm , 'SUM', 'T', tmask, wndm_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0_wp )
CALL crs_dom_ope( rnf , 'MAX', 'T', tmask, rnf_crs , psgn=1.0_wp )
CALL crs_dom_ope( qsr , 'SUM', 'T', tmask, qsr_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0_wp )
......
......@@ -53,7 +53,7 @@ CONTAINS
INTEGER :: dia_ar5_alloc
!!----------------------------------------------------------------------
!
ALLOCATE( thick0(jpi,jpj) , sn0(jpi,jpj,jpk), STAT=dia_ar5_alloc )
ALLOCATE( thick0(A2D(0)) , sn0(A2D(0),jpk), STAT=dia_ar5_alloc )
!
CALL mpp_sum ( 'diaar5', dia_ar5_alloc )
IF( dia_ar5_alloc /= 0 ) CALL ctl_stop( 'STOP', 'dia_ar5_alloc: failed to allocate arrays' )
......@@ -75,9 +75,10 @@ CONTAINS
REAL(wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass, zsst
REAL(wp) :: zaw, zbw, zrw
!
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zarea_ssh , zbotpres ! 2D workspace
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d, zpe ! 2D workspace
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z3d, zrhd, ztpot, zgdept ! 3D workspace (zgdept: needed to use the substitute)
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zarea_ssh ! 2D workspace
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d ! 2D workspace
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z3d
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zrhd, zgdept ! 3D workspace (zgdept: needed to use the substitute)
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn ! 4D workspace
!!--------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_ar5')
......@@ -85,10 +86,12 @@ CONTAINS
IF( kt == nit000 ) CALL dia_ar5_init
IF( l_ar5 ) THEN
ALLOCATE( zarea_ssh(jpi,jpj), zbotpres(jpi,jpj), z2d(jpi,jpj) )
ALLOCATE( zarea_ssh(A2D(0)), z2d(A2D(0)), z3d(A2D(0),jpk) )
ALLOCATE( zrhd(jpi,jpj,jpk) )
ALLOCATE( ztsn(jpi,jpj,jpk,jpts) )
zarea_ssh(:,:) = e1e2t(:,:) * ssh(:,:,Kmm)
zarea_ssh(:,:) = e1e2t(A2D(0)) * ssh(A2D(0),Kmm)
ztsn(:,:,:,:) = 0._wp
zrhd(:,:,:) = 0._wp
ENDIF
!
CALL iom_put( 'e2u' , e2u (:,:) )
......@@ -96,19 +99,19 @@ CONTAINS
CALL iom_put( 'areacello', e1e2t(:,:) )
!
IF( iom_use( 'volcello' ) .OR. iom_use( 'masscello' ) ) THEN
zrhd(:,:,jpk) = 0._wp ! ocean volume ; rhd is used as workspace
DO jk = 1, jpkm1
zrhd(:,:,jk) = e1e2t(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk)
END DO
DO jk = 1, jpk
z3d(:,:,jk) = rho0 * e3t(:,:,jk,Kmm) * tmask(:,:,jk)
END DO
CALL iom_put( 'volcello' , zrhd(:,:,:) ) ! WARNING not consistent with CMIP DR where volcello is at ca. 2000
z3d(:,:,jpk) = 0._wp ! ocean volume ; rhd is used as workspace
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z3d(ji,jj,jk) = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
END_3D
CALL iom_put( 'volcello' , z3d(:,:,:) ) ! WARNING not consistent with CMIP DR where volcello is at ca. 2000
DO_3D( 0, 0, 0, 0, 1, jpk )
z3d(ji,jj,jk) = rho0 * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
END_3D
CALL iom_put( 'masscello' , z3d (:,:,:) ) ! ocean mass
ENDIF
!
IF( iom_use( 'e3tb' ) ) THEN ! bottom layer thickness
DO_2D( 1, 1, 1, 1 )
DO_2D( 0, 0, 0, 0 )
ikb = mbkt(ji,jj)
z2d(ji,jj) = e3t(ji,jj,ikb,Kmm)
END_2D
......@@ -128,61 +131,63 @@ CONTAINS
IF( iom_use( 'botpres' ) .OR. iom_use( 'sshthster' ) .OR. iom_use( 'sshsteric' ) ) THEN
!
ztsn(:,:,:,jp_tem) = ts(:,:,:,jp_tem,Kmm) ! thermosteric ssh
ztsn(:,:,:,jp_sal) = sn0(:,:,:)
ztsn(A2D(0),:,jp_tem) = ts(A2D(0),:,jp_tem,Kmm) ! thermosteric ssh
ztsn(A2D(0),:,jp_sal) = sn0(:,:,:)
ALLOCATE( zgdept(jpi,jpj,jpk) )
DO jk = 1, jpk
zgdept(:,:,jk) = gdept(:,:,jk,Kmm)
END DO
CALL eos( ztsn, zrhd, zgdept) ! now in situ density using initial salinity
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
zgdept(ji,jj,jk) = gdept(ji,jj,jk,Kmm)
END_3D
CALL eos( ztsn, zrhd, zgdept ) ! now in situ density using initial salinity
!
zbotpres(:,:) = 0._wp ! no atmospheric surface pressure, levitating sea-ice
DO jk = 1, jpkm1
zbotpres(:,:) = zbotpres(:,:) + e3t(:,:,jk,Kmm) * zrhd(:,:,jk)
END DO
z2d(:,:) = 0._wp ! no atmospheric surface pressure, levitating sea-ice
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * zrhd(ji,jj,jk)
END_3D
IF( ln_linssh ) THEN
IF( ln_isfcav ) THEN
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
DO_2D( 0, 0, 0, 0 )
iks = mikt(ji,jj)
zbotpres(ji,jj) = zbotpres(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,iks) + riceload(ji,jj)
z2d(ji,jj) = z2d(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,iks) + riceload(ji,jj)
END_2D
ELSE
zbotpres(:,:) = zbotpres(:,:) + ssh(:,:,Kmm) * zrhd(:,:,1)
DO_2D( 0, 0, 0, 0 )
z2d(ji,jj) = z2d(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,1)
END_2D
END IF
!!gm
!!gm riceload should be added in both ln_linssh=T or F, no?
!!gm
END IF
!
zarho = glob_sum( 'diaar5', e1e2t(:,:) * zbotpres(:,:) )
zarho = glob_sum( 'diaar5', e1e2t(A2D(0)) * z2d(:,:) )
zssh_steric = - zarho / area_tot
CALL iom_put( 'sshthster', zssh_steric )
! ! steric sea surface height
zbotpres(:,:) = 0._wp ! no atmospheric surface pressure, levitating sea-ice
DO jk = 1, jpkm1
zbotpres(:,:) = zbotpres(:,:) + e3t(:,:,jk,Kmm) * rhd(:,:,jk)
END DO
z2d(:,:) = 0._wp ! no atmospheric surface pressure, levitating sea-ice
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + e3t(ji,jj,jk,Kmm) * rhd(ji,jj,jk)
END_3D
IF( ln_linssh ) THEN
IF ( ln_isfcav ) THEN
DO ji = 1,jpi
DO jj = 1,jpj
iks = mikt(ji,jj)
zbotpres(ji,jj) = zbotpres(ji,jj) + ssh(ji,jj,Kmm) * rhd(ji,jj,iks) + riceload(ji,jj)
END DO
END DO
DO_2D( 0, 0, 0, 0 )
iks = mikt(ji,jj)
z2d(ji,jj) = z2d(ji,jj) + ssh(ji,jj,Kmm) * rhd(ji,jj,iks) + riceload(ji,jj)
END_2D
ELSE
zbotpres(:,:) = zbotpres(:,:) + ssh(:,:,Kmm) * rhd(:,:,1)
DO_2D( 0, 0, 0, 0 )
z2d(ji,jj) = z2d(ji,jj) + ssh(ji,jj,Kmm) * rhd(ji,jj,1)
END_2D
END IF
END IF
!
zarho = glob_sum( 'diaar5', e1e2t(:,:) * zbotpres(:,:) )
zarho = glob_sum( 'diaar5', e1e2t(A2D(0)) * z2d(:,:) )
zssh_steric = - zarho / area_tot
CALL iom_put( 'sshsteric', zssh_steric )
! ! ocean bottom pressure
zztmp = rho0 * grav * 1.e-4_wp ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa
zbotpres(:,:) = zztmp * ( zbotpres(:,:) + ssh(:,:,Kmm) + thick0(:,:) )
CALL iom_put( 'botpres', zbotpres )
z2d(:,:) = zztmp * ( z2d(:,:) + ssh(A2D(0),Kmm) + thick0(:,:) )
CALL iom_put( 'botpres', z2d )
!
DEALLOCATE( zgdept )
!
......@@ -191,7 +196,7 @@ CONTAINS
IF( iom_use( 'masstot' ) .OR. iom_use( 'temptot' ) .OR. iom_use( 'saltot' ) ) THEN
! ! Mean density anomalie, temperature and salinity
ztsn(:,:,:,:) = 0._wp ! ztsn(:,:,1,jp_tem/sal) is used here as 2D Workspace for temperature & salinity
DO_3D( 1, 1, 1, 1, 1, jpkm1 )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zztmp = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)
ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zztmp * ts(ji,jj,jk,jp_tem,Kmm)
ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zztmp * ts(ji,jj,jk,jp_sal,Kmm)
......@@ -199,16 +204,16 @@ CONTAINS
IF( ln_linssh ) THEN
IF( ln_isfcav ) THEN
DO ji = 1, jpi
DO jj = 1, jpj
iks = mikt(ji,jj)
ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zarea_ssh(ji,jj) * ts(ji,jj,iks,jp_tem,Kmm)
ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zarea_ssh(ji,jj) * ts(ji,jj,iks,jp_sal,Kmm)
END DO
END DO
DO_2D( 0, 0, 0, 0 )
iks = mikt(ji,jj)
ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zarea_ssh(ji,jj) * ts(ji,jj,iks,jp_tem,Kmm)
ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zarea_ssh(ji,jj) * ts(ji,jj,iks,jp_sal,Kmm)
END_2D
ELSE
ztsn(:,:,1,jp_tem) = ztsn(:,:,1,jp_tem) + zarea_ssh(:,:) * ts(:,:,1,jp_tem,Kmm)
ztsn(:,:,1,jp_sal) = ztsn(:,:,1,jp_sal) + zarea_ssh(:,:) * ts(:,:,1,jp_sal,Kmm)
DO_2D( 0, 0, 0, 0 )
ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zarea_ssh(ji,jj) * ts(ji,jj,1,jp_tem,Kmm)
ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zarea_ssh(ji,jj) * ts(ji,jj,1,jp_sal,Kmm)
END_2D
END IF
ENDIF
!
......@@ -226,37 +231,35 @@ CONTAINS
IF( iom_use( 'toce_pot') .OR. iom_use( 'temptot_pot' ) .OR. iom_use( 'sst_pot' ) &
.OR. iom_use( 'ssttot' ) .OR. iom_use( 'tosmint_pot' ) ) THEN
!
ALLOCATE( ztpot(jpi,jpj,jpk) )
ztpot(:,:,jpk) = 0._wp
z3d(:,:,jpk) = 0._wp
DO jk = 1, jpkm1
ztpot(:,:,jk) = eos_pt_from_ct( ts(:,:,jk,jp_tem,Kmm), ts(:,:,jk,jp_sal,Kmm) )
z3d(:,:,jk) = eos_pt_from_ct( ts(A2D(0),jk,jp_tem,Kmm), ts(A2D(0),jk,jp_sal,Kmm) )
END DO
!
CALL iom_put( 'toce_pot', ztpot(:,:,:) ) ! potential temperature (TEOS-10 case)
CALL iom_put( 'sst_pot' , ztpot(:,:,1) ) ! surface temperature
CALL iom_put( 'toce_pot', z3d(:,:,:) ) ! potential temperature (TEOS-10 case)
CALL iom_put( 'sst_pot' , z3d(:,:,1) ) ! surface temperature
!
IF( iom_use( 'temptot_pot' ) ) THEN ! Output potential temperature in case we use TEOS-10
z2d(:,:) = 0._wp
DO jk = 1, jpkm1
z2d(:,:) = z2d(:,:) + e1e2t(:,:) * e3t(:,:,jk,Kmm) * ztpot(:,:,jk)
END DO
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * z3d(ji,jj,jk)
END_3D
ztemp = glob_sum( 'diaar5', z2d(:,:) )
CALL iom_put( 'temptot_pot', ztemp / zvol )
ENDIF
!
IF( iom_use( 'ssttot' ) ) THEN ! Output potential temperature in case we use TEOS-10
zsst = glob_sum( 'diaar5', e1e2t(:,:) * ztpot(:,:,1) )
zsst = glob_sum( 'diaar5', e1e2t(A2D(0)) * z3d(:,:,1) )
CALL iom_put( 'ssttot', zsst / area_tot )
ENDIF
! Vertical integral of temperature
IF( iom_use( 'tosmint_pot') ) THEN
z2d(:,:) = 0._wp
DO_3D( 1, 1, 1, 1, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + rho0 * e3t(ji,jj,jk,Kmm) * ztpot(ji,jj,jk)
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + rho0 * e3t(ji,jj,jk,Kmm) * z3d(ji,jj,jk)
END_3D
CALL iom_put( 'tosmint_pot', z2d )
ENDIF
DEALLOCATE( ztpot )
ENDIF
ELSE
IF( iom_use('ssttot') ) THEN ! Output sst in case we use EOS-80
......@@ -269,33 +272,31 @@ CONTAINS
! Work done against stratification by vertical mixing
! Exclude points where rn2 is negative as convection kicks in here and
! work is not being done against stratification
ALLOCATE( zpe(jpi,jpj) )
zpe(:,:) = 0._wp
z2d(:,:) = 0._wp
IF( ln_zdfddm ) THEN
DO_3D( 1, 1, 1, 1, 2, jpk )
DO_3D( 0, 0, 0, 0, 2, jpk )
IF( rn2(ji,jj,jk) > 0._wp ) THEN
zrw = ( gdept(ji,jj,jk,Kmm) - gdepw(ji,jj,jk,Kmm) ) / e3w(ji,jj,jk,Kmm)
!
zaw = rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem)* zrw
zbw = rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal)* zrw
!
zpe(ji, jj) = zpe(ji,jj) &
z2d(ji, jj) = z2d(ji,jj) &
& - grav * ( avt(ji,jj,jk) * zaw * (ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) &
& - avs(ji,jj,jk) * zbw * (ts(ji,jj,jk-1,jp_sal,Kmm) - ts(ji,jj,jk,jp_sal,Kmm) ) )
ENDIF
END_3D
ELSE
DO_3D( 1, 1, 1, 1, 1, jpk )
zpe(ji,jj) = zpe(ji,jj) + avt(ji,jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * rho0 * e3w(ji,jj,jk,Kmm)
DO_3D( 0, 0, 0, 0, 1, jpk )
z2d(ji,jj) = z2d(ji,jj) + avt(ji,jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * rho0 * e3w(ji,jj,jk,Kmm)
END_3D
ENDIF
CALL iom_put( 'tnpeo', zpe )
DEALLOCATE( zpe )
CALL iom_put( 'tnpeo', z2d )
ENDIF
IF( l_ar5 ) THEN
DEALLOCATE( zarea_ssh , zbotpres, z2d )
DEALLOCATE( ztsn )
DEALLOCATE( zarea_ssh , z2d, z3d )
DEALLOCATE( ztsn )
ENDIF
!
IF( ln_timing ) CALL timing_stop('dia_ar5')
......@@ -316,33 +317,33 @@ CONTAINS
REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in) :: pvflx ! v-flux of advection/diffusion
!
INTEGER :: ji, jj, jk
REAL(wp), DIMENSION(A2D(nn_hls)) :: z2d
REAL(wp), DIMENSION(A2D(0)) :: z2d
!!----------------------------------------------------------------------
z2d(:,:) = puflx(:,:,1)
z2d(:,:) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + puflx(ji,jj,jk)
END_3D
IF( cptr == 'adv' ) THEN
IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rho0_rcp * z2d(:,:) ) ! advective heat transport in i-direction
IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rho0 * z2d(:,:) ) ! advective salt transport in i-direction
IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rho0_rcp * z2d(:,:) ) ! advective heat transport in i-direction
IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rho0 * z2d(:,:) ) ! advective salt transport in i-direction
ELSE IF( cptr == 'ldf' ) THEN
IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rho0_rcp * z2d(:,:) ) ! diffusive heat transport in i-direction
IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rho0 * z2d(:,:) ) ! diffusive salt transport in i-direction
IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rho0_rcp * z2d(:,:) ) ! diffusive heat transport in i-direction
IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rho0 * z2d(:,:) ) ! diffusive salt transport in i-direction
ENDIF
!
z2d(:,:) = pvflx(:,:,1)
z2d(:,:) = 0._wp
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
z2d(ji,jj) = z2d(ji,jj) + pvflx(ji,jj,jk)
END_3D
IF( cptr == 'adv' ) THEN
IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rho0_rcp * z2d(:,:) ) ! advective heat transport in j-direction
IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rho0 * z2d(:,:) ) ! advective salt transport in j-direction
IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rho0_rcp * z2d(:,:) ) ! advective heat transport in j-direction
IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rho0 * z2d(:,:) ) ! advective salt transport in j-direction
ELSE IF( cptr == 'ldf' ) THEN
IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rho0_rcp * z2d(:,:) ) ! diffusive heat transport in j-direction
IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rho0 * z2d(:,:) ) ! diffusive salt transport in j-direction
IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rho0_rcp * z2d(:,:) ) ! diffusive heat transport in j-direction
IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rho0 * z2d(:,:) ) ! diffusive salt transport in j-direction
ENDIF
END SUBROUTINE dia_ar5_hst
......@@ -380,10 +381,10 @@ CONTAINS
area_tot = glob_sum( 'diaar5', e1e2t(:,:) )
ALLOCATE( zvol0(jpi,jpj) )
ALLOCATE( zvol0(A2D(0)) )
zvol0 (:,:) = 0._wp
thick0(:,:) = 0._wp
DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! interpolation of salinity at the last ocean level (i.e. the partial step)
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! interpolation of salinity at the last ocean level (i.e. the partial step)
zztmp = tmask(ji,jj,jk) * e3t_0(ji,jj,jk)
zvol0 (ji,jj) = zvol0 (ji,jj) + zztmp * e1e2t(ji,jj)
thick0(ji,jj) = thick0(ji,jj) + zztmp
......@@ -392,16 +393,16 @@ CONTAINS
DEALLOCATE( zvol0 )
IF( iom_use( 'sshthster' ) ) THEN
ALLOCATE( zsaldta(jpi,jpj,jpk,jpts) )
ALLOCATE( zsaldta(A2D(0),jpk,jpts) )
CALL iom_open ( 'sali_ref_clim_monthly', inum )
CALL iom_get ( inum, jpdom_global, 'vosaline' , zsaldta(:,:,:,1), 1 )
CALL iom_get ( inum, jpdom_global, 'vosaline' , zsaldta(:,:,:,2), 12 )
CALL iom_close( inum )
sn0(:,:,:) = 0.5_wp * ( zsaldta(:,:,:,1) + zsaldta(:,:,:,2) )
sn0(:,:,:) = sn0(:,:,:) * tmask(:,:,:)
sn0(:,:,:) = sn0(:,:,:) * tmask(A2D(0),:)
IF( ln_zps ) THEN ! z-coord. partial steps
DO_2D( 1, 1, 1, 1 ) ! interpolation of salinity at the last ocean level (i.e. the partial step)
DO_2D( 0, 0, 0, 0 ) ! interpolation of salinity at the last ocean level (i.e. the partial step)
ik = mbkt(ji,jj)
IF( ik > 1 ) THEN
zztmp = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
......
......@@ -51,19 +51,19 @@ CONTAINS
INTEGER, INTENT(in) :: kt ! ocean time-step index
INTEGER, INTENT(in) :: Kmm ! ocean time level index
!
INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zCu_max, zCv_max, zCw_max ! local scalars
INTEGER , DIMENSION(3) :: iloc_u , iloc_v , iloc_w , iloc ! workspace
REAL(wp), DIMENSION(jpi,jpj,jpk) :: zCu_cfl, zCv_cfl, zCw_cfl ! workspace
LOGICAL , DIMENSION(jpi,jpj,jpk) :: llmsk
INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zCu_max, zCv_max, zCw_max ! local scalars
INTEGER , DIMENSION(3) :: iloc_u , iloc_v , iloc_w , iloc ! workspace
REAL(wp), DIMENSION(A2D(0),jpk) :: zCu_cfl, zCv_cfl, zCw_cfl ! workspace
LOGICAL , DIMENSION(A2D(0),jpk) :: llmsk
!!----------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('dia_cfl')
!
llmsk( 1:nn_hls,:,:) = .FALSE. ! exclude halos from the checked region
llmsk(Nie0+1: jpi,:,:) = .FALSE.
llmsk(:, 1:nn_hls,:) = .FALSE.
llmsk(:,Nje0+1: jpj,:) = .FALSE.
!llmsk( 1:nn_hls,:,:) = .FALSE. ! exclude halos from the checked region
!llmsk(Nie0+1: jpi,:,:) = .FALSE.
!llmsk(:, 1:nn_hls,:) = .FALSE.
!llmsk(:,Nje0+1: jpj,:) = .FALSE.
!
DO_3D( 0, 0, 0, 0, 1, jpk ) ! calculate Courant numbers
zCu_cfl(ji,jj,jk) = ABS( uu(ji,jj,jk,Kmm) ) * rDt / e1u (ji,jj) ! for i-direction
......
......@@ -414,9 +414,9 @@ CONTAINS
!verify if the point is on the local domain:(1,Nie0)*(1,Nje0)
IF( iiloc >= 1 .AND. iiloc <= Nie0 .AND. &
ijloc >= 1 .AND. ijloc <= Nje0 )THEN
iptloc = iptloc + 1 ! count local points
secs(jsec)%listPoint(iptloc) = POINT_SECTION(mi0(iiglo),mj0(ijglo)) ! store local coordinates
secs(jsec)%direction(iptloc) = directemp(jpt) ! store local direction
iptloc = iptloc + 1 ! count local points
secs(jsec)%listPoint(iptloc) = POINT_SECTION(mi0(iiglo,nn_hls),mj0(ijglo,nn_hls)) ! store local coordinates
secs(jsec)%direction(iptloc) = directemp(jpt) ! store local direction
ENDIF
!
END DO
......
......@@ -5,11 +5,11 @@ MODULE diadetide
!!======================================================================
!! History : ! 2019 (S. Mueller)
!!----------------------------------------------------------------------
USE par_oce , ONLY : wp, jpi, jpj
USE in_out_manager , ONLY : lwp, numout
USE iom , ONLY : iom_put
USE dom_oce , ONLY : rn_Dt, nsec_day
USE phycst , ONLY : rpi
USE par_oce
USE in_out_manager
USE iom
USE dom_oce
USE phycst
USE tide_mod
#if defined key_xios
USE xios
......@@ -24,6 +24,8 @@ MODULE diadetide
PUBLIC :: dia_detide_init, dia_detide
!! * Substitutions
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2019)
!! $Id$
......@@ -90,9 +92,9 @@ CONTAINS
!!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt
REAL(wp), DIMENSION(jpi,jpj) :: zwght_2D
REAL(wp), DIMENSION(A2D(0)) :: zwght_2D
REAL(wp) :: zwght, ztmp
INTEGER :: jn
INTEGER :: ji, jj, jn
! Compute detiding weight at the current time-step; the daily total weight
! is one, and the daily summation of a diagnosed field multiplied by this
......@@ -104,7 +106,10 @@ CONTAINS
zwght = zwght + 1.0_wp / REAL( ndiadetide, KIND=wp )
END IF
END DO
zwght_2D(:,:) = zwght
DO_2D( 0, 0, 0, 0 )
zwght_2D(ji,jj) = zwght
END_2D
CALL iom_put( "diadetide_weight", zwght_2D)
END SUBROUTINE dia_detide
......
......@@ -50,6 +50,7 @@ MODULE diahsb
REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: tmask_ini
!! * Substitutions
# include "do_loop_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
......@@ -82,42 +83,73 @@ CONTAINS
REAL(wp) :: z_frc_trd_v ! - -
REAL(wp) :: z_wn_trd_t , z_wn_trd_s ! - -
REAL(wp) :: z_ssh_hc , z_ssh_sc ! - -
REAL(wp), DIMENSION(jpi,jpj,13) :: ztmp
REAL(wp), DIMENSION(jpi,jpj,jpkm1,4) :: ztmpk
REAL(wp), DIMENSION(17) :: zbg
REAL(wp), DIMENSION(A2D(0),13) :: ztmp
REAL(wp), DIMENSION(A2D(0),jpkm1,4) :: ztmpk
REAL(wp), DIMENSION(17) :: zbg
!!---------------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_hsb')
!
ztmp (:,:,:) = 0._wp ! should be better coded
ztmpk(:,:,:,:) = 0._wp ! should be better coded
!
ts(:,:,:,1,Kmm) = ts(:,:,:,1,Kmm) * tmask(:,:,:) ; ts(:,:,:,1,Kbb) = ts(:,:,:,1,Kbb) * tmask(:,:,:) ;
ts(:,:,:,2,Kmm) = ts(:,:,:,2,Kmm) * tmask(:,:,:) ; ts(:,:,:,2,Kbb) = ts(:,:,:,2,Kbb) * tmask(:,:,:) ;
DO_2D( 0, 0, 0, 0 )
ztmp (ji,jj,:) = 0._wp ! should be better coded
ztmpk(ji,jj,:,:) = 0._wp ! should be better coded
!
ts(ji,jj,:,1,Kmm) = ts(ji,jj,:,1,Kmm) * tmask(ji,jj,:)
ts(ji,jj,:,1,Kbb) = ts(ji,jj,:,1,Kbb) * tmask(ji,jj,:)
!
ts(ji,jj,:,2,Kmm) = ts(ji,jj,:,2,Kmm) * tmask(ji,jj,:)
ts(ji,jj,:,2,Kbb) = ts(ji,jj,:,2,Kbb) * tmask(ji,jj,:)
END_2D
!
! ------------------------- !
! 1 - Trends due to forcing !
! ------------------------- !
! prepare trends
ztmp(:,:,1) = - r1_rho0 * ( emp(:,:) - rnf(:,:) - fwfisf_cav(:,:) - fwfisf_par(:,:) ) * surf(:,:) ! volume
ztmp(:,:,2) = sbc_tsc(:,:,jp_tem) * surf(:,:) ! heat
ztmp(:,:,3) = sbc_tsc(:,:,jp_sal) * surf(:,:) ! salt
IF( ln_rnf ) ztmp(:,:,4) = rnf_tsc(:,:,jp_tem) * surf(:,:) ! runoff temp
IF( ln_rnf_sal ) ztmp(:,:,5) = rnf_tsc(:,:,jp_sal) * surf(:,:) ! runoff salt
IF( ln_isf ) ztmp(:,:,6) = ( risf_cav_tsc(:,:,jp_tem) + risf_par_tsc(:,:,jp_tem) ) * surf(:,:) ! isf temp
IF( ln_traqsr ) ztmp(:,:,7) = r1_rho0_rcp * qsr(:,:) * surf(:,:) ! penetrative solar radiation
IF( ln_trabbc ) ztmp(:,:,8) = qgh_trd0(:,:) * surf(:,:) ! geothermal heat
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,1) = - r1_rho0 * ( emp(ji,jj) & ! volume
& - rnf(ji,jj) &
& - fwfisf_cav(ji,jj) &
& - fwfisf_par(ji,jj) ) * surf(ji,jj)
ztmp(ji,jj,2) = sbc_tsc(ji,jj,jp_tem) * surf(ji,jj) ! heat
ztmp(ji,jj,3) = sbc_tsc(ji,jj,jp_sal) * surf(ji,jj) ! salt
END_2D
IF( ln_rnf ) THEN
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,4) = rnf_tsc(ji,jj,jp_tem) * surf(ji,jj) ! runoff temp
END_2D
END IF
IF( ln_rnf_sal ) THEN
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,5) = rnf_tsc(ji,jj,jp_sal) * surf(ji,jj) ! runoff salt
END_2D
END IF
IF( ln_isf ) THEN
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,6) = ( risf_cav_tsc(ji,jj,jp_tem) &
& + risf_par_tsc(ji,jj,jp_tem) ) * surf(ji,jj) ! isf temp
END_2D
END IF
IF( ln_traqsr ) THEN
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,7) = r1_rho0_rcp * qsr(ji,jj) * surf(ji,jj) ! penetrative solar radiation
END_2D
END IF
IF( ln_trabbc ) THEN
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,8) = qgh_trd0(ji,jj) * surf(ji,jj) ! geothermal heat
END_2D
END IF
!
IF( ln_linssh ) THEN ! Advection flux through fixed surface (z=0)
IF( ln_isfcav ) THEN
DO ji=1,jpi
DO jj=1,jpj
ztmp(ji,jj,9 ) = - surf(ji,jj) * ww(ji,jj,mikt(ji,jj)) * ts(ji,jj,mikt(ji,jj),jp_tem,Kbb)
ztmp(ji,jj,10) = - surf(ji,jj) * ww(ji,jj,mikt(ji,jj)) * ts(ji,jj,mikt(ji,jj),jp_sal,Kbb)
END DO
END DO
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,9 ) = - surf(ji,jj) * ww(ji,jj,mikt(ji,jj)) * ts(ji,jj,mikt(ji,jj),jp_tem,Kbb)
ztmp(ji,jj,10) = - surf(ji,jj) * ww(ji,jj,mikt(ji,jj)) * ts(ji,jj,mikt(ji,jj),jp_sal,Kbb)
END_2D
ELSE
ztmp(:,:,9 ) = - surf(:,:) * ww(:,:,1) * ts(:,:,1,jp_tem,Kbb)
ztmp(:,:,10) = - surf(:,:) * ww(:,:,1) * ts(:,:,1,jp_sal,Kbb)
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,9 ) = - surf(ji,jj) * ww(ji,jj,1) * ts(ji,jj,1,jp_tem,Kbb)
ztmp(ji,jj,10) = - surf(ji,jj) * ww(ji,jj,1) * ts(ji,jj,1,jp_sal,Kbb)
END_2D
END IF
ENDIF
......@@ -152,20 +184,22 @@ CONTAINS
! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl)
!
! ! volume variation (calculated with ssh)
ztmp(:,:,11) = surf(:,:)*ssh(:,:,Kmm) - surf_ini(:,:)*ssh_ini(:,:)
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,11) = surf(ji,jj)*ssh(ji,jj,Kmm) - surf_ini(ji,jj)*ssh_ini(ji,jj)
END_2D
! ! heat & salt content variation (associated with ssh)
IF( ln_linssh ) THEN ! linear free surface case
IF( ln_isfcav ) THEN ! ISF case
DO ji = 1, jpi
DO jj = 1, jpj
ztmp(ji,jj,12) = surf(ji,jj) * ( ts(ji,jj,mikt(ji,jj),jp_tem,Kmm) * ssh(ji,jj,Kmm) - ssh_hc_loc_ini(ji,jj) )
ztmp(ji,jj,13) = surf(ji,jj) * ( ts(ji,jj,mikt(ji,jj),jp_sal,Kmm) * ssh(ji,jj,Kmm) - ssh_sc_loc_ini(ji,jj) )
END DO
END DO
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,12) = surf(ji,jj) * ( ts(ji,jj,mikt(ji,jj),jp_tem,Kmm) * ssh(ji,jj,Kmm) - ssh_hc_loc_ini(ji,jj) )
ztmp(ji,jj,13) = surf(ji,jj) * ( ts(ji,jj,mikt(ji,jj),jp_sal,Kmm) * ssh(ji,jj,Kmm) - ssh_sc_loc_ini(ji,jj) )
END_2D
ELSE ! no under ice-shelf seas
ztmp(:,:,12) = surf(:,:) * ( ts(:,:,1,jp_tem,Kmm) * ssh(:,:,Kmm) - ssh_hc_loc_ini(:,:) )
ztmp(:,:,13) = surf(:,:) * ( ts(:,:,1,jp_sal,Kmm) * ssh(:,:,Kmm) - ssh_sc_loc_ini(:,:) )
DO_2D( 0, 0, 0, 0 )
ztmp(ji,jj,12) = surf(ji,jj) * ( ts(ji,jj,1,jp_tem,Kmm) * ssh(ji,jj,Kmm) - ssh_hc_loc_ini(ji,jj) )
ztmp(ji,jj,13) = surf(ji,jj) * ( ts(ji,jj,1,jp_sal,Kmm) * ssh(ji,jj,Kmm) - ssh_sc_loc_ini(ji,jj) )
END_2D
END IF
ENDIF
......@@ -185,19 +219,27 @@ CONTAINS
! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl)
!
DO jk = 1, jpkm1 ! volume
ztmpk(:,:,jk,1) = surf (:,:) * e3t(:,:,jk,Kmm)*tmask(:,:,jk) &
& - surf_ini(:,:) * e3t_ini(:,:,jk )*tmask_ini(:,:,jk)
DO_2D( 0, 0, 0, 0 )
ztmpk(ji,jj,jk,1) = surf (ji,jj) * e3t(ji,jj,jk,Kmm)*tmask(ji,jj,jk) &
& - surf_ini(ji,jj) * e3t_ini(ji,jj,jk )*tmask_ini(ji,jj,jk)
END_2D
END DO
DO jk = 1, jpkm1 ! heat
ztmpk(:,:,jk,2) = ( surf (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_tem,Kmm) &
& - surf_ini(:,:) * hc_loc_ini(:,:,jk) )
DO_2D( 0, 0, 0, 0 )
ztmpk(ji,jj,jk,2) = ( surf (ji,jj) * e3t(ji,jj,jk,Kmm)*ts(ji,jj,jk,jp_tem,Kmm) &
& - surf_ini(ji,jj) * hc_loc_ini(ji,jj,jk) )
END_2D
END DO
DO jk = 1, jpkm1 ! salt
ztmpk(:,:,jk,3) = ( surf (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_sal,Kmm) &
& - surf_ini(:,:) * sc_loc_ini(:,:,jk) )
DO_2D( 0, 0, 0, 0 )
ztmpk(ji,jj,jk,3) = ( surf (ji,jj) * e3t(ji,jj,jk,Kmm)*ts(ji,jj,jk,jp_sal,Kmm) &
& - surf_ini(ji,jj) * sc_loc_ini(ji,jj,jk) )
END_2D
END DO
DO jk = 1, jpkm1 ! total ocean volume
ztmpk(:,:,jk,4) = surf(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk)
DO_2D( 0, 0, 0, 0 )
ztmpk(ji,jj,jk,4) = surf(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
END_2D
END DO
! global sum
......@@ -315,32 +357,34 @@ CONTAINS
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' dia_hsb_rst : initialise hsb at initial state '
IF(lwp) WRITE(numout,*)
surf_ini(:,:) = e1e2t(:,:) * tmask_i(:,:) ! initial ocean surface
ssh_ini(:,:) = ssh(:,:,Kmm) ! initial ssh
DO jk = 1, jpk
! if ice sheet/oceqn coupling, need to mask ini variables here (mask could change at the next NEMO instance).
e3t_ini (:,:,jk) = e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial vertical scale factors
tmask_ini (:,:,jk) = tmask(:,:,jk) ! initial mask
hc_loc_ini(:,:,jk) = ts(:,:,jk,jp_tem,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial heat content
sc_loc_ini(:,:,jk) = ts(:,:,jk,jp_sal,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial salt content
END DO
DO_2D( 0, 0, 0, 0 )
surf_ini(ji,jj) = e1e2t(ji,jj) * tmask_i(ji,jj) ! initial ocean surface
ssh_ini(ji,jj) = ssh(ji,jj,Kmm) ! initial ssh
END_2D
! if ice sheet/oceqn coupling, need to mask ini variables here (mask could change at the next NEMO instance).
DO_3D( 0, 0, 0, 0, 1, jpk )
e3t_ini (ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial vertical scale factors
tmask_ini (ji,jj,jk) = tmask(ji,jj,jk) ! initial mask
hc_loc_ini(ji,jj,jk) = ts(ji,jj,jk,jp_tem,Kmm) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial heat content
sc_loc_ini(ji,jj,jk) = ts(ji,jj,jk,jp_sal,Kmm) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial salt content
END_3D
frc_v = 0._wp ! volume trend due to forcing
frc_t = 0._wp ! heat content - - - -
frc_s = 0._wp ! salt content - - - -
IF( ln_linssh ) THEN
IF( ln_isfcav ) THEN
DO ji = 1, jpi
DO jj = 1, jpj
ssh_hc_loc_ini(ji,jj) = ts(ji,jj,mikt(ji,jj),jp_tem,Kmm) * ssh(ji,jj,Kmm) ! initial heat content in ssh
ssh_sc_loc_ini(ji,jj) = ts(ji,jj,mikt(ji,jj),jp_sal,Kmm) * ssh(ji,jj,Kmm) ! initial salt content in ssh
END DO
END DO
ELSE
ssh_hc_loc_ini(:,:) = ts(:,:,1,jp_tem,Kmm) * ssh(:,:,Kmm) ! initial heat content in ssh
ssh_sc_loc_ini(:,:) = ts(:,:,1,jp_sal,Kmm) * ssh(:,:,Kmm) ! initial salt content in ssh
DO_2D( 0, 0, 0, 0 )
ssh_hc_loc_ini(ji,jj) = ts(ji,jj,mikt(ji,jj),jp_tem,Kmm) * ssh(ji,jj,Kmm) ! initial heat content in ssh
ssh_sc_loc_ini(ji,jj) = ts(ji,jj,mikt(ji,jj),jp_sal,Kmm) * ssh(ji,jj,Kmm) ! initial salt content in ssh
END_2D
ELSE
DO_2D( 0, 0, 0, 0 )
ssh_hc_loc_ini(ji,jj) = ts(ji,jj,1,jp_tem,Kmm) * ssh(ji,jj,Kmm) ! initial heat content in ssh
ssh_sc_loc_ini(ji,jj) = ts(ji,jj,1,jp_sal,Kmm) * ssh(ji,jj,Kmm) ! initial salt content in ssh
END_2D
END IF
frc_wn_t = 0._wp ! initial heat content misfit due to free surface
frc_wn_s = 0._wp ! initial salt content misfit due to free surface
frc_wn_t = 0._wp ! initial heat content misfit due to free surface
frc_wn_s = 0._wp ! initial salt content misfit due to free surface
ENDIF
ENDIF
!
......@@ -388,6 +432,7 @@ CONTAINS
INTEGER, INTENT(in) :: Kmm ! time level index
!
INTEGER :: ierror, ios ! local integer
INTEGER :: ji, jj ! loop index
!!
NAMELIST/namhsb/ ln_diahsb
!!----------------------------------------------------------------------
......@@ -413,13 +458,13 @@ CONTAINS
! ------------------- !
! 1 - Allocate memory !
! ------------------- !
ALLOCATE( hc_loc_ini(jpi,jpj,jpk), sc_loc_ini(jpi,jpj,jpk), surf_ini(jpi,jpj), &
& e3t_ini(jpi,jpj,jpk), surf(jpi,jpj), ssh_ini(jpi,jpj), tmask_ini(jpi,jpj,jpk),STAT=ierror )
ALLOCATE( hc_loc_ini(A2D(0),jpk), sc_loc_ini(A2D(0),jpk), surf_ini(A2D(0)), &
& e3t_ini(A2D(0),jpk), surf(A2D(0)), ssh_ini(A2D(0)), tmask_ini(A2D(0),jpk),STAT=ierror )
IF( ierror > 0 ) THEN
CALL ctl_stop( 'dia_hsb_init: unable to allocate hc_loc_ini' ) ; RETURN
ENDIF
IF( ln_linssh ) ALLOCATE( ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj),STAT=ierror )
IF( ln_linssh ) ALLOCATE( ssh_hc_loc_ini(A2D(0)), ssh_sc_loc_ini(A2D(0)),STAT=ierror )
IF( ierror > 0 ) THEN
CALL ctl_stop( 'dia_hsb: unable to allocate ssh_hc_loc_ini' ) ; RETURN
ENDIF
......@@ -427,7 +472,10 @@ CONTAINS
! ----------------------------------------------- !
! 2 - Time independant variables and file opening !
! ----------------------------------------------- !
surf(:,:) = e1e2t(:,:) * tmask_i(:,:) ! masked surface grid cell area
DO_2D( 0, 0, 0, 0 )
surf(ji,jj) = e1e2t(ji,jj) * smask0_i(ji,jj) ! masked surface grid cell area
END_2D
surf_tot = glob_sum( 'diahsb', surf(:,:) ) ! total ocean surface area
IF( ln_bdy ) CALL ctl_warn( 'dia_hsb_init: heat/salt budget does not consider open boundary fluxes' )
......
......@@ -86,22 +86,22 @@ CONTAINS
INTEGER, INTENT( in ) :: kt ! ocean time-step index
INTEGER, INTENT( in ) :: Kmm ! ocean time level index
!!
INTEGER :: ji, jj, jk ! dummy loop arguments
REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth
REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth
REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth
REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop
REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace
REAL(wp), DIMENSION(jpi,jpj) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2
REAL(wp), DIMENSION(jpi,jpj) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2
REAL(wp), DIMENSION(jpi,jpj) :: zrho10_3 ! MLD: rho = rho10m + zrho3
REAL(wp), DIMENSION(jpi,jpj) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC)
REAL(wp), DIMENSION(jpi,jpj) :: ztinv ! max of temperature inversion
REAL(wp), DIMENSION(jpi,jpj) :: zdepinv ! depth of temperature inversion
REAL(wp), DIMENSION(jpi,jpj) :: zrho0_3 ! MLD rho = rho(surf) = 0.03
REAL(wp), DIMENSION(jpi,jpj) :: zrho0_1 ! MLD rho = rho(surf) = 0.01
REAL(wp), DIMENSION(jpi,jpj) :: zmaxdzT ! max of dT/dz
REAL(wp), DIMENSION(jpi,jpj) :: zdelr ! delta rho equivalent to deltaT = 0.2
INTEGER :: ji, jj, jk ! dummy loop arguments
REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth
REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth
REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth
REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop
REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace
REAL(wp), DIMENSION(A2D(0)) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2
REAL(wp), DIMENSION(A2D(0)) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2
REAL(wp), DIMENSION(A2D(0)) :: zrho10_3 ! MLD: rho = rho10m + zrho3
REAL(wp), DIMENSION(A2D(0)) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC)
REAL(wp), DIMENSION(A2D(0)) :: ztinv ! max of temperature inversion
REAL(wp), DIMENSION(A2D(0)) :: zdepinv ! depth of temperature inversion
REAL(wp), DIMENSION(A2D(0)) :: zrho0_3 ! MLD rho = rho(surf) = 0.03
REAL(wp), DIMENSION(A2D(0)) :: zrho0_1 ! MLD rho = rho(surf) = 0.01
REAL(wp), DIMENSION(A2D(0)) :: zmaxdzT ! max of dT/dz
REAL(wp), DIMENSION(A2D(0)) :: zdelr ! delta rho equivalent to deltaT = 0.2
!!----------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_hth')
......@@ -131,7 +131,7 @@ CONTAINS
IF( iom_use( 'mlddzt' ) ) zmaxdzT(:,:) = 0._wp
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) &
& .OR. iom_use( 'mldr10_3' ) .OR. iom_use( 'pycndep' ) ) THEN
DO_2D( 1, 1, 1, 1 )
DO_2D( 0, 0, 0, 0 )
zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
hth (ji,jj) = zztmp
zabs2 (ji,jj) = zztmp
......@@ -142,7 +142,7 @@ CONTAINS
ENDIF
IF( iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN
IF( nla10 > 1 ) THEN
DO_2D( 1, 1, 1, 1 )
DO_2D( 0, 0, 0, 0 )
zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
zrho0_3(ji,jj) = zztmp
zrho0_1(ji,jj) = zztmp
......@@ -157,7 +157,7 @@ CONTAINS
! MLD: rho = rho(1) + zrho3 !
! MLD: rho = rho(1) + zrho1 !
! ------------------------------------------------------------- !
DO_3DS( 1, 1, 1, 1, jpkm1, 2, -1 ) ! loop from bottom to 2
DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! loop from bottom to 2
!
zzdep = gdepw(ji,jj,jk,Kmm)
zztmp = ( ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) &
......@@ -189,7 +189,7 @@ CONTAINS
!
! Preliminary computation
! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC)
DO_2D( 1, 1, 1, 1 )
DO_2D( 0, 0, 0, 0 )
IF( tmask(ji,jj,nla10) == 1. ) THEN
zu = 1779.50 + 11.250 * ts(ji,jj,nla10,jp_tem,Kmm) - 3.80 * ts(ji,jj,nla10,jp_sal,Kmm) &
& - 0.0745 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm) &
......@@ -213,7 +213,7 @@ CONTAINS
! temperature inversion: max( 0, max of tn - tn(10m) ) !
! depth of temperature inversion !
! ------------------------------------------------------------- !
DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! loop from bottom to nlb10
DO_3DS( 0, 0, 0, 0, jpkm1, nlb10, -1 ) ! loop from bottom to nlb10
!
zzdep = gdepw(ji,jj,jk,Kmm) * tmask(ji,jj,1)
!
......@@ -305,13 +305,16 @@ CONTAINS
!
INTEGER :: ji, jj, jk, iid
REAL(wp) :: zztmp, zzdep
INTEGER, DIMENSION(jpi,jpj) :: iktem
INTEGER, DIMENSION(A2D(0)) :: iktem
! --------------------------------------- !
! search deepest level above ptem !
! --------------------------------------- !
iktem(:,:) = 1
DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom
DO_2D( 0, 0, 0, 0 )
iktem(ji,jj) = 1
END_2D
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom
zztmp = ts(ji,jj,jk,jp_tem,Kmm)
IF( zztmp >= ptem ) iktem(ji,jj) = jk
END_3D
......@@ -319,7 +322,7 @@ CONTAINS
! ------------------------------- !
! Depth of ptem isotherm !
! ------------------------------- !
DO_2D( 1, 1, 1, 1 )
DO_2D( 0, 0, 0, 0 )
!
zzdep = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) ! depth of the ocean bottom
!
......@@ -346,18 +349,29 @@ CONTAINS
REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc
!
INTEGER :: ji, jj, jk, ik
REAL(wp), DIMENSION(jpi,jpj) :: zthick
INTEGER , DIMENSION(jpi,jpj) :: ilevel
REAL(wp), DIMENSION(A2D(0)) :: zthick
INTEGER , DIMENSION(A2D(0)) :: ilevel
! surface boundary condition
IF( .NOT. ln_linssh ) THEN ; zthick(:,:) = 0._wp ; phtc(:,:) = 0._wp
ELSE ; zthick(:,:) = ssh(:,:,Kmm) ; phtc(:,:) = pt(:,:,1) * ssh(:,:,Kmm) * tmask(:,:,1)
IF( .NOT. ln_linssh ) THEN
DO_2D( 0, 0, 0, 0 )
zthick(ji,jj) = 0._wp
phtc (ji,jj) = 0._wp
END_2D
ELSE
DO_2D( 0, 0, 0, 0 )
zthick(ji,jj) = ssh(ji,jj,Kmm)
phtc (ji,jj) = pt(ji,jj,1) * ssh(ji,jj,Kmm) * tmask(ji,jj,1)
END_2D
ENDIF
!
ilevel(:,:) = 1
DO_3D( 1, 1, 1, 1, 1, jpkm1 )
DO_2D( 0, 0, 0, 0 )
ilevel(ji,jj) = 1
END_2D
!
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
IF( ( gdepw(ji,jj,jk+1,Kmm) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN
ilevel(ji,jj) = jk+1
zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm)
......@@ -365,7 +379,7 @@ CONTAINS
ENDIF
END_3D
!
DO_2D( 1, 1, 1, 1 )
DO_2D( 0, 0, 0, 0 )
ik = ilevel(ji,jj)
IF( tmask(ji,jj,ik) == 1 ) THEN
zthick(ji,jj) = MIN ( gdepw(ji,jj,ik+1,Kmm), pdep ) - zthick(ji,jj) ! remaining thickness to reach dephw pdep
......
......@@ -6,7 +6,7 @@ MODULE diamlr
!! History : 4.0 ! 2019 (S. Mueller) Original code
!!----------------------------------------------------------------------
USE par_oce , ONLY : wp, jpi, jpj
USE par_oce , ONLY : wp, jpi, jpj, ntsi, ntei, ntsj, ntej
USE phycst , ONLY : rpi
USE dom_oce , ONLY : adatrj
USE tide_mod
......@@ -407,8 +407,9 @@ CONTAINS
!! ** Purpose : update time used in multiple-linear-regression analysis
!!
!!----------------------------------------------------------------------
REAL(wp), DIMENSION(jpi,jpj) :: zadatrj2d
REAL(wp), DIMENSION(A2D(0)) :: zadatrj2d
!!----------------------------------------------------------------------
INTEGER :: ji, jj
IF( ln_timing ) CALL timing_start('dia_mlr')
......@@ -417,7 +418,9 @@ CONTAINS
!
! A 2-dimensional field of constant value is sent, and subsequently used directly
! or transformed to a scalar or a constant 3-dimensional field as required.
zadatrj2d(:,:) = adatrj*86400.0_wp
DO_2D( 0, 0, 0, 0 )
zadatrj2d(ji,jj) = adatrj*86400.0_wp
END_2D
IF ( iom_use('diamlr_time') ) CALL iom_put('diamlr_time', zadatrj2d)
!
IF( ln_timing ) CALL timing_stop('dia_mlr')
......
......@@ -78,7 +78,7 @@ CONTAINS
!!----------------------------------------------------------------------
INTEGER , INTENT(in) :: kt ! ocean time-step index
INTEGER , INTENT(in) :: Kmm ! time level index
REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport
REAL(wp), DIMENSION(A2D(0),jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport
!!----------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('dia_ptr')
......@@ -110,13 +110,13 @@ CONTAINS
!!----------------------------------------------------------------------
!! ** Purpose : Calculate diagnostics and send to XIOS
!!----------------------------------------------------------------------
INTEGER , INTENT(in) :: kt ! ocean time-step index
INTEGER , INTENT(in) :: Kmm ! time level index
REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport
INTEGER , INTENT(in) :: kt ! ocean time-step index
INTEGER , INTENT(in) :: Kmm ! time level index
REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL :: pvtr ! j-effective transport (used only by PRESENT)
!
INTEGER :: ji, jj, jk, jn ! dummy loop indices
REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace
REAL(wp), DIMENSION(jpj) :: zvsum, ztsum, zssum ! 1D workspace
REAL(wp), DIMENSION(Nis0:Nie0,Njs0:Nje0) :: z2d ! 2D workspace
REAL(wp), DIMENSION( Njs0:Nje0) :: zvsum, ztsum, zssum ! 1D workspace
!
!overturning calculation
REAL(wp), DIMENSION(:,:,: ), ALLOCATABLE :: sjk, r1_sjk, v_msf ! i-mean i-k-surface and its inverse
......@@ -126,19 +126,19 @@ CONTAINS
REAL(wp), DIMENSION(:,:,: ), ALLOCATABLE :: z3dtr
!!----------------------------------------------------------------------
!
ALLOCATE( z3dtr(jpi,jpj,nbasin) )
ALLOCATE( z3dtr(Nis0:Nie0,Njs0:Nje0,nbasin) )
IF( PRESENT( pvtr ) ) THEN
IF( iom_use( 'zomsf' ) ) THEN ! effective MSF
ALLOCATE( z4d1(jpi,jpj,jpk,nbasin) )
ALLOCATE( z4d1(Nis0:Nie0,Njs0:Nje0,jpk,nbasin) )
!
DO jn = 1, nbasin ! by sub-basins
z4d1(1,:,:,jn) = pvtr_int(:,:,jp_vtr,jn) ! zonal cumulative effective transport excluding closed seas
z4d1(Nis0,:,:,jn) = pvtr_int(:,:,jp_vtr,jn) ! zonal cumulative effective transport excluding closed seas
DO jk = jpkm1, 1, -1
z4d1(1,:,jk,jn) = z4d1(1,:,jk+1,jn) - z4d1(1,:,jk,jn) ! effective j-Stream-Function (MSF)
z4d1(Nis0,:,jk,jn) = z4d1(Nis0,:,jk+1,jn) - z4d1(Nis0,:,jk,jn) ! effective j-Stream-Function (MSF)
END DO
DO ji = 2, jpi
z4d1(ji,:,:,jn) = z4d1(1,:,:,jn)
DO ji = Nis0+1, Nie0
z4d1(ji,:,:,jn) = z4d1(Nis0,:,:,jn)
ENDDO
END DO
CALL iom_put( 'zomsf', z4d1 * rc_sv )
......@@ -146,8 +146,8 @@ CONTAINS
DEALLOCATE( z4d1 )
ENDIF
IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
ALLOCATE( sjk(jpj,jpk,nbasin), r1_sjk(jpj,jpk,nbasin), v_msf(jpj,jpk,nbasin), &
& zt_jk(jpj,jpk,nbasin), zs_jk(jpj,jpk,nbasin) )
ALLOCATE( sjk( Njs0:Nje0,jpk,nbasin), r1_sjk(Njs0:Nje0,jpk,nbasin), v_msf(Njs0:Nje0,jpk,nbasin), &
& zt_jk(Njs0:Nje0,jpk,nbasin), zs_jk( Njs0:Nje0,jpk,nbasin) )
!
DO jn = 1, nbasin
sjk(:,:,jn) = pvtr_int(:,:,jp_msk,jn)
......@@ -162,16 +162,16 @@ CONTAINS
!
ENDDO
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sophtove', z3dtr )
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sopstove', z3dtr )
......@@ -181,7 +181,7 @@ CONTAINS
IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
! Calculate barotropic heat and salt transport here
ALLOCATE( sjk(jpj,1,nbasin), r1_sjk(jpj,1,nbasin) )
ALLOCATE( sjk(A1Dj(0),1,nbasin), r1_sjk(A1Dj(0),1,nbasin) )
!
DO jn = 1, nbasin
sjk(:,1,jn) = SUM( pvtr_int(:,:,jp_msk,jn), 2 )
......@@ -196,16 +196,16 @@ CONTAINS
!
ENDDO
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sophtbtr', z3dtr )
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sopstbtr', z3dtr )
......@@ -218,28 +218,28 @@ CONTAINS
pvtr_int(:,:,:,:) = 0._wp
ELSE
IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' ) ) THEN ! i-mean i-k-surface
ALLOCATE( z4d1(jpi,jpj,jpk,nbasin), z4d2(jpi,jpj,jpk,nbasin) )
ALLOCATE( z4d1(Nis0:Nie0,Njs0:Nje0,jpk,nbasin), z4d2(Nis0:Nie0,Njs0:Nje0,jpk,nbasin) )
!
DO jn = 1, nbasin
z4d1(1,:,:,jn) = pzon_int(:,:,jp_msk,jn)
DO ji = 2, jpi
z4d1(ji,:,:,jn) = z4d1(1,:,:,jn)
z4d1(Nis0,:,:,jn) = pzon_int(:,:,jp_msk,jn)
DO ji = Nis0+1, Nie0
z4d1(ji,:,:,jn) = z4d1(Nis0,:,:,jn)
ENDDO
ENDDO
CALL iom_put( 'zosrf', z4d1 )
!
DO jn = 1, nbasin
z4d2(1,:,:,jn) = pzon_int(:,:,jp_tem,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 )
DO ji = 2, jpi
z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
z4d2(Nis0,:,:,jn) = pzon_int(:,:,jp_tem,jn) / MAX( z4d1(Nis0,:,:,jn), 10.e-15 )
DO ji = Nis0+1, Nie0
z4d2(ji,:,:,jn) = z4d2(Nis0,:,:,jn)
ENDDO
ENDDO
CALL iom_put( 'zotem', z4d2 )
!
DO jn = 1, nbasin
z4d2(1,:,:,jn) = pzon_int(:,:,jp_sal,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 )
DO ji = 2, jpi
z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
z4d2(Nis0,:,:,jn) = pzon_int(:,:,jp_sal,jn) / MAX( z4d1(Nis0,:,:,jn), 10.e-15 )
DO ji = Nis0+1, Nie0
z4d2(ji,:,:,jn) = z4d2(Nis0,:,:,jn)
ENDDO
ENDDO
CALL iom_put( 'zosal', z4d2 )
......@@ -251,16 +251,16 @@ CONTAINS
IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN
!
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sophtadv', z3dtr )
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sopstadv', z3dtr )
......@@ -269,16 +269,16 @@ CONTAINS
IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN
!
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sophtldf', z3dtr )
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sopstldf', z3dtr )
......@@ -287,16 +287,16 @@ CONTAINS
IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN
!
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sophteiv', z3dtr )
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sopsteiv', z3dtr )
......@@ -304,16 +304,16 @@ CONTAINS
!
IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sophtvtr', z3dtr )
DO jn = 1, nbasin
z3dtr(1,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = 2, jpi
z3dtr(ji,:,jn) = z3dtr(1,:,jn)
z3dtr(Nis0,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg)
DO ji = Nis0+1, Nie0
z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn)
ENDDO
ENDDO
CALL iom_put( 'sopstvtr', z3dtr )
......@@ -349,8 +349,8 @@ CONTAINS
!! ** Action : pvtr_int - terms for volume streamfunction, heat/salt transport barotropic/overturning terms
!! pzon_int - terms for i mean temperature/salinity
!!----------------------------------------------------------------------
INTEGER , INTENT(in) :: Kmm ! time level index
REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(in), OPTIONAL :: pvtr ! j-effective transport
INTEGER , INTENT(in) :: Kmm ! time level index
REAL(wp), DIMENSION(A2D(0),jpk), INTENT(in), OPTIONAL :: pvtr ! j-effective transport
REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zmask ! 3D workspace
REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE :: zts ! 4D workspace
REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: sjk, v_msf ! Zonal sum: i-k surface area, j-effective transport
......@@ -362,7 +362,7 @@ CONTAINS
IF( PRESENT( pvtr ) ) THEN
! i sum of effective j transport excluding closed seas
IF( iom_use( 'zomsf' ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
ALLOCATE( v_msf(A1Dj(nn_hls),jpk,nbasin) )
ALLOCATE( v_msf(A1Dj(0),jpk,nbasin) )
DO jn = 1, nbasin
v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )
......@@ -374,16 +374,16 @@ CONTAINS
ENDIF
! i sum of j surface area, j surface area - temperature/salinity product on V grid
IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR. &
IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR. &
& iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
ALLOCATE( zmask(A2D(nn_hls),jpk), zts(A2D(nn_hls),jpk,jpts), &
& sjk(A1Dj(nn_hls),jpk,nbasin), &
& zt_jk(A1Dj(nn_hls),jpk,nbasin), zs_jk(A1Dj(nn_hls),jpk,nbasin) )
ALLOCATE( zmask( A2D(0),jpk ), zts( A2D(0),jpk,jpts ), &
& sjk( A1Dj(0),jpk,nbasin), &
& zt_jk(A1Dj(0),jpk,nbasin), zs_jk(A1Dj(0),jpk,nbasin) )
zmask(:,:,:) = 0._wp
zts(:,:,:,:) = 0._wp
DO_3D( 1, 1, 1, 0, 1, jpkm1 )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
zmask(ji,jj,jk) = vmask(ji,jj,jk) * zvfc
zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc !Tracers averaged onto V grid
......@@ -405,14 +405,14 @@ CONTAINS
ELSE
! i sum of j surface area - temperature/salinity product on T grid
IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' ) ) THEN
ALLOCATE( zmask(A2D(nn_hls),jpk), zts(A2D(nn_hls),jpk,jpts), &
& sjk(A1Dj(nn_hls),jpk,nbasin), &
& zt_jk(A1Dj(nn_hls),jpk,nbasin), zs_jk(A1Dj(nn_hls),jpk,nbasin) )
ALLOCATE( zmask( A2D(0),jpk ), zts( A2D(0),jpk,jpts ), &
& sjk( A1Dj(0),jpk,nbasin), &
& zt_jk(A1Dj(0),jpk,nbasin), zs_jk(A1Dj(0),jpk,nbasin) )
zmask(:,:,:) = 0._wp
zts(:,:,:,:) = 0._wp
DO_3D( 1, 1, 1, 1, 1, jpkm1 )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zsfc = e1t(ji,jj) * e3t(ji,jj,jk,Kmm)
zmask(ji,jj,jk) = tmask(ji,jj,jk) * zsfc
zts(ji,jj,jk,jp_tem) = ts(ji,jj,jk,jp_tem,Kmm) * zsfc
......@@ -434,11 +434,12 @@ CONTAINS
! i-k sum of j surface area - temperature/salinity product on V grid
IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN
! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations
ALLOCATE( zts(A2D(nn_hls),jpk,jpts) )
zts(:,:,:,:) = 0._wp
DO_3D( 1, 1, 1, 0, 1, jpkm1 )
DO_3D( 0, 0, 0, 0, 1, jpkm1 )
zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc !Tracers averaged onto V grid
zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
......@@ -538,11 +539,12 @@ CONTAINS
!! Wrapper for heat and salt transport calculations to calculate them for each basin
!! Called from all advection and/or diffusion routines
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: ktra ! tracer index
CHARACTER(len=3) , INTENT(in) :: cptr ! transport type 'adv'/'ldf'/'eiv'
REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in) :: pvflx ! 3D input array of advection/diffusion
REAL(wp), DIMENSION(A1Dj(nn_hls),nbasin) :: zsj !
INTEGER :: jn !
INTEGER , INTENT(in) :: ktra ! tracer index
CHARACTER(len=3) , INTENT(in) :: cptr ! transport type 'adv'/'ldf'/'eiv'
! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations
REAL(wp), DIMENSION(A2D(0),jpk) , INTENT(in) :: pvflx ! 3D input array of advection/diffusion
REAL(wp), DIMENSION(A1Dj(0),nbasin) :: zsj !
INTEGER :: jn !
DO jn = 1, nbasin
zsj(:,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
......@@ -576,13 +578,13 @@ CONTAINS
!!
!! ** Action : phstr
!!----------------------------------------------------------------------
REAL(wp), DIMENSION(jpj,nbasin) , INTENT(inout) :: phstr !
REAL(wp), DIMENSION(A1Dj(nn_hls),nbasin), INTENT(in) :: pva !
REAL(wp), DIMENSION(Njs0:Nje0,nbasin), INTENT(inout) :: phstr !
REAL(wp), DIMENSION(A1Dj(0) ,nbasin), INTENT(in ) :: pva !
INTEGER :: jj
#if ! defined key_mpi_off
INTEGER, DIMENSION(1) :: ish1d
INTEGER, DIMENSION(2) :: ish2d
REAL(wp), DIMENSION(jpj*nbasin) :: zwork
INTEGER, DIMENSION(1) :: ish1d
INTEGER, DIMENSION(2) :: ish2d
REAL(wp), DIMENSION(:), ALLOCATABLE :: zwork
#endif
DO jj = ntsj, ntej
......@@ -591,11 +593,13 @@ CONTAINS
#if ! defined key_mpi_off
IF( .NOT. l_istiled .OR. ntile == nijtile ) THEN
ish1d(1) = jpj*nbasin
ish2d(1) = jpj ; ish2d(2) = nbasin
ALLOCATE( zwork(Nj_0*nbasin) )
ish1d(1) = Nj_0*nbasin
ish2d(1) = Nj_0 ; ish2d(2) = nbasin
zwork(:) = RESHAPE( phstr(:,:), ish1d )
CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl )
phstr(:,:) = RESHAPE( zwork, ish2d )
DEALLOCATE( zwork )
ENDIF
#endif
END SUBROUTINE ptr_sum_2d
......@@ -612,13 +616,13 @@ CONTAINS
!!
!! ** Action : phstr
!!----------------------------------------------------------------------
REAL(wp), DIMENSION(jpj,jpk,nbasin) , INTENT(inout) :: phstr !
REAL(wp), DIMENSION(A1Dj(nn_hls),jpk,nbasin), INTENT(in) :: pva !
INTEGER :: jj, jk
REAL(wp), DIMENSION(Njs0:Nje0,jpk,nbasin), INTENT(inout) :: phstr !
REAL(wp), DIMENSION(A1Dj(0) ,jpk,nbasin), INTENT(in ) :: pva !
INTEGER :: jj, jk
#if ! defined key_mpi_off
INTEGER, DIMENSION(1) :: ish1d
INTEGER, DIMENSION(3) :: ish3d
REAL(wp), DIMENSION(jpj*jpk*nbasin) :: zwork
INTEGER, DIMENSION(1) :: ish1d
INTEGER, DIMENSION(3) :: ish3d
REAL(wp), DIMENSION(:), ALLOCATABLE :: zwork
#endif
DO jk = 1, jpk
......@@ -629,11 +633,13 @@ CONTAINS
#if ! defined key_mpi_off
IF( .NOT. l_istiled .OR. ntile == nijtile ) THEN
ish1d(1) = jpj*jpk*nbasin
ish3d(1) = jpj ; ish3d(2) = jpk ; ish3d(3) = nbasin
ALLOCATE( zwork(Nj_0*jpk*nbasin) )
ish1d(1) = Nj_0*jpk*nbasin
ish3d(1) = Nj_0 ; ish3d(2) = jpk ; ish3d(3) = nbasin
zwork(:) = RESHAPE( phstr(:,:,:), ish1d )
CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl )
phstr(:,:,:) = RESHAPE( zwork, ish3d )
DEALLOCATE( zwork )
ENDIF
#endif
END SUBROUTINE ptr_sum_3d
......@@ -651,13 +657,13 @@ CONTAINS
! nbasin has been initialized in iom_init to define the axis "basin"
!
IF( .NOT. ALLOCATED( btmsk ) ) THEN
ALLOCATE( btmsk(jpi,jpj,nbasin) , btmsk34(jpi,jpj,nbasin), &
& hstr_adv(jpj,jpts,nbasin), hstr_eiv(jpj,jpts,nbasin), &
& hstr_ove(jpj,jpts,nbasin), hstr_btr(jpj,jpts,nbasin), &
& hstr_ldf(jpj,jpts,nbasin), hstr_vtr(jpj,jpts,nbasin), STAT=ierr(1) )
ALLOCATE( btmsk(jpi,jpj,nbasin) , btmsk34(jpi,jpj,nbasin), &
& hstr_adv(Njs0:Nje0,jpts,nbasin), hstr_eiv(Njs0:Nje0,jpts,nbasin), &
& hstr_ove(Njs0:Nje0,jpts,nbasin), hstr_btr(Njs0:Nje0,jpts,nbasin), &
& hstr_ldf(Njs0:Nje0,jpts,nbasin), hstr_vtr(Njs0:Nje0,jpts,nbasin), STAT=ierr(1) )
!
ALLOCATE( pvtr_int(jpj,jpk,jpts+2,nbasin), &
& pzon_int(jpj,jpk,jpts+1,nbasin), STAT=ierr(2) )
ALLOCATE( pvtr_int(Njs0:Nje0,jpk,jpts+2,nbasin), &
& pzon_int(Njs0:Nje0,jpk,jpts+1,nbasin), STAT=ierr(2) )
!
dia_ptr_alloc = MAXVAL( ierr )
CALL mpp_sum( 'diaptr', dia_ptr_alloc )
......@@ -677,11 +683,12 @@ CONTAINS
!!
!! ** Action : - p_fval: i-k-mean poleward flux of pvflx
!!----------------------------------------------------------------------
REAL(wp), INTENT(in), DIMENSION(A2D(nn_hls),jpk) :: pvflx ! mask flux array at V-point
REAL(wp), INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask
! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations
REAL(wp), INTENT(in), DIMENSION(A2D(0),jpk) :: pvflx ! mask flux array at V-point
REAL(wp), INTENT(in), DIMENSION(jpi,jpj ) :: pmsk ! Optional 2D basin mask
!
INTEGER :: ji, jj, jk ! dummy loop arguments
REAL(wp), DIMENSION(A1Dj(nn_hls)) :: p_fval ! function value
INTEGER :: ji, jj, jk ! dummy loop arguments
REAL(wp), DIMENSION(A1Dj(0)) :: p_fval ! function value
!!--------------------------------------------------------------------
!
p_fval(:) = 0._wp
......@@ -702,11 +709,12 @@ CONTAINS
!!
!! ** Action : - p_fval: i-k-mean poleward flux of pvflx
!!----------------------------------------------------------------------
REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls)) :: pvflx ! mask flux array at V-point
REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask
! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations
REAL(wp) , INTENT(in), DIMENSION(A2D(0) ) :: pvflx ! mask flux array at V-point
REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask
!
INTEGER :: ji,jj ! dummy loop arguments
REAL(wp), DIMENSION(A1Dj(nn_hls)) :: p_fval ! function value
INTEGER :: ji, jj ! dummy loop arguments
REAL(wp), DIMENSION(A1Dj(0)) :: p_fval ! function value
!!--------------------------------------------------------------------
!
p_fval(:) = 0._wp
......@@ -725,14 +733,13 @@ CONTAINS
!!
!! ** Action : - p_fval: j-cumulated sum of pva
!!----------------------------------------------------------------------
REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pva ! mask flux array at V-point
REAL(wp) , INTENT(in), DIMENSION(A2D(0)) :: pva ! mask flux array at V-point
!
INTEGER :: ji,jj,jc ! dummy loop arguments
INTEGER :: ijpj ! ???
REAL(wp), DIMENSION(jpi,jpj) :: p_fval ! function value
INTEGER :: ji,jj,jc ! dummy loop arguments
INTEGER :: ijpj ! ???
REAL(wp), DIMENSION(A2D(0)) :: p_fval ! function value
!!--------------------------------------------------------------------
!
ijpj = jpj ! ???
p_fval(:,:) = 0._wp
DO jc = 1, jpnj ! looping over all processors in j axis
DO_2D( 0, 0, 0, 0 )
......@@ -756,11 +763,11 @@ CONTAINS
!!----------------------------------------------------------------------
!!
IMPLICIT none
REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls),jpk) :: pta ! mask flux array at V-point
REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask
REAL(wp) , INTENT(in), DIMENSION(A2D(0) ,jpk) :: pta ! mask flux array at V-point
REAL(wp) , INTENT(in), DIMENSION(jpi,jpj ) :: pmsk ! Optional 2D basin mask
!!
INTEGER :: ji, jj, jk ! dummy loop arguments
REAL(wp), DIMENSION(A1Dj(nn_hls),jpk) :: p_fval ! return function value
REAL(wp), DIMENSION(A1Dj(0),jpk) :: p_fval ! return function value
!!--------------------------------------------------------------------
!
p_fval(:,:) = 0._wp
......