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.

trcrst.F90 in trunk/NEMO/TOP_SRC – NEMO

Context Navigation

source: trunk/NEMO/TOP_SRC/trcrst.F90 @ 959

Last change on this file since 959 was 945, checked in by cetlod, 16 years ago
Update modules for new version of TOP model, see ticket 144
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 15.9 KB

Rev	Line
[268]	1	MODULE trcrst
[335]	2	!!======================================================================
[945]	3	!! * MODULE trcrst *
	4	!! TOP : create, write, read the restart files for passive tracers
[335]	5	!!======================================================================
[945]	6	!! History : 1.0 ! 2007-02 (C. Ethe) adaptation from the ocean
[274]	7	!!----------------------------------------------------------------------
[945]	8	#if defined key_top
[335]	9	!!----------------------------------------------------------------------
[945]	10	!! 'key_top' TOP models
	11	!!----------------------------------------------------------------------
	12	!! trc_rst_opn : open restart file
	13	!! trc_rst_read : read restart file
	14	!! trc_rst_wri : write restart file
	15	!!----------------------------------------------------------------------
[335]	16	USE oce_trc
	17	USE trc
	18	USE sms
[945]	19	USE trcsms_cfc ! CFC variables
[335]	20	USE trctrp_lec
[433]	21	USE lib_mpp
[616]	22	USE iom
[335]	23
	24	IMPLICIT NONE
	25	PRIVATE
	26
[945]	27	PUBLIC trc_rst_opn ! called by ???
	28	PUBLIC trc_rst_read ! called by ???
	29	PUBLIC trc_rst_wri ! called by ???
[335]	30
[616]	31	LOGICAL, PUBLIC :: lrst_trc !: logical to control the trc restart write
	32	INTEGER, PUBLIC :: numrtr, numrtw !: logical unit for trc restart (read and write)
[350]	33
[945]	34	REAL(wp) :: &
	35	alkmean = 2426. , & ! mean value of alkalinity ( Glodap ; for Goyet 2391. )
	36	po4mean = 2.165 , & ! mean value of phosphates
	37	no3mean = 30.90 , & ! mean value of nitrate
	38	siomean = 91.51 ! mean value of silicate
[616]	39
[350]	40	!! * Substitutions
[945]	41	# include "top_substitute.h90"
	42	!!----------------------------------------------------------------------
	43	!! NEMO/TOP 1.0 , LOCEAN-IPSL (2005)
	44	!! $Id$
	45	!! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
	46	!!----------------------------------------------------------------------
[335]	47
[268]	48	CONTAINS
[335]	49
[616]	50	SUBROUTINE trc_rst_opn( kt )
	51	!!----------------------------------------------------------------------
	52	!! * trc_rst_opn *
	53	!!
	54	!! ** purpose : output of sea-trc variable in a netcdf file
	55	!!----------------------------------------------------------------------
	56	INTEGER, INTENT(in) :: kt ! number of iteration
	57	!
	58	CHARACTER(LEN=20) :: clkt ! ocean time-step define as a character
	59	CHARACTER(LEN=50) :: clname ! trc output restart file name
	60	!!----------------------------------------------------------------------
	61	!
	62	IF( kt == nit000 ) THEN
	63	lrst_trc = .FALSE.
[945]	64	# if defined key_off_tra
[616]	65	nitrst = nitend ! in online version, already done in rst_opn routine defined in restart.F90 module
[945]	66	# endif
[616]	67	ENDIF
	68
	69	IF( kt == nitrst - ndttrc .OR. nitend - nit000 + 1 < 2 * ndttrc ) THEN
	70	! beware if model runs less than 2*ndttrc time step
	71	! beware of the format used to write kt (default is i8.8, that should be large enough)
[945]	72	IF( nitrst > 1.0e9 ) THEN ; WRITE(clkt,* ) nitrst
	73	ELSE ; WRITE(clkt,'(i8.8)') nitrst
[616]	74	ENDIF
	75	! create the file
	76	IF(lwp) WRITE(numout,*)
	77	clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_restart_trc"
	78	IF(lwp) WRITE(numout,*) ' open trc restart.output NetCDF file: '//clname
	79	CALL iom_open( clname, numrtw, ldwrt = .TRUE., kiolib = jprstlib )
	80	lrst_trc = .TRUE.
	81	ENDIF
	82	!
	83	END SUBROUTINE trc_rst_opn
	84
	85
	86	SUBROUTINE trc_rst_read
[945]	87	!!----------------------------------------------------------------------
	88	!! * trc_rst_opn *
[335]	89	!!
[945]	90	!! ** purpose : read passive tracer fields in restart files
	91	!!----------------------------------------------------------------------
	92	INTEGER :: jn
	93	INTEGER :: iarak0
	94	REAL(wp) :: zkt, zarak0
	95	# if defined key_pisces
	96	REAL(wp) :: ztrasum
	97	INTEGER :: ji, jj, jk
	98	REAL(wp) :: caralk, bicarb, co3
	99	# endif
	100	!!----------------------------------------------------------------------
[268]	101
[945]	102	IF(lwp) WRITE(numout,*)
	103	IF(lwp) WRITE(numout,*) 'trc_rst_read : read the TOP restart file'
	104	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
[350]	105
[945]	106	IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN ; iarak0 = 1
	107	ELSE ; iarak0 = 0
[268]	108	ENDIF
	109
[945]	110	IF(lwp) WRITE(numout,*)
	111	IF(lwp) WRITE(numout,*) ' the present run starts at the time step nit000 = ', nit000
	112	IF(lwp .AND. iarak0 == 1 ) WRITE(numout,*) ' and needs previous fields for Arakawa sheme ??? '
[268]	113
	114
	115	! Time domain : restart
	116	! -------------------------
	117	IF(lwp) WRITE(numout,*)
	118	IF(lwp) WRITE(numout,) ' ** passive tracer restart option'
	119	SELECT CASE ( nrsttr )
	120	CASE ( 0 )
[945]	121	IF(lwp) WRITE(numout,*) ' nrsttr = 0 no control of nit000'
[268]	122	CASE ( 1 )
[945]	123	IF(lwp) WRITE(numout,*) ' nrsttr = 1 we control the date of nit000'
[268]	124	CASE ( 2 )
[945]	125	IF(lwp) WRITE(numout,*) ' nrsttr = 2 the date adatrj is read in restart file'
[268]	126	CASE DEFAULT
	127	IF(lwp) WRITE(numout,*) ' ===>>>> nrsttr not equal 0, 1 or 2 : no control of the date'
[945]	128	IF(lwp) WRITE(numout,*) ' ======= ========='
[268]	129	END SELECT
	130
[945]	131	CALL iom_open( 'restart.trc', numrtr, kiolib = jprstlib )
[268]	132
[616]	133	CALL iom_get( numrtr, 'kt' , zkt )
	134	CALL iom_get( numrtr, 'arak0', zarak0 )
[268]	135
[494]	136	IF(lwp) WRITE(numout,*)
	137	IF(lwp) WRITE(numout,*) ' Info on the restart file read : '
[945]	138	IF(lwp) WRITE(numout,*) ' time-step : ', NINT( zkt )
	139	IF(lwp) WRITE(numout,*) ' arakawa option : ', NINT( zarak0 )
[494]	140
	141
[945]	142	IF( nittrc000 - NINT( zkt ) /= 1 .AND. nrsttr /= 0 ) & ! control of date
	143	& CALL ctl_stop( ' ===>>>> : problem with nit000 for the restart', &
	144	& ' verify the restart file or rerun with nrstdt = 0 (namelist)' )
[268]	145
[945]	146	IF( iarak0 /= NINT( zarak0 ) ) & ! Control of the scheme
	147	& CALL ctl_stop( ' ===>>>> : problem with advection scheme', &
	148	& ' it must be the same type for both restart and previous run', &
	149	& ' centered or euler ' )
[268]	150
	151
[945]	152	! READ prognostic variables and computes diagnostic variable
[494]	153	DO jn = 1, jptra
[945]	154	CALL iom_get( numrtr, jpdom_local, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
	155	CALL iom_get( numrtr, jpdom_local, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
[335]	156	END DO
[945]	157	# if defined key_lobster
[616]	158	CALL iom_get( numrtr, jpdom_local, 'SEDB'//ctrcnm(jpdet), sedpocb(:,:) )
	159	CALL iom_get( numrtr, jpdom_local, 'SEDN'//ctrcnm(jpdet), sedpocn(:,:) )
[945]	160	# elif defined key_pisces
[616]	161	CALL iom_get( numrtr, jpdom_local, 'Silicalim', xksi(:,:) )
[945]	162	CALL iom_get( numrtr, jpdom_local, 'Silicamax', xksimax(:,:) )
	163	# elif defined key_cfc
[561]	164	DO jn = 1, jptra
[945]	165	CALL iom_get( numrtr, jpdom_local, 'qint'//ctrcnm(jn), qint(:,:,jn) )
	166	CALL iom_get( numrtr, jpdom_local, 'qtr'//ctrcnm(jn) , qtr( :,:,jn) )
[616]	167	END DO
[945]	168	# endif
[268]	169
[945]	170	# if defined key_pisces
	171	! ! --------------------------- !
	172	IF( cp_cfg == "orca" .AND. .NOT. lk_trc_c1d ) THEN ! ORCA condiguration (not 1D) !
	173	! ! --------------------------- !
	174	! set total alkalinity, phosphate, NO3 & silicate
	175	! total alkalinity
	176	! -----------------------------------------------
	177	ztrasum = 0.e0
	178	DO jk = 1, jpk
	179	DO jj = 1, jpj
	180	DO ji = 1, jpi
	181	ztrasum = ztrasum + trn(ji,jj,jk,jptal) * tmask(ji,jj,jk) * tmask_i(ji,jj) &
	182	# if defined key_off_degrad
	183	& * facvol(ji,jj,jk) &
	184	# endif
	185	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	186	END DO
[350]	187	END DO
	188	END DO
[945]	189	IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain
[350]	190
[433]	191
[945]	192	ztrasum = ztrasum / areatot * 1.e6
	193	IF(lwp) WRITE(numout,*) 'TALK moyen ', ztrasum
	194	trn(:,:,:,jptal) = trn(:,:,:,jptal) * alkmean / ztrasum
[350]	195
[945]	196	! phosphate
	197	! ---------
	198	ztrasum = 0.e0
	199	DO jk = 1, jpk
	200	DO jj = 1, jpj
	201	DO ji = 1, jpi
	202	ztrasum = ztrasum + trn(ji,jj,jk,jppo4) * tmask(ji,jj,jk) * tmask_i(ji,jj) &
	203	# if defined key_off_degrad
	204	& * facvol(ji,jj,jk) &
	205	# endif
	206	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	207	END DO
[350]	208	END DO
	209	END DO
[945]	210	IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain
[350]	211
[945]	212	ztrasum = ztrasum / areatot * 1.e6 / 122.
	213	IF(lwp) WRITE(numout,*) 'PO4 moyen ', ztrasum
	214	trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * po4mean / ztrasum
[433]	215
[945]	216	! NO3
	217	! ---
	218	ztrasum = 0.e0
	219	DO jk = 1, jpk
	220	DO jj = 1, jpj
	221	DO ji = 1, jpi
	222	ztrasum = ztrasum + trn(ji,jj,jk,jpno3) * tmask(ji,jj,jk) * tmask_i(ji,jj) &
	223	# if defined key_off_degrad
	224	& * facvol(ji,jj,jk) &
	225	# endif
	226	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	227	END DO
[350]	228	END DO
	229	END DO
[945]	230	IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain
[350]	231
[945]	232	ztrasum = ztrasum / areatot * 1.e6 / 7.6
	233	IF(lwp) WRITE(numout,*) 'NO3 moyen ', ztrasum
	234	trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * no3mean / ztrasum
[433]	235
[945]	236	! Silicate
	237	! --------
	238	ztrasum = 0.e0
	239	DO jk = 1, jpk
	240	DO jj = 1, jpj
	241	DO ji = 1, jpi
	242	ztrasum = ztrasum + trn(ji,jj,jk,jpsil) * tmask(ji,jj,jk) * tmask_i(ji,jj) &
	243	# if defined key_off_degrad
	244	& * facvol(ji,jj,jk) &
	245	# endif
	246	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
	247	END DO
[350]	248	END DO
	249	END DO
[945]	250	IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain
[350]	251
[945]	252	IF(lwp) WRITE(numout,) 'SiO3 moyen ', ztrasum / areatot 1.e6
	253	ztrasum = ztrasum / areatot * 1.e6
	254	trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * siomean / ztrasum )
	255	!
	256	ENDIF
[433]	257
[945]	258	!#if defined key_kriest
[616]	259	! !! Initialize number of particles from a standart restart file
	260	! !! The name of big organic particles jpgoc has been only change
	261	! !! and replace by jpnum but the values here are concentration
	262	! trn(:,:,:,jppoc) = trn(:,:,:,jppoc) + trn(:,:,:,jpnum)
	263	! trn(:,:,:,jpnum) = trn(:,:,:,jppoc) / ( 6. * xkr_massp )
	264	!#endif
[945]	265	!! Set hi (???) from total alcalinity, borat (???), akb3 (???) and ak23 (???)
	266	!! ---------------------------------------------------------------------
	267	DO jk = 1, jpk
	268	DO jj = 1, jpj
[268]	269	DO ji = 1,jpi
[945]	270	caralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.e-8 / ( rtrn + akb3(ji,jj,jk) ) )
	271	co3 = ( caralk - trn(ji,jj,jk,jpdic) ) * tmask(ji,jj,jk) &
	272	& + 0.5e-3 * ( 1.- tmask(ji,jj,jk) )
	273	bicarb = 2.* trn(ji,jj,jk,jpdic) - caralk
	274	hi(ji,jj,jk) = ( ak23(ji,jj,jk) * bicarb / co3 ) * tmask(ji,jj,jk) &
	275	& + 1.0e-9 * ( 1.- tmask(ji,jj,jk) )
	276	END DO
	277	END DO
	278	END DO
	279	# endif
[268]	280
[616]	281	CALL iom_close( numrtr )
[945]	282	!
	283	END SUBROUTINE trc_rst_read
[494]	284
	285
[945]	286	SUBROUTINE trc_rst_wri( kt )
	287	!!----------------------------------------------------------------------
	288	!! * trc_rst_wri *
[335]	289	!!
[945]	290	!! ** purpose : write passive tracer fields in restart files
	291	!!----------------------------------------------------------------------
	292	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[335]	293	!!
[945]	294	INTEGER :: ji, jj, jk, jn
[616]	295	REAL(wp) :: zdiag_var, zdiag_varmin, zdiag_varmax, zdiag_tot
	296	REAL(wp) :: zder
[945]	297	!!----------------------------------------------------------------------
[268]	298
[945]	299	IF( MOD( kt, nstock ) == 0 .OR. kt == nitend ) THEN
[268]	300
[945]	301	! 0. initialisations
	302	! ------------------
	303	IF(lwp) WRITE(numout,*)
	304	IF(lwp) WRITE(numout,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp
	305	IF(lwp) WRITE(numout,*) '~~~~~~~'
[268]	306
	307
[616]	308	CALL iom_rstput( kt, nitrst, numrtw, 'kt' , REAL( kt, wp ) )
	309
[335]	310	IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN
[616]	311	CALL iom_rstput( kt, nitrst, numrtw, 'arak0', 1. )
[335]	312	ELSE
[616]	313	CALL iom_rstput( kt, nitrst, numrtw, 'arak0', 0. )
[335]	314	ENDIF
[268]	315
	316
[616]	317	! prognostic variables
	318	! --------------------
[945]	319	DO jn = 1, jptra
[616]	320	CALL iom_rstput( kt, nitrst, numrtw, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) )
	321	CALL iom_rstput( kt, nitrst, numrtw, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) )
[335]	322	END DO
[268]	323
[945]	324	#if defined key_lobster
[616]	325	CALL iom_rstput( kt, nitrst, numrtw, 'SEDB'//ctrcnm(jpdet), sedpocb(:,:) )
	326	CALL iom_rstput( kt, nitrst, numrtw, 'SEDN'//ctrcnm(jpdet), sedpocn(:,:) )
[945]	327	#elif defined key_pisces
[616]	328	CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) )
[945]	329	CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) )
[268]	330
[616]	331	#elif defined key_cfc
[945]	332	DO jn = 1, jptra
[616]	333	CALL iom_rstput( kt, nitrst, numrtw, 'qint'//ctrcnm(jn), qint(:,:,jn) )
[945]	334	CALL iom_rstput( kt, nitrst, numrtw, 'qtr'//ctrcnm(jn) , qtr (:,:,jn) )
[616]	335	END DO
	336	#endif
[268]	337
[945]	338	IF(lwp) WRITE(numout,*) '----TRACER STAT----'
[616]	339
[945]	340	zdiag_tot = 0.e0
	341	DO jn = 1, jptra
	342	zdiag_var = 0.e0
	343	zdiag_varmin = 0.e0
	344	zdiag_varmax = 0.e0
	345	DO ji = 1, jpi
	346	DO jj = 1, jpj
	347	DO jk = 1,jpk
	348	zdiag_var = zdiag_var + trn(ji,jj,jk,jn) * tmask(ji,jj,jk) * tmask_i(ji,jj) &
[616]	349	#if defined key_off_degrad
	350	& * facvol(ji,jj,jk) &
	351	#endif
	352	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
[335]	353	END DO
	354	END DO
[268]	355	END DO
	356
[945]	357	zdiag_varmin = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
	358	zdiag_varmax = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) )
[268]	359
[616]	360	IF( lk_mpp ) THEN
[945]	361	CALL mpp_min( zdiag_varmin ) ! min over the global domain
	362	CALL mpp_max( zdiag_varmax ) ! max over the global domain
	363	CALL mpp_sum( zdiag_var ) ! sum over the global domain
[433]	364	END IF
[268]	365
[945]	366	zdiag_tot = zdiag_tot + zdiag_var
	367	zdiag_var = zdiag_var / areatot
[433]	368
[945]	369	IF(lwp) WRITE(numout,*) ' MEAN NO ', jn, ctrcnm(jn), ' = ', zdiag_var, &
	370	& ' MIN = ', zdiag_varmin, ' MAX = ', zdiag_varmax
[268]	371	END DO
	372
[945]	373	zder = ( ( zdiag_tot - trai ) / ( trai + 1.e-12 ) ) * 100._wp
	374	IF(lwp) WRITE(numout,*) ' Integral of all tracers over the full domain = ', zdiag_tot
	375	IF(lwp) WRITE(numout,*) ' Drift of the sum of all tracers =', zder, ' %'
[433]	376
[616]	377	CALL iom_close(numrtw)
[945]	378	!
[268]	379	ENDIF
[945]	380	!
[616]	381	END SUBROUTINE trc_rst_wri
[268]	382
	383	#else
[945]	384	!!----------------------------------------------------------------------
	385	!! Dummy module : No passive tracer
	386	!!----------------------------------------------------------------------
[335]	387	CONTAINS
[945]	388	SUBROUTINE trc_rst_read ! Empty routines
[616]	389	END SUBROUTINE trc_rst_read
[945]	390	SUBROUTINE trc_rst_wri( kt )
[335]	391	INTEGER, INTENT ( in ) :: kt
[616]	392	WRITE(,) 'trc_rst_wri: You should not have seen this print! error?', kt
[945]	393	END SUBROUTINE trc_rst_wri
[268]	394	#endif
[945]	395
	396	!!======================================================================
[335]	397	END MODULE trcrst

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

Download in other formats: