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
Expand all Hide whitespace changes
Inline Side-by-side
Showing
with 1621 additions and 743 deletions
src/ICE/iceupdate.F90
View file @ fc409945
...@@ -314,7 +314,7 @@ CONTAINS ...@@ -314,7 +314,7 @@ CONTAINS
!! !!
!! ** Action : * at each ice time step (every nn_fsbc time step): !! ** Action : * at each ice time step (every nn_fsbc time step):
!! - compute the modulus of ice-ocean relative velocity !! - 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 | !! tmod_io = rhoco * | U_ice-U_oce |
!! - update the modulus of stress at ocean surface !! - update the modulus of stress at ocean surface
!! taum = (1-a) * taum + a * tmod_io * | U_ice-U_oce | !! taum = (1-a) * taum + a * tmod_io * | U_ice-U_oce |
...@@ -325,19 +325,19 @@ CONTAINS ...@@ -325,19 +325,19 @@ CONTAINS
!! !!
!! NB: - ice-ocean rotation angle no more allowed !! NB: - ice-ocean rotation angle no more allowed
!! - here we make an approximation: taum is only computed every ice time step !! - 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 !! This avoids mutiple average to pass from U,V grids to T grids
!! to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... !! 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 !! - taum : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes
!!--------------------------------------------------------------------- !!---------------------------------------------------------------------
INTEGER , INTENT(in) :: kt ! ocean time-step index INTEGER , INTENT(in) :: kt ! ocean time-step index
REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents
! !
INTEGER :: ji, jj ! dummy loop indices INTEGER :: ji, jj ! dummy loop indices
REAL(wp) :: zat_u, zutau_ice, zu_t, zmodt ! local scalar REAL(wp) :: zutau_ice, zu_t, zmodt ! local scalar
REAL(wp) :: zat_v, zvtau_ice, zv_t, zrhoco ! - - REAL(wp) :: zvtau_ice, zv_t, zrhoco ! - -
REAL(wp) :: zflagi ! - - REAL(wp) :: zflagi ! - -
!!--------------------------------------------------------------------- !!---------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('iceupdate') IF( ln_timing ) CALL timing_start('iceupdate')
...@@ -352,8 +352,8 @@ CONTAINS ...@@ -352,8 +352,8 @@ CONTAINS
IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step) 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( 0, 0, 0, 0 ) !* update the modulus of stress at ocean surface (T-point)
! ! 2*(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) 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) 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 ! ! |U_ice-U_oce|^2
zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t ) zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t )
! ! update the ocean stress modulus ! ! update the ocean stress modulus
...@@ -377,19 +377,14 @@ CONTAINS ...@@ -377,19 +377,14 @@ CONTAINS
ENDIF ENDIF
! !
DO_2D( 0, 0, 0, 0 ) !* update the stress WITHOUT an ice-ocean rotation angle 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 ! ! 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) ) 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) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,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 ! ! stresses at the ocean surface
utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_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 - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice vtau(ji,jj) = ( 1._wp - at_i(ji,jj) ) * vtau_oce(ji,jj) + at_i(ji,jj) * zvtau_ice
END_2D END_2D
CALL lbc_lnk( 'iceupdate', utau, 'U', -1.0_wp, vtau, 'V', -1.0_wp ) ! lateral boundary condition CALL lbc_lnk( 'iceupdate', utau, 'T', -1.0_wp, vtau, 'T', -1.0_wp ) ! lateral boundary condition
! !
IF( ln_timing ) CALL timing_stop('iceupdate') IF( ln_timing ) CALL timing_stop('iceupdate')
! !
......
src/OCE/CRS/crsfld.F90
View file @ fc409945
...@@ -211,8 +211,8 @@ CONTAINS ...@@ -211,8 +211,8 @@ CONTAINS
! sbc fields ! 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( 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( 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', 'V', vmask, vtau_crs , p_e12=e1v , p_surf_crs=e1v_crs , psgn=1.0_wp ) 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( 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( 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 ) CALL crs_dom_ope( qsr , 'SUM', 'T', tmask, qsr_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0_wp )
......
src/OCE/DIA/diadetide.F90
View file @ fc409945
...@@ -5,11 +5,11 @@ MODULE diadetide ...@@ -5,11 +5,11 @@ MODULE diadetide
!!====================================================================== !!======================================================================
!! History : ! 2019 (S. Mueller) !! History : ! 2019 (S. Mueller)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
USE par_oce , ONLY : wp, jpi, jpj USE par_oce
USE in_out_manager , ONLY : lwp, numout USE in_out_manager
USE iom , ONLY : iom_put USE iom
USE dom_oce , ONLY : rn_Dt, nsec_day USE dom_oce
USE phycst , ONLY : rpi USE phycst
USE tide_mod USE tide_mod
#if defined key_xios #if defined key_xios
USE xios USE xios
...@@ -24,6 +24,8 @@ MODULE diadetide ...@@ -24,6 +24,8 @@ MODULE diadetide
PUBLIC :: dia_detide_init, dia_detide PUBLIC :: dia_detide_init, dia_detide
!! * Substitutions
# include "do_loop_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2019) !! NEMO/OCE 4.0 , NEMO Consortium (2019)
!! $Id$ !! $Id$
...@@ -90,9 +92,9 @@ CONTAINS ...@@ -90,9 +92,9 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt INTEGER, INTENT(in) :: kt
REAL(wp), DIMENSION(jpi,jpj) :: zwght_2D REAL(wp), DIMENSION(A2D(0)) :: zwght_2D
REAL(wp) :: zwght, ztmp REAL(wp) :: zwght, ztmp
INTEGER :: jn INTEGER :: ji, jj, jn
! Compute detiding weight at the current time-step; the daily total weight ! 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 ! is one, and the daily summation of a diagnosed field multiplied by this
...@@ -104,7 +106,10 @@ CONTAINS ...@@ -104,7 +106,10 @@ CONTAINS
zwght = zwght + 1.0_wp / REAL( ndiadetide, KIND=wp ) zwght = zwght + 1.0_wp / REAL( ndiadetide, KIND=wp )
END IF END IF
END DO 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) CALL iom_put( "diadetide_weight", zwght_2D)
END SUBROUTINE dia_detide END SUBROUTINE dia_detide
......
src/OCE/DIA/diahsb.F90
View file @ fc409945
...@@ -50,6 +50,7 @@ MODULE diahsb ...@@ -50,6 +50,7 @@ MODULE diahsb
REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: tmask_ini REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: tmask_ini
!! * Substitutions !! * Substitutions
# include "do_loop_substitute.h90"
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
...@@ -82,30 +83,61 @@ CONTAINS ...@@ -82,30 +83,61 @@ CONTAINS
REAL(wp) :: z_frc_trd_v ! - - REAL(wp) :: z_frc_trd_v ! - -
REAL(wp) :: z_wn_trd_t , z_wn_trd_s ! - - REAL(wp) :: z_wn_trd_t , z_wn_trd_s ! - -
REAL(wp) :: z_ssh_hc , z_ssh_sc ! - - REAL(wp) :: z_ssh_hc , z_ssh_sc ! - -
REAL(wp), DIMENSION(jpi,jpj,13) :: ztmp REAL(wp), DIMENSION(A2D(0),13) :: ztmp
REAL(wp), DIMENSION(jpi,jpj,jpkm1,4) :: ztmpk REAL(wp), DIMENSION(A2D(0),jpkm1,4) :: ztmpk
REAL(wp), DIMENSION(17) :: zbg REAL(wp), DIMENSION(17) :: zbg
!!--------------------------------------------------------------------------- !!---------------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_hsb') IF( ln_timing ) CALL timing_start('dia_hsb')
! !
ztmp (:,:,:) = 0._wp ! should be better coded DO_2D( 0, 0, 0, 0 )
ztmpk(:,:,:,:) = 0._wp ! should be better coded ztmp (ji,jj,:) = 0._wp ! should be better coded
! ztmpk(ji,jj,:,:) = 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(:,:,:) ; 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 ! ! 1 - Trends due to forcing !
! ------------------------- ! ! ------------------------- !
! prepare trends ! prepare trends
ztmp(:,:,1) = - r1_rho0 * ( emp(:,:) - rnf(:,:) - fwfisf_cav(:,:) - fwfisf_par(:,:) ) * surf(:,:) ! volume DO_2D( 0, 0, 0, 0 )
ztmp(:,:,2) = sbc_tsc(:,:,jp_tem) * surf(:,:) ! heat ztmp(ji,jj,1) = - r1_rho0 * ( emp(ji,jj) & ! volume
ztmp(:,:,3) = sbc_tsc(:,:,jp_sal) * surf(:,:) ! salt & - rnf(ji,jj) &
IF( ln_rnf ) ztmp(:,:,4) = rnf_tsc(:,:,jp_tem) * surf(:,:) ! runoff temp & - fwfisf_cav(ji,jj) &
IF( ln_rnf_sal ) ztmp(:,:,5) = rnf_tsc(:,:,jp_sal) * surf(:,:) ! runoff salt & - fwfisf_par(ji,jj) ) * surf(ji,jj)
IF( ln_isf ) ztmp(:,:,6) = ( risf_cav_tsc(:,:,jp_tem) + risf_par_tsc(:,:,jp_tem) ) * surf(:,:) ! isf temp ztmp(ji,jj,2) = sbc_tsc(ji,jj,jp_tem) * surf(ji,jj) ! heat
IF( ln_traqsr ) ztmp(:,:,7) = r1_rho0_rcp * qsr(:,:) * surf(:,:) ! penetrative solar radiation ztmp(ji,jj,3) = sbc_tsc(ji,jj,jp_sal) * surf(ji,jj) ! salt
IF( ln_trabbc ) ztmp(:,:,8) = qgh_trd0(:,:) * surf(:,:) ! geothermal heat 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_linssh ) THEN ! Advection flux through fixed surface (z=0)
IF( ln_isfcav ) THEN IF( ln_isfcav ) THEN
...@@ -116,8 +148,10 @@ CONTAINS ...@@ -116,8 +148,10 @@ CONTAINS
END DO END DO
END DO END DO
ELSE ELSE
ztmp(:,:,9 ) = - surf(:,:) * ww(:,:,1) * ts(:,:,1,jp_tem,Kbb) DO_2D( 0, 0, 0, 0 )
ztmp(:,:,10) = - surf(:,:) * ww(:,:,1) * ts(:,:,1,jp_sal,Kbb) 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 END IF
ENDIF ENDIF
...@@ -152,7 +186,9 @@ CONTAINS ...@@ -152,7 +186,9 @@ CONTAINS
! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl) ! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl)
! !
! ! volume variation (calculated with ssh) ! ! 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) ! ! heat & salt content variation (associated with ssh)
IF( ln_linssh ) THEN ! linear free surface case IF( ln_linssh ) THEN ! linear free surface case
...@@ -164,8 +200,10 @@ CONTAINS ...@@ -164,8 +200,10 @@ CONTAINS
END DO END DO
END DO END DO
ELSE ! no under ice-shelf seas ELSE ! no under ice-shelf seas
ztmp(:,:,12) = surf(:,:) * ( ts(:,:,1,jp_tem,Kmm) * ssh(:,:,Kmm) - ssh_hc_loc_ini(:,:) ) DO_2D( 0, 0, 0, 0 )
ztmp(:,:,13) = surf(:,:) * ( ts(:,:,1,jp_sal,Kmm) * ssh(:,:,Kmm) - ssh_sc_loc_ini(:,:) ) 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 END IF
ENDIF ENDIF
...@@ -185,19 +223,27 @@ CONTAINS ...@@ -185,19 +223,27 @@ CONTAINS
! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl) ! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl)
! !
DO jk = 1, jpkm1 ! volume DO jk = 1, jpkm1 ! volume
ztmpk(:,:,jk,1) = surf (:,:) * e3t(:,:,jk,Kmm)*tmask(:,:,jk) & DO_2D( 0, 0, 0, 0 )
& - surf_ini(:,:) * e3t_ini(:,:,jk )*tmask_ini(:,:,jk) 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 END DO
DO jk = 1, jpkm1 ! heat DO jk = 1, jpkm1 ! heat
ztmpk(:,:,jk,2) = ( surf (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_tem,Kmm) & DO_2D( 0, 0, 0, 0 )
& - surf_ini(:,:) * hc_loc_ini(:,:,jk) ) 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 END DO
DO jk = 1, jpkm1 ! salt DO jk = 1, jpkm1 ! salt
ztmpk(:,:,jk,3) = ( surf (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_sal,Kmm) & DO_2D( 0, 0, 0, 0 )
& - surf_ini(:,:) * sc_loc_ini(:,:,jk) ) 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 END DO
DO jk = 1, jpkm1 ! total ocean volume 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 END DO
! global sum ! global sum
...@@ -315,14 +361,18 @@ CONTAINS ...@@ -315,14 +361,18 @@ CONTAINS
IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' dia_hsb_rst : initialise hsb at initial state ' IF(lwp) WRITE(numout,*) ' dia_hsb_rst : initialise hsb at initial state '
IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*)
surf_ini(:,:) = e1e2t(:,:) * tmask_i(:,:) ! initial ocean surface DO_2D( 0, 0, 0, 0 )
ssh_ini(:,:) = ssh(:,:,Kmm) ! initial ssh 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
DO jk = 1, jpk DO jk = 1, jpk
! if ice sheet/oceqn coupling, need to mask ini variables here (mask could change at the next NEMO instance). ! 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 DO_2D( 0, 0, 0, 0 )
tmask_ini (:,:,jk) = tmask(:,:,jk) ! initial mask e3t_ini (ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial vertical scale factors
hc_loc_ini(:,:,jk) = ts(:,:,jk,jp_tem,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial heat content tmask_ini (ji,jj,jk) = tmask(ji,jj,jk) ! initial mask
sc_loc_ini(:,:,jk) = ts(:,:,jk,jp_sal,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial salt content 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_2D
END DO END DO
frc_v = 0._wp ! volume trend due to forcing frc_v = 0._wp ! volume trend due to forcing
frc_t = 0._wp ! heat content - - - - frc_t = 0._wp ! heat content - - - -
...@@ -334,13 +384,15 @@ CONTAINS ...@@ -334,13 +384,15 @@ CONTAINS
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_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 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
END DO END DO
ELSE ELSE
ssh_hc_loc_ini(:,:) = ts(:,:,1,jp_tem,Kmm) * ssh(:,:,Kmm) ! initial heat content in ssh DO_2D( 0, 0, 0, 0 )
ssh_sc_loc_ini(:,:) = ts(:,:,1,jp_sal,Kmm) * ssh(:,:,Kmm) ! initial salt content in ssh 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 END IF
frc_wn_t = 0._wp ! initial heat 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 frc_wn_s = 0._wp ! initial salt content misfit due to free surface
ENDIF ENDIF
ENDIF ENDIF
! !
...@@ -388,6 +440,7 @@ CONTAINS ...@@ -388,6 +440,7 @@ CONTAINS
INTEGER, INTENT(in) :: Kmm ! time level index INTEGER, INTENT(in) :: Kmm ! time level index
! !
INTEGER :: ierror, ios ! local integer INTEGER :: ierror, ios ! local integer
INTEGER :: ji, jj ! loop index
!! !!
NAMELIST/namhsb/ ln_diahsb NAMELIST/namhsb/ ln_diahsb
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
...@@ -427,7 +480,10 @@ CONTAINS ...@@ -427,7 +480,10 @@ CONTAINS
! ----------------------------------------------- ! ! ----------------------------------------------- !
! 2 - Time independant variables and file opening ! ! 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) * tmask_i(ji,jj) ! masked surface grid cell area
END_2D
surf_tot = glob_sum( 'diahsb', surf(:,:) ) ! total ocean surface area 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' ) IF( ln_bdy ) CALL ctl_warn( 'dia_hsb_init: heat/salt budget does not consider open boundary fluxes' )
......
src/OCE/DIA/diahth.F90
View file @ fc409945
...@@ -86,22 +86,22 @@ CONTAINS ...@@ -86,22 +86,22 @@ CONTAINS
INTEGER, INTENT( in ) :: kt ! ocean time-step index INTEGER, INTENT( in ) :: kt ! ocean time-step index
INTEGER, INTENT( in ) :: Kmm ! ocean time level index INTEGER, INTENT( in ) :: Kmm ! ocean time level index
!! !!
INTEGER :: ji, jj, jk ! dummy loop arguments INTEGER :: ji, jj, jk ! dummy loop arguments
REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth 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) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth
REAL(wp) :: ztem2 = 0.2_wp ! temperature 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) :: zztmp, zzdep ! temporary scalars inside do loop
REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace
REAL(wp), DIMENSION(jpi,jpj) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 REAL(wp), DIMENSION(A2D(0)) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2
REAL(wp), DIMENSION(jpi,jpj) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 REAL(wp), DIMENSION(A2D(0)) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2
REAL(wp), DIMENSION(jpi,jpj) :: zrho10_3 ! MLD: rho = rho10m + zrho3 REAL(wp), DIMENSION(A2D(0)) :: 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(A2D(0)) :: 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(A2D(0)) :: ztinv ! max of temperature inversion
REAL(wp), DIMENSION(jpi,jpj) :: zdepinv ! depth of temperature inversion REAL(wp), DIMENSION(A2D(0)) :: zdepinv ! depth of temperature inversion
REAL(wp), DIMENSION(jpi,jpj) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 REAL(wp), DIMENSION(A2D(0)) :: zrho0_3 ! MLD rho = rho(surf) = 0.03
REAL(wp), DIMENSION(jpi,jpj) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 REAL(wp), DIMENSION(A2D(0)) :: zrho0_1 ! MLD rho = rho(surf) = 0.01
REAL(wp), DIMENSION(jpi,jpj) :: zmaxdzT ! max of dT/dz REAL(wp), DIMENSION(A2D(0)) :: zmaxdzT ! max of dT/dz
REAL(wp), DIMENSION(jpi,jpj) :: zdelr ! delta rho equivalent to deltaT = 0.2 REAL(wp), DIMENSION(A2D(0)) :: zdelr ! delta rho equivalent to deltaT = 0.2
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_hth') IF( ln_timing ) CALL timing_start('dia_hth')
...@@ -131,7 +131,7 @@ CONTAINS ...@@ -131,7 +131,7 @@ CONTAINS
IF( iom_use( 'mlddzt' ) ) zmaxdzT(:,:) = 0._wp IF( iom_use( 'mlddzt' ) ) zmaxdzT(:,:) = 0._wp
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) & IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) &
& .OR. iom_use( 'mldr10_3' ) .OR. iom_use( 'pycndep' ) ) THEN & .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) zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
hth (ji,jj) = zztmp hth (ji,jj) = zztmp
zabs2 (ji,jj) = zztmp zabs2 (ji,jj) = zztmp
...@@ -142,7 +142,7 @@ CONTAINS ...@@ -142,7 +142,7 @@ CONTAINS
ENDIF ENDIF
IF( iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN IF( iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN
IF( nla10 > 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) zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
zrho0_3(ji,jj) = zztmp zrho0_3(ji,jj) = zztmp
zrho0_1(ji,jj) = zztmp zrho0_1(ji,jj) = zztmp
...@@ -157,7 +157,7 @@ CONTAINS ...@@ -157,7 +157,7 @@ CONTAINS
! MLD: rho = rho(1) + zrho3 ! ! MLD: rho = rho(1) + zrho3 !
! MLD: rho = rho(1) + zrho1 ! ! 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) zzdep = gdepw(ji,jj,jk,Kmm)
zztmp = ( ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) & zztmp = ( ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) &
...@@ -189,7 +189,7 @@ CONTAINS ...@@ -189,7 +189,7 @@ CONTAINS
! !
! Preliminary computation ! Preliminary computation
! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) ! 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 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) & 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) & & - 0.0745 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm) &
...@@ -213,7 +213,7 @@ CONTAINS ...@@ -213,7 +213,7 @@ CONTAINS
! temperature inversion: max( 0, max of tn - tn(10m) ) ! ! temperature inversion: max( 0, max of tn - tn(10m) ) !
! depth of temperature inversion ! ! 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) zzdep = gdepw(ji,jj,jk,Kmm) * tmask(ji,jj,1)
! !
...@@ -305,13 +305,16 @@ CONTAINS ...@@ -305,13 +305,16 @@ CONTAINS
! !
INTEGER :: ji, jj, jk, iid INTEGER :: ji, jj, jk, iid
REAL(wp) :: zztmp, zzdep REAL(wp) :: zztmp, zzdep
INTEGER, DIMENSION(jpi,jpj) :: iktem INTEGER, DIMENSION(A2D(0)) :: iktem
! --------------------------------------- ! ! --------------------------------------- !
! search deepest level above ptem ! ! search deepest level above ptem !
! --------------------------------------- ! ! --------------------------------------- !
iktem(:,:) = 1 DO_2D( 0, 0, 0, 0 )
DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom 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) zztmp = ts(ji,jj,jk,jp_tem,Kmm)
IF( zztmp >= ptem ) iktem(ji,jj) = jk IF( zztmp >= ptem ) iktem(ji,jj) = jk
END_3D END_3D
...@@ -319,7 +322,7 @@ CONTAINS ...@@ -319,7 +322,7 @@ CONTAINS
! ------------------------------- ! ! ------------------------------- !
! Depth of ptem isotherm ! ! 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 zzdep = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) ! depth of the ocean bottom
! !
...@@ -346,18 +349,29 @@ CONTAINS ...@@ -346,18 +349,29 @@ CONTAINS
REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc
! !
INTEGER :: ji, jj, jk, ik INTEGER :: ji, jj, jk, ik
REAL(wp), DIMENSION(jpi,jpj) :: zthick REAL(wp), DIMENSION(A2D(0)) :: zthick
INTEGER , DIMENSION(jpi,jpj) :: ilevel INTEGER , DIMENSION(A2D(0)) :: ilevel
! surface boundary condition ! surface boundary condition
IF( .NOT. ln_linssh ) THEN ; zthick(:,:) = 0._wp ; phtc(:,:) = 0._wp IF( .NOT. ln_linssh ) THEN
ELSE ; zthick(:,:) = ssh(:,:,Kmm) ; phtc(:,:) = pt(:,:,1) * ssh(:,:,Kmm) * tmask(:,:,1) 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 ENDIF
! !
ilevel(:,:) = 1 DO_2D( 0, 0, 0, 0 )
DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 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 IF( ( gdepw(ji,jj,jk+1,Kmm) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN
ilevel(ji,jj) = jk+1 ilevel(ji,jj) = jk+1
zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm) zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm)
...@@ -365,7 +379,7 @@ CONTAINS ...@@ -365,7 +379,7 @@ CONTAINS
ENDIF ENDIF
END_3D END_3D
! !
DO_2D( 1, 1, 1, 1 ) DO_2D( 0, 0, 0, 0 )
ik = ilevel(ji,jj) ik = ilevel(ji,jj)
IF( tmask(ji,jj,ik) == 1 ) THEN 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 zthick(ji,jj) = MIN ( gdepw(ji,jj,ik+1,Kmm), pdep ) - zthick(ji,jj) ! remaining thickness to reach dephw pdep
......
src/OCE/DIA/diamlr.F90
View file @ fc409945
...@@ -6,7 +6,7 @@ MODULE diamlr ...@@ -6,7 +6,7 @@ MODULE diamlr
!! History : 4.0 ! 2019 (S. Mueller) Original code !! 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 phycst , ONLY : rpi
USE dom_oce , ONLY : adatrj USE dom_oce , ONLY : adatrj
USE tide_mod USE tide_mod
...@@ -407,8 +407,9 @@ CONTAINS ...@@ -407,8 +407,9 @@ CONTAINS
!! ** Purpose : update time used in multiple-linear-regression analysis !! ** 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') IF( ln_timing ) CALL timing_start('dia_mlr')
...@@ -417,7 +418,9 @@ CONTAINS ...@@ -417,7 +418,9 @@ CONTAINS
! !
! A 2-dimensional field of constant value is sent, and subsequently used directly ! 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. ! 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 ( iom_use('diamlr_time') ) CALL iom_put('diamlr_time', zadatrj2d)
! !
IF( ln_timing ) CALL timing_stop('dia_mlr') IF( ln_timing ) CALL timing_stop('dia_mlr')
......
src/OCE/DIA/diaptr.F90
View file @ fc409945
This diff is collapsed.
src/OCE/DIA/diawri.F90
View file @ fc409945
...@@ -135,8 +135,8 @@ CONTAINS ...@@ -135,8 +135,8 @@ CONTAINS
ENDIF ENDIF
! initialize arrays ! initialize arrays
z2d(:,:) = 0._wp z2d(A2D(0)) = 0._wp
z3d(:,:,:) = 0._wp z3d(A2D(0),:) = 0._wp
! Output of initial vertical scale factor ! Output of initial vertical scale factor
CALL iom_put("e3t_0", e3t_0(:,:,:) ) CALL iom_put("e3t_0", e3t_0(:,:,:) )
...@@ -868,7 +868,11 @@ CONTAINS ...@@ -868,7 +868,11 @@ CONTAINS
CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "J/m2" , & ! htc3 CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "J/m2" , & ! htc3
& jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout )
#endif #endif
CALL histdef( nid_T, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
& jpi, jpj, nh_T, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
CALL histdef( nid_T, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
& jpi, jpj, nh_T, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
!
CALL histend( nid_T, snc4chunks=snc4set ) CALL histend( nid_T, snc4chunks=snc4set )
! !!! nid_U : 3D ! !!! nid_U : 3D
...@@ -878,10 +882,7 @@ CONTAINS ...@@ -878,10 +882,7 @@ CONTAINS
CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current" , "m/s" , & ! usd CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current" , "m/s" , & ! usd
& jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) & jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )
ENDIF ENDIF
! !!! nid_U : 2D !
CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau
& jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
CALL histend( nid_U, snc4chunks=snc4set ) CALL histend( nid_U, snc4chunks=snc4set )
! !!! nid_V : 3D ! !!! nid_V : 3D
...@@ -891,10 +892,7 @@ CONTAINS ...@@ -891,10 +892,7 @@ CONTAINS
CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current" , "m/s" , & ! vsd CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current" , "m/s" , & ! vsd
& jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) & jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )
ENDIF ENDIF
! !!! nid_V : 2D !
CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau
& jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout )
CALL histend( nid_V, snc4chunks=snc4set ) CALL histend( nid_V, snc4chunks=snc4set )
! !!! nid_W : 3D ! !!! nid_W : 3D
...@@ -1066,12 +1064,12 @@ CONTAINS ...@@ -1066,12 +1064,12 @@ CONTAINS
CALL histwrite( nid_T, "so28chgt", it, hd28 , ndim_hT, ndex_hT ) ! depth of the 28 isotherm CALL histwrite( nid_T, "so28chgt", it, hd28 , ndim_hT, ndex_hT ) ! depth of the 28 isotherm
CALL histwrite( nid_T, "sohtc300", it, htc3 , ndim_hT, ndex_hT ) ! first 300m heaat content CALL histwrite( nid_T, "sohtc300", it, htc3 , ndim_hT, ndex_hT ) ! first 300m heaat content
#endif #endif
CALL histwrite( nid_T, "sozotaux", it, utau , ndim_hT, ndex_hT ) ! i-wind stress
CALL histwrite( nid_T, "sometauy", it, vtau , ndim_hT, ndex_hT ) ! j-wind stress
CALL histwrite( nid_U, "vozocrtx", it, uu(:,:,:,Kmm) , ndim_U , ndex_U ) ! i-current CALL histwrite( nid_U, "vozocrtx", it, uu(:,:,:,Kmm) , ndim_U , ndex_U ) ! i-current
CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress
CALL histwrite( nid_V, "vomecrty", it, vv(:,:,:,Kmm) , ndim_V , ndex_V ) ! j-current CALL histwrite( nid_V, "vomecrty", it, vv(:,:,:,Kmm) , ndim_V , ndex_V ) ! j-current
CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress
IF( ln_zad_Aimp ) THEN IF( ln_zad_Aimp ) THEN
DO_3D( 0, 0, 0, 0, 1, jpk ) DO_3D( 0, 0, 0, 0, 1, jpk )
......
src/OCE/DYN/dynadv.F90
View file @ fc409945
...@@ -7,6 +7,7 @@ MODULE dynadv ...@@ -7,6 +7,7 @@ MODULE dynadv
!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase !! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
!! 3.6 ! 2015-05 (N. Ducousso, G. Madec) add Hollingsworth scheme as an option !! 3.6 ! 2015-05 (N. Ducousso, G. Madec) add Hollingsworth scheme as an option
!! 4.0 ! 2017-07 (G. Madec) add a linear dynamics option !! 4.0 ! 2017-07 (G. Madec) add a linear dynamics option
!! 4.5 ! 2022-06 (S. Techene, G, Madec) refactorization to reduce local memory usage
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
...@@ -50,7 +51,7 @@ MODULE dynadv ...@@ -50,7 +51,7 @@ MODULE dynadv
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
SUBROUTINE dyn_adv( kt, Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) SUBROUTINE dyn_adv( kt, Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw )
!!--------------------------------------------------------------------- !!---------------------------------------------------------------------
!! *** ROUTINE dyn_adv *** !! *** ROUTINE dyn_adv ***
!! !!
...@@ -64,7 +65,6 @@ CONTAINS ...@@ -64,7 +65,6 @@ CONTAINS
!! (see dynvor module). !! (see dynvor module).
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt , Kbb, Kmm, Krhs ! ocean time step and level indices INTEGER , INTENT(in ) :: kt , Kbb, Kmm, Krhs ! ocean time step and level indices
INTEGER , OPTIONAL , INTENT(in ) :: no_zad ! no vertical advection compotation
REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity
REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum Eq. REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum Eq.
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
...@@ -73,12 +73,23 @@ CONTAINS ...@@ -73,12 +73,23 @@ CONTAINS
! !
SELECT CASE( n_dynadv ) !== compute advection trend and add it to general trend ==! SELECT CASE( n_dynadv ) !== compute advection trend and add it to general trend ==!
CASE( np_VEC_c2 ) != vector form =! CASE( np_VEC_c2 ) != vector form =!
CALL dyn_keg ( kt, nn_dynkeg , Kmm, puu, pvv, Krhs ) ! horizontal gradient of kinetic energy ! !* horizontal gradient of kinetic energy
CALL dyn_zad ( kt , Kmm, puu, pvv, Krhs ) ! vertical advection IF (nn_hls==1) THEN ! halo 1 case
CALL dyn_keg_hls1( kt, nn_dynkeg , Kmm, puu, pvv, Krhs ) ! lbc needed with Hollingsworth scheme
ELSE ! halo 2 case
CALL dyn_keg ( kt, nn_dynkeg , Kmm, puu, pvv, Krhs )
ENDIF
CALL dyn_zad ( kt , Kmm, puu, pvv, Krhs ) !* vertical advection
!
CASE( np_FLX_c2 ) != flux form =! CASE( np_FLX_c2 ) != flux form =!
CALL dyn_adv_cen2( kt , Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) ! 2nd order centered scheme CALL dyn_adv_cen2( kt , Kmm, puu, pvv, Krhs, pau, pav, paw ) !* 2nd order centered scheme
CASE( np_FLX_ubs ) !
CALL dyn_adv_ubs ( kt , Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) ! 3rd order UBS scheme (UP3) CASE( np_FLX_ubs ) !* 3rd order UBS scheme (UP3)
IF (nn_hls==1) THEN ! halo 1 case
CALL dyn_adv_ubs_hls1( kt , Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw )
ELSE ! halo 2 case
CALL dyn_adv_ubs ( kt , Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw )
ENDIF
END SELECT END SELECT
! !
IF( ln_timing ) CALL timing_stop( 'dyn_adv' ) IF( ln_timing ) CALL timing_stop( 'dyn_adv' )
......
src/OCE/DYN/dynadv_cen2.F90
View file @ fc409945
...@@ -6,6 +6,7 @@ MODULE dynadv_cen2 ...@@ -6,6 +6,7 @@ MODULE dynadv_cen2
!!====================================================================== !!======================================================================
!! History : 2.0 ! 2006-08 (G. Madec, S. Theetten) Original code !! History : 2.0 ! 2006-08 (G. Madec, S. Theetten) Original code
!! 3.2 ! 2009-07 (R. Benshila) Suppression of rigid-lid option !! 3.2 ! 2009-07 (R. Benshila) Suppression of rigid-lid option
!! 4.5 ! 2022-06 (S. Techene, G, Madec) refactorization to reduce local memory usage
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
...@@ -35,7 +36,7 @@ MODULE dynadv_cen2 ...@@ -35,7 +36,7 @@ MODULE dynadv_cen2
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
SUBROUTINE dyn_adv_cen2( kt, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) SUBROUTINE dyn_adv_cen2( kt, Kmm, puu, pvv, Krhs, pau, pav, paw )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE dyn_adv_cen2 *** !! *** ROUTINE dyn_adv_cen2 ***
!! !!
...@@ -51,15 +52,17 @@ CONTAINS ...@@ -51,15 +52,17 @@ CONTAINS
!! ** Action : (puu,pvv)(:,:,:,Krhs) updated with the advective trend !! ** Action : (puu,pvv)(:,:,:,Krhs) updated with the advective trend
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt , Kmm, Krhs ! ocean time-step and level indices INTEGER , INTENT(in ) :: kt , Kmm, Krhs ! ocean time-step and level indices
INTEGER , OPTIONAL , INTENT(in ) :: no_zad ! no vertical advection computation
REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation
REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity
! !
INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zzu, zzv ! local scalars REAL(wp) :: zzu, zzfu_kp1 ! local scalars
REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zfu_t, zfu_f, zfu_uw, zfu REAL(wp) :: zzv, zzfv_kp1 ! - -
REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zfv_t, zfv_f, zfv_vw, zfv, zfw REAL(wp), DIMENSION(A2D(1)) :: zfu_t, zfu_f, zfu
REAL(wp), DIMENSION(A2D(1)) :: zfv_t, zfv_f, zfv
REAL(wp), DIMENSION(A2D(1)) :: zfu_uw, zfv_vw, zfw
REAL(wp), DIMENSION(:,:,:) , POINTER :: zpt_u, zpt_v, zpt_w REAL(wp), DIMENSION(:,:,:) , POINTER :: zpt_u, zpt_v, zpt_w
REAL(wp), DIMENSION(:,:,:) , ALLOCATABLE :: zu_trd, zv_trd
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile
...@@ -71,8 +74,9 @@ CONTAINS ...@@ -71,8 +74,9 @@ CONTAINS
ENDIF ENDIF
! !
IF( l_trddyn ) THEN ! trends: store the input trends IF( l_trddyn ) THEN ! trends: store the input trends
zfu_uw(:,:,:) = puu(:,:,:,Krhs) ALLOCATE( zu_trd(A2D(0),jpkm1), zv_trd(A2D(0),jpkm1) )
zfv_vw(:,:,:) = pvv(:,:,:,Krhs) zu_trd(A2D(0),:) = puu(A2D(0),:,Krhs)
zv_trd(A2D(0),:) = pvv(A2D(0),:,Krhs)
ENDIF ENDIF
! !
IF( PRESENT( pau ) ) THEN ! RK3: advective velocity (pau,pav,paw) /= advected velocity (puu,pvv,ww) IF( PRESENT( pau ) ) THEN ! RK3: advective velocity (pau,pav,paw) /= advected velocity (puu,pvv,ww)
...@@ -89,84 +93,86 @@ CONTAINS ...@@ -89,84 +93,86 @@ CONTAINS
! !
DO jk = 1, jpkm1 ! horizontal transport DO jk = 1, jpkm1 ! horizontal transport
DO_2D( 1, 1, 1, 1 ) DO_2D( 1, 1, 1, 1 )
zfu(ji,jj,jk) = 0.25_wp * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * zpt_u(ji,jj,jk) zfu(ji,jj) = 0.25_wp * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * zpt_u(ji,jj,jk)
zfv(ji,jj,jk) = 0.25_wp * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * zpt_v(ji,jj,jk) zfv(ji,jj) = 0.25_wp * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * zpt_v(ji,jj,jk)
END_2D END_2D
DO_2D( 1, 0, 1, 0 ) ! horizontal momentum fluxes (at T- and F-point) DO_2D( 1, 0, 1, 0 ) ! horizontal momentum fluxes (at T- and F-point)
zfu_t(ji+1,jj ,jk) = ( zfu(ji,jj,jk) + zfu(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji+1,jj ,jk,Kmm) ) zfu_t(ji+1,jj ) = ( zfu(ji,jj) + zfu(ji+1,jj) ) * ( puu(ji,jj,jk,Kmm) + puu(ji+1,jj ,jk,Kmm) )
zfv_f(ji ,jj ,jk) = ( zfv(ji,jj,jk) + zfv(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji ,jj+1,jk,Kmm) ) zfv_f(ji ,jj ) = ( zfv(ji,jj) + zfv(ji+1,jj) ) * ( puu(ji,jj,jk,Kmm) + puu(ji ,jj+1,jk,Kmm) )
zfu_f(ji ,jj ,jk) = ( zfu(ji,jj,jk) + zfu(ji,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji+1,jj ,jk,Kmm) ) zfu_f(ji ,jj ) = ( zfu(ji,jj) + zfu(ji,jj+1) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji+1,jj ,jk,Kmm) )
zfv_t(ji ,jj+1,jk) = ( zfv(ji,jj,jk) + zfv(ji,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji ,jj+1,jk,Kmm) ) zfv_t(ji ,jj+1) = ( zfv(ji,jj) + zfv(ji,jj+1) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji ,jj+1,jk,Kmm) )
END_2D END_2D
DO_2D( 0, 0, 0, 0 ) ! divergence of horizontal momentum fluxes DO_2D( 0, 0, 0, 0 ) ! divergence of horizontal momentum fluxes
puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_t(ji+1,jj,jk) - zfu_t(ji,jj ,jk) & puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_t(ji+1,jj) - zfu_t(ji,jj ) &
& + zfv_f(ji ,jj,jk) - zfv_f(ji,jj-1,jk) ) * r1_e1e2u(ji,jj) & & + zfv_f(ji ,jj) - zfv_f(ji,jj-1) ) * r1_e1e2u(ji,jj) &
& / e3u(ji,jj,jk,Kmm) & / e3u(ji,jj,jk,Kmm)
pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfu_f(ji,jj ,jk) - zfu_f(ji-1,jj,jk) & pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfu_f(ji,jj ) - zfu_f(ji-1,jj) &
& + zfv_t(ji,jj+1,jk) - zfv_t(ji ,jj,jk) ) * r1_e1e2v(ji,jj) & & + zfv_t(ji,jj+1) - zfv_t(ji ,jj) ) * r1_e1e2v(ji,jj) &
& / e3v(ji,jj,jk,Kmm) & / e3v(ji,jj,jk,Kmm)
END_2D END_2D
END DO END DO
! !
IF( l_trddyn ) THEN ! trends: send trend to trddyn for diagnostic IF( l_trddyn ) THEN ! trends: send trend to trddyn for diagnostic
zfu_uw(:,:,:) = puu(:,:,:,Krhs) - zfu_uw(:,:,:) zu_trd(A2D(0),:) = puu(A2D(0),:,Krhs) - zu_trd(A2D(0),:)
zfv_vw(:,:,:) = pvv(:,:,:,Krhs) - zfv_vw(:,:,:) zv_trd(A2D(0),:) = pvv(A2D(0),:,Krhs) - zv_trd(A2D(0),:)
CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt, Kmm ) CALL trd_dyn( zu_trd, zv_trd, jpdyn_keg, kt, Kmm )
zfu_t(:,:,:) = puu(:,:,:,Krhs) zu_trd(A2D(0),:) = puu(A2D(0),:,Krhs)
zfv_t(:,:,:) = pvv(:,:,:,Krhs) zv_trd(A2D(0),:) = pvv(A2D(0),:,Krhs)
ENDIF ENDIF
! !
IF( PRESENT( no_zad ) ) THEN !== No vertical advection ==! (except if linear free surface) ! !== Vertical advection ==!
! == !
IF( ln_linssh ) THEN ! linear free surface: advection through the surface z=0 ! ! surface vertical fluxes
DO_2D( 0, 0, 0, 0 ) !
zzu = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) IF( ln_linssh ) THEN ! linear free surface: advection through the surface z=0
zzv = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) DO_2D( 0, 0, 0, 0 )
puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) - zzu * r1_e1e2u(ji,jj) & zfu_uw(ji,jj) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm)
& / e3u(ji,jj,1,Kmm) zfv_vw(ji,jj) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm)
pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) - zzv * r1_e1e2v(ji,jj) & END_2D
& / e3v(ji,jj,1,Kmm) ELSE ! non linear free: surface advective fluxes set to zero
END_2D DO_2D( 0, 0, 0, 0 )
ENDIF zfu_uw(ji,jj) = 0._wp
! zfv_vw(ji,jj) = 0._wp
ELSE !== Vertical advection ==! END_2D
ENDIF
!
DO jk = 1, jpk-2 ! divergence of advective fluxes
! !
DO_2D( 0, 0, 0, 0 ) ! surface/bottom advective fluxes set to zero DO_2D( 0, 1, 0, 1 ) ! 1/4 * Vertical transport at level k+1
zfu_uw(ji,jj,jpk) = 0._wp ; zfv_vw(ji,jj,jpk) = 0._wp zfw(ji,jj) = 0.25_wp * e1e2t(ji,jj) * zpt_w(ji,jj,jk+1)
zfu_uw(ji,jj, 1 ) = 0._wp ; zfv_vw(ji,jj, 1 ) = 0._wp
END_2D END_2D
IF( ln_linssh ) THEN ! linear free surface: advection through the surface z=0 DO_2D( 0, 0, 0, 0 )
DO_2D( 0, 0, 0, 0 ) ! ! vertical flux at level k+1
zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) zzfu_kp1 = ( zfw(ji,jj) + zfw(ji+1,jj ) ) * ( puu(ji,jj,jk+1,Kmm) + puu(ji,jj,jk,Kmm) )
zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) zzfv_kp1 = ( zfw(ji,jj) + zfw(ji ,jj+1) ) * ( pvv(ji,jj,jk+1,Kmm) + pvv(ji,jj,jk,Kmm) )
END_2D ! ! divergence of vertical momentum flux
ENDIF puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj) - zzfu_kp1 ) * r1_e1e2u(ji,jj) &
DO jk = 2, jpkm1 ! interior advective fluxes
DO_2D( 0, 1, 0, 1 ) ! 1/4 * Vertical transport
zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * zpt_w(ji,jj,jk)
END_2D
DO_2D( 0, 0, 0, 0 )
zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji+1,jj ,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) )
zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji ,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) )
END_2D
END DO
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! divergence of vertical momentum flux divergence
puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) &
& / e3u(ji,jj,jk,Kmm) & / e3u(ji,jj,jk,Kmm)
pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) & pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj) - zzfv_kp1 ) * r1_e1e2v(ji,jj) &
& / e3v(ji,jj,jk,Kmm) & / e3v(ji,jj,jk,Kmm)
END_3D ! ! store vertical flux for next level calculation
! zfu_uw(ji,jj) = zzfu_kp1
IF( l_trddyn ) THEN ! trends: send trend to trddyn for diagnostic zfv_vw(ji,jj) = zzfv_kp1
zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:) END_2D
zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:) END DO
CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm ) !
ENDIF jk = jpkm1
! ! Control print DO_2D( 0, 0, 0, 0 )
IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' cen2 adv - Ua: ', mask1=umask, & puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - zfu_uw(ji,jj) * r1_e1e2u(ji,jj) &
& tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) & / e3u(ji,jj,jk,Kmm)
! pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - zfv_vw(ji,jj) * r1_e1e2v(ji,jj) &
& / e3v(ji,jj,jk,Kmm)
END_2D
!
IF( l_trddyn ) THEN ! trends: send trend to trddyn for diagnostic
zu_trd(A2D(0),:) = puu(A2D(0),:,Krhs) - zu_trd(A2D(0),:)
zv_trd(A2D(0),:) = pvv(A2D(0),:,Krhs) - zv_trd(A2D(0),:)
CALL trd_dyn( zu_trd, zv_trd, jpdyn_zad, kt, Kmm )
DEALLOCATE( zu_trd, zv_trd )
ENDIF ENDIF
! ! Control print
IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' cen2 adv - Ua: ', mask1=umask, &
& tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
! !
END SUBROUTINE dyn_adv_cen2 END SUBROUTINE dyn_adv_cen2
......
src/OCE/DYN/dynadv_ubs.F90
View file @ fc409945
This diff is collapsed.
src/OCE/DYN/dynatf.F90
View file @ fc409945
...@@ -331,7 +331,7 @@ CONTAINS ...@@ -331,7 +331,7 @@ CONTAINS
ALLOCATE(zutau(jpi,jpj)) ALLOCATE(zutau(jpi,jpj))
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 0, 0, 0 )
jk = miku(ji,jj) jk = miku(ji,jj)
zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * puu(ji,jj,jk,Kaa) zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( puu(ji-1,jj,jk,Kaa) + puu(ji,jj,jk,Kaa) )
END_2D END_2D
CALL iom_put( "utau", zutau(:,:) ) CALL iom_put( "utau", zutau(:,:) )
DEALLOCATE(zutau) DEALLOCATE(zutau)
...@@ -345,7 +345,7 @@ CONTAINS ...@@ -345,7 +345,7 @@ CONTAINS
ALLOCATE(zvtau(jpi,jpj)) ALLOCATE(zvtau(jpi,jpj))
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 0, 0, 0 )
jk = mikv(ji,jj) jk = mikv(ji,jj)
zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * pvv(ji,jj,jk,Kaa) zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( pvv(ji,jj-1,jk,Kaa) + pvv(ji,jj,jk,Kaa) )
END_2D END_2D
CALL iom_put( "vtau", zvtau(:,:) ) CALL iom_put( "vtau", zvtau(:,:) )
DEALLOCATE(zvtau) DEALLOCATE(zvtau)
......
src/OCE/DYN/dynatf_qco.F90
View file @ fc409945
...@@ -248,7 +248,7 @@ CONTAINS ...@@ -248,7 +248,7 @@ CONTAINS
ALLOCATE(zutau(jpi,jpj)) ALLOCATE(zutau(jpi,jpj))
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 0, 0, 0 )
jk = miku(ji,jj) jk = miku(ji,jj)
zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * puu(ji,jj,jk,Kaa) zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( puu(ji-1,jj,jk,Kaa) + puu(ji,jj,jk,Kaa) )
END_2D END_2D
CALL iom_put( "utau", zutau(:,:) ) CALL iom_put( "utau", zutau(:,:) )
DEALLOCATE(zutau) DEALLOCATE(zutau)
...@@ -262,7 +262,7 @@ CONTAINS ...@@ -262,7 +262,7 @@ CONTAINS
ALLOCATE(zvtau(jpi,jpj)) ALLOCATE(zvtau(jpi,jpj))
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 0, 0, 0 )
jk = mikv(ji,jj) jk = mikv(ji,jj)
zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * pvv(ji,jj,jk,Kaa) zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( pvv(ji,jj-1,jk,Kaa) + pvv(ji,jj,jk,Kaa) )
END_2D END_2D
CALL iom_put( "vtau", zvtau(:,:) ) CALL iom_put( "vtau", zvtau(:,:) )
DEALLOCATE(zvtau) DEALLOCATE(zvtau)
......
src/OCE/DYN/dynkeg.F90
View file @ fc409945
...@@ -6,7 +6,8 @@ MODULE dynkeg ...@@ -6,7 +6,8 @@ MODULE dynkeg
!! History : 1.0 ! 1987-09 (P. Andrich, M.-A. Foujols) Original code !! History : 1.0 ! 1987-09 (P. Andrich, M.-A. Foujols) Original code
!! 7.0 ! 1997-05 (G. Madec) Split dynber into dynkeg and dynhpg !! 7.0 ! 1997-05 (G. Madec) Split dynber into dynkeg and dynhpg
!! NEMO 1.0 ! 2002-07 (G. Madec) F90: Free form and module !! NEMO 1.0 ! 2002-07 (G. Madec) F90: Free form and module
!! 3.6 ! 2015-05 (N. Ducousso, G. Madec) add Hollingsworth scheme as an option !! 3.6 ! 2015-05 (N. Ducousso, G. Madec) add Hollingsworth scheme as an option
!! 4.5 ! 2022-06 (S. Techene, G, Madec) refactorization to reduce local memory usage
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
...@@ -27,7 +28,8 @@ MODULE dynkeg ...@@ -27,7 +28,8 @@ MODULE dynkeg
IMPLICIT NONE IMPLICIT NONE
PRIVATE PRIVATE
PUBLIC dyn_keg ! routine called by step module PUBLIC dyn_keg ! routine called by step module
PUBLIC dyn_keg_hls1 ! routine called by step module
INTEGER, PARAMETER, PUBLIC :: nkeg_C2 = 0 !: 2nd order centered scheme (standard scheme) INTEGER, PARAMETER, PUBLIC :: nkeg_C2 = 0 !: 2nd order centered scheme (standard scheme)
INTEGER, PARAMETER, PUBLIC :: nkeg_HW = 1 !: Hollingsworth et al., QJRMS, 1983 INTEGER, PARAMETER, PUBLIC :: nkeg_HW = 1 !: Hollingsworth et al., QJRMS, 1983
...@@ -44,6 +46,123 @@ MODULE dynkeg ...@@ -44,6 +46,123 @@ MODULE dynkeg
CONTAINS CONTAINS
SUBROUTINE dyn_keg( kt, kscheme, Kmm, puu, pvv, Krhs ) SUBROUTINE dyn_keg( kt, kscheme, Kmm, puu, pvv, Krhs )
!!----------------------------------------------------------------------
!! *** ROUTINE dyn_keg ***
!!
!! ** Purpose : Compute the now momentum trend due to the horizontal
!! gradient of the horizontal kinetic energy and add it to the
!! general momentum trend.
!!
!! ** Method : * kscheme = nkeg_C2 : 2nd order centered scheme that
!! conserve kinetic energy. Compute the now horizontal kinetic energy
!! zhke = 1/2 [ mi-1( un^2 ) + mj-1( vn^2 ) ]
!! * kscheme = nkeg_HW : Hollingsworth correction following
!! Arakawa (2001). The now horizontal kinetic energy is given by:
!! zhke = 1/6 [ mi-1( 2 * un^2 + ((u(j+1)+u(j-1))/2)^2 )
!! + mj-1( 2 * vn^2 + ((v(i+1)+v(i-1))/2)^2 ) ]
!!
!! Take its horizontal gradient and add it to the general momentum
!! trend.
!! u(rhs) = u(rhs) - 1/e1u di[ zhke ]
!! v(rhs) = v(rhs) - 1/e2v dj[ zhke ]
!!
!! ** Action : - Update the (puu(:,:,:,Krhs), pvv(:,:,:,Krhs)) with the hor. ke gradient trend
!! - send this trends to trd_dyn (l_trddyn=T) for post-processing
!!
!! ** References : Arakawa, A., International Geophysics 2001.
!! Hollingsworth et al., Quart. J. Roy. Meteor. Soc., 1983.
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: kscheme ! =0/1 type of KEG scheme
INTEGER , INTENT(in ) :: Kmm, Krhs ! ocean time level indices
REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation
!
INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zu, zv ! local scalars
REAL(wp), DIMENSION(:,: ) , ALLOCATABLE :: zhke
REAL(wp), DIMENSION(:,:,:) , ALLOCATABLE :: zu_trd, zv_trd
!!----------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('dyn_keg')
!
IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile
IF( kt == nit000 ) THEN
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) 'dyn_keg : kinetic energy gradient trend, scheme number=', kscheme
IF(lwp) WRITE(numout,*) '~~~~~~~'
ENDIF
ENDIF
!
IF( l_trddyn ) THEN ! Save the input trends
ALLOCATE( zu_trd(A2D(0),jpk), zv_trd(A2D(0),jpk) )
zu_trd(A2D(0),:) = puu(A2D(0),:,Krhs)
zv_trd(A2D(0),:) = pvv(A2D(0),:,Krhs)
ENDIF
!
SELECT CASE ( kscheme )
!
CASE ( nkeg_C2 ) !== Standard scheme ==!
ALLOCATE( zhke(A2D(1)) )
DO jk = 1, jpkm1
DO_2D( 0, 1, 0, 1 ) !* Horizontal kinetic energy at T-point
zu = puu(ji-1,jj ,jk,Kmm) * puu(ji-1,jj ,jk,Kmm) &
& + puu(ji ,jj ,jk,Kmm) * puu(ji ,jj ,jk,Kmm)
zv = pvv(ji ,jj-1,jk,Kmm) * pvv(ji ,jj-1,jk,Kmm) &
& + pvv(ji ,jj ,jk,Kmm) * pvv(ji ,jj ,jk,Kmm)
zhke(ji,jj) = 0.25_wp * ( zv + zu )
END_2D
!
DO_2D( 0, 0, 0, 0 ) !* grad( KE ) added to the general momentum trends
puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zhke(ji+1,jj ) - zhke(ji,jj) ) * r1_e1u(ji,jj)
pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zhke(ji ,jj+1) - zhke(ji,jj) ) * r1_e2v(ji,jj)
END_2D
END DO
DEALLOCATE( zhke )
!
CASE ( nkeg_HW ) !* Hollingsworth scheme
ALLOCATE( zhke(A2D(1)) )
DO jk = 1, jpkm1
DO_2D( 0, 1, 0, 1 )
! round brackets added to fix the order of floating point operations
! needed to ensure halo 1 - halo 2 compatibility
zu = ( puu(ji-1,jj ,jk,Kmm) * puu(ji-1,jj ,jk,Kmm) &
& + puu(ji ,jj ,jk,Kmm) * puu(ji ,jj ,jk,Kmm) ) * 8._wp &
& + ( ( puu(ji-1,jj-1,jk,Kmm) + puu(ji-1,jj+1,jk,Kmm) ) * ( puu(ji-1,jj-1,jk,Kmm) + puu(ji-1,jj+1,jk,Kmm) ) &
& + ( puu(ji ,jj-1,jk,Kmm) + puu(ji ,jj+1,jk,Kmm) ) * ( puu(ji ,jj-1,jk,Kmm) + puu(ji ,jj+1,jk,Kmm) ) &
& ) ! bracket for halo 1 - halo 2 compatibility
zv = ( pvv(ji ,jj-1,jk,Kmm) * pvv(ji ,jj-1,jk,Kmm) &
& + pvv(ji ,jj ,jk,Kmm) * pvv(ji ,jj ,jk,Kmm) ) * 8._wp &
& + ( ( pvv(ji-1,jj-1,jk,Kmm) + pvv(ji+1,jj-1,jk,Kmm) ) * ( pvv(ji-1,jj-1,jk,Kmm) + pvv(ji+1,jj-1,jk,Kmm) ) &
& + ( pvv(ji-1,jj ,jk,Kmm) + pvv(ji+1,jj ,jk,Kmm) ) * ( pvv(ji-1,jj ,jk,Kmm) + pvv(ji+1,jj ,jk,Kmm) ) &
& ) ! bracket for halo 1 - halo 2 compatibility
zhke(ji,jj) = r1_48 * ( zv + zu )
END_2D
!
DO_2D( 0, 0, 0, 0 ) !* grad( KE ) added to the general momentum trends
puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zhke(ji+1,jj ) - zhke(ji,jj) ) * r1_e1u(ji,jj)
pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zhke(ji ,jj+1) - zhke(ji,jj) ) * r1_e2v(ji,jj)
END_2D
END DO
DEALLOCATE( zhke )
!
END SELECT
!
IF( l_trddyn ) THEN ! save the Kinetic Energy trends for diagnostic
zu_trd(A2D(0),:) = puu(A2D(0),:,Krhs) - zu_trd(A2D(0),:)
zv_trd(A2D(0),:) = pvv(A2D(0),:,Krhs) - zv_trd(A2D(0),:)
CALL trd_dyn( zu_trd, zv_trd, jpdyn_keg, kt, Kmm )
DEALLOCATE( zu_trd, zv_trd )
ENDIF
!
IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' keg - Ua: ', mask1=umask, &
& tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
!
IF( ln_timing ) CALL timing_stop('dyn_keg')
!
END SUBROUTINE dyn_keg
SUBROUTINE dyn_keg_hls1( kt, kscheme, Kmm, puu, pvv, Krhs )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE dyn_keg *** !! *** ROUTINE dyn_keg ***
!! !!
...@@ -86,7 +205,7 @@ CONTAINS ...@@ -86,7 +205,7 @@ CONTAINS
IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile
IF( kt == nit000 ) THEN IF( kt == nit000 ) THEN
IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) 'dyn_keg : kinetic energy gradient trend, scheme number=', kscheme IF(lwp) WRITE(numout,*) 'dyn_keg_hls1 : kinetic energy gradient trend, scheme number=', kscheme
IF(lwp) WRITE(numout,*) '~~~~~~~' IF(lwp) WRITE(numout,*) '~~~~~~~'
ENDIF ENDIF
ENDIF ENDIF
...@@ -147,7 +266,7 @@ CONTAINS ...@@ -147,7 +266,7 @@ CONTAINS
! !
IF( ln_timing ) CALL timing_stop('dyn_keg') IF( ln_timing ) CALL timing_stop('dyn_keg')
! !
END SUBROUTINE dyn_keg END SUBROUTINE dyn_keg_hls1
!!====================================================================== !!======================================================================
END MODULE dynkeg END MODULE dynkeg
src/OCE/DYN/dynspg_ts.F90
View file @ fc409945
...@@ -334,14 +334,14 @@ CONTAINS ...@@ -334,14 +334,14 @@ CONTAINS
! ! ------------------ ! ! ! ------------------ !
IF( ln_bt_fw ) THEN IF( ln_bt_fw ) THEN
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 0, 0, 0 )
zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utauU(ji,jj) * r1_hu(ji,jj,Kmm)
zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm) zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtauV(ji,jj) * r1_hv(ji,jj,Kmm)
END_2D END_2D
ELSE ELSE
zztmp = r1_rho0 * r1_2 zztmp = r1_rho0 * r1_2
DO_2D( 0, 0, 0, 0 ) DO_2D( 0, 0, 0, 0 )
zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utauU(ji,jj) ) * r1_hu(ji,jj,Kmm)
zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv(ji,jj,Kmm) zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtauV(ji,jj) ) * r1_hv(ji,jj,Kmm)
END_2D END_2D
ENDIF ENDIF
! !
......
src/OCE/DYN/dynvor.F90
View file @ fc409945
This diff is collapsed.
src/OCE/DYN/dynzdf.F90
View file @ fc409945
This diff is collapsed.
src/OCE/IOM/iom.F90
View file @ fc409945
...@@ -1327,14 +1327,21 @@ CONTAINS ...@@ -1327,14 +1327,21 @@ CONTAINS
IF( irankpv == 1 ) ishape(1:1) = SHAPE(pv_r1d) IF( irankpv == 1 ) ishape(1:1) = SHAPE(pv_r1d)
IF( irankpv == 2 ) ishape(1:2) = SHAPE(pv_r2d) IF( irankpv == 2 ) ishape(1:2) = SHAPE(pv_r2d)
IF( irankpv == 3 ) ishape(1:3) = SHAPE(pv_r3d) IF( irankpv == 3 ) ishape(1:3) = SHAPE(pv_r3d)
ix1 = 1 ; ix2 = icnt(1) ; iy1 = 1 ; iy2 = icnt(2) ! index of the array to be read
ctmp1 = 'd' ctmp1 = 'd'
ELSE ELSE
IF( irankpv == 2 ) THEN IF( irankpv == 2 ) ishape(1:2) = SHAPE(pv_r2d)
ishape(1:2) = SHAPE(pv_r2d(Nis0:Nie0,Njs0:Nje0 )) ; ctmp1 = 'd(Nis0:Nie0,Njs0:Nje0)' IF( irankpv == 3 ) ishape(1:3) = SHAPE(pv_r3d)
ENDIF IF( ishape(1) == jpi .AND. ishape(2) == jpj ) THEN
IF( irankpv == 3 ) THEN ishape(1:2) = (/ Ni_0, Nj_0 /)
ishape(1:3) = SHAPE(pv_r3d(Nis0:Nie0,Njs0:Nje0,:)) ; ctmp1 = 'd(Nis0:Nie0,Njs0:Nje0,:)' ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0 ! index of the array to be read
ctmp1 = 'd(Nis0:Nie0,Njs0:Nje0'
ELSE
ix1 = 1 ; ix2 = icnt(1) ; iy1 = 1 ; iy2 = icnt(2) ! index of the array to be read
ctmp1 = 'd(:,:'
ENDIF ENDIF
IF( irankpv == 3 ) ctmp1 = TRIM(ctmp1)//',:'
ctmp1 = TRIM(ctmp1)//')'
ENDIF ENDIF
DO jl = 1, irankpv DO jl = 1, irankpv
WRITE( ctmp2, FMT="(', ', i1,'): ', i5,' /= icnt(', i1,'):', i5)" ) jl, ishape(jl), jl, icnt(jl) WRITE( ctmp2, FMT="(', ', i1,'): ', i5,' /= icnt(', i1,'):', i5)" ) jl, ishape(jl), jl, icnt(jl)
...@@ -1347,11 +1354,6 @@ CONTAINS ...@@ -1347,11 +1354,6 @@ CONTAINS
!- !-
IF( idvar > 0 .AND. istop == nstop ) THEN ! no additional errors until this point... IF( idvar > 0 .AND. istop == nstop ) THEN ! no additional errors until this point...
! !
! find the right index of the array to be read
IF( idom /= jpdom_unknown ) THEN ; ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0
ELSE ; ix1 = 1 ; ix2 = icnt(1) ; iy1 = 1 ; iy2 = icnt(2)
ENDIF
CALL iom_nf90_get( kiomid, idvar, inbdim, istart, icnt, ix1, ix2, iy1, iy2, pv_r1d, pv_r2d, pv_r3d ) CALL iom_nf90_get( kiomid, idvar, inbdim, istart, icnt, ix1, ix2, iy1, iy2, pv_r1d, pv_r2d, pv_r3d )
IF( istop == nstop ) THEN ! no additional errors until this point... IF( istop == nstop ) THEN ! no additional errors until this point...
...@@ -1362,9 +1364,11 @@ CONTAINS ...@@ -1362,9 +1364,11 @@ CONTAINS
zsgn = 1._wp zsgn = 1._wp
IF( PRESENT(psgn ) ) zsgn = psgn IF( PRESENT(psgn ) ) zsgn = psgn
!--- overlap areas and extra hallows (mpp) !--- overlap areas and extra hallows (mpp)
IF( PRESENT(pv_r2d) .AND. idom /= jpdom_unknown .AND. cl_type /= 'Z' ) THEN llok = idom /= jpdom_unknown .AND. cl_type /= 'Z' &
& .AND. ix1 == Nis0 .AND. ix2 == Nie0 .AND. iy1 == Njs0 .AND. iy2 == Nje0
IF( PRESENT(pv_r2d) .AND. llok ) THEN
CALL lbc_lnk( 'iom', pv_r2d, cl_type, zsgn, kfillmode = kfill ) CALL lbc_lnk( 'iom', pv_r2d, cl_type, zsgn, kfillmode = kfill )
ELSEIF( PRESENT(pv_r3d) .AND. idom /= jpdom_unknown .AND. cl_type /= 'Z' ) THEN ELSEIF( PRESENT(pv_r3d) .AND. llok ) THEN
CALL lbc_lnk( 'iom', pv_r3d, cl_type, zsgn, kfillmode = kfill ) CALL lbc_lnk( 'iom', pv_r3d, cl_type, zsgn, kfillmode = kfill )
ENDIF ENDIF
! !
...@@ -2336,14 +2340,14 @@ CONTAINS ...@@ -2336,14 +2340,14 @@ CONTAINS
idb(jn) = -nn_hls ! Tile data offset (halo size) idb(jn) = -nn_hls ! Tile data offset (halo size)
END DO END DO
! Tile_[ij]begin are defined with respect to the processor data domain, so data_[ij]begin is added
CALL iom_set_domain_attr("grid_"//cdgrd, ntiles=nijtile, & CALL iom_set_domain_attr("grid_"//cdgrd, ntiles=nijtile, &
& tile_ibegin=ntsi_a(1:nijtile) + idb(:) - 1, tile_jbegin=ntsj_a(1:nijtile) + idb(:) - 1, & & tile_ibegin=ntsi_a(1:nijtile) - nn_hls - 1, tile_jbegin=ntsj_a(1:nijtile) - nn_hls - 1, &
& tile_ni=ini(:), tile_nj=inj(:), & & tile_ni=ini(:), tile_nj=inj(:), &
& tile_data_ibegin=idb(:), tile_data_jbegin=idb(:), & & tile_data_ibegin=idb(:), tile_data_jbegin=idb(:), &
& tile_data_ni=ini(:) - 2 * idb(:), tile_data_nj=inj(:) - 2 * idb(:)) & tile_data_ni=ini(:) - 2 * idb(:), tile_data_nj=inj(:) - 2 * idb(:))
idb(:) = 0
CALL iom_set_domain_attr("grid_"//cdgrd//"_inner", ntiles=nijtile, & CALL iom_set_domain_attr("grid_"//cdgrd//"_inner", ntiles=nijtile, &
& tile_ibegin=ntsi_a(1:nijtile) + idb(:) - 1, tile_jbegin=ntsj_a(1:nijtile) + idb(:) - 1, & & tile_ibegin=ntsi_a(1:nijtile) - nn_hls - 1, tile_jbegin=ntsj_a(1:nijtile) - nn_hls - 1, &
& tile_ni=ini(:), tile_nj=inj(:), & & tile_ni=ini(:), tile_nj=inj(:), &
& tile_data_ibegin=idb(:), tile_data_jbegin=idb(:), & & tile_data_ibegin=idb(:), tile_data_jbegin=idb(:), &
& tile_data_ni=ini(:) - 2 * idb(:), tile_data_nj=inj(:) - 2 * idb(:)) & tile_data_ni=ini(:) - 2 * idb(:), tile_data_nj=inj(:) - 2 * idb(:))
...@@ -2453,7 +2457,7 @@ CONTAINS ...@@ -2453,7 +2457,7 @@ CONTAINS
! CALL dom_ngb( -168.53_wp, 65.03_wp, ix, iy, 'T' ) ! i-line that passes through Bering Strait: Reference latitude (used in plots) ! CALL dom_ngb( -168.53_wp, 65.03_wp, ix, iy, 'T' ) ! i-line that passes through Bering Strait: Reference latitude (used in plots)
CALL dom_ngb( 180.0_wp, 90.0_wp, ix, iy, 'T' ) ! i-line that passes near the North Pole : Reference latitude (used in plots) CALL dom_ngb( 180.0_wp, 90.0_wp, ix, iy, 'T' ) ! i-line that passes near the North Pole : Reference latitude (used in plots)
CALL iom_set_domain_attr("gznl", ni_glo=Ni0glo, nj_glo=Nj0glo, ibegin=mig0(Nis0)-1, jbegin=mjg0(Njs0)-1, ni=Ni_0, nj=Nj_0) CALL iom_set_domain_attr("gznl", ni_glo=Ni0glo, nj_glo=Nj0glo, ibegin=mig0(Nis0)-1, jbegin=mjg0(Njs0)-1, ni=Ni_0, nj=Nj_0)
CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin = -nn_hls, data_ni = jpi, data_jbegin = -nn_hls, data_nj = jpj) CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin=0, data_ni=Ni_0, data_jbegin=0, data_nj=Nj_0)
CALL iom_set_domain_attr("gznl", lonvalue = real(zlon, dp), & CALL iom_set_domain_attr("gznl", lonvalue = real(zlon, dp), &
& latvalue = real(RESHAPE(plat(Nis0:Nie0, Njs0:Nje0),(/ Ni_0*Nj_0 /)),dp)) & latvalue = real(RESHAPE(plat(Nis0:Nie0, Njs0:Nje0),(/ Ni_0*Nj_0 /)),dp))
CALL iom_set_zoom_domain_attr("ptr", ibegin=ix-1, jbegin=0, ni=1, nj=Nj0glo) CALL iom_set_zoom_domain_attr("ptr", ibegin=ix-1, jbegin=0, ni=1, nj=Nj0glo)
......
src/OCE/IOM/iom_nf90.F90
View file @ fc409945
...@@ -40,6 +40,8 @@ MODULE iom_nf90 ...@@ -40,6 +40,8 @@ MODULE iom_nf90
MODULE PROCEDURE iom_nf90_rp0123d_dp MODULE PROCEDURE iom_nf90_rp0123d_dp
END INTERFACE END INTERFACE
!! * Substitutions
# include "do_loop_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: iom_nf90.F90 14433 2021-02-11 08:06:49Z smasson $ !! $Id: iom_nf90.F90 14433 2021-02-11 08:06:49Z smasson $
...@@ -544,7 +546,7 @@ CONTAINS ...@@ -544,7 +546,7 @@ CONTAINS
INTEGER :: idvar ! variable id INTEGER :: idvar ! variable id
INTEGER :: jd ! dimension loop counter INTEGER :: jd ! dimension loop counter
INTEGER :: ix1, ix2, iy1, iy2 ! subdomain indexes INTEGER :: ix1, ix2, iy1, iy2 ! subdomain indexes
INTEGER, DIMENSION(4) :: idimsz ! dimensions size INTEGER, DIMENSION(3) :: ishape ! dimensions size
INTEGER, DIMENSION(4) :: idimid ! dimensions id INTEGER, DIMENSION(4) :: idimid ! dimensions id
CHARACTER(LEN=256) :: clinfo ! info character CHARACTER(LEN=256) :: clinfo ! info character
INTEGER :: if90id ! nf90 file identifier INTEGER :: if90id ! nf90 file identifier
...@@ -627,11 +629,9 @@ CONTAINS ...@@ -627,11 +629,9 @@ CONTAINS
itype = NF90_DOUBLE itype = NF90_DOUBLE
ENDIF ENDIF
IF( PRESENT(pv_r0d) ) THEN IF( PRESENT(pv_r0d) ) THEN
CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, & CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, iom_file(kiomid)%nvid(idvar) ), clinfo )
& iom_file(kiomid)%nvid(idvar) ), clinfo )
ELSE ELSE
CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, idimid(1:idims), & CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, idimid(1:idims), iom_file(kiomid)%nvid(idvar) ), clinfo )
& iom_file(kiomid)%nvid(idvar) ), clinfo )
ENDIF ENDIF
lchunk = .false. lchunk = .false.
IF( snc4set%luse .AND. idims == 4 ) lchunk = .true. IF( snc4set%luse .AND. idims == 4 ) lchunk = .true.
...@@ -673,23 +673,13 @@ CONTAINS ...@@ -673,23 +673,13 @@ CONTAINS
ENDIF ENDIF
! on what kind of domain must the data be written? ! on what kind of domain must the data be written?
IF( PRESENT(pv_r2d) .OR. PRESENT(pv_r3d) ) THEN IF( PRESENT(pv_r2d) .OR. PRESENT(pv_r3d) ) THEN
idimsz(1:2) = iom_file(kiomid)%dimsz(1:2,idvar)
IF( idimsz(1) == Ni_0 .AND. idimsz(2) == Nj_0 ) THEN
ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0
ELSEIF( idimsz(1) == jpi .AND. idimsz(2) == jpj ) THEN
ix1 = 1 ; ix2 = jpi ; iy1 = 1 ; iy2 = jpj
ELSEIF( idimsz(1) == jpi .AND. idimsz(2) == jpj ) THEN
ix1 = 1 ; ix2 = jpi ; iy1 = 1 ; iy2 = jpj
ELSE
CALL ctl_stop( 'iom_nf90_rp0123d: should have been an impossible case...' )
ENDIF
! write dimension variables if it is not already done ! write dimension variables if it is not already done
! ============= ! =============
! trick: is defined to 0 => dimension variable are defined but not yet written ! trick: is defined to 0 => dimension variable are defined but not yet written
IF( iom_file(kiomid)%dimsz(1, 4) == 0 ) THEN ! time_counter = 0 IF( iom_file(kiomid)%dimsz(1, 4) == 0 ) THEN ! time_counter = 0
CALL iom_nf90_check( NF90_PUT_VAR( if90id, 1, glamt(ix1:ix2, iy1:iy2) ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, 1, glamt(A2D(0)) ), clinfo )
CALL iom_nf90_check( NF90_PUT_VAR( if90id, 2, gphit(ix1:ix2, iy1:iy2) ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, 2, gphit(A2D(0)) ), clinfo )
SELECT CASE (iom_file(kiomid)%comp) SELECT CASE (iom_file(kiomid)%comp)
CASE ('OCE') CASE ('OCE')
CALL iom_nf90_check( NF90_PUT_VAR( if90id, 3, gdept_1d ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, 3, gdept_1d ), clinfo )
...@@ -704,6 +694,17 @@ CONTAINS ...@@ -704,6 +694,17 @@ CONTAINS
iom_file(kiomid)%dimsz(1, 4) = 1 ! so we don't enter this IF case any more... iom_file(kiomid)%dimsz(1, 4) = 1 ! so we don't enter this IF case any more...
IF(lwp) WRITE(numout,*) TRIM(clinfo)//' write dimension variables done' IF(lwp) WRITE(numout,*) TRIM(clinfo)//' write dimension variables done'
ENDIF ENDIF
IF( PRESENT(pv_r2d) ) ishape(1:2) = SHAPE(pv_r2d)
IF( PRESENT(pv_r3d) ) ishape(1:3) = SHAPE(pv_r3d)
IF( ishape(1) == Ni_0 .AND. ishape(2) == Nj_0 ) THEN
ix1 = 1 ; ix2 = Ni_0 ; iy1 = 1 ; iy2 = Nj_0
ELSEIF( ishape(1) == jpi .AND. ishape(2) == jpj ) THEN
ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0
ELSE
CALL ctl_stop( 'iom_nf90_rp0123d: should have been an impossible case...' )
ENDIF
ENDIF ENDIF
! write the data ! write the data
...@@ -712,7 +713,7 @@ CONTAINS ...@@ -712,7 +713,7 @@ CONTAINS
CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r0d ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r0d ), clinfo )
ELSEIF( PRESENT(pv_r1d) ) THEN ELSEIF( PRESENT(pv_r1d) ) THEN
CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r1d(:) ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r1d(:) ), clinfo )
ELSEIF( PRESENT(pv_r2d) ) THEN ELSEIF( PRESENT(pv_r2d) ) THEN
CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r2d(ix1:ix2,iy1:iy2) ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r2d(ix1:ix2,iy1:iy2) ), clinfo )
ELSEIF( PRESENT(pv_r3d) ) THEN ELSEIF( PRESENT(pv_r3d) ) THEN
CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r3d(ix1:ix2,iy1:iy2,:) ), clinfo ) CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r3d(ix1:ix2,iy1:iy2,:) ), clinfo )
......
src/OCE/IOM/prtctl.F90
View file @ fc409945
...@@ -137,9 +137,14 @@ CONTAINS ...@@ -137,9 +137,14 @@ CONTAINS
INTEGER :: jn, jl, kdir INTEGER :: jn, jl, kdir
INTEGER :: iis, iie, jjs, jje INTEGER :: iis, iie, jjs, jje
INTEGER :: itra, inum INTEGER :: itra, inum
INTEGER, DIMENSION(4) :: ishape
REAL(2*wp) :: zsum1, zsum2, zvctl1, zvctl2 REAL(2*wp) :: zsum1, zsum2, zvctl1, zvctl2
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( ( ktab2d_1 * ktab3d_1 * ktab4d_1 * ktab2d_2 * ktab3d_2 ) /= 0 ) THEN
CALL ctl_stop( 'prt_ctl is not working with tiles' )
ENDIF
! Arrays, scalars initialization ! Arrays, scalars initialization
cl1 = '' cl1 = ''
cl2 = '' cl2 = ''
...@@ -157,12 +162,19 @@ CONTAINS ...@@ -157,12 +162,19 @@ CONTAINS
! Loop over each sub-domain, i.e. the total number of processors ijsplt ! Loop over each sub-domain, i.e. the total number of processors ijsplt
DO jl = 1, SIZE(nall_ictls) DO jl = 1, SIZE(nall_ictls)
! define shoter names... IF( PRESENT(tab2d_1) ) ishape(1:2) = SHAPE(tab2d_1)
iis = MAX( nall_ictls(jl), ntsi ) IF( PRESENT(tab3d_1) ) ishape(1:3) = SHAPE(tab3d_1)
iie = MIN( nall_ictle(jl), ntei ) IF( PRESENT(tab4d_1) ) ishape(1:4) = SHAPE(tab4d_1)
jjs = MAX( nall_jctls(jl), ntsj ) IF( ishape(1) == jpi .AND. ishape(2) == jpj ) THEN
jje = MIN( nall_jctle(jl), ntej ) iis = Nis0 ; iie = Nie0 ; jjs = Njs0 ; jje = Nje0
ELSE
iis = 1 ; iie = ishape(1) ; jjs = 1 ; jje = ishape(2)
ENDIF
iis = MAX( nall_ictls(jl), iis )
iie = MIN( nall_ictle(jl), iie )
jjs = MAX( nall_jctls(jl), jjs )
jje = MIN( nall_jctle(jl), jje )
IF( PRESENT(clinfo) ) THEN ; inum = numprt_top(jl) IF( PRESENT(clinfo) ) THEN ; inum = numprt_top(jl)
ELSE ; inum = numprt_oce(jl) ELSE ; inum = numprt_oce(jl)
...@@ -188,32 +200,32 @@ CONTAINS ...@@ -188,32 +200,32 @@ CONTAINS
! 2D arrays ! 2D arrays
IF( PRESENT(tab2d_1) ) THEN IF( PRESENT(tab2d_1) ) THEN
IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) * mask1(iis:iie,jjs:jje,1) ) IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) * mask1(A2D(0),1) )
ELSE ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) ) ELSE ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) )
ENDIF ENDIF
ENDIF ENDIF
IF( PRESENT(tab2d_2) ) THEN IF( PRESENT(tab2d_2) ) THEN
IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) * mask2(iis:iie,jjs:jje,1) ) IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) * mask2(A2D(0),1) )
ELSE ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) ) ELSE ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) )
ENDIF ENDIF
ENDIF ENDIF
! 3D arrays ! 3D arrays
IF( PRESENT(tab3d_1) ) THEN IF( PRESENT(tab3d_1) ) THEN
IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) * mask1(iis:iie,jjs:jje,1:kdir) ) IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) * mask1(A2D(0),1:kdir) )
ELSE ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) ) ELSE ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) )
ENDIF ENDIF
ENDIF ENDIF
IF( PRESENT(tab3d_2) ) THEN IF( PRESENT(tab3d_2) ) THEN
IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) * mask2(iis:iie,jjs:jje,1:kdir) ) IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) * mask2(A2D(0),1:kdir) )
ELSE ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) ) ELSE ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) )
ENDIF ENDIF
ENDIF ENDIF
! 4D arrays ! 4D arrays
IF( PRESENT(tab4d_1) ) THEN IF( PRESENT(tab4d_1) ) THEN
IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) * mask1(iis:iie,jjs:jje,1:kdir) ) IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) * mask1(A2D(0),1:kdir) )
ELSE ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) ) ELSE ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) )
ENDIF ENDIF
ENDIF ENDIF
......