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.

dtadyn.F90 in NEMO/branches/2019/dev_r12114_ticket_2263/src/OFF – NEMO

Context Navigation

source: NEMO/branches/2019/dev_r12114_ticket_2263/src/OFF/dtadyn.F90 @ 12172

Last change on this file since 12172 was 12116, checked in by smasson, 5 years ago
dev_r12114_ticket_2263: new version of fldread calendar, see #2263
Property svn:keywords set to `Id`
File size: 41.0 KB

Rev	Line
[325]	1	MODULE dtadyn
	2	!!======================================================================
	3	!! * MODULE dtadyn *
[2528]	4	!! Off-line : interpolation of the physical fields
	5	!!======================================================================
	6	!! History : OPA ! 1992-01 (M. Imbard) Original code
	7	!! 8.0 ! 1998-04 (L.Bopp MA Foujols) slopes for isopyc.
	8	!! - ! 1998-05 (L. Bopp) read output of coupled run
	9	!! 8.2 ! 2001-01 (M. Levy et M. Benjelloul) add netcdf FORMAT
	10	!! NEMO 1.0 ! 2005-03 (O. Aumont and A. El Moussaoui) F90
	11	!! - ! 2005-12 (C. Ethe) Adapted for DEGINT
	12	!! 3.0 ! 2007-06 (C. Ethe) use of iom module
	13	!! 3.3 ! 2010-11 (C. Ethe) Full reorganization of the off-line: phasing with the on-line
[3294]	14	!! 3.4 ! 2011-05 (C. Ethe) Use of fldread
[2528]	15	!!----------------------------------------------------------------------
[325]	16
	17	!!----------------------------------------------------------------------
[3294]	18	!! dta_dyn_init : initialization, namelist read, and SAVEs control
[325]	19	!! dta_dyn : Interpolation of the fields
	20	!!----------------------------------------------------------------------
	21	USE oce ! ocean dynamics and tracers variables
[2528]	22	USE c1d ! 1D configuration: lk_c1d
	23	USE dom_oce ! ocean domain: variables
[7646]	24	USE domvvl ! variable volume
[2528]	25	USE zdf_oce ! ocean vertical physics: variables
	26	USE sbc_oce ! surface module: variables
[3294]	27	USE trc_oce ! share ocean/biogeo variables
[325]	28	USE phycst ! physical constants
[2528]	29	USE trabbl ! active tracer: bottom boundary layer
	30	USE ldfslp ! lateral diffusion: iso-neutral slopes
[7646]	31	USE sbcrnf ! river runoffs
	32	USE ldftra ! ocean tracer lateral physics
[2528]	33	USE zdfmxl ! vertical physics: mixed layer depth
	34	USE eosbn2 ! equation of state - Brunt Vaisala frequency
[325]	35	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[2528]	36	USE zpshde ! z-coord. with partial steps: horizontal derivatives
	37	USE in_out_manager ! I/O manager
	38	USE iom ! I/O library
[325]	39	USE lib_mpp ! distributed memory computing library
[3294]	40	USE prtctl ! print control
	41	USE fldread ! read input fields
	42	USE timing ! Timing
[7646]	43	USE trc, ONLY : ln_rsttr, numrtr, numrtw, lrst_trc
[325]	44
	45	IMPLICIT NONE
	46	PRIVATE
	47
[10222]	48	PUBLIC dta_dyn_init ! called by opa.F90
	49	PUBLIC dta_dyn ! called by step.F90
	50	PUBLIC dta_dyn_sed_init ! called by opa.F90
	51	PUBLIC dta_dyn_sed ! called by step.F90
	52	PUBLIC dta_dyn_swp ! called by step.F90
[325]	53
[7646]	54	CHARACTER(len=100) :: cn_dir !: Root directory for location of ssr files
	55	LOGICAL :: ln_dynrnf !: read runoff data in file (T) or set to zero (F)
	56	LOGICAL :: ln_dynrnf_depth !: read runoff data in file (T) or set to zero (F)
	57	REAL(wp) :: fwbcorr
[325]	58
[7646]	59
	60	INTEGER , PARAMETER :: jpfld = 20 ! maximum number of fields to read
[3294]	61	INTEGER , SAVE :: jf_tem ! index of temperature
	62	INTEGER , SAVE :: jf_sal ! index of salinity
[7646]	63	INTEGER , SAVE :: jf_uwd ! index of u-transport
	64	INTEGER , SAVE :: jf_vwd ! index of v-transport
	65	INTEGER , SAVE :: jf_wwd ! index of v-transport
[3294]	66	INTEGER , SAVE :: jf_avt ! index of Kz
	67	INTEGER , SAVE :: jf_mld ! index of mixed layer deptht
	68	INTEGER , SAVE :: jf_emp ! index of water flux
[7646]	69	INTEGER , SAVE :: jf_empb ! index of water flux
[3294]	70	INTEGER , SAVE :: jf_qsr ! index of solar radiation
	71	INTEGER , SAVE :: jf_wnd ! index of wind speed
	72	INTEGER , SAVE :: jf_ice ! index of sea ice cover
[4570]	73	INTEGER , SAVE :: jf_rnf ! index of river runoff
[7646]	74	INTEGER , SAVE :: jf_fmf ! index of downward salt flux
[3294]	75	INTEGER , SAVE :: jf_ubl ! index of u-bbl coef
	76	INTEGER , SAVE :: jf_vbl ! index of v-bbl coef
[7646]	77	INTEGER , SAVE :: jf_div ! index of e3t
[325]	78
[7646]	79
	80	TYPE(FLD), ALLOCATABLE, SAVE, DIMENSION(:) :: sf_dyn ! structure of input fields (file informations, fields read)
[3294]	81	! !
	82	REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: uslpdta ! zonal isopycnal slopes
	83	REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: vslpdta ! meridional isopycnal slopes
	84	REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wslpidta ! zonal diapycnal slopes
	85	REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wslpjdta ! meridional diapycnal slopes
[1735]	86
[7646]	87	INTEGER, SAVE :: nprevrec, nsecdyn
[325]	88
[343]	89	!!----------------------------------------------------------------------
[9598]	90	!! NEMO/OFF 4.0 , NEMO Consortium (2018)
[2528]	91	!! $Id$
[10068]	92	!! Software governed by the CeCILL license (see ./LICENSE)
[343]	93	!!----------------------------------------------------------------------
[325]	94	CONTAINS
	95
[1501]	96	SUBROUTINE dta_dyn( kt )
[325]	97	!!----------------------------------------------------------------------
	98	!! * ROUTINE dta_dyn *
	99	!!
[3294]	100	!! ** Purpose : Prepares dynamics and physics fields from a NEMO run
	101	!! for an off-line simulation of passive tracers
[325]	102	!!
[3294]	103	!! ** Method : calculates the position of data
	104	!! - computes slopes if needed
	105	!! - interpolates data if needed
[2528]	106	!!----------------------------------------------------------------------
	107	INTEGER, INTENT(in) :: kt ! ocean time-step index
[9212]	108	!
[7646]	109	INTEGER :: ji, jj, jk
[9212]	110	REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zemp
	111	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zhdivtr
[325]	112	!!----------------------------------------------------------------------
[3294]	113	!
[9124]	114	IF( ln_timing ) CALL timing_start( 'dta_dyn')
[3294]	115	!
[7646]	116	nsecdyn = nsec_year + nsec1jan000 ! number of seconds between Jan. 1st 00h of nit000 year and the middle of time step
[3294]	117	!
[7646]	118	IF( kt == nit000 ) THEN ; nprevrec = 0
	119	ELSE ; nprevrec = sf_dyn(jf_tem)%nrec_a(2)
[325]	120	ENDIF
[7646]	121	CALL fld_read( kt, 1, sf_dyn ) != read data at kt time step ==!
	122	!
	123	IF( l_ldfslp .AND. .NOT.lk_c1d ) CALL dta_dyn_slp( kt ) ! Computation of slopes
	124	!
	125	tsn(:,:,:,jp_tem) = sf_dyn(jf_tem)%fnow(:,:,:) * tmask(:,:,:) ! temperature
	126	tsn(:,:,:,jp_sal) = sf_dyn(jf_sal)%fnow(:,:,:) * tmask(:,:,:) ! salinity
	127	wndm(:,:) = sf_dyn(jf_wnd)%fnow(:,:,1) * tmask(:,:,1) ! wind speed - needed for gas exchange
	128	fmmflx(:,:) = sf_dyn(jf_fmf)%fnow(:,:,1) * tmask(:,:,1) ! downward salt flux (v3.5+)
	129	fr_i(:,:) = sf_dyn(jf_ice)%fnow(:,:,1) * tmask(:,:,1) ! Sea-ice fraction
	130	qsr (:,:) = sf_dyn(jf_qsr)%fnow(:,:,1) * tmask(:,:,1) ! solar radiation
	131	emp (:,:) = sf_dyn(jf_emp)%fnow(:,:,1) * tmask(:,:,1) ! E-P
	132	IF( ln_dynrnf ) THEN
	133	rnf (:,:) = sf_dyn(jf_rnf)%fnow(:,:,1) * tmask(:,:,1) ! E-P
	134	IF( ln_dynrnf_depth .AND. .NOT. ln_linssh ) CALL dta_dyn_hrnf
[3294]	135	ENDIF
[325]	136	!
[7646]	137	un(:,:,:) = sf_dyn(jf_uwd)%fnow(:,:,:) * umask(:,:,:) ! effective u-transport
	138	vn(:,:,:) = sf_dyn(jf_vwd)%fnow(:,:,:) * vmask(:,:,:) ! effective v-transport
	139	wn(:,:,:) = sf_dyn(jf_wwd)%fnow(:,:,:) * tmask(:,:,:) ! effective v-transport
	140	!
	141	IF( .NOT.ln_linssh ) THEN
[9212]	142	ALLOCATE( zemp(jpi,jpj) , zhdivtr(jpi,jpj,jpk) )
	143	zhdivtr(:,:,:) = sf_dyn(jf_div)%fnow(:,:,:) * tmask(:,:,:) ! effective u-transport
	144	emp_b (:,:) = sf_dyn(jf_empb)%fnow(:,:,1) * tmask(:,:,1) ! E-P
	145	zemp (:,:) = ( 0.5_wp * ( emp(:,:) + emp_b(:,:) ) + rnf(:,:) + fwbcorr ) * tmask(:,:,1)
[7646]	146	CALL dta_dyn_ssh( kt, zhdivtr, sshb, zemp, ssha, e3t_a(:,:,:) ) != ssh, vertical scale factor & vertical transport
[9212]	147	DEALLOCATE( zemp , zhdivtr )
[7646]	148	! Write in the tracer restart file
[9212]	149	! *********************************
[7646]	150	IF( lrst_trc ) THEN
[3827]	151	IF(lwp) WRITE(numout,*)
[9212]	152	IF(lwp) WRITE(numout,*) 'dta_dyn_ssh : ssh field written in tracer restart file at it= ', kt,' date= ', ndastp
	153	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[7646]	154	CALL iom_rstput( kt, nitrst, numrtw, 'sshn', ssha )
	155	CALL iom_rstput( kt, nitrst, numrtw, 'sshb', sshn )
[1501]	156	ENDIF
[3294]	157	ENDIF
[325]	158	!
[5131]	159	CALL eos ( tsn, rhd, rhop, gdept_0(:,:,:) ) ! In any case, we need rhop
	160	CALL eos_rab( tsn, rab_n ) ! now local thermal/haline expension ratio at T-points
	161	CALL bn2 ( tsn, rab_n, rn2 ) ! before Brunt-Vaisala frequency need for zdfmxl
	162
	163	rn2b(:,:,:) = rn2(:,:,:) ! need for zdfmxl
[7646]	164	CALL zdf_mxl( kt ) ! In any case, we need mxl
[2528]	165	!
[9019]	166	hmld(:,:) = sf_dyn(jf_mld)%fnow(:,:,1) * tmask(:,:,1) ! mixed layer depht
	167	avt(:,:,:) = sf_dyn(jf_avt)%fnow(:,:,:) * tmask(:,:,:) ! vertical diffusive coefficient
[10213]	168	avs(:,:,:) = avt(:,:,:)
[7646]	169	!
[9019]	170	IF( ln_trabbl .AND. .NOT.lk_c1d ) THEN ! diffusive Bottom boundary layer param
	171	ahu_bbl(:,:) = sf_dyn(jf_ubl)%fnow(:,:,1) * umask(:,:,1) ! bbl diffusive coef
	172	ahv_bbl(:,:) = sf_dyn(jf_vbl)%fnow(:,:,1) * vmask(:,:,1)
	173	ENDIF
[2762]	174	!
[7646]	175	!
	176	CALL eos( tsn, rhd, rhop, gdept_0(:,:,:) ) ! In any case, we need rhop
	177	!
[3294]	178	IF(ln_ctl) THEN ! print control
[9440]	179	CALL prt_ctl(tab3d_1=tsn(:,:,:,jp_tem), clinfo1=' tn - : ', mask1=tmask, kdim=jpk )
	180	CALL prt_ctl(tab3d_1=tsn(:,:,:,jp_sal), clinfo1=' sn - : ', mask1=tmask, kdim=jpk )
	181	CALL prt_ctl(tab3d_1=un , clinfo1=' un - : ', mask1=umask, kdim=jpk )
	182	CALL prt_ctl(tab3d_1=vn , clinfo1=' vn - : ', mask1=vmask, kdim=jpk )
	183	CALL prt_ctl(tab3d_1=wn , clinfo1=' wn - : ', mask1=tmask, kdim=jpk )
	184	CALL prt_ctl(tab3d_1=avt , clinfo1=' kz - : ', mask1=tmask, kdim=jpk )
[7646]	185	CALL prt_ctl(tab3d_1=uslp , clinfo1=' slp - u : ', tab3d_2=vslp, clinfo2=' v : ', kdim=jpk)
	186	CALL prt_ctl(tab3d_1=wslpi , clinfo1=' slp - wi: ', tab3d_2=wslpj, clinfo2=' wj: ', kdim=jpk)
[2528]	187	ENDIF
	188	!
[9124]	189	IF( ln_timing ) CALL timing_stop( 'dta_dyn')
[3294]	190	!
[325]	191	END SUBROUTINE dta_dyn
	192
[2528]	193
[3294]	194	SUBROUTINE dta_dyn_init
[325]	195	!!----------------------------------------------------------------------
[3294]	196	!! * ROUTINE dta_dyn_init *
[325]	197	!!
[3294]	198	!! ** Purpose : Initialisation of the dynamical data
	199	!! ** Method : - read the data namdta_dyn namelist
[325]	200	!!----------------------------------------------------------------------
[3294]	201	INTEGER :: ierr, ierr0, ierr1, ierr2, ierr3 ! return error code
	202	INTEGER :: ifpr ! dummy loop indice
	203	INTEGER :: jfld ! dummy loop arguments
	204	INTEGER :: inum, idv, idimv ! local integer
[4147]	205	INTEGER :: ios ! Local integer output status for namelist read
[7646]	206	INTEGER :: ji, jj, jk
	207	REAL(wp) :: zcoef
	208	INTEGER :: nkrnf_max
	209	REAL(wp) :: hrnf_max
[3294]	210	!!
[7646]	211	CHARACTER(len=100) :: cn_dir ! Root directory for location of core files
	212	TYPE(FLD_N), DIMENSION(jpfld) :: slf_d ! array of namelist informations on the fields to read
	213	TYPE(FLD_N) :: sn_uwd, sn_vwd, sn_wwd, sn_empb, sn_emp ! informations about the fields to be read
	214	TYPE(FLD_N) :: sn_tem , sn_sal , sn_avt ! " "
	215	TYPE(FLD_N) :: sn_mld, sn_qsr, sn_wnd , sn_ice , sn_fmf ! " "
	216	TYPE(FLD_N) :: sn_ubl, sn_vbl, sn_rnf ! " "
	217	TYPE(FLD_N) :: sn_div ! informations about the fields to be read
[9212]	218	!!
[7646]	219	NAMELIST/namdta_dyn/cn_dir, ln_dynrnf, ln_dynrnf_depth, fwbcorr, &
[9212]	220	& sn_uwd, sn_vwd, sn_wwd, sn_emp, &
	221	& sn_avt, sn_tem, sn_sal, sn_mld , sn_qsr , &
	222	& sn_wnd, sn_ice, sn_fmf, &
	223	& sn_ubl, sn_vbl, sn_rnf, &
[7646]	224	& sn_empb, sn_div
[9212]	225	!!----------------------------------------------------------------------
[7646]	226	!
[4147]	227	REWIND( numnam_ref ) ! Namelist namdta_dyn in reference namelist : Offline: init. of dynamical data
	228	READ ( numnam_ref, namdta_dyn, IOSTAT = ios, ERR = 901)
[11536]	229	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdta_dyn in reference namelist' )
[4147]	230	REWIND( numnam_cfg ) ! Namelist namdta_dyn in configuration namelist : Offline: init. of dynamical data
	231	READ ( numnam_cfg, namdta_dyn, IOSTAT = ios, ERR = 902 )
[11536]	232	902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namdta_dyn in configuration namelist' )
[4624]	233	IF(lwm) WRITE ( numond, namdta_dyn )
[3294]	234	! ! store namelist information in an array
	235	! ! Control print
[325]	236	IF(lwp) THEN
	237	WRITE(numout,*)
[3294]	238	WRITE(numout,*) 'dta_dyn : offline dynamics '
	239	WRITE(numout,*) '~~~~~~~ '
	240	WRITE(numout,*) ' Namelist namdta_dyn'
[7646]	241	WRITE(numout,*) ' runoffs option enabled (T) or not (F) ln_dynrnf = ', ln_dynrnf
	242	WRITE(numout,*) ' runoffs is spread in vertical ln_dynrnf_depth = ', ln_dynrnf_depth
	243	WRITE(numout,*) ' annual global mean of empmr for ssh correction fwbcorr = ', fwbcorr
[325]	244	WRITE(numout,*)
	245	ENDIF
[3294]	246	!
[7646]	247	jf_uwd = 1 ; jf_vwd = 2 ; jf_wwd = 3 ; jf_emp = 4 ; jf_avt = 5
	248	jf_tem = 6 ; jf_sal = 7 ; jf_mld = 8 ; jf_qsr = 9
	249	jf_wnd = 10 ; jf_ice = 11 ; jf_fmf = 12 ; jfld = jf_fmf
[3294]	250	!
[7646]	251	slf_d(jf_uwd) = sn_uwd ; slf_d(jf_vwd) = sn_vwd ; slf_d(jf_wwd) = sn_wwd
	252	slf_d(jf_emp) = sn_emp ; slf_d(jf_avt) = sn_avt
	253	slf_d(jf_tem) = sn_tem ; slf_d(jf_sal) = sn_sal ; slf_d(jf_mld) = sn_mld
	254	slf_d(jf_qsr) = sn_qsr ; slf_d(jf_wnd) = sn_wnd ; slf_d(jf_ice) = sn_ice
	255	slf_d(jf_fmf) = sn_fmf
[3294]	256	!
[7646]	257	IF( .NOT.ln_linssh ) THEN
[9212]	258	jf_div = jfld + 1 ; jf_empb = jfld + 2 ; jfld = jf_empb
	259	slf_d(jf_div) = sn_div ; slf_d(jf_empb) = sn_empb
[7646]	260	ENDIF
	261	!
[9019]	262	IF( ln_trabbl ) THEN
[9212]	263	jf_ubl = jfld + 1 ; jf_vbl = jfld + 2 ; jfld = jf_vbl
	264	slf_d(jf_ubl) = sn_ubl ; slf_d(jf_vbl) = sn_vbl
[7646]	265	ENDIF
	266	!
[5385]	267	IF( ln_dynrnf ) THEN
[9212]	268	jf_rnf = jfld + 1 ; jfld = jf_rnf
	269	slf_d(jf_rnf) = sn_rnf
[4570]	270	ELSE
[7646]	271	rnf(:,:) = 0._wp
[4570]	272	ENDIF
	273
[3294]	274	ALLOCATE( sf_dyn(jfld), STAT=ierr ) ! set sf structure
[7646]	275	IF( ierr > 0 ) THEN
[3294]	276	CALL ctl_stop( 'dta_dyn: unable to allocate sf structure' ) ; RETURN
	277	ENDIF
[5768]	278	! ! fill sf with slf_i and control print
	279	CALL fld_fill( sf_dyn, slf_d, cn_dir, 'dta_dyn_init', 'Data in file', 'namdta_dyn' )
[7646]	280	!
[3294]	281	! Open file for each variable to get his number of dimension
	282	DO ifpr = 1, jfld
[12116]	283	CALL fld_def( sf_dyn(ifpr) )
	284	CALL iom_open( sf_dyn(ifpr)%clname, sf_dyn(ifpr)%num )
[5768]	285	idv = iom_varid( sf_dyn(ifpr)%num , slf_d(ifpr)%clvar ) ! id of the variable sdjf%clvar
	286	idimv = iom_file ( sf_dyn(ifpr)%num )%ndims(idv) ! number of dimension for variable sdjf%clvar
[12116]	287	CALL iom_close( sf_dyn(ifpr)%num ) ! close file if already open
[5768]	288	ierr1=0
[3294]	289	IF( idimv == 3 ) THEN ! 2D variable
	290	ALLOCATE( sf_dyn(ifpr)%fnow(jpi,jpj,1) , STAT=ierr0 )
	291	IF( slf_d(ifpr)%ln_tint ) ALLOCATE( sf_dyn(ifpr)%fdta(jpi,jpj,1,2) , STAT=ierr1 )
	292	ELSE ! 3D variable
	293	ALLOCATE( sf_dyn(ifpr)%fnow(jpi,jpj,jpk) , STAT=ierr0 )
	294	IF( slf_d(ifpr)%ln_tint ) ALLOCATE( sf_dyn(ifpr)%fdta(jpi,jpj,jpk,2), STAT=ierr1 )
[2528]	295	ENDIF
[3294]	296	IF( ierr0 + ierr1 > 0 ) THEN
	297	CALL ctl_stop( 'dta_dyn_init : unable to allocate sf_dyn array structure' ) ; RETURN
	298	ENDIF
	299	END DO
[325]	300	!
[5836]	301	IF( l_ldfslp .AND. .NOT.lk_c1d ) THEN ! slopes
[3294]	302	IF( sf_dyn(jf_tem)%ln_tint ) THEN ! time interpolation
	303	ALLOCATE( uslpdta (jpi,jpj,jpk,2), vslpdta (jpi,jpj,jpk,2), &
	304	& wslpidta(jpi,jpj,jpk,2), wslpjdta(jpi,jpj,jpk,2), STAT=ierr2 )
[7646]	305	!
	306	IF( ierr2 > 0 ) THEN
	307	CALL ctl_stop( 'dta_dyn_init : unable to allocate slope arrays' ) ; RETURN
	308	ENDIF
[3294]	309	ENDIF
[2528]	310	ENDIF
[7646]	311	!
	312	IF( .NOT.ln_linssh ) THEN
	313	IF( .NOT. sf_dyn(jf_uwd)%ln_clim .AND. ln_rsttr .AND. & ! Restart: read in restart file
	314	iom_varid( numrtr, 'sshn', ldstop = .FALSE. ) > 0 ) THEN
	315	IF(lwp) WRITE(numout,*) ' sshn forcing fields read in the restart file for initialisation'
	316	CALL iom_get( numrtr, jpdom_autoglo, 'sshn', sshn(:,:) )
	317	CALL iom_get( numrtr, jpdom_autoglo, 'sshb', sshb(:,:) )
	318	ELSE
	319	IF(lwp) WRITE(numout,*) ' sshn forcing fields read in the restart file for initialisation'
	320	CALL iom_open( 'restart', inum )
	321	CALL iom_get( inum, jpdom_autoglo, 'sshn', sshn(:,:) )
	322	CALL iom_get( inum, jpdom_autoglo, 'sshb', sshb(:,:) )
	323	CALL iom_close( inum ) ! close file
	324	ENDIF
	325	!
	326	DO jk = 1, jpkm1
	327	e3t_n(:,:,jk) = e3t_0(:,:,jk) * ( 1._wp + sshn(:,:) * tmask(:,:,1) / ( ht_0(:,:) + 1.0 - tmask(:,:,1) ) )
	328	ENDDO
	329	e3t_a(:,:,jpk) = e3t_0(:,:,jpk)
	330
	331	! Horizontal scale factor interpolations
	332	! --------------------------------------
	333	CALL dom_vvl_interpol( e3t_n(:,:,:), e3u_n(:,:,:), 'U' )
	334	CALL dom_vvl_interpol( e3t_n(:,:,:), e3v_n(:,:,:), 'V' )
	335
	336	! Vertical scale factor interpolations
	337	! ------------------------------------
	338	CALL dom_vvl_interpol( e3t_n(:,:,:), e3w_n(:,:,:), 'W' )
	339
	340	e3t_b(:,:,:) = e3t_n(:,:,:)
	341	e3u_b(:,:,:) = e3u_n(:,:,:)
	342	e3v_b(:,:,:) = e3v_n(:,:,:)
	343
	344	! t- and w- points depth
	345	! ----------------------
	346	gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1)
	347	gdepw_n(:,:,1) = 0.0_wp
	348
	349	DO jk = 2, jpk
	350	DO jj = 1,jpj
	351	DO ji = 1,jpi
	352	! zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) ! 0 everywhere
	353	! tmask = wmask, ie everywhere expect at jk = mikt
	354	! 1 for jk =
	355	! mikt
	356	zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
	357	gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1)
	358	gdept_n(ji,jj,jk) = zcoef * ( gdepw_n(ji,jj,jk ) + 0.5 * e3w_n(ji,jj,jk)) &
	359	& + (1-zcoef) * ( gdept_n(ji,jj,jk-1) + e3w_n(ji,jj,jk))
	360	END DO
	361	END DO
	362	END DO
	363
	364	gdept_b(:,:,:) = gdept_n(:,:,:)
	365	gdepw_b(:,:,:) = gdepw_n(:,:,:)
	366	!
[495]	367	ENDIF
[2715]	368	!
[7646]	369	IF( ln_dynrnf .AND. ln_dynrnf_depth ) THEN ! read depht over which runoffs are distributed
	370	IF(lwp) WRITE(numout,*)
	371	IF(lwp) WRITE(numout,*) ' read in the file depht over which runoffs are distributed'
	372	CALL iom_open ( "runoffs", inum ) ! open file
	373	CALL iom_get ( inum, jpdom_data, 'rodepth', h_rnf ) ! read the river mouth array
	374	CALL iom_close( inum ) ! close file
	375	!
	376	nk_rnf(:,:) = 0 ! set the number of level over which river runoffs are applied
	377	DO jj = 1, jpj
	378	DO ji = 1, jpi
	379	IF( h_rnf(ji,jj) > 0._wp ) THEN
	380	jk = 2
	381	DO WHILE ( jk /= mbkt(ji,jj) .AND. gdept_0(ji,jj,jk) < h_rnf(ji,jj) ) ; jk = jk + 1
	382	END DO
	383	nk_rnf(ji,jj) = jk
	384	ELSEIF( h_rnf(ji,jj) == -1._wp ) THEN ; nk_rnf(ji,jj) = 1
	385	ELSEIF( h_rnf(ji,jj) == -999._wp ) THEN ; nk_rnf(ji,jj) = mbkt(ji,jj)
	386	ELSE
	387	CALL ctl_stop( 'sbc_rnf_init: runoff depth not positive, and not -999 or -1, rnf value in file fort.999' )
	388	WRITE(999,*) 'ji, jj, h_rnf(ji,jj) :', ji, jj, h_rnf(ji,jj)
	389	ENDIF
	390	END DO
	391	END DO
	392	DO jj = 1, jpj ! set the associated depth
	393	DO ji = 1, jpi
	394	h_rnf(ji,jj) = 0._wp
	395	DO jk = 1, nk_rnf(ji,jj)
	396	h_rnf(ji,jj) = h_rnf(ji,jj) + e3t_n(ji,jj,jk)
	397	END DO
	398	END DO
	399	END DO
	400	ELSE ! runoffs applied at the surface
	401	nk_rnf(:,:) = 1
	402	h_rnf (:,:) = e3t_n(:,:,1)
	403	ENDIF
	404	nkrnf_max = MAXVAL( nk_rnf(:,:) )
	405	hrnf_max = MAXVAL( h_rnf(:,:) )
	406	IF( lk_mpp ) THEN
[10425]	407	CALL mpp_max( 'dtadyn', nkrnf_max ) ! max over the global domain
	408	CALL mpp_max( 'dtadyn', hrnf_max ) ! max over the global domain
[7646]	409	ENDIF
	410	IF(lwp) WRITE(numout,*) ' '
	411	IF(lwp) WRITE(numout,*) ' max depht of runoff : ', hrnf_max,' max level : ', nkrnf_max
	412	IF(lwp) WRITE(numout,*) ' '
	413	!
[2528]	414	CALL dta_dyn( nit000 )
	415	!
[1501]	416	END SUBROUTINE dta_dyn_init
	417
[10222]	418	SUBROUTINE dta_dyn_sed( kt )
	419	!!----------------------------------------------------------------------
	420	!! * ROUTINE dta_dyn *
	421	!!
	422	!! ** Purpose : Prepares dynamics and physics fields from a NEMO run
	423	!! for an off-line simulation of passive tracers
	424	!!
	425	!! ** Method : calculates the position of data
	426	!! - computes slopes if needed
	427	!! - interpolates data if needed
	428	!!----------------------------------------------------------------------
	429	INTEGER, INTENT(in) :: kt ! ocean time-step index
	430	!
	431	!!----------------------------------------------------------------------
	432	!
	433	IF( ln_timing ) CALL timing_start( 'dta_dyn_sed')
	434	!
	435	nsecdyn = nsec_year + nsec1jan000 ! number of seconds between Jan. 1st 00h of nit000 year and the middle of time step
	436	!
	437	IF( kt == nit000 ) THEN ; nprevrec = 0
	438	ELSE ; nprevrec = sf_dyn(jf_tem)%nrec_a(2)
	439	ENDIF
	440	CALL fld_read( kt, 1, sf_dyn ) != read data at kt time step ==!
	441	!
	442	tsn(:,:,:,jp_tem) = sf_dyn(jf_tem)%fnow(:,:,:) * tmask(:,:,:) ! temperature
	443	tsn(:,:,:,jp_sal) = sf_dyn(jf_sal)%fnow(:,:,:) * tmask(:,:,:) ! salinity
	444	!
	445	CALL eos ( tsn, rhd, rhop, gdept_0(:,:,:) ) ! In any case, we need rhop
[9212]	446
[10222]	447	IF(ln_ctl) THEN ! print control
	448	CALL prt_ctl(tab3d_1=tsn(:,:,:,jp_tem), clinfo1=' tn - : ', mask1=tmask, kdim=jpk )
	449	CALL prt_ctl(tab3d_1=tsn(:,:,:,jp_sal), clinfo1=' sn - : ', mask1=tmask, kdim=jpk )
	450	ENDIF
	451	!
	452	IF( ln_timing ) CALL timing_stop( 'dta_dyn_sed')
	453	!
	454	END SUBROUTINE dta_dyn_sed
	455
	456
	457	SUBROUTINE dta_dyn_sed_init
	458	!!----------------------------------------------------------------------
	459	!! * ROUTINE dta_dyn_init *
	460	!!
	461	!! ** Purpose : Initialisation of the dynamical data
	462	!! ** Method : - read the data namdta_dyn namelist
	463	!!----------------------------------------------------------------------
	464	INTEGER :: ierr, ierr0, ierr1, ierr2, ierr3 ! return error code
	465	INTEGER :: ifpr ! dummy loop indice
	466	INTEGER :: jfld ! dummy loop arguments
	467	INTEGER :: inum, idv, idimv ! local integer
	468	INTEGER :: ios ! Local integer output status for namelist read
	469	!!
	470	CHARACTER(len=100) :: cn_dir ! Root directory for location of core files
	471	TYPE(FLD_N), DIMENSION(2) :: slf_d ! array of namelist informations on the fields to read
	472	TYPE(FLD_N) :: sn_tem , sn_sal ! " "
	473	!!
	474	NAMELIST/namdta_dyn/cn_dir, ln_dynrnf, ln_dynrnf_depth, fwbcorr, &
	475	& sn_tem, sn_sal
	476	!!----------------------------------------------------------------------
	477	!
	478	REWIND( numnam_ref ) ! Namelist namdta_dyn in reference namelist : Offline: init. of dynamical data
	479	READ ( numnam_ref, namdta_dyn, IOSTAT = ios, ERR = 901)
[11536]	480	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdta_dyn in reference namelist' )
[10222]	481	REWIND( numnam_cfg ) ! Namelist namdta_dyn in configuration namelist : Offline: init. of dynamical data
	482	READ ( numnam_cfg, namdta_dyn, IOSTAT = ios, ERR = 902 )
[11536]	483	902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namdta_dyn in configuration namelist' )
[10222]	484	IF(lwm) WRITE ( numond, namdta_dyn )
	485	! ! store namelist information in an array
	486	! ! Control print
	487	IF(lwp) THEN
	488	WRITE(numout,*)
	489	WRITE(numout,*) 'dta_dyn : offline dynamics '
	490	WRITE(numout,*) '~~~~~~~ '
	491	WRITE(numout,*) ' Namelist namdta_dyn'
	492	WRITE(numout,*) ' runoffs option enabled (T) or not (F) ln_dynrnf = ', ln_dynrnf
	493	WRITE(numout,*) ' runoffs is spread in vertical ln_dynrnf_depth = ', ln_dynrnf_depth
	494	WRITE(numout,*) ' annual global mean of empmr for ssh correction fwbcorr = ', fwbcorr
	495	WRITE(numout,*)
	496	ENDIF
	497	!
	498	jf_tem = 1 ; jf_sal = 2 ; jfld = jf_sal
	499	!
	500	slf_d(jf_tem) = sn_tem ; slf_d(jf_sal) = sn_sal
	501	!
	502	ALLOCATE( sf_dyn(jfld), STAT=ierr ) ! set sf structure
	503	IF( ierr > 0 ) THEN
	504	CALL ctl_stop( 'dta_dyn: unable to allocate sf structure' ) ; RETURN
	505	ENDIF
	506	! ! fill sf with slf_i and control print
	507	CALL fld_fill( sf_dyn, slf_d, cn_dir, 'dta_dyn_init', 'Data in file', 'namdta_dyn' )
	508	!
	509	! Open file for each variable to get his number of dimension
	510	DO ifpr = 1, jfld
[12116]	511	CALL fld_def( sf_dyn(ifpr) )
	512	CALL iom_open( sf_dyn(ifpr)%clname, sf_dyn(ifpr)%num )
[10222]	513	idv = iom_varid( sf_dyn(ifpr)%num , slf_d(ifpr)%clvar ) ! id of the variable sdjf%clvar
	514	idimv = iom_file ( sf_dyn(ifpr)%num )%ndims(idv) ! number of dimension for variable sdjf%clvar
[12116]	515	CALL iom_close( sf_dyn(ifpr)%num ) ! close file if already open
[10222]	516	ierr1=0
	517	IF( idimv == 3 ) THEN ! 2D variable
	518	ALLOCATE( sf_dyn(ifpr)%fnow(jpi,jpj,1) , STAT=ierr0 )
	519	IF( slf_d(ifpr)%ln_tint ) ALLOCATE( sf_dyn(ifpr)%fdta(jpi,jpj,1,2) , STAT=ierr1 )
	520	ELSE ! 3D variable
	521	ALLOCATE( sf_dyn(ifpr)%fnow(jpi,jpj,jpk) , STAT=ierr0 )
	522	IF( slf_d(ifpr)%ln_tint ) ALLOCATE( sf_dyn(ifpr)%fdta(jpi,jpj,jpk,2), STAT=ierr1 )
	523	ENDIF
	524	IF( ierr0 + ierr1 > 0 ) THEN
	525	CALL ctl_stop( 'dta_dyn_init : unable to allocate sf_dyn array structure' ) ; RETURN
	526	ENDIF
	527	END DO
	528	!
	529	CALL dta_dyn_sed( nit000 )
	530	!
	531	END SUBROUTINE dta_dyn_sed_init
	532
[7646]	533	SUBROUTINE dta_dyn_swp( kt )
	534	!!---------------------------------------------------------------------
	535	!! * ROUTINE dta_dyn_swp *
	536	!!
[9212]	537	!! ** Purpose : Swap and the data and compute the vertical scale factor
	538	!! at U/V/W pointand the depht
[7646]	539	!!---------------------------------------------------------------------
	540	INTEGER, INTENT(in) :: kt ! time step
[9212]	541	!
[7646]	542	INTEGER :: ji, jj, jk
	543	REAL(wp) :: zcoef
	544	!!---------------------------------------------------------------------
[6140]	545
[7646]	546	IF( kt == nit000 ) THEN
	547	IF(lwp) WRITE(numout,*)
	548	IF(lwp) WRITE(numout,*) 'ssh_swp : Asselin time filter and swap of sea surface height'
	549	IF(lwp) WRITE(numout,*) '~~~~~~~ '
	550	ENDIF
	551
	552	sshb(:,:) = sshn(:,:) + atfp * ( sshb(:,:) - 2 * sshn(:,:) + ssha(:,:)) ! before <-- now filtered
	553	sshn(:,:) = ssha(:,:)
	554
	555	e3t_n(:,:,:) = e3t_a(:,:,:)
	556
	557	! Reconstruction of all vertical scale factors at now and before time steps
	558	! =============================================================================
	559
	560	! Horizontal scale factor interpolations
	561	! --------------------------------------
	562	CALL dom_vvl_interpol( e3t_n(:,:,:), e3u_n(:,:,:), 'U' )
	563	CALL dom_vvl_interpol( e3t_n(:,:,:), e3v_n(:,:,:), 'V' )
	564
	565	! Vertical scale factor interpolations
	566	! ------------------------------------
	567	CALL dom_vvl_interpol( e3t_n(:,:,:), e3w_n (:,:,:), 'W' )
	568
	569	e3t_b(:,:,:) = e3t_n(:,:,:)
	570	e3u_b(:,:,:) = e3u_n(:,:,:)
	571	e3v_b(:,:,:) = e3v_n(:,:,:)
	572
	573	! t- and w- points depth
	574	! ----------------------
	575	gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1)
	576	gdepw_n(:,:,1) = 0.0_wp
[9212]	577	!
[7646]	578	DO jk = 2, jpk
	579	DO jj = 1,jpj
	580	DO ji = 1,jpi
[9212]	581	zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
	582	gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1)
	583	gdept_n(ji,jj,jk) = zcoef * ( gdepw_n(ji,jj,jk ) + 0.5 * e3w_n(ji,jj,jk)) &
	584	& + (1-zcoef) * ( gdept_n(ji,jj,jk-1) + e3w_n(ji,jj,jk))
	585	END DO
	586	END DO
	587	END DO
	588	!
[7646]	589	gdept_b(:,:,:) = gdept_n(:,:,:)
	590	gdepw_b(:,:,:) = gdepw_n(:,:,:)
	591	!
	592	END SUBROUTINE dta_dyn_swp
[9212]	593
[7646]	594
	595	SUBROUTINE dta_dyn_ssh( kt, phdivtr, psshb, pemp, pssha, pe3ta )
[1501]	596	!!----------------------------------------------------------------------
[7646]	597	!! * ROUTINE dta_dyn_wzv *
	598	!!
	599	!! ** Purpose : compute the after ssh (ssha) and the now vertical velocity
[1501]	600	!!
[7646]	601	!! ** Method : Using the incompressibility hypothesis,
	602	!! - the ssh increment is computed by integrating the horizontal divergence
	603	!! and multiply by the time step.
[1501]	604	!!
[7646]	605	!! - compute the after scale factor : repartition of ssh INCREMENT proportionnaly
	606	!! to the level thickness ( z-star case )
	607	!!
	608	!! - the vertical velocity is computed by integrating the horizontal divergence
	609	!! from the bottom to the surface minus the scale factor evolution.
	610	!! The boundary conditions are w=0 at the bottom (no flux)
	611	!!
	612	!! ** action : ssha / e3t_a / wn
	613	!!
	614	!! Reference : Leclair, M., and G. Madec, 2009, Ocean Modelling.
[2528]	615	!!----------------------------------------------------------------------
[7646]	616	INTEGER, INTENT(in ) :: kt ! time-step
	617	REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: phdivtr ! horizontal divergence transport
	618	REAL(wp), DIMENSION(jpi,jpj) , OPTIONAL, INTENT(in ) :: psshb ! now ssh
	619	REAL(wp), DIMENSION(jpi,jpj) , OPTIONAL, INTENT(in ) :: pemp ! evaporation minus precipitation
	620	REAL(wp), DIMENSION(jpi,jpj) , OPTIONAL, INTENT(inout) :: pssha ! after ssh
	621	REAL(wp), DIMENSION(jpi,jpj,jpk), OPTIONAL, INTENT(out) :: pe3ta ! after vertical scale factor
[9212]	622	!
[7646]	623	INTEGER :: jk
	624	REAL(wp), DIMENSION(jpi,jpj) :: zhdiv
	625	REAL(wp) :: z2dt
	626	!!----------------------------------------------------------------------
[3294]	627	!
[7646]	628	z2dt = 2._wp * rdt
	629	!
	630	zhdiv(:,:) = 0._wp
	631	DO jk = 1, jpkm1
	632	zhdiv(:,:) = zhdiv(:,:) + phdivtr(:,:,jk) * tmask(:,:,jk)
	633	END DO
	634	! ! Sea surface elevation time-stepping
	635	pssha(:,:) = ( psshb(:,:) - z2dt * ( r1_rau0 * pemp(:,:) + zhdiv(:,:) ) ) * ssmask(:,:)
	636	! !
	637	! ! After acale factors at t-points ( z_star coordinate )
	638	DO jk = 1, jpkm1
	639	pe3ta(:,:,jk) = e3t_0(:,:,jk) * ( 1._wp + pssha(:,:) * tmask(:,:,1) / ( ht_0(:,:) + 1.0 - tmask(:,:,1) ) )
	640	END DO
	641	!
	642	END SUBROUTINE dta_dyn_ssh
	643
	644
	645	SUBROUTINE dta_dyn_hrnf
	646	!!----------------------------------------------------------------------
	647	!! * ROUTINE sbc_rnf *
[1501]	648	!!
[7646]	649	!! ** Purpose : update the horizontal divergence with the runoff inflow
	650	!!
	651	!! ** Method :
	652	!! CAUTION : rnf is positive (inflow) decreasing the
	653	!! divergence and expressed in m/s
	654	!!
	655	!! ** Action : phdivn decreased by the runoff inflow
[2528]	656	!!----------------------------------------------------------------------
[7646]	657	!!
	658	INTEGER :: ji, jj, jk ! dummy loop indices
	659	!!----------------------------------------------------------------------
[2528]	660	!
[7646]	661	DO jj = 1, jpj ! update the depth over which runoffs are distributed
	662	DO ji = 1, jpi
	663	h_rnf(ji,jj) = 0._wp
	664	DO jk = 1, nk_rnf(ji,jj) ! recalculates h_rnf to be the depth in metres
	665	h_rnf(ji,jj) = h_rnf(ji,jj) + e3t_n(ji,jj,jk) ! to the bottom of the relevant grid box
[1501]	666	END DO
[7646]	667	END DO
[2528]	668	END DO
[5836]	669	!
[7646]	670	END SUBROUTINE dta_dyn_hrnf
	671
	672
	673
	674	SUBROUTINE dta_dyn_slp( kt )
	675	!!---------------------------------------------------------------------
	676	!! * ROUTINE dta_dyn_slp *
	677	!!
	678	!! ** Purpose : Computation of slope
	679	!!
	680	!!---------------------------------------------------------------------
	681	INTEGER, INTENT(in) :: kt ! time step
	682	!
	683	INTEGER :: ji, jj ! dummy loop indices
	684	REAL(wp) :: ztinta ! ratio applied to after records when doing time interpolation
	685	REAL(wp) :: ztintb ! ratio applied to before records when doing time interpolation
	686	INTEGER :: iswap
[9212]	687	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zuslp, zvslp, zwslpi, zwslpj
	688	REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: zts
[7646]	689	!!---------------------------------------------------------------------
	690	!
	691	IF( sf_dyn(jf_tem)%ln_tint ) THEN ! Computes slopes (here avt is used as workspace)
	692	IF( kt == nit000 ) THEN
	693	IF(lwp) WRITE(numout,*) ' Compute new slopes at kt = ', kt
	694	zts(:,:,:,jp_tem) = sf_dyn(jf_tem)%fdta(:,:,:,1) * tmask(:,:,:) ! temperature
	695	zts(:,:,:,jp_sal) = sf_dyn(jf_sal)%fdta(:,:,:,1) * tmask(:,:,:) ! salinity
	696	avt(:,:,:) = sf_dyn(jf_avt)%fdta(:,:,:,1) * tmask(:,:,:) ! vertical diffusive coef.
	697	CALL compute_slopes( kt, zts, zuslp, zvslp, zwslpi, zwslpj )
	698	uslpdta (:,:,:,1) = zuslp (:,:,:)
	699	vslpdta (:,:,:,1) = zvslp (:,:,:)
	700	wslpidta(:,:,:,1) = zwslpi(:,:,:)
	701	wslpjdta(:,:,:,1) = zwslpj(:,:,:)
	702	!
	703	zts(:,:,:,jp_tem) = sf_dyn(jf_tem)%fdta(:,:,:,2) * tmask(:,:,:) ! temperature
	704	zts(:,:,:,jp_sal) = sf_dyn(jf_sal)%fdta(:,:,:,2) * tmask(:,:,:) ! salinity
	705	avt(:,:,:) = sf_dyn(jf_avt)%fdta(:,:,:,2) * tmask(:,:,:) ! vertical diffusive coef.
	706	CALL compute_slopes( kt, zts, zuslp, zvslp, zwslpi, zwslpj )
	707	uslpdta (:,:,:,2) = zuslp (:,:,:)
	708	vslpdta (:,:,:,2) = zvslp (:,:,:)
	709	wslpidta(:,:,:,2) = zwslpi(:,:,:)
	710	wslpjdta(:,:,:,2) = zwslpj(:,:,:)
	711	ELSE
	712	!
	713	iswap = 0
	714	IF( sf_dyn(jf_tem)%nrec_a(2) - nprevrec /= 0 ) iswap = 1
	715	IF( nsecdyn > sf_dyn(jf_tem)%nrec_b(2) .AND. iswap == 1 ) THEN ! read/update the after data
	716	IF(lwp) WRITE(numout,*) ' Compute new slopes at kt = ', kt
	717	uslpdta (:,:,:,1) = uslpdta (:,:,:,2) ! swap the data
	718	vslpdta (:,:,:,1) = vslpdta (:,:,:,2)
	719	wslpidta(:,:,:,1) = wslpidta(:,:,:,2)
	720	wslpjdta(:,:,:,1) = wslpjdta(:,:,:,2)
	721	!
	722	zts(:,:,:,jp_tem) = sf_dyn(jf_tem)%fdta(:,:,:,2) * tmask(:,:,:) ! temperature
	723	zts(:,:,:,jp_sal) = sf_dyn(jf_sal)%fdta(:,:,:,2) * tmask(:,:,:) ! salinity
	724	avt(:,:,:) = sf_dyn(jf_avt)%fdta(:,:,:,2) * tmask(:,:,:) ! vertical diffusive coef.
	725	CALL compute_slopes( kt, zts, zuslp, zvslp, zwslpi, zwslpj )
	726	!
	727	uslpdta (:,:,:,2) = zuslp (:,:,:)
	728	vslpdta (:,:,:,2) = zvslp (:,:,:)
	729	wslpidta(:,:,:,2) = zwslpi(:,:,:)
	730	wslpjdta(:,:,:,2) = zwslpj(:,:,:)
	731	ENDIF
	732	ENDIF
	733	ENDIF
	734	!
	735	IF( sf_dyn(jf_tem)%ln_tint ) THEN
	736	ztinta = REAL( nsecdyn - sf_dyn(jf_tem)%nrec_b(2), wp ) &
	737	& / REAL( sf_dyn(jf_tem)%nrec_a(2) - sf_dyn(jf_tem)%nrec_b(2), wp )
	738	ztintb = 1. - ztinta
	739	IF( l_ldfslp .AND. .NOT.lk_c1d ) THEN ! Computes slopes (here avt is used as workspace)
	740	uslp (:,:,:) = ztintb * uslpdta (:,:,:,1) + ztinta * uslpdta (:,:,:,2)
	741	vslp (:,:,:) = ztintb * vslpdta (:,:,:,1) + ztinta * vslpdta (:,:,:,2)
	742	wslpi(:,:,:) = ztintb * wslpidta(:,:,:,1) + ztinta * wslpidta(:,:,:,2)
	743	wslpj(:,:,:) = ztintb * wslpjdta(:,:,:,1) + ztinta * wslpjdta(:,:,:,2)
	744	ENDIF
	745	ELSE
	746	zts(:,:,:,jp_tem) = sf_dyn(jf_tem)%fnow(:,:,:) * tmask(:,:,:) ! temperature
	747	zts(:,:,:,jp_sal) = sf_dyn(jf_sal)%fnow(:,:,:) * tmask(:,:,:) ! salinity
	748	avt(:,:,:) = sf_dyn(jf_avt)%fnow(:,:,:) * tmask(:,:,:) ! vertical diffusive coef.
	749	CALL compute_slopes( kt, zts, zuslp, zvslp, zwslpi, zwslpj )
	750	!
	751	IF( l_ldfslp .AND. .NOT.lk_c1d ) THEN ! Computes slopes (here avt is used as workspace)
	752	uslp (:,:,:) = zuslp (:,:,:)
	753	vslp (:,:,:) = zvslp (:,:,:)
	754	wslpi(:,:,:) = zwslpi(:,:,:)
	755	wslpj(:,:,:) = zwslpj(:,:,:)
	756	ENDIF
	757	ENDIF
	758	!
	759	END SUBROUTINE dta_dyn_slp
[1501]	760
[9212]	761
[7646]	762	SUBROUTINE compute_slopes( kt, pts, puslp, pvslp, pwslpi, pwslpj )
[1501]	763	!!---------------------------------------------------------------------
[3294]	764	!! * ROUTINE dta_dyn_slp *
[1501]	765	!!
[9212]	766	!! ** Purpose : Computation of slope
[1501]	767	!!---------------------------------------------------------------------
[3294]	768	INTEGER , INTENT(in ) :: kt ! time step
	769	REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! temperature/salinity
	770	REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(out) :: puslp ! zonal isopycnal slopes
	771	REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(out) :: pvslp ! meridional isopycnal slopes
	772	REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(out) :: pwslpi ! zonal diapycnal slopes
	773	REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(out) :: pwslpj ! meridional diapycnal slopes
[1501]	774	!!---------------------------------------------------------------------
[9212]	775	!
[7646]	776	IF( l_ldfslp .AND. .NOT.lk_c1d ) THEN ! Computes slopes (here avt is used as workspace)
[5836]	777	CALL eos ( pts, rhd, rhop, gdept_0(:,:,:) )
	778	CALL eos_rab( pts, rab_n ) ! now local thermal/haline expension ratio at T-points
	779	CALL bn2 ( pts, rab_n, rn2 ) ! now Brunt-Vaisala
[5131]	780
[6140]	781	! Partial steps: before Horizontal DErivative
	782	IF( ln_zps .AND. .NOT. ln_isfcav) &
	783	& CALL zps_hde ( kt, jpts, pts, gtsu, gtsv, & ! Partial steps: before horizontal gradient
	784	& rhd, gru , grv ) ! of t, s, rd at the last ocean level
	785	IF( ln_zps .AND. ln_isfcav) &
[7646]	786	& CALL zps_hde_isf( kt, jpts, pts, gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
	787	& rhd, gru , grv , grui, grvi ) ! of t, s, rd at the first ocean level
[4990]	788
[5836]	789	rn2b(:,:,:) = rn2(:,:,:) ! need for zdfmxl
	790	CALL zdf_mxl( kt ) ! mixed layer depth
	791	CALL ldf_slp( kt, rhd, rn2 ) ! slopes
[7646]	792	puslp (:,:,:) = uslp (:,:,:)
	793	pvslp (:,:,:) = vslp (:,:,:)
	794	pwslpi(:,:,:) = wslpi(:,:,:)
	795	pwslpj(:,:,:) = wslpj(:,:,:)
[5836]	796	ELSE
	797	puslp (:,:,:) = 0. ! to avoid warning when compiling
	798	pvslp (:,:,:) = 0.
	799	pwslpi(:,:,:) = 0.
	800	pwslpj(:,:,:) = 0.
	801	ENDIF
[2528]	802	!
[7646]	803	END SUBROUTINE compute_slopes
[9212]	804
[2528]	805	!!======================================================================
[325]	806	END MODULE dtadyn

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

Download in other formats: