Context Navigation

trasbc.F90 @ 7731

Last change on this file since 7731 was 7731, checked in by dford, 7 years ago
Merge in revisions 6625:7726 of dev_r5518_v3.4_asm_nemovar_community, so this branch will be identical to revison 7726 of dev_r5518_v3.6_asm_nemovar_community.
File size: 16.9 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), POINTER, DIMENSION(:,:,:) :: ztrdt, ztrds
126	!!----------------------------------------------------------------------
127	!
128	IF( nn_timing == 1 ) CALL timing_start('tra_sbc')
129	!
130	IF( kt == nit000 ) THEN
131	IF(lwp) WRITE(numout,*)
132	IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition'
133	IF(lwp) WRITE(numout,*) '~~~~~~~ '
134	ENDIF
135
136	IF( l_trdtra ) THEN !* Save ta and sa trends
137	CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds )
138	ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
139	ztrds(:,:,:) = tsa(:,:,:,jp_sal)
140	ENDIF
141
142	!!gm IF( .NOT.ln_traqsr ) qsr(:,:) = 0.e0 ! no solar radiation penetration
143	IF( .NOT.ln_traqsr ) THEN ! no solar radiation penetration
144	qns(:,:) = qns(:,:) + qsr(:,:) ! total heat flux in qns
145	qsr(:,:) = 0.e0 ! qsr set to zero
146	ENDIF
147
148	!----------------------------------------
149	! EMP, SFX and QNS effects
150	!----------------------------------------
151	! Set before sbc tracer content fields
152	! ************************************
153	IF( kt == nit000 ) THEN ! Set the forcing field at nit000 - 1
154	! ! -----------------------------------
155	IF( ln_rstart .AND. & ! Restart: read in restart file
156	& iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN
157	IF(lwp) WRITE(numout,*) ' nit000-1 surface tracer content forcing fields red in the restart file'
158	zfact = 0.5_wp
159	CALL iom_get( numror, jpdom_autoglo, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem) ) ! before heat content sbc trend
160	CALL iom_get( numror, jpdom_autoglo, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal) ) ! before salt content sbc trend
161	ELSE ! No restart or restart not found: Euler forward time stepping
162	zfact = 1._wp
163	sbc_tsc_b(:,:,:) = 0._wp
164	ENDIF
165	ELSE ! Swap of forcing fields
166	! ! ----------------------
167	zfact = 0.5_wp
168	sbc_tsc_b(:,:,:) = sbc_tsc(:,:,:)
169	ENDIF
170	! Compute now sbc tracer content fields
171	! *************************************
172
173	! Concentration dilution effect on (t,s) due to
174	! evaporation, precipitation and qns, but not river runoff
175
176	IF( lk_vvl ) THEN ! Variable Volume case ==>> heat content of mass flux is in qns
177	DO jj = 1, jpj
178	DO ji = 1, jpi
179	sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj) ! non solar heat flux
180	sbc_tsc(ji,jj,jp_sal) = r1_rau0 * sfx(ji,jj) ! salt flux due to freezing/melting
181	END DO
182	END DO
183	ELSE ! Constant Volume case ==>> Concentration dilution effect
184	DO jj = 2, jpj
185	DO ji = fs_2, fs_jpim1 ! vector opt.
186	! temperature : heat flux
187	sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj) & ! non solar heat flux
188	& + r1_rau0 * emp(ji,jj) * tsn(ji,jj,1,jp_tem) ! concent./dilut. effect
189	! salinity : salt flux + concent./dilut. effect (both in sfx)
190	sbc_tsc(ji,jj,jp_sal) = r1_rau0 * ( sfx(ji,jj) & ! salt flux (freezing/melting)
191	& + emp(ji,jj) * tsn(ji,jj,1,jp_sal) ) ! concent./dilut. effect
192	END DO
193	END DO
194	IF( iom_use('emp_x_sst') ) CALL iom_put( "emp_x_sst", emp (:,:) * tsn(:,:,1,jp_tem) ) ! c/d term on sst
195	IF( iom_use('emp_x_sss') ) CALL iom_put( "emp_x_sss", emp (:,:) * tsn(:,:,1,jp_sal) ) ! c/d term on sss
196	ENDIF
197	! Concentration dilution effect on (t,s) due to evapouration, precipitation and qns, but not river runoff
198	DO jn = 1, jpts
199	DO jj = 2, jpj
200	DO ji = fs_2, fs_jpim1 ! vector opt.
201	z1_e3t = zfact / fse3t(ji,jj,1)
202	tsa(ji,jj,1,jn) = tsa(ji,jj,1,jn) + ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) * z1_e3t
203	END DO
204	END DO
205	END DO
206	! Write in the ocean restart file
207	! *******************************
208	IF( lrst_oce ) THEN
209	IF(lwp) WRITE(numout,*)
210	IF(lwp) WRITE(numout,*) 'sbc : ocean surface tracer content forcing fields written in ocean restart file ', &
211	& 'at it= ', kt,' date= ', ndastp
212	IF(lwp) WRITE(numout,*) '~~~~'
213	CALL iom_rstput( kt, nitrst, numrow, 'sbc_hc_b', sbc_tsc(:,:,jp_tem) )
214	CALL iom_rstput( kt, nitrst, numrow, 'sbc_sc_b', sbc_tsc(:,:,jp_sal) )
215	ENDIF
216	!
217	!
218	!----------------------------------------
219	! Ice Shelf effects (ISF)
220	! tbl treated as in Losh (2008) JGR
221	!----------------------------------------
222	!
223	IF( nn_isf > 0 ) THEN
224	zfact = 0.5e0
225	DO jj = 2, jpj
226	DO ji = fs_2, fs_jpim1
227
228	ikt = misfkt(ji,jj)
229	ikb = misfkb(ji,jj)
230
231	! level fully include in the ice shelf boundary layer
232	! if isfdiv, we have to remove heat flux due to inflow at 0oC (as in rnf when you add rnf at sst)
233	! sign - because fwf sign of evapo (rnf sign of precip)
234	DO jk = ikt, ikb - 1
235	! compute tfreez for the temperature correction (we add water at freezing temperature)
236	! compute trend
237	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &
238	& + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)) * r1_hisf_tbl(ji,jj)
239	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &
240	& + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj)
241	END DO
242
243	! level partially include in ice shelf boundary layer
244	! compute tfreez for the temperature correction (we add water at freezing temperature)
245	! compute trend
246	tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem) &
247	& + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)
248	tsa(ji,jj,ikb,jp_sal) = tsa(ji,jj,ikb,jp_sal) &
249	& + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)
250	END DO
251	END DO
252	IF( lrst_oce ) THEN
253	IF(lwp) WRITE(numout,*)
254	IF(lwp) WRITE(numout,*) 'sbc : isf surface tracer content forcing fields written in ocean restart file ', &
255	& 'at it= ', kt,' date= ', ndastp
256	IF(lwp) WRITE(numout,*) '~~~~'
257	CALL iom_rstput( kt, nitrst, numrow, 'fwf_isf_b', fwfisf(:,:) )
258	CALL iom_rstput( kt, nitrst, numrow, 'isf_hc_b' , risf_tsc(:,:,jp_tem) )
259	CALL iom_rstput( kt, nitrst, numrow, 'isf_sc_b' , risf_tsc(:,:,jp_sal) )
260	ENDIF
261	END IF
262	!
263	!----------------------------------------
264	! River Runoff effects
265	!----------------------------------------
266	!
267	IF( ln_rnf ) THEN ! input of heat and salt due to river runoff
268	zfact = 0.5_wp
269	DO jj = 2, jpj
270	DO ji = fs_2, fs_jpim1
271	IF( rnf(ji,jj) /= 0._wp ) THEN
272	zdep = zfact / h_rnf(ji,jj)
273	DO jk = 1, nk_rnf(ji,jj)
274	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &
275	& + ( rnf_tsc_b(ji,jj,jp_tem) + rnf_tsc(ji,jj,jp_tem) ) * zdep
276	IF( ln_rnf_sal ) tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &
277	& + ( rnf_tsc_b(ji,jj,jp_sal) + rnf_tsc(ji,jj,jp_sal) ) * zdep
278	END DO
279	ENDIF
280	END DO
281	END DO
282	ENDIF
283
284	#if defined key_asminc
285	! WARNING: THIS MAY WELL NOT BE REQUIRED - WE DON'T WANT TO CHANGE T&S BUT THIS MAY COMPENSATE ANOTHER TERM...
286	! Rate of change in e3t for each level is ssh_iau*e3t_0/ht_0
287	! Contribution to tsa should be rate of change in level / per m of ocean? (hence the division by fse3t_n)
288	IF( ln_sshinc ) THEN ! input of heat and salt due to assimilation
289	DO jj = 2, jpj
290	DO ji = fs_2, fs_jpim1
291	zdep = ssh_iau(ji,jj) / ( ht_0(ji,jj) + 1.0 - ssmask(ji, jj) )
292	DO jk = 1, jpkm1
293	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &
294	& + tsn(ji,jj,jk,jp_tem) * zdep * ( e3t_0(ji,jj,jk) / fse3t_n(ji,jj,jk) )
295	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &
296	& + tsn(ji,jj,jk,jp_sal) * zdep * ( e3t_0(ji,jj,jk) / fse3t_n(ji,jj,jk) )
297	END DO
298	END DO
299	END DO
300	ENDIF
301	#endif
302
303	IF( l_trdtra ) THEN ! send trends for further diagnostics
304	ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
305	ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
306	CALL trd_tra( kt, 'TRA', jp_tem, jptra_nsr, ztrdt )
307	CALL trd_tra( kt, 'TRA', jp_sal, jptra_nsr, ztrds )
308	CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds )
309	ENDIF
310	!
311	IF(ln_ctl) CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' sbc - Ta: ', mask1=tmask, &
312	& tab3d_2=tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )
313	!
314	IF( nn_timing == 1 ) CALL timing_stop('tra_sbc')
315	!
316	END SUBROUTINE tra_sbc
317
318	!!======================================================================
319	END MODULE trasbc

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

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

Context Navigation

source: branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90 @ 7731

Download in other formats: