New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

trcnxt.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP – NEMO

Context Navigation

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcnxt.F90 @ 10963

Last change on this file since 10963 was 10963, checked in by acc, 5 years ago
2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps : Convert TOP routines in top-level TOP directory and all knock on effects of these conversions. SETTE tested.
Property svn:keywords set to `Id`
File size: 12.5 KB

Rev	Line
[941]	1	MODULE trcnxt
	2	!!======================================================================
	3	!! * MODULE trcnxt *
	4	!! Ocean passive tracers: time stepping on passives tracers
	5	!!======================================================================
[2528]	6	!! History : 7.0 ! 1991-11 (G. Madec) Original code
	7	!! ! 1993-03 (M. Guyon) symetrical conditions
	8	!! ! 1995-02 (M. Levy) passive tracers
	9	!! ! 1996-02 (G. Madec & M. Imbard) opa release 8.0
	10	!! 8.0 ! 1996-04 (A. Weaver) Euler forward step
	11	!! 8.2 ! 1999-02 (G. Madec, N. Grima) semi-implicit pressure grad.
	12	!! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module
	13	!! ! 2002-08 (G. Madec) F90: Free form and module
	14	!! ! 2002-11 (C. Talandier, A-M Treguier) Open boundaries
	15	!! ! 2004-03 (C. Ethe) passive tracers
	16	!! ! 2007-02 (C. Deltel) Diagnose ML trends for passive tracers
	17	!! 2.0 ! 2006-02 (L. Debreu, C. Mazauric) Agrif implementation
	18	!! 3.0 ! 2008-06 (G. Madec) time stepping always done in trazdf
	19	!! 3.1 ! 2009-02 (G. Madec, R. Benshila) re-introduce the vvl option
	20	!! 3.3 ! 2010-06 (C. Ethe, G. Madec) Merge TRA-TRC
[1175]	21	!!----------------------------------------------------------------------
[941]	22	#if defined key_top
	23	!!----------------------------------------------------------------------
	24	!! 'key_top' TOP models
	25	!!----------------------------------------------------------------------
	26	!! trc_nxt : time stepping on passive tracers
	27	!!----------------------------------------------------------------------
	28	USE oce_trc ! ocean dynamics and tracers variables
[2528]	29	USE trc ! ocean passive tracers variables
[4990]	30	USE trd_oce
[2528]	31	USE trdtra
	32	USE tranxt
[7646]	33	USE bdy_oce , ONLY: ln_bdy
[6140]	34	USE trcbdy ! BDY open boundaries
[2528]	35	# if defined key_agrif
[941]	36	USE agrif_top_interp
[2528]	37	# endif
[9019]	38	!
	39	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	40	USE prtctl_trc ! Print control for debbuging
[941]	41
	42	IMPLICIT NONE
	43	PRIVATE
	44
[9019]	45	PUBLIC trc_nxt ! routine called by step.F90
[2528]	46
[7881]	47	REAL(wp) :: rfact1, rfact2
[7872]	48
[941]	49	!!----------------------------------------------------------------------
[10067]	50	!! NEMO/TOP 4.0 , NEMO Consortium (2018)
[7753]	51	!! $Id$
[10068]	52	!! Software governed by the CeCILL license (see ./LICENSE)
[941]	53	!!----------------------------------------------------------------------
	54	CONTAINS
	55
[10954]	56	SUBROUTINE trc_nxt( kt, Kbb, Kmm, Krhs )
[941]	57	!!----------------------------------------------------------------------
	58	!! * ROUTINE trcnxt *
	59	!!
	60	!! ** Purpose : Compute the passive tracers fields at the
	61	!! next time-step from their temporal trends and swap the fields.
	62	!!
	63	!! ** Method : Apply lateral boundary conditions on (ua,va) through
	64	!! call to lbc_lnk routine
	65	!! default:
	66	!! arrays swap
	67	!! (trn) = (tra) ; (tra) = (0,0)
	68	!! (trb) = (trn)
	69	!!
	70	!! For Arakawa or TVD Scheme :
	71	!! A Asselin time filter applied on now tracers (trn) to avoid
	72	!! the divergence of two consecutive time-steps and tr arrays
	73	!! to prepare the next time_step:
	74	!! (trb) = (trn) + atfp [ (trb) + (tra) - 2 (trn) ]
	75	!! (trn) = (tra) ; (tra) = (0,0)
	76	!!
	77	!!
	78	!! ** Action : - update trb, trn
	79	!!----------------------------------------------------------------------
[1271]	80	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[10954]	81	INTEGER, INTENT( in ) :: Kbb, Kmm, Krhs ! time level indices
[2715]	82	!
[7753]	83	INTEGER :: jk, jn ! dummy loop indices
[941]	84	REAL(wp) :: zfact ! temporary scalar
	85	CHARACTER (len=22) :: charout
[9019]	86	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztrdt ! 4D workspace
[941]	87	!!----------------------------------------------------------------------
[3294]	88	!
[9124]	89	IF( ln_timing ) CALL timing_start('trc_nxt')
[3294]	90	!
	91	IF( kt == nittrc000 .AND. lwp ) THEN
[941]	92	WRITE(numout,*)
	93	WRITE(numout,*) 'trc_nxt : time stepping on passive tracers'
	94	ENDIF
[6140]	95	!
[5656]	96	#if defined key_agrif
	97	CALL Agrif_trc ! AGRIF zoom boundaries
	98	#endif
[9081]	99	! Update after tracer on domain lateral boundaries
[10425]	100	CALL lbc_lnk( 'trcnxt', tra(:,:,:,:), 'T', 1. )
[941]	101
[10963]	102	IF( ln_bdy ) CALL trc_bdy( kt, Kbb, Kmm, Krhs )
[941]	103
[6140]	104	IF( l_trdtrc ) THEN ! trends: store now fields before the Asselin filter application
[9019]	105	ALLOCATE( ztrdt(jpi,jpj,jpk,jptra) )
[10096]	106	ztrdt(:,:,:,:) = 0._wp
	107	IF( ln_traldf_iso ) THEN ! diagnose the "pure" Kz diffusive trend
	108	DO jn = 1, jptra
[10946]	109	CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_zdfp, ztrdt(:,:,:,jn) )
[10096]	110	ENDDO
	111	ENDIF
	112
	113	! total trend for the non-time-filtered variables.
	114	zfact = 1.0 / rdttrc
	115	! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from tsn terms
	116	IF( ln_linssh ) THEN ! linear sea surface height only
	117	DO jn = 1, jptra
	118	DO jk = 1, jpkm1
	119	ztrdt(:,:,jk,jn) = ( tra(:,:,jk,jn)e3t_a(:,:,jk) / e3t_n(:,:,jk) - trn(:,:,jk,jn)) zfact
	120	END DO
	121	END DO
	122	ELSE
	123	DO jn = 1, jptra
	124	DO jk = 1, jpkm1
	125	ztrdt(:,:,jk,jn) = ( tra(:,:,jk,jn) - trn(:,:,jk,jn) ) * zfact
	126	END DO
	127	END DO
	128	ENDIF
	129	!
	130	DO jn = 1, jptra
[10946]	131	CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_tot, ztrdt(:,:,:,jn) )
[10096]	132	ENDDO
	133	!
	134	IF( ln_linssh ) THEN ! linear sea surface height only
	135	! Store now fields before applying the Asselin filter
	136	! in order to calculate Asselin filter trend later.
	137	ztrdt(:,:,:,:) = trn(:,:,:,:)
	138	ENDIF
	139
[2528]	140	ENDIF
[6140]	141	! ! Leap-Frog + Asselin filter time stepping
[8399]	142	IF( (neuler == 0 .AND. kt == nittrc000) .OR. ln_top_euler ) THEN ! Euler time-stepping (only swap)
[2528]	143	DO jn = 1, jptra
	144	DO jk = 1, jpkm1
[7753]	145	trn(:,:,jk,jn) = tra(:,:,jk,jn)
[8399]	146	trb(:,:,jk,jn) = trn(:,:,jk,jn)
[2528]	147	END DO
	148	END DO
[10096]	149	IF (l_trdtrc .AND. .NOT. ln_linssh ) THEN ! Zero Asselin filter contribution must be explicitly written out since for vvl
	150	! ! Asselin filter is output by tra_nxt_vvl that is not called on this time step
	151	ztrdt(:,:,:,:) = 0._wp
	152	DO jn = 1, jptra
[10946]	153	CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
[10096]	154	ENDDO
	155	END IF
	156	!
[7872]	157	ELSE
[7881]	158	IF( .NOT. l_offline ) THEN ! Leap-Frog + Asselin filter time stepping
[10954]	159	IF( ln_linssh ) THEN ; CALL tra_nxt_fix( kt, Kmm, nittrc000, 'TRC', trb, trn, tra, jptra ) ! linear ssh
	160	ELSE ; CALL tra_nxt_vvl( kt, Kbb, Kmm, Krhs, nittrc000, rdttrc, 'TRC', trb, trn, tra, &
	161	& sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh
[7872]	162	ENDIF
	163	ELSE
	164	CALL trc_nxt_off( kt ) ! offline
[2528]	165	ENDIF
[6140]	166	!
[10425]	167	CALL lbc_lnk_multi( 'trcnxt', trb(:,:,:,:), 'T', 1._wp, trn(:,:,:,:), 'T', 1._wp, tra(:,:,:,:), 'T', 1._wp )
[2528]	168	ENDIF
[6140]	169	!
[10096]	170	IF( l_trdtrc .AND. ln_linssh ) THEN ! trend of the Asselin filter (tb filtered - tb)/dt )
[2528]	171	DO jn = 1, jptra
	172	DO jk = 1, jpkm1
[6140]	173	zfact = 1._wp / r2dttrc
[7753]	174	ztrdt(:,:,jk,jn) = ( trb(:,:,jk,jn) - ztrdt(:,:,jk,jn) ) * zfact
[1175]	175	END DO
[10946]	176	CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) )
[2528]	177	END DO
	178	END IF
[10096]	179	IF( l_trdtrc ) DEALLOCATE( ztrdt )
[2528]	180	!
[941]	181	IF(ln_ctl) THEN ! print mean trends (used for debugging)
	182	WRITE(charout, FMT="('nxt')")
	183	CALL prt_ctl_trc_info(charout)
	184	CALL prt_ctl_trc(tab4d=trn, mask=tmask, clinfo=ctrcnm)
	185	ENDIF
[2528]	186	!
[9124]	187	IF( ln_timing ) CALL timing_stop('trc_nxt')
[3294]	188	!
[941]	189	END SUBROUTINE trc_nxt
	190
[7872]	191
	192	SUBROUTINE trc_nxt_off( kt )
	193	!!----------------------------------------------------------------------
	194	!! * ROUTINE tra_nxt_vvl *
	195	!!
	196	!! ** Purpose : Time varying volume: apply the Asselin time filter
	197	!! and swap the tracer fields.
	198	!!
	199	!! ** Method : - Apply a thickness weighted Asselin time filter on now fields.
	200	!! - save in (ta,sa) a thickness weighted average over the three
	201	!! time levels which will be used to compute rdn and thus the semi-
	202	!! implicit hydrostatic pressure gradient (ln_dynhpg_imp = T)
	203	!! - swap tracer fields to prepare the next time_step.
	204	!! This can be summurized for tempearture as:
	205	!! ztm = ( e3t_ntn + rbcp[ e3t_btb - 2 e3t_ntn + e3t_a*ta ] ) ln_dynhpg_imp = T
	206	!! /( e3t_n + rbcp*[ e3t_b - 2 e3t_n + e3t_a ] )
	207	!! ztm = 0 otherwise
	208	!! tb = ( e3t_ntn + atfp[ e3t_btb - 2 e3t_ntn + e3t_a*ta ] )
	209	!! /( e3t_n + atfp*[ e3t_b - 2 e3t_n + e3t_a ] )
	210	!! tn = ta
	211	!! ta = zt (NB: reset to 0 after eos_bn2 call)
	212	!!
	213	!! ** Action : - (tb,sb) and (tn,sn) ready for the next time step
	214	!! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T)
	215	!!----------------------------------------------------------------------
	216	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	217	!!
	218	INTEGER :: ji, jj, jk, jn ! dummy loop indices
	219	REAL(wp) :: ztc_a , ztc_n , ztc_b , ztc_f , ztc_d ! local scalar
	220	REAL(wp) :: ze3t_b, ze3t_n, ze3t_a, ze3t_f, ze3t_d ! - -
	221	!!----------------------------------------------------------------------
	222	!
	223	IF( kt == nittrc000 ) THEN
	224	IF(lwp) WRITE(numout,*)
	225	IF(lwp) WRITE(numout,*) 'trc_nxt_off : time stepping'
	226	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
	227	IF( .NOT. ln_linssh ) THEN
	228	rfact1 = atfp * rdttrc
	229	rfact2 = rfact1 / rau0
	230	ENDIF
	231	!
	232	ENDIF
	233	!
	234	DO jn = 1, jptra
	235	DO jk = 1, jpkm1
	236	DO jj = 1, jpj
	237	DO ji = 1, jpi
	238	ze3t_b = e3t_b(ji,jj,jk)
	239	ze3t_n = e3t_n(ji,jj,jk)
	240	ze3t_a = e3t_a(ji,jj,jk)
	241	! ! tracer content at Before, now and after
	242	ztc_b = trb(ji,jj,jk,jn) * ze3t_b
	243	ztc_n = trn(ji,jj,jk,jn) * ze3t_n
	244	ztc_a = tra(ji,jj,jk,jn) * ze3t_a
	245	!
	246	ze3t_d = ze3t_a - 2. * ze3t_n + ze3t_b
	247	ztc_d = ztc_a - 2. * ztc_n + ztc_b
	248	!
	249	ze3t_f = ze3t_n + atfp * ze3t_d
	250	ztc_f = ztc_n + atfp * ztc_d
	251	!
	252	IF( .NOT. ln_linssh .AND. jk == mikt(ji,jj) ) THEN ! first level
	253	ze3t_f = ze3t_f - rfact2 * ( emp_b(ji,jj) - emp(ji,jj) )
	254	ztc_f = ztc_f - rfact1 * ( sbc_trc(ji,jj,jn) - sbc_trc_b(ji,jj,jn) )
	255	ENDIF
	256
	257	ze3t_f = 1.e0 / ze3t_f
	258	trb(ji,jj,jk,jn) = ztc_f * ze3t_f ! ptb <-- ptn filtered
	259	trn(ji,jj,jk,jn) = tra(ji,jj,jk,jn) ! ptn <-- pta
	260	!
	261	END DO
	262	END DO
	263	END DO
	264	!
	265	END DO
	266	!
	267	END SUBROUTINE trc_nxt_off
	268
[941]	269	#else
	270	!!----------------------------------------------------------------------
	271	!! Default option Empty module
	272	!!----------------------------------------------------------------------
	273	CONTAINS
	274	SUBROUTINE trc_nxt( kt )
	275	INTEGER, INTENT(in) :: kt
	276	WRITE(,) 'trc_nxt: You should not have seen this print! error?', kt
	277	END SUBROUTINE trc_nxt
	278	#endif
	279	!!======================================================================
	280	END MODULE trcnxt

Note: See TracBrowser for help on using the repository browser.

Download in other formats: