Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

diaar5.F90 @ 7910

Last change on this file since 7910 was 7910, checked in by timgraham, 7 years ago
All wrk_alloc removed
Property svn:keywords set to `Id`
File size: 15.9 KB

Rev	Line
[1756]	1	MODULE diaar5
	2	!!======================================================================
	3	!! * MODULE diaar5 *
	4	!! AR5 diagnostics
	5	!!======================================================================
[2528]	6	!! History : 3.2 ! 2009-11 (S. Masson) Original code
	7	!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase + merge TRC-TRA
[1756]	8	!!----------------------------------------------------------------------
[2528]	9	!! dia_ar5 : AR5 diagnostics
	10	!! dia_ar5_init : initialisation of AR5 diagnostics
[1756]	11	!!----------------------------------------------------------------------
	12	USE oce ! ocean dynamics and active tracers
	13	USE dom_oce ! ocean space and time domain
[2528]	14	USE eosbn2 ! equation of state (eos_bn2 routine)
[1756]	15	USE lib_mpp ! distribued memory computing library
	16	USE iom ! I/O manager library
[3294]	17	USE timing ! preformance summary
[5253]	18	USE fldread ! type FLD_N
	19	USE phycst ! physical constant
	20	USE in_out_manager ! I/O manager
[7646]	21	USE zdfddm
	22	USE zdf_oce
[1756]	23
	24	IMPLICIT NONE
	25	PRIVATE
	26
[2528]	27	PUBLIC dia_ar5 ! routine called in step.F90 module
[2715]	28	PUBLIC dia_ar5_alloc ! routine called in nemogcm.F90 module
[7646]	29	PUBLIC dia_ar5_hst ! heat/salt transport
[1756]	30
[2528]	31	REAL(wp) :: vol0 ! ocean volume (interior domain)
	32	REAL(wp) :: area_tot ! total ocean surface (interior domain)
[2715]	33	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,: ) :: area ! cell surface (interior domain)
	34	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,: ) :: thick0 ! ocean thickness (interior domain)
	35	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sn0 ! initial salinity
[7646]	36
	37	LOGICAL :: l_ar5
[1756]	38
[7646]	39	!! * Substitutions
	40	# include "zdfddm_substitute.h90"
	41	# include "vectopt_loop_substitute.h90"
[1756]	42	!!----------------------------------------------------------------------
[2528]	43	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
	44	!! $Id$
	45	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1756]	46	!!----------------------------------------------------------------------
	47	CONTAINS
	48
[2715]	49	FUNCTION dia_ar5_alloc()
	50	!!----------------------------------------------------------------------
	51	!! * ROUTINE dia_ar5_alloc *
	52	!!----------------------------------------------------------------------
	53	INTEGER :: dia_ar5_alloc
	54	!!----------------------------------------------------------------------
	55	!
	56	ALLOCATE( area(jpi,jpj), thick0(jpi,jpj) , sn0(jpi,jpj,jpk) , STAT=dia_ar5_alloc )
	57	!
	58	IF( lk_mpp ) CALL mpp_sum ( dia_ar5_alloc )
	59	IF( dia_ar5_alloc /= 0 ) CALL ctl_warn('dia_ar5_alloc: failed to allocate arrays')
	60	!
	61	END FUNCTION dia_ar5_alloc
	62
	63
[1756]	64	SUBROUTINE dia_ar5( kt )
	65	!!----------------------------------------------------------------------
	66	!! * ROUTINE dia_ar5 *
	67	!!
[2528]	68	!! ** Purpose : compute and output some AR5 diagnostics
[1756]	69	!!----------------------------------------------------------------------
[2715]	70	!
[1756]	71	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[2715]	72	!
[1756]	73	INTEGER :: ji, jj, jk ! dummy loop arguments
	74	REAL(wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass
[7646]	75	REAL(wp) :: zaw, zbw, zrw
[3294]	76	!
[7910]	77	REAL(wp), DIMENSION(jpi,jpj) :: zarea_ssh , zbotpres ! 2D workspace
	78	REAL(wp), DIMENSION(jpi,jpj) :: zpe ! 2D workspace
	79	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrhd , zrhop ! 3D workspace
	80	REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: ztsn ! 4D workspace
[1756]	81	!!--------------------------------------------------------------------
[3294]	82	IF( nn_timing == 1 ) CALL timing_start('dia_ar5')
	83
[7646]	84	IF( kt == nit000 ) CALL dia_ar5_init
[1756]	85
[7646]	86	IF( l_ar5 ) THEN
[7753]	87	zarea_ssh(:,:) = area(:,:) * sshn(:,:)
[7646]	88	ENDIF
	89	!
	90	IF( iom_use( 'voltot' ) .OR. iom_use( 'sshtot' ) .OR. iom_use( 'sshdyn' ) ) THEN
	91	! ! total volume of liquid seawater
	92	zvolssh = SUM( zarea_ssh(:,:) )
	93	IF( lk_mpp ) CALL mpp_sum( zvolssh )
	94	zvol = vol0 + zvolssh
[1756]	95
[7646]	96	CALL iom_put( 'voltot', zvol )
	97	CALL iom_put( 'sshtot', zvolssh / area_tot )
	98	CALL iom_put( 'sshdyn', sshn(:,:) - (zvolssh / area_tot) )
	99	!
	100	ENDIF
[1756]	101
[7646]	102	IF( iom_use( 'botpres' ) .OR. iom_use( 'sshthster' ) .OR. iom_use( 'sshsteric' ) ) THEN
	103	!
[7753]	104	ztsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) ! thermosteric ssh
	105	ztsn(:,:,:,jp_sal) = sn0(:,:,:)
[7646]	106	CALL eos( ztsn, zrhd, gdept_n(:,:,:) ) ! now in situ density using initial salinity
	107	!
[7753]	108	zbotpres(:,:) = 0._wp ! no atmospheric surface pressure, levitating sea-ice
[7646]	109	DO jk = 1, jpkm1
[7753]	110	zbotpres(:,:) = zbotpres(:,:) + e3t_n(:,:,jk) * zrhd(:,:,jk)
[7646]	111	END DO
	112	IF( ln_linssh ) THEN
	113	IF( ln_isfcav ) THEN
	114	DO ji = 1, jpi
	115	DO jj = 1, jpj
	116	zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj)
	117	END DO
[5120]	118	END DO
[7646]	119	ELSE
[7753]	120	zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1)
[7646]	121	END IF
[6140]	122	!!gm
	123	!!gm riceload should be added in both ln_linssh=T or F, no?
	124	!!gm
[7646]	125	END IF
	126	!
[7753]	127	zarho = SUM( area(:,:) * zbotpres(:,:) )
[7646]	128	IF( lk_mpp ) CALL mpp_sum( zarho )
	129	zssh_steric = - zarho / area_tot
	130	CALL iom_put( 'sshthster', zssh_steric )
[1756]	131
[7646]	132	! ! steric sea surface height
	133	CALL eos( tsn, zrhd, zrhop, gdept_n(:,:,:) ) ! now in situ and potential density
[7753]	134	zrhop(:,:,jpk) = 0._wp
[7646]	135	CALL iom_put( 'rhop', zrhop )
	136	!
[7753]	137	zbotpres(:,:) = 0._wp ! no atmospheric surface pressure, levitating sea-ice
[7646]	138	DO jk = 1, jpkm1
[7753]	139	zbotpres(:,:) = zbotpres(:,:) + e3t_n(:,:,jk) * zrhd(:,:,jk)
[7646]	140	END DO
	141	IF( ln_linssh ) THEN
	142	IF ( ln_isfcav ) THEN
	143	DO ji = 1,jpi
	144	DO jj = 1,jpj
	145	zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj)
	146	END DO
[5120]	147	END DO
[7646]	148	ELSE
[7753]	149	zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1)
[7646]	150	END IF
[5120]	151	END IF
[7646]	152	!
[7753]	153	zarho = SUM( area(:,:) * zbotpres(:,:) )
[7646]	154	IF( lk_mpp ) CALL mpp_sum( zarho )
	155	zssh_steric = - zarho / area_tot
	156	CALL iom_put( 'sshsteric', zssh_steric )
[1756]	157
[7646]	158	! ! ocean bottom pressure
	159	zztmp = rau0 * grav * 1.e-4_wp ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa
[7753]	160	zbotpres(:,:) = zztmp * ( zbotpres(:,:) + sshn(:,:) + thick0(:,:) )
[7646]	161	CALL iom_put( 'botpres', zbotpres )
	162	!
	163	ENDIF
[1756]	164
[7646]	165	IF( iom_use( 'masstot' ) .OR. iom_use( 'temptot' ) .OR. iom_use( 'saltot' ) ) THEN
	166	! ! Mean density anomalie, temperature and salinity
	167	ztemp = 0._wp
	168	zsal = 0._wp
	169	DO jk = 1, jpkm1
	170	DO jj = 1, jpj
	171	DO ji = 1, jpi
	172	zztmp = area(ji,jj) * e3t_n(ji,jj,jk)
	173	ztemp = ztemp + zztmp * tsn(ji,jj,jk,jp_tem)
	174	zsal = zsal + zztmp * tsn(ji,jj,jk,jp_sal)
	175	END DO
[1756]	176	END DO
	177	END DO
[7646]	178	IF( ln_linssh ) THEN
	179	IF( ln_isfcav ) THEN
	180	DO ji = 1, jpi
	181	DO jj = 1, jpj
	182	ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem)
	183	zsal = zsal + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal)
	184	END DO
[5120]	185	END DO
[7646]	186	ELSE
	187	ztemp = ztemp + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_tem) )
	188	zsal = zsal + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_sal) )
	189	END IF
	190	ENDIF
	191	IF( lk_mpp ) THEN
	192	CALL mpp_sum( ztemp )
	193	CALL mpp_sum( zsal )
[5120]	194	END IF
[7646]	195	!
	196	zmass = rau0 * ( zarho + zvol ) ! total mass of liquid seawater
	197	ztemp = ztemp / zvol ! potential temperature in liquid seawater
	198	zsal = zsal / zvol ! Salinity of liquid seawater
	199	!
	200	CALL iom_put( 'masstot', zmass )
	201	CALL iom_put( 'temptot', ztemp )
	202	CALL iom_put( 'saltot' , zsal )
	203	!
[1756]	204	ENDIF
[7646]	205
	206	IF( iom_use( 'tnpeo' )) THEN
	207	! Work done against stratification by vertical mixing
	208	! Exclude points where rn2 is negative as convection kicks in here and
	209	! work is not being done against stratification
[7753]	210	zpe(:,:) = 0._wp
[7646]	211	IF( lk_zdfddm ) THEN
	212	DO ji=1,jpi
	213	DO jj=1,jpj
	214	DO jk=1,jpk
	215	zrw = ( gdepw_n(ji,jj,jk ) - gdept_n(ji,jj,jk) ) &
	216	& / ( gdept_n(ji,jj,jk-1) - gdept_n(ji,jj,jk) )
	217	!
	218	zaw = rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem)* zrw
	219	zbw = rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal)* zrw
	220	!
	221	zpe(ji, jj) = zpe(ji, jj) - MIN(0._wp, rn2(ji,jj,jk)) * &
	222	& grav * (avt(ji,jj,jk) * zaw * (tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) &
	223	& - fsavs(ji,jj,jk) * zbw * (tsn(ji,jj,jk-1,jp_sal) - tsn(ji,jj,jk,jp_sal) ) )
	224
	225	ENDDO
	226	ENDDO
	227	ENDDO
	228	ELSE
	229	DO ji = 1, jpi
	230	DO jj = 1, jpj
	231	DO jk = 1, jpk
	232	zpe(ji,jj) = zpe(ji,jj) + avt(ji, jj, jk) * MIN(0._wp,rn2(ji, jj, jk)) * rau0 * e3w_n(ji, jj, jk)
	233	ENDDO
	234	ENDDO
	235	ENDDO
	236	ENDIF
	237	CALL lbc_lnk( zpe, 'T', 1._wp)
	238	CALL iom_put( 'tnpeo', zpe )
	239	ENDIF
[1756]	240	!
[7646]	241	IF( l_ar5 ) THEN
	242	ENDIF
[1756]	243	!
[3294]	244	IF( nn_timing == 1 ) CALL timing_stop('dia_ar5')
	245	!
[1756]	246	END SUBROUTINE dia_ar5
	247
[7646]	248	SUBROUTINE dia_ar5_hst( ktra, cptr, pua, pva )
	249	!!----------------------------------------------------------------------
	250	!! * ROUTINE dia_ar5_htr *
	251	!!----------------------------------------------------------------------
	252	!! Wrapper for heat transport calculations
	253	!! Called from all advection and/or diffusion routines
	254	!!----------------------------------------------------------------------
	255	INTEGER , INTENT(in ) :: ktra ! tracer index
	256	CHARACTER(len=3) , INTENT(in) :: cptr ! transport type 'adv'/'ldf'
	257	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: pua ! 3D input array of advection/diffusion
	258	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: pva ! 3D input array of advection/diffusion
	259	!
	260	INTEGER :: ji, jj, jk
[7910]	261	REAL(wp), DIMENSION(jpi,jpj) :: z2d
[1756]	262
[7646]	263
	264
	265	z2d(:,:) = pua(:,:,1)
	266	DO jk = 1, jpkm1
	267	DO jj = 2, jpjm1
	268	DO ji = fs_2, fs_jpim1 ! vector opt.
	269	z2d(ji,jj) = z2d(ji,jj) + pua(ji,jj,jk)
	270	END DO
	271	END DO
	272	END DO
	273	CALL lbc_lnk( z2d, 'U', -1. )
	274	IF( cptr == 'adv' ) THEN
	275	IF( ktra == jp_tem ) CALL iom_put( "uadv_heattr" , rau0_rcp * z2d ) ! advective heat transport in i-direction
	276	IF( ktra == jp_sal ) CALL iom_put( "uadv_salttr" , rau0 * z2d ) ! advective salt transport in i-direction
	277	ENDIF
	278	IF( cptr == 'ldf' ) THEN
	279	IF( ktra == jp_tem ) CALL iom_put( "udiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in i-direction
	280	IF( ktra == jp_sal ) CALL iom_put( "udiff_salttr" , rau0 * z2d ) ! diffusive salt transport in i-direction
	281	ENDIF
	282	!
	283	z2d(:,:) = pva(:,:,1)
	284	DO jk = 1, jpkm1
	285	DO jj = 2, jpjm1
	286	DO ji = fs_2, fs_jpim1 ! vector opt.
	287	z2d(ji,jj) = z2d(ji,jj) + pva(ji,jj,jk)
	288	END DO
	289	END DO
	290	END DO
	291	CALL lbc_lnk( z2d, 'V', -1. )
	292	IF( cptr == 'adv' ) THEN
	293	IF( ktra == jp_tem ) CALL iom_put( "vadv_heattr" , rau0_rcp * z2d ) ! advective heat transport in j-direction
	294	IF( ktra == jp_sal ) CALL iom_put( "vadv_salttr" , rau0 * z2d ) ! advective salt transport in j-direction
	295	ENDIF
	296	IF( cptr == 'ldf' ) THEN
	297	IF( ktra == jp_tem ) CALL iom_put( "vdiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in j-direction
	298	IF( ktra == jp_sal ) CALL iom_put( "vdiff_salttr" , rau0 * z2d ) ! diffusive salt transport in j-direction
	299	ENDIF
	300
	301
	302	END SUBROUTINE dia_ar5_hst
	303
	304
[1756]	305	SUBROUTINE dia_ar5_init
	306	!!----------------------------------------------------------------------
	307	!! * ROUTINE dia_ar5_init *
	308	!!
[2528]	309	!! ** Purpose : initialization for AR5 diagnostic computation
[1756]	310	!!----------------------------------------------------------------------
	311	INTEGER :: inum
	312	INTEGER :: ik
	313	INTEGER :: ji, jj, jk ! dummy loop indices
[7753]	314	REAL(wp) :: zztmp
[7910]	315	REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: zsaldta ! Jan/Dec levitus salinity
[5253]	316	!
[1756]	317	!!----------------------------------------------------------------------
	318	!
[3294]	319	IF( nn_timing == 1 ) CALL timing_start('dia_ar5_init')
	320	!
[7646]	321	l_ar5 = .FALSE.
	322	IF( iom_use( 'voltot' ) .OR. iom_use( 'sshtot' ) .OR. iom_use( 'sshdyn' ) .OR. &
	323	& iom_use( 'masstot' ) .OR. iom_use( 'temptot' ) .OR. iom_use( 'saltot' ) .OR. &
	324	& iom_use( 'botpres' ) .OR. iom_use( 'sshthster' ) .OR. iom_use( 'sshsteric' ) ) L_ar5 = .TRUE.
	325
	326	IF( l_ar5 ) THEN
	327	!
	328	! ! allocate dia_ar5 arrays
	329	IF( dia_ar5_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' )
[2715]	330
[7753]	331	area(:,:) = e1e2t(:,:) * tmask_i(:,:)
[1756]	332
[7646]	333	area_tot = SUM( area(:,:) ) ; IF( lk_mpp ) CALL mpp_sum( area_tot )
[1756]	334
[7646]	335	vol0 = 0._wp
[7753]	336	thick0(:,:) = 0._wp
[7646]	337	DO jk = 1, jpkm1
[7753]	338	vol0 = vol0 + SUM( area (:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) )
	339	thick0(:,:) = thick0(:,:) + tmask_i(:,:) * tmask(:,:,jk) * e3t_0(:,:,jk)
[7646]	340	END DO
	341	IF( lk_mpp ) CALL mpp_sum( vol0 )
[5253]	342
[7646]	343	CALL iom_open ( 'sali_ref_clim_monthly', inum )
	344	CALL iom_get ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,1), 1 )
	345	CALL iom_get ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,2), 12 )
	346	CALL iom_close( inum )
[6665]	347
[7753]	348	sn0(:,:,:) = 0.5_wp * ( zsaldta(:,:,:,1) + zsaldta(:,:,:,2) )
	349	sn0(:,:,:) = sn0(:,:,:) * tmask(:,:,:)
[7646]	350	IF( ln_zps ) THEN ! z-coord. partial steps
	351	DO jj = 1, jpj ! interpolation of salinity at the last ocean level (i.e. the partial step)
	352	DO ji = 1, jpi
	353	ik = mbkt(ji,jj)
	354	IF( ik > 1 ) THEN
	355	zztmp = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
	356	sn0(ji,jj,ik) = ( 1._wp - zztmp ) * sn0(ji,jj,ik) + zztmp * sn0(ji,jj,ik-1)
	357	ENDIF
	358	END DO
[1756]	359	END DO
[7646]	360	ENDIF
	361	!
	362	!
[1756]	363	ENDIF
	364	!
[3294]	365	IF( nn_timing == 1 ) CALL timing_stop('dia_ar5_init')
	366	!
[1756]	367	END SUBROUTINE dia_ar5_init
	368
	369	!!======================================================================
	370	END MODULE diaar5

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

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

Context Navigation

source: branches/2017/dev_r7881_no_wrk_alloc/NEMOGCM/NEMO/OPA_SRC/DIA/diaar5.F90 @ 7910

Download in other formats: