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.

limsbc.F90 in branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

Context Navigation

source: branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90 @ 3517

Last change on this file since 3517 was 3517, checked in by gm, 11 years ago
gm: Branch: dev_r3385_NOCS04_HAMF; #665. update sbccpl ; change LIM3 from equivalent salt flux to salt flux and mass flux
Property svn:keywords set to `Id`
File size: 25.4 KB

Line
1	MODULE limsbc
2	!!======================================================================
3	!! * MODULE limsbc *
4	!! computation of the flux at the sea ice/ocean interface
5	!!======================================================================
6	!! History : - ! 2006-07 (M. Vancoppelle) LIM3 original code
7	!! 3.0 ! 2008-03 (C. Tallandier) surface module
8	!! - ! 2008-04 (C. Tallandier) split in 2 + new ice-ocean coupling
9	!! 3.3 ! 2010-05 (G. Madec) decrease ocean & ice reference salinities in the Baltic sea
10	!! ! + simplification of the ice-ocean stress calculation
11	!! 3.4 ! 2011-02 (G. Madec) dynamical allocation
12	!! 3.5 ! 2012-10 (A. Coward, G. Madec) salt fluxes ; ice+snow mass
13	!!----------------------------------------------------------------------
14	#if defined key_lim3
15	!!----------------------------------------------------------------------
16	!! 'key_lim3' LIM 3.0 sea-ice model
17	!!----------------------------------------------------------------------
18	!! lim_sbc_alloc : allocate the limsbc arrays
19	!! lim_sbc_init : initialisation
20	!! lim_sbc_flx : updates mass, heat and salt fluxes at the ocean surface
21	!! lim_sbc_tau : update i- and j-stresses, and its modulus at the ocean surface
22	!!----------------------------------------------------------------------
23	USE par_oce ! ocean parameters
24	USE par_ice ! ice parameters
25	USE dom_oce ! ocean domain
26	USE sbc_ice ! Surface boundary condition: sea-ice fields
27	USE sbc_oce ! Surface boundary condition: ocean fields
28	USE phycst ! physical constants
29	USE albedo ! albedo parameters
30	USE ice ! LIM sea-ice variables
31	USE lbclnk ! ocean lateral boundary condition
32	USE in_out_manager ! I/O manager
33	USE lib_mpp ! MPP library
34	USE wrk_nemo ! work arrays
35	USE prtctl ! Print control
36	USE cpl_oasis3, ONLY : lk_cpl
37	USE oce, ONLY : sshn, sshb, snwice_mass, snwice_mass_b, snwice_fmass, sshu_b, sshv_b, sshu_n, sshv_n, sshf_n
38	USE dom_ice, ONLY : tms
39
40	IMPLICIT NONE
41	PRIVATE
42
43	PUBLIC lim_sbc_init ! called by ice_init
44	PUBLIC lim_sbc_flx ! called by sbc_ice_lim
45	PUBLIC lim_sbc_tau ! called by sbc_ice_lim
46
47	REAL(wp) :: epsi16 = 1.e-16_wp ! constant values
48	REAL(wp) :: rzero = 0._wp
49	REAL(wp) :: rone = 1._wp
50
51	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau_oce, vtau_oce ! air-ocean surface i- & j-stress [N/m2]
52	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: tmod_io ! modulus of the ice-ocean velocity [m/s]
53	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: soce_0 , sice_0 ! cst SSS and ice salinity (levitating sea-ice)
54
55	!! * Substitutions
56	# include "vectopt_loop_substitute.h90"
57	!!----------------------------------------------------------------------
58	!! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
59	!! $Id$
60	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
61	!!----------------------------------------------------------------------
62	CONTAINS
63
64	INTEGER FUNCTION lim_sbc_alloc()
65	!!-------------------------------------------------------------------
66	!! * ROUTINE lim_sbc_alloc *
67	!!-------------------------------------------------------------------
68	ALLOCATE( soce_0(jpi,jpj) , utau_oce(jpi,jpj) , &
69	& sice_0(jpi,jpj) , vtau_oce(jpi,jpj) , tmod_io(jpi,jpj), STAT=lim_sbc_alloc)
70	!
71	IF( lk_mpp ) CALL mpp_sum( lim_sbc_alloc )
72	IF( lim_sbc_alloc /= 0 ) CALL ctl_warn('lim_sbc_alloc: failed to allocate arrays')
73	END FUNCTION lim_sbc_alloc
74
75
76	SUBROUTINE lim_sbc_flx( kt )
77	!!-------------------------------------------------------------------
78	!! * ROUTINE lim_sbc_flx *
79	!!
80	!! ** Purpose : Update the surface ocean boundary condition for heat
81	!! salt and mass over areas where sea-ice is non-zero
82	!!
83	!! ** Action : - computes the heat and freshwater/salt fluxes
84	!! at the ice-ocean interface.
85	!! - Update the ocean sbc
86	!!
87	!! ** Outputs : - qsr : sea heat flux: solar
88	!! - qns : sea heat flux: non solar
89	!! - emp : freshwater budget: volume flux
90	!! - sfx : salt flux
91	!! - fr_i : ice fraction
92	!! - tn_ice : sea-ice surface temperature
93	!! - alb_ice : sea-ice alberdo (lk_cpl=T)
94	!!
95	!! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90.
96	!! Tartinville et al. 2001 Ocean Modelling, 3, 95-108.
97	!!---------------------------------------------------------------------
98	INTEGER, INTENT(in) :: kt ! number of iteration
99	!
100	INTEGER :: ji, jj ! dummy loop indices
101	INTEGER :: ierr ! local integer
102	INTEGER :: ifvt, i1mfr, idfr ! some switches
103	INTEGER :: iflt, ial, iadv, ifral, ifrdv
104	REAL(wp) :: zinda, zemp, zemp_snow, zfmm ! local scalars
105	REAL(wp), POINTER, DIMENSION(:,:) :: zfcm1 , zfcm2 ! solar/non solar heat fluxes
106	REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace
107	!!---------------------------------------------------------------------
108
109	CALL wrk_alloc( jpi, jpj, zfcm1 , zfcm2 )
110	IF( lk_cpl ) CALL wrk_alloc( jpi, jpj, jpl, zalb, zalbp )
111
112	!------------------------------------------!
113	! heat flux at the ocean surface !
114	!------------------------------------------!
115	! pfrld is the lead fraction at the previous time step (actually between TRP and THD)
116	! changed to old_frld and old ht_i
117
118	DO jj = 1, jpj
119	DO ji = 1, jpi
120	zinda = 1.0 - MAX( rzero , SIGN( rone , - ( 1.0 - pfrld(ji,jj) ) ) )
121	ifvt = zinda * MAX( rzero , SIGN( rone, -phicif (ji,jj) ) ) !subscripts are bad here
122	i1mfr = 1.0 - MAX( rzero , SIGN( rone , - ( at_i(ji,jj) ) ) )
123	idfr = 1.0 - MAX( rzero , SIGN( rone , ( 1.0 - at_i(ji,jj) ) - pfrld(ji,jj) ) )
124	iflt = zinda * (1 - i1mfr) * (1 - ifvt )
125	ial = ifvt * i1mfr + ( 1 - ifvt ) * idfr
126	iadv = ( 1 - i1mfr ) * zinda
127	ifral = ( 1 - i1mfr * ( 1 - ial ) )
128	ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv
129
130	! switch --- 1.0 ---------------- 0.0 --------------------
131	! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
132	! zinda \| if pfrld = 1 \| if pfrld < 1 \|
133	! -> ifvt\| if pfrld old_ht_i
134	! i1mfr \| if frld = 1 \| if frld < 1 \|
135	! idfr \| if frld <= pfrld \| if frld > pfrld \|
136	! iflt \|
137	! ial \|
138	! iadv \|
139	! ifral
140	! ifrdv
141
142	! computation the solar flux at ocean surface
143	zfcm1(ji,jj) = pfrld(ji,jj) * qsr(ji,jj) + ( 1. - pfrld(ji,jj) ) * fstric(ji,jj)
144	! fstric Solar flux transmitted trough the ice
145	! qsr Net short wave heat flux on free ocean
146	! new line
147	fscmbq(ji,jj) = ( 1.0 - pfrld(ji,jj) ) * fstric(ji,jj)
148
149	! computation the non solar heat flux at ocean surface
150	zfcm2(ji,jj) = - zfcm1(ji,jj) &
151	& + iflt * ( fscmbq(ji,jj) ) & ! total abl -> fscmbq is given to the ocean
152	! fscmbq and ffltbif are obsolete
153	! & + iflt * ffltbif(ji,jj) !!! only if one category is used
154	& + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice &
155	& + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice &
156	& + fhmec(ji,jj) & ! new contribution due to snow melt when ridging!!
157	& + fheat_rpo(ji,jj) & ! contribution from ridge formation
158	& + fheat_res(ji,jj)
159	! fscmbq Part of the solar radiation transmitted through the ice and going to the ocean computed in limthd_zdf.F90
160	! ffltbif Total heat content of the ice (brine pockets+ice) / delta_t
161	! qcmif Energy needed to bring the ocean surface layer until its freezing (ok)
162	! qldif heat balance of the lead (or of the open ocean)
163	! qfvbq i think this is wrong!
164	! ---> Array used to store energy in case of total lateral ablation
165	! qfvbq latent heat uptake/release after accretion/ablation
166	! qdtcn Energy from the turbulent oceanic heat flux heat flux coming in the lead
167
168	IF( num_sal == 2 ) zfcm2(ji,jj) = zfcm2(ji,jj) + fhbri(ji,jj) ! add contribution due to brine drainage
169
170	! bottom radiative component is sent to the computation of the oceanic heat flux
171	fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj)
172
173	! used to compute the oceanic heat flux at the next time step
174	qsr(ji,jj) = zfcm1(ji,jj) ! solar heat flux
175	qns(ji,jj) = zfcm2(ji,jj) - fdtcn(ji,jj) ! non solar heat flux
176	! ! fdtcn : turbulent oceanic heat flux
177
178	!!gm this IF prevents the vertorisation of the whole loop
179	IF ( ( ji == jiindx ) .AND. ( jj == jjindx) ) THEN
180	WRITE(numout,*) ' lim_sbc : heat fluxes '
181	WRITE(numout,*) ' qsr : ', qsr(jiindx,jjindx)
182	WRITE(numout,*) ' zfcm1 : ', zfcm1(jiindx,jjindx)
183	WRITE(numout,*) ' pfrld : ', pfrld(jiindx,jjindx)
184	WRITE(numout,*) ' fstric : ', fstric (jiindx,jjindx)
185	WRITE(numout,*)
186	WRITE(numout,*) ' qns : ', qns(jiindx,jjindx)
187	WRITE(numout,*) ' zfcm2 : ', zfcm2(jiindx,jjindx)
188	WRITE(numout,*) ' zfcm1 : ', zfcm1(jiindx,jjindx)
189	WRITE(numout,*) ' ifral : ', ifral
190	WRITE(numout,*) ' ial : ', ial
191	WRITE(numout,*) ' qcmif : ', qcmif(jiindx,jjindx)
192	WRITE(numout,*) ' qldif : ', qldif(jiindx,jjindx)
193	WRITE(numout,) ' qcmif / dt: ', qcmif(jiindx,jjindx) r1_rdtice
194	WRITE(numout,) ' qldif / dt: ', qldif(jiindx,jjindx) r1_rdtice
195	WRITE(numout,*) ' ifrdv : ', ifrdv
196	WRITE(numout,*) ' qfvbq : ', qfvbq(jiindx,jjindx)
197	WRITE(numout,*) ' qdtcn : ', qdtcn(jiindx,jjindx)
198	WRITE(numout,) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) r1_rdtice
199	WRITE(numout,) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) r1_rdtice
200	WRITE(numout,*) ' '
201	WRITE(numout,*) ' fdtcn : ', fdtcn(jiindx,jjindx)
202	WRITE(numout,*) ' fhmec : ', fhmec(jiindx,jjindx)
203	WRITE(numout,*) ' fheat_rpo : ', fheat_rpo(jiindx,jjindx)
204	WRITE(numout,*) ' fhbri : ', fhbri(jiindx,jjindx)
205	WRITE(numout,*) ' fheat_res : ', fheat_res(jiindx,jjindx)
206	ENDIF
207	!!gm end
208	END DO
209	END DO
210
211	!------------------------------------------!
212	! mass flux at the ocean surface !
213	!------------------------------------------!
214
215	!!gm optimisation: this loop have to be merged with the previous one
216	DO jj = 1, jpj
217	DO ji = 1, jpi
218	! case of realistic freshwater flux (Tartinville et al., 2001) (presently ACTIVATED)
219	! -------------------------------------------------------------------------------------
220	! The idea of this approach is that the system that we consider is the ICE-OCEAN system
221	! Thus FW flux = External ( E-P+snow melt)
222	! Salt flux = Exchanges in the ice-ocean system then converted into FW
223	! Associated to Ice formation AND Ice melting
224	! Even if i see Ice melting as a FW and SALT flux
225	!
226
227	! computing freshwater exchanges at the ice/ocean interface
228	zemp = emp(ji,jj) * ( 1.0 - at_i(ji,jj) ) & ! evaporation over oceanic fraction
229	& - tprecip(ji,jj) * at_i(ji,jj) & ! all precipitation reach the ocean
230	& + sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) ) & ! except solid precip intercepted by sea-ice
231	& - fmmec(ji,jj) ! snow falling when ridging
232
233	! mass flux at the ocean/ice interface (sea ice fraction)
234	zemp_snw = rdm_snw(ji,jj) * r1_rdtice ! snow melting = pure water that enters the ocean
235	zfmm = rdm_ice(ji,jj) * r1_rdtice ! Freezing minus mesting
236
237	emp(ji,jj) = zemp + zemp_snw + zfmm ! mass flux + F/M mass flux (always ice/ocean mass exchange)
238
239	! correcting brine salt fluxes (zinda = 1 if pfrld=1 , =0 otherwise)
240	zinda = 1.0 - MAX( rzero , SIGN( rone , - ( 1.0 - pfrld(ji,jj) ) ) )
241	fsbri(ji,jj) = zinda * fsbri(ji,jj)
242	END DO
243	END DO
244
245	!------------------------------------------!
246	! salt flux at the ocean surface !
247	!------------------------------------------!
248
249	IF( num_sal == 2 ) THEN ! variable ice salinity: brine drainage included in the salt flux
250	sfx (:,:) = fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:) + fsbri(:,:)
251	ELSE ! constant ice salinity:
252	sfx (:,:) = fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:)
253	ENDIF
254	!-----------------------------------------------!
255	! mass of snow and ice per unit area !
256	!-----------------------------------------------!
257	IF( nn_ice_embd /= 0 ) THEN ! embedded sea-ice (mass required)
258	snwice_mass_b(:,:) = snwice_mass(:,:) ! save mass from the previous ice time step
259	! ! new mass per unit area
260	snwice_mass (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:) )
261	! ! time evolution of snow+ice mass
262	snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) * r1_rdtice
263	ENDIF
264
265	!-----------------------------------------------!
266	! Storing the transmitted variables !
267	!-----------------------------------------------!
268	fr_i (:,:) = at_i(:,:) ! Sea-ice fraction
269	tn_ice(:,:,:) = t_su(:,:,:) ! Ice surface temperature
270
271	!------------------------------------------------!
272	! Computation of snow/ice and ocean albedo !
273	!------------------------------------------------!
274	IF( lk_cpl ) THEN ! coupled case
275	CALL albedo_ice( t_su, ht_i, ht_s, zalbp, zalb ) ! snow/ice albedo
276	!
277	alb_ice(:,:,:) = 0.5_wp * zalbp(:,:,:) + 0.5_wp * zalb (:,:,:) ! Ice albedo (mean clear and overcast skys)
278	ENDIF
279
280	IF(ln_ctl) THEN
281	CALL prt_ctl( tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ' )
282	CALL prt_ctl( tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=sfx , clinfo2=' sfx : ' )
283	CALL prt_ctl( tab2d_1=fr_i , clinfo1=' lim_sbc: fr_i : ' )
284	CALL prt_ctl( tab3d_1=tn_ice, clinfo1=' lim_sbc: tn_ice : ', kdim=jpl )
285	ENDIF
286	!
287	CALL wrk_dealloc( jpi, jpj, zfcm1 , zfcm2 )
288	IF( lk_cpl ) CALL wrk_dealloc( jpi, jpj, jpl, zalb, zalbp )
289	!
290	END SUBROUTINE lim_sbc_flx
291
292
293	SUBROUTINE lim_sbc_tau( kt , pu_oce, pv_oce )
294	!!-------------------------------------------------------------------
295	!! * ROUTINE lim_sbc_tau *
296	!!
297	!! ** Purpose : Update the ocean surface stresses due to the ice
298	!!
299	!! ** Action : * at each ice time step (every nn_fsbc time step):
300	!! - compute the modulus of ice-ocean relative velocity
301	!! (rhoCd) at T-point (C-grid) or I-point (B-grid)
302	!! tmod_io = rhoco * \| U_ice-U_oce \|
303	!! - update the modulus of stress at ocean surface
304	!! taum = frld * taum + (1-frld) * tmod_io * \| U_ice-U_oce \|
305	!! * at each ocean time step (every kt):
306	!! compute linearized ice-ocean stresses as
307	!! Utau = tmod_io * \| U_ice - pU_oce \|
308	!! using instantaneous current ocean velocity (usually before)
309	!!
310	!! NB: - ice-ocean rotation angle no more allowed
311	!! - here we make an approximation: taum is only computed every ice time step
312	!! This avoids mutiple average to pass from T -> U,V grids and next from U,V grids
313	!! to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton...
314	!!
315	!! ** Outputs : - utau, vtau : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes
316	!! - taum : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes
317	!!---------------------------------------------------------------------
318	INTEGER , INTENT(in) :: kt ! ocean time-step index
319	REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents
320	!!
321	INTEGER :: ji, jj ! dummy loop indices
322	REAL(wp) :: zat_u, zutau_ice, zu_t, zmodt ! local scalar
323	REAL(wp) :: zat_v, zvtau_ice, zv_t ! - -
324	!!---------------------------------------------------------------------
325	!
326	IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step)
327	!CDIR NOVERRCHK
328	DO jj = 2, jpjm1 !* update the modulus of stress at ocean surface (T-point)
329	!CDIR NOVERRCHK
330	DO ji = fs_2, fs_jpim1
331	! ! 2*(U_ice-U_oce) at T-point
332	zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj)
333	zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1)
334	! ! \|U_ice-U_oce\|^2
335	zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t )
336	! ! update the ocean stress modulus
337	taum(ji,jj) = ( 1._wp - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zmodt
338	tmod_io(ji,jj) = rhoco * SQRT( zmodt ) ! rhoco * \|U_ice-U_oce\| at T-point
339	END DO
340	END DO
341	CALL lbc_lnk( taum, 'T', 1. ) ; CALL lbc_lnk( tmod_io, 'T', 1. )
342	!
343	utau_oce(:,:) = utau(:,:) !* save the air-ocean stresses at ice time-step
344	vtau_oce(:,:) = vtau(:,:)
345	!
346	ENDIF
347	!
348	! !== every ocean time-step ==!
349	!
350	DO jj = 2, jpjm1 !* update the stress WITHOUT a ice-ocean rotation angle
351	DO ji = fs_2, fs_jpim1 ! Vect. Opt.
352	zat_u = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5_wp ! ice area at u and V-points
353	zat_v = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5_wp
354	! ! linearized quadratic drag formulation
355	zutau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) )
356	zvtau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) )
357	! ! stresses at the ocean surface
358	utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_ice
359	vtau(ji,jj) = ( 1._wp - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice
360	END DO
361	END DO
362	CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition
363	!
364	IF(ln_ctl) CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau : ', mask1=umask, &
365	& tab2d_2=vtau, clinfo2=' vtau : ' , mask2=vmask )
366	!
367	END SUBROUTINE lim_sbc_tau
368
369
370	SUBROUTINE lim_sbc_init
371	!!-------------------------------------------------------------------
372	!! * ROUTINE lim_sbc_init *
373	!!
374	!! ** Purpose : Preparation of the file ice_evolu for the output of
375	!! the temporal evolution of key variables
376	!!
377	!! ** input : Namelist namicedia
378	!!-------------------------------------------------------------------
379	!
380	INTEGER :: ji, jj ! dummy loop indices
381	REAL(wp) :: zcoefu, zcoefv, zcoeff ! local scalar
382	IF(lwp) WRITE(numout,*)
383	IF(lwp) WRITE(numout,*) 'lim_sbc_init : LIM-3 sea-ice - surface boundary condition'
384	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ '
385
386	! ! allocate lim_sbc array
387	IF( lim_sbc_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'lim_sbc_init : unable to allocate standard arrays' )
388	!
389	soce_0(:,:) = soce ! constant SSS and ice salinity used in levitating sea-ice case
390	sice_0(:,:) = sice
391	!
392	IF( cp_cfg == "orca" ) THEN ! decrease ocean & ice reference salinities in the Baltic sea
393	WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND. &
394	& 54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp )
395	soce_0(:,:) = 4._wp
396	sice_0(:,:) = 2._wp
397	END WHERE
398	ENDIF
399	! ! embedded sea ice
400	IF( nn_ice_embd /= 0 ) THEN ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass
401	snwice_mass (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:) )
402	snwice_mass_b(:,:) = snwice_mass(:,:)
403	ELSE
404	snwice_mass (:,:) = 0.0_wp ! no mass exchanges
405	snwice_mass_b(:,:) = 0.0_wp ! no mass exchanges
406	ENDIF
407	IF( nn_ice_embd == 2 .AND. & ! full embedment (case 2) & no restart
408	& .NOT. ln_rstart ) THEN ! deplete the initial ssh below sea-ice area
409	sshn(:,:) = sshn(:,:) - snwice_mass(:,:) * r1_rau0
410	sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0
411	!
412	! Note: Changed the initial values of sshb and sshn=> need to recompute ssh[u,v,f]_[b,n]
413	! which were previously set in domvvl
414	IF ( lk_vvl ) THEN ! Is this necessary? embd 2 should be restricted to vvl only???
415	DO jj = 1, jpjm1
416	DO ji = 1, jpim1 ! caution: use of Vector Opt. not possible
417	zcoefu = 0.5 * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) )
418	zcoefv = 0.5 * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) )
419	zcoeff = 0.25 * umask(ji,jj,1) * umask(ji,jj+1,1)
420	sshu_b(ji,jj) = zcoefu * ( e1t(ji ,jj) * e2t(ji ,jj) * sshb(ji ,jj) &
421	& + e1t(ji+1,jj) * e2t(ji+1,jj) * sshb(ji+1,jj) )
422	sshv_b(ji,jj) = zcoefv * ( e1t(ji,jj ) * e2t(ji,jj ) * sshb(ji,jj ) &
423	& + e1t(ji,jj+1) * e2t(ji,jj+1) * sshb(ji,jj+1) )
424	sshu_n(ji,jj) = zcoefu * ( e1t(ji ,jj) * e2t(ji ,jj) * sshn(ji ,jj) &
425	& + e1t(ji+1,jj) * e2t(ji+1,jj) * sshn(ji+1,jj) )
426	sshv_n(ji,jj) = zcoefv * ( e1t(ji,jj ) * e2t(ji,jj ) * sshn(ji,jj ) &
427	& + e1t(ji,jj+1) * e2t(ji,jj+1) * sshn(ji,jj+1) )
428	END DO
429	END DO
430	CALL lbc_lnk( sshu_b, 'U', 1. ) ; CALL lbc_lnk( sshu_n, 'U', 1. )
431	CALL lbc_lnk( sshv_b, 'V', 1. ) ; CALL lbc_lnk( sshv_n, 'V', 1. )
432	DO jj = 1, jpjm1
433	DO ji = 1, jpim1 ! NO Vector Opt.
434	sshf_n(ji,jj) = 0.5 * umask(ji,jj,1) * umask(ji,jj+1,1) &
435	& / ( e1f(ji,jj ) * e2f(ji,jj ) ) &
436	& * ( e1u(ji,jj ) * e2u(ji,jj ) * sshu_n(ji,jj ) &
437	& + e1u(ji,jj+1) * e2u(ji,jj+1) * sshu_n(ji,jj+1) )
438	END DO
439	END DO
440	CALL lbc_lnk( sshf_n, 'F', 1. )
441	ENDIF
442	ENDIF
443	!
444	END SUBROUTINE lim_sbc_init
445
446	#else
447	!!----------------------------------------------------------------------
448	!! Default option : Dummy module NO LIM 3.0 sea-ice model
449	!!----------------------------------------------------------------------
450	CONTAINS
451	SUBROUTINE lim_sbc ! Dummy routine
452	END SUBROUTINE lim_sbc
453	#endif
454
455	!!======================================================================
456	END MODULE limsbc

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

Download in other formats: