Skip to content
Snippets Groups Projects
Commit 50ecaa32 authored by Guillaume Samson's avatar Guillaume Samson :snowman2:
Browse files

Merge branch '124-fct-is-not-monotonic-with-non-linear-free-surface' into 'main' 

Resolve "FCT is not monotonic with non-linear free surface"

Closes #124

See merge request nemo/nemo!191
parents 873812c0 7dc2bfb1
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
src/OCE/TRA/traadv.F90
+ 6
− 6
View file @ 50ecaa32
......@@ -78,7 +78,7 @@ MODULE traadv
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE tra_adv( kt, Kbb, Kmm, pts, Krhs, pau, pav, paw )
SUBROUTINE tra_adv( kt, Kbb, Kmm, Kaa, pts, Krhs, pau, pav, paw )
!!----------------------------------------------------------------------
!! *** ROUTINE tra_adv ***
!!
......@@ -86,10 +86,10 @@ CONTAINS
!!
!! ** Method : - Update ts(Krhs) with the advective trend following nadv
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! time level indices
REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity
REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt) , INTENT(inout) :: pts ! active tracers and RHS of tracer equation
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Kaa, Krhs ! time level indices
REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity
REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt) , INTENT(inout) :: pts ! active tracers and RHS of tracer equation
!
INTEGER :: ji, jj, jk ! dummy loop index
REAL(wp), DIMENSION(:,:,:), POINTER :: zptu, zptv, zptw
......@@ -193,7 +193,7 @@ CONTAINS
CASE ( np_CEN ) ! Centered scheme : 2nd / 4th order
CALL tra_adv_cen ( kt, nit000, 'TRA', zuu, zvv, zww, Kmm, pts, jpts, Krhs, nn_cen_h, nn_cen_v )
CASE ( np_FCT ) ! FCT scheme : 2nd / 4th order
CALL tra_adv_fct ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_fct_h, nn_fct_v )
CALL tra_adv_fct ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, Kaa, pts, jpts, Krhs, nn_fct_h, nn_fct_v )
CASE ( np_MUS ) ! MUSCL
CALL tra_adv_mus( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, ln_mus_ups )
CASE ( np_UBS ) ! UBS
......
src/OCE/TRA/traadv_fct.F90
+ 24
− 24
View file @ 50ecaa32
......@@ -54,7 +54,7 @@ MODULE traadv_fct
CONTAINS
SUBROUTINE tra_adv_fct( kt, kit000, cdtype, p2dt, pU, pV, pW, &
& Kbb, Kmm, pt, kjpt, Krhs, kn_fct_h, kn_fct_v )
& Kbb, Kmm, Kaa, pt, kjpt, Krhs, kn_fct_h, kn_fct_v )
!!----------------------------------------------------------------------
!! *** ROUTINE tra_adv_fct ***
!!
......@@ -70,8 +70,8 @@ CONTAINS
!! - send trends to trdtra module for further diagnostics (l_trdtra=T)
!! - poleward advective heat and salt transport (ln_diaptr=T)
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! ocean time level indices
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Kaa, Krhs ! ocean time level indices
INTEGER , INTENT(in ) :: kit000 ! first time step index
CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
INTEGER , INTENT(in ) :: kjpt ! number of tracers
......@@ -93,7 +93,7 @@ CONTAINS
!!----------------------------------------------------------------------
!
#if defined key_loop_fusion
CALL tra_adv_fct_lf ( kt, nit000, cdtype, p2dt, pU, pV, pW, Kbb, Kmm, pt, kjpt, Krhs, kn_fct_h, kn_fct_v )
CALL tra_adv_fct_lf ( kt, nit000, cdtype, p2dt, pU, pV, pW, Kbb, Kmm, Kaa, pt, kjpt, Krhs, kn_fct_h, kn_fct_v )
#else
IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile
IF( kt == kit000 ) THEN
......@@ -143,9 +143,9 @@ CONTAINS
ALLOCATE(zwdia(A2D(nn_hls),jpk), zwinf(A2D(nn_hls),jpk), zwsup(A2D(nn_hls),jpk))
DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )
zwdia(ji,jj,jk) = 1._wp + p2dt * ( MAX( wi(ji,jj,jk) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) &
& / e3t(ji,jj,jk,Krhs)
zwinf(ji,jj,jk) = p2dt * MIN( wi(ji,jj,jk ) , 0._wp ) / e3t(ji,jj,jk,Krhs)
zwsup(ji,jj,jk) = -p2dt * MAX( wi(ji,jj,jk+1) , 0._wp ) / e3t(ji,jj,jk,Krhs)
& / e3t(ji,jj,jk,Kaa)
zwinf(ji,jj,jk) = p2dt * MIN( wi(ji,jj,jk ) , 0._wp ) / e3t(ji,jj,jk,Kaa)
zwsup(ji,jj,jk) = -p2dt * MAX( wi(ji,jj,jk+1) , 0._wp ) / e3t(ji,jj,jk,Kaa)
END_3D
END IF
!
......@@ -189,7 +189,7 @@ CONTAINS
pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra &
& / e3t(ji,jj,jk,Kmm ) * tmask(ji,jj,jk)
zwi(ji,jj,jk) = ( e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * ztra ) &
& / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk)
& / e3t(ji,jj,jk,Kaa) * tmask(ji,jj,jk)
END_3D
IF ( ll_zAimp ) THEN
......@@ -309,7 +309,7 @@ CONTAINS
ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) &
& + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) &
& + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj)
ztw(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk)
ztw(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Kaa) * tmask(ji,jj,jk)
END_3D
!
CALL tridia_solver( zwdia, zwsup, zwinf, ztw, ztw , 0 )
......@@ -323,7 +323,7 @@ CONTAINS
!
! !== monotonicity algorithm ==!
!
CALL nonosc( Kmm, pt(:,:,:,jn,Kbb), zwx, zwy, zwz, zwi, p2dt )
CALL nonosc( Kaa, pt(:,:,:,jn,Kbb), zwx, zwy, zwz, zwi, p2dt )
!
! !== final trend with corrected fluxes ==!
!
......@@ -332,7 +332,7 @@ CONTAINS
& + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) &
& + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj)
pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra / e3t(ji,jj,jk,Kmm)
zwi(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk)
zwi(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Kaa) * tmask(ji,jj,jk)
END_3D
!
IF ( ll_zAimp ) THEN
......@@ -384,7 +384,7 @@ CONTAINS
END SUBROUTINE tra_adv_fct
SUBROUTINE nonosc( Kmm, pbef, paa, pbb, pcc, paft, p2dt )
SUBROUTINE nonosc( Kaa, pbef, paa, pbb, pcc, paft, p2dt )
!!---------------------------------------------------------------------
!! *** ROUTINE nonosc ***
!!
......@@ -397,7 +397,7 @@ CONTAINS
!! drange (1995) multi-dimensional forward-in-time and upstream-
!! in-space based differencing for fluid
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: Kmm ! time level index
INTEGER , INTENT(in ) :: Kaa ! time level index
REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pbef ! before field
REAL(wp), DIMENSION(A2D(nn_hls) ,jpk), INTENT(in ) :: paft ! after field
......@@ -451,7 +451,7 @@ CONTAINS
& + MAX( 0., pcc(ji ,jj ,jk ) ) - MIN( 0., pcc(ji ,jj ,jk+1) )
! up & down beta terms
zbt = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) / p2dt
zbt = e1e2t(ji,jj) * e3t(ji,jj,jk,Kaa) / p2dt
zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt
zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt
END_2D
......@@ -689,7 +689,7 @@ CONTAINS
out = 0.5 * ( zfp * pt(ji,jj,jk,jn,Kbb) + zfm * pt(ji,jj+1,jk,jn,Kbb) )
SUBROUTINE tra_adv_fct_lf( kt, kit000, cdtype, p2dt, pU, pV, pW, &
& Kbb, Kmm, pt, kjpt, Krhs, kn_fct_h, kn_fct_v )
& Kbb, Kmm, Kaa, pt, kjpt, Krhs, kn_fct_h, kn_fct_v )
!!----------------------------------------------------------------------
!! *** ROUTINE tra_adv_fct ***
!!
......@@ -705,8 +705,8 @@ CONTAINS
!! - send trends to trdtra module for further diagnostics (l_trdtra=T)
!! - poleward advective heat and salt transport (ln_diaptr=T)
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! ocean time level indices
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Kaa, Krhs ! ocean time level indices
INTEGER , INTENT(in ) :: kit000 ! first time step index
CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
INTEGER , INTENT(in ) :: kjpt ! number of tracers
......@@ -766,9 +766,9 @@ CONTAINS
ALLOCATE(zwdia(jpi,jpj,jpk), zwinf(jpi,jpj,jpk),zwsup(jpi,jpj,jpk))
DO_3D( 1, 1, 1, 1, 1, jpkm1 )
zwdia(ji,jj,jk) = 1._wp + p2dt * ( MAX( wi(ji,jj,jk) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) &
& / e3t(ji,jj,jk,Krhs)
zwinf(ji,jj,jk) = p2dt * MIN( wi(ji,jj,jk ) , 0._wp ) / e3t(ji,jj,jk,Krhs)
zwsup(ji,jj,jk) = -p2dt * MAX( wi(ji,jj,jk+1) , 0._wp ) / e3t(ji,jj,jk,Krhs)
& / e3t(ji,jj,jk,Kaa)
zwinf(ji,jj,jk) = p2dt * MIN( wi(ji,jj,jk ) , 0._wp ) / e3t(ji,jj,jk,Kaa)
zwsup(ji,jj,jk) = -p2dt * MAX( wi(ji,jj,jk+1) , 0._wp ) / e3t(ji,jj,jk,Kaa)
END_3D
END IF
!
......@@ -813,7 +813,7 @@ CONTAINS
pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra &
& / e3t(ji,jj,jk,Kmm ) * tmask(ji,jj,jk)
zwi(ji,jj,jk) = ( e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * ztra ) &
& / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk)
& / e3t(ji,jj,jk,Kaa ) * tmask(ji,jj,jk)
END_2D
END DO
......@@ -929,7 +929,7 @@ CONTAINS
ztra = - ( zwx_3d(ji,jj,jk) - zwx_3d(ji-1,jj ,jk ) &
& + zwy_3d(ji,jj,jk) - zwy_3d(ji ,jj-1,jk ) &
& + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj)
ztw(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk)
ztw(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Kaa) * tmask(ji,jj,jk)
END_3D
!
CALL tridia_solver( zwdia, zwsup, zwinf, ztw, ztw , 0 )
......@@ -943,7 +943,7 @@ CONTAINS
!
! !== monotonicity algorithm ==!
!
CALL nonosc( Kmm, pt(:,:,:,jn,Kbb), zwx_3d, zwy_3d, zwz, zwi, p2dt )
CALL nonosc( Kaa, pt(:,:,:,jn,Kbb), zwx_3d, zwy_3d, zwz, zwi, p2dt )
!
! !== final trend with corrected fluxes ==!
!
......@@ -952,7 +952,7 @@ CONTAINS
& + zwy_3d(ji,jj,jk) - zwy_3d(ji ,jj-1,jk ) &
& + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj)
pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra / e3t(ji,jj,jk,Kmm)
zwi(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk)
zwi(ji,jj,jk) = zwi(ji,jj,jk) + p2dt * ztra / e3t(ji,jj,jk,Kaa) * tmask(ji,jj,jk)
END_3D
!
IF ( ll_zAimp ) THEN
......
src/OCE/step.F90
+ 1
− 1
View file @ 50ecaa32
......@@ -336,7 +336,7 @@ CONTAINS
DO jtile = 1, nijtile
IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = jtile )
CALL tra_adv ( kstp, Nbb, Nnn, ts, Nrhs ) ! hor. + vert. advection ==> RHS
CALL tra_adv ( kstp, Nbb, Nnn, Naa, ts, Nrhs ) ! hor. + vert. advection ==> RHS
IF( ln_zdfmfc ) CALL tra_mfc ( kstp, Nbb, ts, Nrhs ) ! Mass Flux Convection
IF( ln_zdfosm ) THEN
CALL tra_osm ( kstp, Nnn, ts, Nrhs ) ! OSMOSIS non-local tracer fluxes ==> RHS
......
src/OCE/stpmlf.F90
+ 1
− 1
View file @ 50ecaa32
......@@ -351,7 +351,7 @@ CONTAINS
DO jtile = 1, nijtile
IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = jtile )
CALL tra_adv ( kstp, Nbb, Nnn, ts, Nrhs ) ! hor. + vert. advection ==> RHS
CALL tra_adv ( kstp, Nbb, Nnn, Naa, ts, Nrhs ) ! hor. + vert. advection ==> RHS
IF( ln_zdfmfc ) CALL tra_mfc ( kstp, Nbb, ts, Nrhs ) ! Mass Flux Convection
IF( ln_zdfosm ) THEN
CALL tra_osm ( kstp, Nnn, ts, Nrhs ) ! OSMOSIS non-local tracer fluxes ==> RHS
......
src/OCE/stprk3_stg.F90
+ 3
− 3
View file @ 50ecaa32
......@@ -271,7 +271,7 @@ CONTAINS
!
CALL trc_sbc_RK3( kstp, Kmm, tr, Krhs, kstg ) ! surface boundary condition
!
CALL trc_adv ( kstp, Kbb, Kmm, tr, Krhs, zaU, zaV, ww ) ! horizontal & vertical advection
CALL trc_adv ( kstp, Kbb, Kmm, Kaa, tr, Krhs, zaU, zaV, ww ) ! horizontal & vertical advection
!
! !== time integration ==! ∆t = rn_Dt/3 (stg1) or rn_Dt/2 (stg2)
DO jn = 1, jptra
......@@ -299,7 +299,7 @@ CONTAINS
!
CALL trc_sbc_RK3( kstp, Kmm, tr, Krhs, kstg ) ! surface boundary condition
!
CALL trc_adv ( kstp, Kbb, Kmm, tr, Krhs, zaU, zaV, ww ) ! horizontal & vertical advection
CALL trc_adv ( kstp, Kbb, Kmm, Kaa, tr, Krhs, zaU, zaV, ww ) ! horizontal & vertical advection
!
CALL trc_sms ( kstp, Kbb, Kbb, Krhs ) ! tracers: sinks and sources
CALL trc_trp ( kstp, Kbb, Kmm, Krhs, Kaa ) ! transport of passive tracers (without advection)
......@@ -320,7 +320,7 @@ CONTAINS
CALL eos( ts(:,:,:,:,Kmm), rhd, rhop, gdept_0 ) ! now in potential density for tra_mle computation
!===>>> CAUTION here may be without GM velocity but stokes drift required ! 0 barotropic divergence for GM != 0 barotropic divergence for SD
!!st consistence 2D / 3D - flux de masse
CALL tra_adv( kstp, Kbb, Kmm, ts, Krhs, zaU, zaV, ww ) ! hor. + vert. advection ==> RHS
CALL tra_adv( kstp, Kbb, Kmm, Kaa, ts, Krhs, zaU, zaV, ww ) ! hor. + vert. advection ==> RHS
!===>>>>>> stg1&2: Verify the necessity of these trends (we may need it as there are in the RHS of dynspg_ts ?)
!!gm ====>>>> needed for heat and salt fluxes associated with mass/volume flux
......
src/TOP/TRP/trcadv.F90
+ 3
− 3
View file @ 50ecaa32
......@@ -72,7 +72,7 @@ MODULE trcadv
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE trc_adv( kt, Kbb, Kmm, ptr, Krhs, pau, pav, paw )
SUBROUTINE trc_adv( kt, Kbb, Kmm, Kaa, ptr, Krhs, pau, pav, paw )
!!----------------------------------------------------------------------
!! *** ROUTINE trc_adv ***
!!
......@@ -81,7 +81,7 @@ CONTAINS
!! ** Method : - Update tr(Krhs) with the advective trend following nadv
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! time level indices
INTEGER , INTENT(in ) :: Kbb, Kmm, Kaa, Krhs ! time level indices
REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity
REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt) , INTENT(inout) :: ptr ! passive tracers and RHS of tracer equation
!
......@@ -158,7 +158,7 @@ CONTAINS
CASE ( np_CEN ) ! Centered : 2nd / 4th order
CALL tra_adv_cen ( kt, nittrc000,'TRC', zuu, zvv, zww, Kmm, ptr, jptra, Krhs, nn_cen_h, nn_cen_v )
CASE ( np_FCT ) ! FCT : 2nd / 4th order
CALL tra_adv_fct( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, ptr, jptra, Krhs, nn_fct_h, nn_fct_v )
CALL tra_adv_fct( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, Kaa, ptr, jptra, Krhs, nn_fct_h, nn_fct_v )
CASE ( np_MUS ) ! MUSCL
CALL tra_adv_mus( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, ptr, jptra, Krhs, ln_mus_ups )
CASE ( np_UBS ) ! UBS
......
src/TOP/TRP/trctrp.F90
+ 1
− 1
View file @ 50ecaa32
......@@ -89,7 +89,7 @@ CONTAINS
#endif
#if ! defined key_RK3
! ! MLF only: add the advection trend to the RHS
CALL trc_adv ( kt, Kbb, Kmm, tr, Krhs ) ! horizontal & vertical advection
CALL trc_adv ( kt, Kbb, Kmm, Kaa, tr, Krhs ) ! horizontal & vertical advection
#endif
CALL trc_ldf ( kt, Kbb, Kmm, tr, Krhs ) ! lateral mixing
CALL trc_zdf ( kt, Kbb, Kmm, Krhs, tr, Kaa ) ! vert. mixing & after tracer ==> after
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment