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.

traadv_ubs.F90 in NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/OCE/TRA – NEMO

Context Navigation

source: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/OCE/TRA/traadv_ubs.F90 @ 12546

Last change on this file since 12546 was 12546, checked in by orioltp, 4 years ago
Adding precision specification in hardcoded reals and other modifications to allow compilation without forcing reals without precision specification to a certain value through compiler flags
Property svn:keywords set to `Id`
File size: 18.4 KB

Rev	Line
[503]	1	MODULE traadv_ubs
	2	!!==============================================================================
	3	!! * MODULE traadv_ubs *
	4	!! Ocean active tracers: horizontal & vertical advective trend
	5	!!==============================================================================
[2528]	6	!! History : 1.0 ! 2006-08 (L. Debreu, R. Benshila) Original code
	7	!! 3.3 ! 2010-05 (C. Ethe, G. Madec) merge TRC-TRA + switch from velocity to transport
[503]	8	!!----------------------------------------------------------------------
	9
	10	!!----------------------------------------------------------------------
	11	!! tra_adv_ubs : update the tracer trend with the horizontal
	12	!! advection trends using a third order biaised scheme
	13	!!----------------------------------------------------------------------
[3625]	14	USE oce ! ocean dynamics and active tracers
	15	USE dom_oce ! ocean space and time domain
[4990]	16	USE trc_oce ! share passive tracers/Ocean variables
	17	USE trd_oce ! trends: ocean variables
[5836]	18	USE traadv_fct ! acces to routine interp_4th_cpt
[4990]	19	USE trdtra ! trends manager: tracers
	20	USE diaptr ! poleward transport diagnostics
[7646]	21	USE diaar5 ! AR5 diagnostics
[4990]	22	!
[9019]	23	USE iom ! I/O library
[9124]	24	USE in_out_manager ! I/O manager
[9019]	25	USE lib_mpp ! massively parallel library
[3625]	26	USE lbclnk ! ocean lateral boundary condition (or mpp link)
	27	USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
[503]	28
	29	IMPLICIT NONE
	30	PRIVATE
	31
	32	PUBLIC tra_adv_ubs ! routine called by traadv module
	33
[7646]	34	LOGICAL :: l_trd ! flag to compute trends
	35	LOGICAL :: l_ptr ! flag to compute poleward transport
	36	LOGICAL :: l_hst ! flag to compute heat transport
[503]	37
[7646]	38
[503]	39	!! * Substitutions
[12377]	40	# include "do_loop_substitute.h90"
[503]	41	!!----------------------------------------------------------------------
[9598]	42	!! NEMO/OCE 4.0 , NEMO Consortium (2018)
[1152]	43	!! $Id$
[10068]	44	!! Software governed by the CeCILL license (see ./LICENSE)
[503]	45	!!----------------------------------------------------------------------
	46	CONTAINS
	47
[12377]	48	SUBROUTINE tra_adv_ubs( kt, kit000, cdtype, p2dt, pU, pV, pW, &
	49	& Kbb, Kmm, pt, kjpt, Krhs, kn_ubs_v )
[503]	50	!!----------------------------------------------------------------------
	51	!! * ROUTINE tra_adv_ubs *
	52	!!
	53	!! ** Purpose : Compute the now trend due to the advection of tracers
	54	!! and add it to the general trend of passive tracer equations.
	55	!!
[5836]	56	!! ** Method : The 3rd order Upstream Biased Scheme (UBS) is based on an
[3787]	57	!! upstream-biased parabolic interpolation (Shchepetkin and McWilliams 2005)
[519]	58	!! It is only used in the horizontal direction.
	59	!! For example the i-component of the advective fluxes are given by :
[3787]	60	!! ! e2u e3u un ( mi(Tn) - zltu(i ) ) if un(i) >= 0
[4990]	61	!! ztu = ! or
[3787]	62	!! ! e2u e3u un ( mi(Tn) - zltu(i+1) ) if un(i) < 0
[519]	63	!! where zltu is the second derivative of the before temperature field:
	64	!! zltu = 1/e3t di[ e2u e3u / e1u di[Tb] ]
[5836]	65	!! This results in a dissipatively dominant (i.e. hyper-diffusive)
[519]	66	!! truncation error. The overall performance of the advection scheme
	67	!! is similar to that reported in (Farrow and Stevens, 1995).
[5836]	68	!! For stability reasons, the first term of the fluxes which corresponds
[519]	69	!! to a second order centered scheme is evaluated using the now velocity
	70	!! (centered in time) while the second term which is the diffusive part
	71	!! of the scheme, is evaluated using the before velocity (forward in time).
	72	!! Note that UBS is not positive. Do not use it on passive tracers.
[5836]	73	!! On the vertical, the advection is evaluated using a FCT scheme,
	74	!! as the UBS have been found to be too diffusive.
[6140]	75	!! kn_ubs_v argument controles whether the FCT is based on
	76	!! a 2nd order centrered scheme (kn_ubs_v=2) or on a 4th order compact
	77	!! scheme (kn_ubs_v=4).
[503]	78	!!
[12377]	79	!! ** Action : - update pt(:,:,:,:,Krhs) with the now advective tracer trends
[6140]	80	!! - send trends to trdtra module for further diagnostcs (l_trdtra=T)
[12377]	81	!! - poleward advective heat and salt transport (ln_diaptr=T)
[519]	82	!!
	83	!! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404.
	84	!! Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731Ð1741.
[503]	85	!!----------------------------------------------------------------------
[12377]	86	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	87	INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! ocean time level indices
	88	INTEGER , INTENT(in ) :: kit000 ! first time step index
	89	CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
	90	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	91	INTEGER , INTENT(in ) :: kn_ubs_v ! number of tracers
	92	REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step
	93	REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume transport components
	94	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation
[2715]	95	!
	96	INTEGER :: ji, jj, jk, jn ! dummy loop indices
[6140]	97	REAL(wp) :: ztra, zbtr, zcoef ! local scalars
[2715]	98	REAL(wp) :: zfp_ui, zfm_ui, zcenut, ztak, zfp_wk, zfm_wk ! - -
	99	REAL(wp) :: zfp_vj, zfm_vj, zcenvt, zeeu, zeev, z_hdivn ! - -
[9019]	100	REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztu, ztv, zltu, zltv, zti, ztw ! 3D workspace
[503]	101	!!----------------------------------------------------------------------
[3294]	102	!
	103	IF( kt == kit000 ) THEN
[503]	104	IF(lwp) WRITE(numout,*)
[2528]	105	IF(lwp) WRITE(numout,*) 'tra_adv_ubs : horizontal UBS advection scheme on ', cdtype
[503]	106	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
	107	ENDIF
[2528]	108	!
[4499]	109	l_trd = .FALSE.
[7646]	110	l_hst = .FALSE.
	111	l_ptr = .FALSE.
	112	IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.
[12377]	113	IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.
[7646]	114	IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
	115	& iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.
[4499]	116	!
[6140]	117	ztw (:,:, 1 ) = 0._wp ! surface & bottom value : set to zero for all tracers
	118	zltu(:,:,jpk) = 0._wp ; zltv(:,:,jpk) = 0._wp
[5836]	119	ztw (:,:,jpk) = 0._wp ; zti (:,:,jpk) = 0._wp
	120	!
[2528]	121	! ! ===========
	122	DO jn = 1, kjpt ! tracer loop
	123	! ! ===========
	124	!
[5836]	125	DO jk = 1, jpkm1 !== horizontal laplacian of before tracer ==!
[12377]	126	DO_2D_10_10
	127	zeeu = e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) * umask(ji,jj,jk)
	128	zeev = e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) * vmask(ji,jj,jk)
	129	ztu(ji,jj,jk) = zeeu * ( pt(ji+1,jj ,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) )
	130	ztv(ji,jj,jk) = zeev * ( pt(ji ,jj+1,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) )
	131	END_2D
	132	DO_2D_00_00
	133	zcoef = 1._wp / ( 6._wp * e3t(ji,jj,jk,Kmm) )
	134	zltu(ji,jj,jk) = ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zcoef
	135	zltv(ji,jj,jk) = ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zcoef
	136	END_2D
[2528]	137	!
[5836]	138	END DO
[12546]	139	CALL lbc_lnk( 'traadv_ubs', zltu, 'T', 1.0_wp ) ; CALL lbc_lnk( 'traadv_ubs', zltv, 'T', 1.0_wp ) ! Lateral boundary cond. (unchanged sgn)
[2528]	140	!
[12377]	141	DO_3D_10_10( 1, jpkm1 )
	142	zfp_ui = pU(ji,jj,jk) + ABS( pU(ji,jj,jk) ) ! upstream transport (x2)
	143	zfm_ui = pU(ji,jj,jk) - ABS( pU(ji,jj,jk) )
	144	zfp_vj = pV(ji,jj,jk) + ABS( pV(ji,jj,jk) )
	145	zfm_vj = pV(ji,jj,jk) - ABS( pV(ji,jj,jk) )
	146	! ! 2nd order centered advective fluxes (x2)
	147	zcenut = pU(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) )
	148	zcenvt = pV(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) )
	149	! ! UBS advective fluxes
	150	ztu(ji,jj,jk) = 0.5 * ( zcenut - zfp_ui * zltu(ji,jj,jk) - zfm_ui * zltu(ji+1,jj,jk) )
	151	ztv(ji,jj,jk) = 0.5 * ( zcenvt - zfp_vj * zltv(ji,jj,jk) - zfm_vj * zltv(ji,jj+1,jk) )
	152	END_3D
[5836]	153	!
[12377]	154	zltu(:,:,:) = pt(:,:,:,jn,Krhs) ! store the initial trends before its update
[5836]	155	!
	156	DO jk = 1, jpkm1 !== add the horizontal advective trend ==!
[12377]	157	DO_2D_00_00
	158	pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) &
	159	& - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk) &
	160	& + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
	161	END_2D
[2528]	162	!
[5836]	163	END DO
	164	!
[12377]	165	zltu(:,:,:) = pt(:,:,:,jn,Krhs) - zltu(:,:,:) ! Horizontal advective trend used in vertical 2nd order FCT case
[5836]	166	! ! and/or in trend diagnostic (l_trd=T)
[4990]	167	!
	168	IF( l_trd ) THEN ! trend diagnostics
[12377]	169	CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztu, pU, pt(:,:,:,jn,Kmm) )
	170	CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztv, pV, pt(:,:,:,jn,Kmm) )
[2528]	171	END IF
[7646]	172	!
	173	! ! "Poleward" heat and salt transports (contribution of upstream fluxes)
	174	IF( l_ptr ) CALL dia_ptr_hst( jn, 'adv', ztv(:,:,:) )
	175	! ! heati/salt transport
	176	IF( l_hst ) CALL dia_ar5_hst( jn, 'adv', ztu(:,:,:), ztv(:,:,:) )
[5836]	177	!
[7646]	178	!
[5836]	179	! !== vertical advective trend ==!
	180	!
	181	SELECT CASE( kn_ubs_v ) ! select the vertical advection scheme
	182	!
	183	CASE( 2 ) ! 2nd order FCT
	184	!
[12377]	185	IF( l_trd ) zltv(:,:,:) = pt(:,:,:,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag.
[5836]	186	!
	187	! !* upstream advection with initial mass fluxes & intermediate update ==!
[12377]	188	DO_3D_11_11( 2, jpkm1 )
	189	zfp_wk = pW(ji,jj,jk) + ABS( pW(ji,jj,jk) )
	190	zfm_wk = pW(ji,jj,jk) - ABS( pW(ji,jj,jk) )
	191	ztw(ji,jj,jk) = 0.5_wp * ( zfp_wk * pt(ji,jj,jk,jn,Kbb) + zfm_wk * pt(ji,jj,jk-1,jn,Kbb) ) * wmask(ji,jj,jk)
	192	END_3D
[6140]	193	IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as ztw has been w-masked)
[5836]	194	IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface
[12377]	195	DO_2D_11_11
	196	ztw(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) ! linear free surface
	197	END_2D
[5836]	198	ELSE ! no cavities: only at the ocean surface
[12377]	199	ztw(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb)
[5836]	200	ENDIF
	201	ENDIF
	202	!
[12377]	203	DO_3D_00_00( 1, jpkm1 )
	204	ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
	205	pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztak
	206	zti(ji,jj,jk) = ( pt(ji,jj,jk,jn,Kbb) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk)
	207	END_3D
[12546]	208	CALL lbc_lnk( 'traadv_ubs', zti, 'T', 1.0_wp ) ! Lateral boundary conditions on zti, zsi (unchanged sign)
[5836]	209	!
	210	! !* anti-diffusive flux : high order minus low order
[12377]	211	DO_3D_11_11( 2, jpkm1 )
	212	ztw(ji,jj,jk) = ( 0.5_wp * pW(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) ) &
	213	& - ztw(ji,jj,jk) ) * wmask(ji,jj,jk)
	214	END_3D
[5836]	215	! ! top ocean value: high order == upstream ==>> zwz=0
[6140]	216	IF( ln_linssh ) ztw(:,:, 1 ) = 0._wp ! only ocean surface as interior zwz values have been w-masked
[5836]	217	!
[12377]	218	CALL nonosc_z( Kmm, pt(:,:,:,jn,Kbb), ztw, zti, p2dt ) ! monotonicity algorithm
[5836]	219	!
	220	CASE( 4 ) ! 4th order COMPACT
[12377]	221	CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw ) ! 4th order compact interpolation of T at w-point
	222	DO_3D_00_00( 2, jpkm1 )
	223	ztw(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk)
	224	END_3D
	225	IF( ln_linssh ) ztw(:,:, 1 ) = pW(:,:,1) * pt(:,:,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work
[5836]	226	!
	227	END SELECT
[2528]	228	!
[12377]	229	DO_3D_00_00( 1, jpkm1 )
	230	pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
	231	END_3D
[5836]	232	!
	233	IF( l_trd ) THEN ! vertical advective trend diagnostics
[12377]	234	DO_3D_00_00( 1, jpkm1 )
	235	zltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk) &
	236	& + pt(ji,jj,jk,jn,Kmm) * ( pW(ji,jj,jk) - pW(ji,jj,jk+1) ) &
	237	& * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
	238	END_3D
	239	CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zltv )
[2528]	240	ENDIF
	241	!
[4990]	242	END DO
[503]	243	!
[2528]	244	END SUBROUTINE tra_adv_ubs
[503]	245
	246
[12377]	247	SUBROUTINE nonosc_z( Kmm, pbef, pcc, paft, p2dt )
[503]	248	!!---------------------------------------------------------------------
	249	!! * ROUTINE nonosc_z *
	250	!!
	251	!! ** Purpose : compute monotonic tracer fluxes from the upstream
	252	!! scheme and the before field by a nonoscillatory algorithm
	253	!!
	254	!! ** Method : ... ???
	255	!! warning : pbef and paft must be masked, but the boundaries
	256	!! conditions on the fluxes are not necessary zalezak (1979)
	257	!! drange (1995) multi-dimensional forward-in-time and upstream-
	258	!! in-space based differencing for fluid
	259	!!----------------------------------------------------------------------
[12377]	260	INTEGER , INTENT(in ) :: Kmm ! time level index
[6140]	261	REAL(wp), INTENT(in ) :: p2dt ! tracer time-step
[2528]	262	REAL(wp), DIMENSION (jpi,jpj,jpk) :: pbef ! before field
[503]	263	REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: paft ! after field
	264	REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: pcc ! monotonic flux in the k direction
[2715]	265	!
	266	INTEGER :: ji, jj, jk ! dummy loop indices
	267	INTEGER :: ikm1 ! local integer
[6140]	268	REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn ! local scalars
[9019]	269	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zbetup, zbetdo ! 3D workspace
[503]	270	!!----------------------------------------------------------------------
[3294]	271	!
[2715]	272	zbig = 1.e+40_wp
	273	zrtrn = 1.e-15_wp
	274	zbetup(:,:,:) = 0._wp ; zbetdo(:,:,:) = 0._wp
[5836]	275	!
[503]	276	! Search local extrema
	277	! --------------------
[5836]	278	! ! large negative value (-zbig) inside land
[503]	279	pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) )
	280	paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) )
[5836]	281	!
	282	DO jk = 1, jpkm1 ! search maximum in neighbourhood
[503]	283	ikm1 = MAX(jk-1,1)
[12377]	284	DO_2D_00_00
	285	zbetup(ji,jj,jk) = MAX( pbef(ji ,jj ,jk ), paft(ji ,jj ,jk ), &
	286	& pbef(ji ,jj ,ikm1), pbef(ji ,jj ,jk+1), &
	287	& paft(ji ,jj ,ikm1), paft(ji ,jj ,jk+1) )
	288	END_2D
[503]	289	END DO
[5836]	290	! ! large positive value (+zbig) inside land
[503]	291	pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )
	292	paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )
[5836]	293	!
	294	DO jk = 1, jpkm1 ! search minimum in neighbourhood
[503]	295	ikm1 = MAX(jk-1,1)
[12377]	296	DO_2D_00_00
	297	zbetdo(ji,jj,jk) = MIN( pbef(ji ,jj ,jk ), paft(ji ,jj ,jk ), &
	298	& pbef(ji ,jj ,ikm1), pbef(ji ,jj ,jk+1), &
	299	& paft(ji ,jj ,ikm1), paft(ji ,jj ,jk+1) )
	300	END_2D
[503]	301	END DO
[5836]	302	! ! restore masked values to zero
[503]	303	pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:)
	304	paft(:,:,:) = paft(:,:,:) * tmask(:,:,:)
[5836]	305	!
	306	! Positive and negative part of fluxes and beta terms
	307	! ---------------------------------------------------
[12377]	308	DO_3D_00_00( 1, jpkm1 )
	309	! positive & negative part of the flux
	310	zpos = MAX( 0., pcc(ji ,jj ,jk+1) ) - MIN( 0., pcc(ji ,jj ,jk ) )
	311	zneg = MAX( 0., pcc(ji ,jj ,jk ) ) - MIN( 0., pcc(ji ,jj ,jk+1) )
	312	! up & down beta terms
	313	zbt = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) / p2dt
	314	zbetup(ji,jj,jk) = ( zbetup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt
	315	zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zbetdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt
	316	END_3D
[5836]	317	!
[503]	318	! monotonic flux in the k direction, i.e. pcc
	319	! -------------------------------------------
[12377]	320	DO_3D_00_00( 2, jpkm1 )
	321	za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj,jk-1) )
	322	zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj,jk-1) )
[12546]	323	zc = 0.5 * ( 1.e0 + SIGN( 1.0_wp, pcc(ji,jj,jk) ) )
[12377]	324	pcc(ji,jj,jk) = pcc(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb )
	325	END_3D
[503]	326	!
	327	END SUBROUTINE nonosc_z
	328
	329	!!======================================================================
	330	END MODULE traadv_ubs

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

Download in other formats: