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.

trasbc.F90 in branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90 @ 9321

Last change on this file since 9321 was 9321, checked in by davestorkey, 6 years ago
UKMO/dev_r5518_GO6_package branch: allow timing of I/O and coupling only. See GMED ticket 374.
File size: 17.5 KB

Line
1	MODULE trasbc
2	!!==============================================================================
3	!! * MODULE trasbc *
4	!! Ocean active tracers: surface boundary condition
5	!!==============================================================================
6	!! History : OPA ! 1998-10 (G. Madec, G. Roullet, M. Imbard) Original code
7	!! 8.2 ! 2001-02 (D. Ludicone) sea ice and free surface
8	!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
9	!! 3.3 ! 2010-04 (M. Leclair, G. Madec) Forcing averaged over 2 time steps
10	!! - ! 2010-09 (C. Ethe, G. Madec) Merge TRA-TRC
11	!! 3.6 ! 2014-11 (P. Mathiot) isf melting forcing
12	!!----------------------------------------------------------------------
13
14	!!----------------------------------------------------------------------
15	!! tra_sbc : update the tracer trend at ocean surface
16	!!----------------------------------------------------------------------
17	USE oce ! ocean dynamics and active tracers
18	USE sbc_oce ! surface boundary condition: ocean
19	USE dom_oce ! ocean space domain variables
20	USE phycst ! physical constant
21	USE sbcmod ! ln_rnf
22	USE sbcrnf ! River runoff
23	USE sbcisf ! Ice shelf
24	USE traqsr ! solar radiation penetration
25	USE trd_oce ! trends: ocean variables
26	USE trdtra ! trends manager: tracers
27	!
28	USE in_out_manager ! I/O manager
29	USE prtctl ! Print control
30	USE iom
31	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
32	USE wrk_nemo ! Memory Allocation
33	USE timing ! Timing
34	USE eosbn2
35	#if defined key_asminc
36	USE asminc ! Assimilation increment
37	#endif
38
39	IMPLICIT NONE
40	PRIVATE
41
42	PUBLIC tra_sbc ! routine called by step.F90
43
44	!! * Substitutions
45	# include "domzgr_substitute.h90"
46	# include "vectopt_loop_substitute.h90"
47	!!----------------------------------------------------------------------
48	!! NEMO/OPA 3.7 , NEMO Consortium (2014)
49	!! $Id$
50	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
51	!!----------------------------------------------------------------------
52	CONTAINS
53
54	SUBROUTINE tra_sbc ( kt )
55	!!----------------------------------------------------------------------
56	!! * ROUTINE tra_sbc *
57	!!
58	!! ** Purpose : Compute the tracer surface boundary condition trend of
59	!! (flux through the interface, concentration/dilution effect)
60	!! and add it to the general trend of tracer equations.
61	!!
62	!! ** Method :
63	!! Following Roullet and Madec (2000), the air-sea flux can be divided
64	!! into three effects: (1) Fext, external forcing;
65	!! (2) Fwi, concentration/dilution effect due to water exchanged
66	!! at the surface by evaporation, precipitations and runoff (E-P-R);
67	!! (3) Fwe, tracer carried with the water that is exchanged.
68	!! - salinity : salt flux only due to freezing/melting
69	!! sa = sa + sfx / rau0 / e3t for k=1
70	!!
71	!! Fext, flux through the air-sea interface for temperature and salt:
72	!! - temperature : heat flux q (w/m2). If penetrative solar
73	!! radiation q is only the non solar part of the heat flux, the
74	!! solar part is added in traqsr.F routine.
75	!! ta = ta + q /(rau0 rcp e3t) for k=1
76	!! - salinity : no salt flux
77	!!
78	!! The formulation for Fwb and Fwi vary according to the free
79	!! surface formulation (linear or variable volume).
80	!! * Linear free surface
81	!! The surface freshwater flux modifies the ocean volume
82	!! and thus the concentration of a tracer and the temperature.
83	!! First order of the effect of surface freshwater exchange
84	!! for salinity, it can be neglected on temperature (especially
85	!! as the temperature of precipitations and runoffs is usually
86	!! unknown).
87	!! - temperature : we assume that the temperature of both
88	!! precipitations and runoffs is equal to the SST, thus there
89	!! is no additional flux since in this case, the concentration
90	!! dilution effect is balanced by the net heat flux associated
91	!! to the freshwater exchange (Fwe+Fwi=0):
92	!! (Tp P - Te E) + SST (P-E) = 0 when Tp=Te=SST
93	!! - salinity : evaporation, precipitation and runoff
94	!! water has a zero salinity but there is a salt flux due to
95	!! freezing/melting, thus:
96	!! sa = sa + emp * sn / rau0 / e3t for k=1
97	!! + sfx / rau0 / e3t
98	!! where emp, the surface freshwater budget (evaporation minus
99	!! precipitation minus runoff) given in kg/m2/s is divided
100	!! by rau0 (density of sea water) to obtain m/s.
101	!! Note: even though Fwe does not appear explicitly for
102	!! temperature in this routine, the heat carried by the water
103	!! exchanged through the surface is part of the total heat flux
104	!! forcing and must be taken into account in the global heat
105	!! balance).
106	!! * nonlinear free surface (variable volume, lk_vvl)
107	!! contrary to the linear free surface case, Fwi is properly
108	!! taken into account by using the true layer thicknesses to
109	!! calculate tracer content and advection. There is no need to
110	!! deal with it in this routine.
111	!! - temperature: Fwe=SST (P-E+R) is added to Fext.
112	!! - salinity: Fwe = 0, there is no surface flux of salt.
113	!!
114	!! ** Action : - Update the 1st level of (ta,sa) with the trend associated
115	!! with the tracer surface boundary condition
116	!! - send trends to trdtra module (l_trdtra=T)
117	!!----------------------------------------------------------------------
118	INTEGER, INTENT(in) :: kt ! ocean time-step index
119	!!
120	INTEGER :: ji, jj, jk, jn ! dummy loop indices
121	INTEGER :: ikt, ikb
122	INTEGER :: nk_isf
123	REAL(wp) :: zfact, z1_e3t, zdep
124	REAL(wp) :: zalpha, zhk
125	REAL(wp) :: zt_frz, zpress
126	REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdt, ztrds
127	!!----------------------------------------------------------------------
128	!
129	IF( nn_timing == 1 ) CALL timing_start('tra_sbc')
130	!
131	IF( kt == nit000 ) THEN
132	IF(lwp) WRITE(numout,*)
133	IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition'
134	IF(lwp) WRITE(numout,*) '~~~~~~~ '
135	ENDIF
136
137	IF( l_trdtra ) THEN !* Save ta and sa trends
138	CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds )
139	ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
140	ztrds(:,:,:) = tsa(:,:,:,jp_sal)
141	ENDIF
142
143	!!gm IF( .NOT.ln_traqsr ) qsr(:,:) = 0.e0 ! no solar radiation penetration
144	IF( .NOT.ln_traqsr ) THEN ! no solar radiation penetration
145	qns(:,:) = qns(:,:) + qsr(:,:) ! total heat flux in qns
146	qsr(:,:) = 0.e0 ! qsr set to zero
147	ENDIF
148
149	!----------------------------------------
150	! EMP, SFX and QNS effects
151	!----------------------------------------
152	! Set before sbc tracer content fields
153	! ************************************
154	IF( kt == nit000 ) THEN ! Set the forcing field at nit000 - 1
155	! ! -----------------------------------
156	IF( ln_rstart .AND. & ! Restart: read in restart file
157	& iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN
158	IF(lwp) WRITE(numout,*) ' nit000-1 surface tracer content forcing fields red in the restart file'
159	zfact = 0.5_wp
160	CALL iom_get( numror, jpdom_autoglo, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem) ) ! before heat content sbc trend
161	CALL iom_get( numror, jpdom_autoglo, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal) ) ! before salt content sbc trend
162	ELSE ! No restart or restart not found: Euler forward time stepping
163	zfact = 1._wp
164	sbc_tsc(:,:,:) = 0._wp
165	sbc_tsc_b(:,:,:) = 0._wp
166	ENDIF
167	ELSE ! Swap of forcing fields
168	! ! ----------------------
169	zfact = 0.5_wp
170	sbc_tsc_b(:,:,:) = sbc_tsc(:,:,:)
171	ENDIF
172	! Compute now sbc tracer content fields
173	! *************************************
174
175	! Concentration dilution effect on (t,s) due to
176	! evaporation, precipitation and qns, but not river runoff
177
178	IF( lk_vvl ) THEN ! Variable Volume case ==>> heat content of mass flux is in qns
179	DO jj = 1, jpj
180	DO ji = 1, jpi
181	sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj) ! non solar heat flux
182	sbc_tsc(ji,jj,jp_sal) = r1_rau0 * sfx(ji,jj) ! salt flux due to freezing/melting
183	END DO
184	END DO
185	ELSE ! Constant Volume case ==>> Concentration dilution effect
186	DO jj = 2, jpj
187	DO ji = fs_2, fs_jpim1 ! vector opt.
188	! temperature : heat flux
189	sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj) & ! non solar heat flux
190	& + r1_rau0 * emp(ji,jj) * tsn(ji,jj,1,jp_tem) ! concent./dilut. effect
191	! salinity : salt flux + concent./dilut. effect (both in sfx)
192	sbc_tsc(ji,jj,jp_sal) = r1_rau0 * ( sfx(ji,jj) & ! salt flux (freezing/melting)
193	& + emp(ji,jj) * tsn(ji,jj,1,jp_sal) ) ! concent./dilut. effect
194	END DO
195	END DO
196	IF( iom_use('emp_x_sst') ) CALL iom_put( "emp_x_sst", emp (:,:) * tsn(:,:,1,jp_tem) ) ! c/d term on sst
197	IF( iom_use('emp_x_sss') ) CALL iom_put( "emp_x_sss", emp (:,:) * tsn(:,:,1,jp_sal) ) ! c/d term on sss
198	ENDIF
199	! Concentration dilution effect on (t,s) due to evapouration, precipitation and qns, but not river runoff
200	DO jn = 1, jpts
201	DO jj = 2, jpj
202	DO ji = fs_2, fs_jpim1 ! vector opt.
203	z1_e3t = zfact / fse3t(ji,jj,1)
204	tsa(ji,jj,1,jn) = tsa(ji,jj,1,jn) + ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) * z1_e3t
205	END DO
206	END DO
207	END DO
208	! Write in the ocean restart file
209	! *******************************
210	IF( lrst_oce ) THEN
211	IF(lwp) WRITE(numout,*)
212	IF(lwp) WRITE(numout,*) 'sbc : ocean surface tracer content forcing fields written in ocean restart file ', &
213	& 'at it= ', kt,' date= ', ndastp
214	IF(lwp) WRITE(numout,*) '~~~~'
215	IF(nn_timing == 2) CALL timing_start('iom_rstput')
216	CALL iom_rstput( kt, nitrst, numrow, 'sbc_hc_b', sbc_tsc(:,:,jp_tem) )
217	CALL iom_rstput( kt, nitrst, numrow, 'sbc_sc_b', sbc_tsc(:,:,jp_sal) )
218	IF(nn_timing == 2) CALL timing_stop('iom_rstput')
219	ENDIF
220	!
221	!
222	!----------------------------------------
223	! Ice Shelf effects (ISF)
224	! tbl treated as in Losh (2008) JGR
225	!----------------------------------------
226	!
227	IF( nn_isf > 0 ) THEN
228	zfact = 0.5e0
229	DO jj = 2, jpj
230	DO ji = fs_2, fs_jpim1
231
232	ikt = misfkt(ji,jj)
233	ikb = misfkb(ji,jj)
234
235	! level fully include in the ice shelf boundary layer
236	! if isfdiv, we have to remove heat flux due to inflow at 0oC (as in rnf when you add rnf at sst)
237	! sign - because fwf sign of evapo (rnf sign of precip)
238	DO jk = ikt, ikb - 1
239	! compute tfreez for the temperature correction (we add water at freezing temperature)
240	! compute trend
241	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &
242	& + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)) * r1_hisf_tbl(ji,jj)
243	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &
244	& + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj)
245	END DO
246
247	! level partially include in ice shelf boundary layer
248	! compute tfreez for the temperature correction (we add water at freezing temperature)
249	! compute trend
250	tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem) &
251	& + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)
252	tsa(ji,jj,ikb,jp_sal) = tsa(ji,jj,ikb,jp_sal) &
253	& + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)
254	END DO
255	END DO
256	IF( lrst_oce ) THEN
257	IF(lwp) WRITE(numout,*)
258	IF(lwp) WRITE(numout,*) 'sbc : isf surface tracer content forcing fields written in ocean restart file ', &
259	& 'at it= ', kt,' date= ', ndastp
260	IF(lwp) WRITE(numout,*) '~~~~'
261	IF(nn_timing == 2) CALL timing_start('iom_rstput')
262	CALL iom_rstput( kt, nitrst, numrow, 'fwf_isf_b', fwfisf(:,:) )
263	CALL iom_rstput( kt, nitrst, numrow, 'isf_hc_b' , risf_tsc(:,:,jp_tem) )
264	CALL iom_rstput( kt, nitrst, numrow, 'isf_sc_b' , risf_tsc(:,:,jp_sal) )
265	IF(nn_timing == 2) CALL timing_stop('iom_rstput')
266	ENDIF
267	END IF
268	!
269	!----------------------------------------
270	! River Runoff effects
271	!----------------------------------------
272	!
273	IF( ln_rnf ) THEN ! input of heat and salt due to river runoff
274	zfact = 0.5_wp
275	DO jj = 2, jpj
276	DO ji = fs_2, fs_jpim1
277	IF( rnf(ji,jj) /= 0._wp ) THEN
278	zdep = zfact / h_rnf(ji,jj)
279	DO jk = 1, nk_rnf(ji,jj)
280	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &
281	& + ( rnf_tsc_b(ji,jj,jp_tem) + rnf_tsc(ji,jj,jp_tem) ) * zdep
282	IF( ln_rnf_sal ) tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &
283	& + ( rnf_tsc_b(ji,jj,jp_sal) + rnf_tsc(ji,jj,jp_sal) ) * zdep
284	END DO
285	ENDIF
286	END DO
287	END DO
288	ENDIF
289
290	IF( iom_use('rnf_x_sst') ) CALL iom_put( "rnf_x_sst", rnf*tsn(:,:,1,jp_tem) ) ! runoff term on sst
291	IF( iom_use('rnf_x_sss') ) CALL iom_put( "rnf_x_sss", rnf*tsn(:,:,1,jp_sal) ) ! runoff term on sss
292
293	#if defined key_asminc
294	! WARNING: THIS MAY WELL NOT BE REQUIRED - WE DON'T WANT TO CHANGE T&S BUT THIS MAY COMPENSATE ANOTHER TERM...
295	! Rate of change in e3t for each level is ssh_iau*e3t_0/ht_0
296	! Contribution to tsa should be rate of change in level / per m of ocean? (hence the division by fse3t_n)
297	IF( ln_sshinc ) THEN ! input of heat and salt due to assimilation
298	DO jj = 2, jpj
299	DO ji = fs_2, fs_jpim1
300	zdep = ssh_iau(ji,jj) / ( ht_0(ji,jj) + 1.0 - ssmask(ji, jj) )
301	DO jk = 1, jpkm1
302	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &
303	& + tsn(ji,jj,jk,jp_tem) * zdep * ( e3t_0(ji,jj,jk) / fse3t_n(ji,jj,jk) )
304	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &
305	& + tsn(ji,jj,jk,jp_sal) * zdep * ( e3t_0(ji,jj,jk) / fse3t_n(ji,jj,jk) )
306	END DO
307	END DO
308	END DO
309	ENDIF
310	#endif
311
312	IF( l_trdtra ) THEN ! send trends for further diagnostics
313	ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
314	ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
315	CALL trd_tra( kt, 'TRA', jp_tem, jptra_nsr, ztrdt )
316	CALL trd_tra( kt, 'TRA', jp_sal, jptra_nsr, ztrds )
317	CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds )
318	ENDIF
319	!
320	IF(ln_ctl) CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' sbc - Ta: ', mask1=tmask, &
321	& tab3d_2=tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )
322	!
323	IF( nn_timing == 1 ) CALL timing_stop('tra_sbc')
324	!
325	END SUBROUTINE tra_sbc
326
327	!!======================================================================
328	END MODULE trasbc

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

Download in other formats: