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.

sbccpl.F90 in branches/UKMO/dev_r8183_ICEMODEL_svn_removed/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/UKMO/dev_r8183_ICEMODEL_svn_removed/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 8825

Last change on this file since 8825 was 8825, checked in by alexwestmohc, 6 years ago
Adding some temporary code to initialise the meltpond variables to zero prior to sending these through the coupler
File size: 162.9 KB

Line
1	MODULE sbccpl
2	!!======================================================================
3	!! * MODULE sbccpl *
4	!! Surface Boundary Condition : momentum, heat and freshwater fluxes in coupled mode
5	!!======================================================================
6	!! History : 2.0 ! 2007-06 (R. Redler, N. Keenlyside, W. Park) Original code split into flxmod & taumod
7	!! 3.0 ! 2008-02 (G. Madec, C Talandier) surface module
8	!! 3.1 ! 2009_02 (G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface
9	!! 3.4 ! 2011_11 (C. Harris) more flexibility + multi-category fields
10	!!----------------------------------------------------------------------
11	!!----------------------------------------------------------------------
12	!! namsbc_cpl : coupled formulation namlist
13	!! sbc_cpl_init : initialisation of the coupled exchanges
14	!! sbc_cpl_rcv : receive fields from the atmosphere over the ocean (ocean only)
15	!! receive stress from the atmosphere over the ocean (ocean-ice case)
16	!! sbc_cpl_ice_tau : receive stress from the atmosphere over ice
17	!! sbc_cpl_ice_flx : receive fluxes from the atmosphere over ice
18	!! sbc_cpl_snd : send fields to the atmosphere
19	!!----------------------------------------------------------------------
20	USE dom_oce ! ocean space and time domain
21	USE sbc_oce ! Surface boundary condition: ocean fields
22	USE trc_oce ! share SMS/Ocean variables
23	USE sbc_ice ! Surface boundary condition: ice fields
24	USE sbcapr ! Stochastic param. : ???
25	USE sbcdcy ! surface boundary condition: diurnal cycle
26	USE sbcwave ! surface boundary condition: waves
27	USE phycst ! physical constants
28	#if defined key_lim3
29	USE ice ! ice variables
30	#endif
31	USE cpl_oasis3 ! OASIS3 coupling
32	USE geo2ocean !
33	USE oce , ONLY : tsn, un, vn, sshn, ub, vb, sshb, fraqsr_1lev
34	USE ocealb !
35	USE eosbn2 !
36	USE sbcrnf, ONLY : l_rnfcpl
37	USE sbcisf , ONLY : l_isfcpl
38	#if defined key_cice
39	USE ice_domain_size, only: ncat
40	#endif
41	#if defined key_lim3
42	USE icethd_dh ! for CALL ice_thd_snwblow
43	#endif
44	!
45	USE in_out_manager ! I/O manager
46	USE iom ! NetCDF library
47	USE lib_mpp ! distribued memory computing library
48	USE wrk_nemo ! work arrays
49	USE timing ! Timing
50	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
51
52	IMPLICIT NONE
53	PRIVATE
54
55	PUBLIC sbc_cpl_init ! routine called by sbcmod.F90
56	PUBLIC sbc_cpl_rcv ! routine called by icestp.F90
57	PUBLIC sbc_cpl_snd ! routine called by step.F90
58	PUBLIC sbc_cpl_ice_tau ! routine called by icestp.F90
59	PUBLIC sbc_cpl_ice_flx ! routine called by icestp.F90
60	PUBLIC sbc_cpl_alloc ! routine called in sbcice_cice.F90
61
62	INTEGER, PARAMETER :: jpr_otx1 = 1 ! 3 atmosphere-ocean stress components on grid 1
63	INTEGER, PARAMETER :: jpr_oty1 = 2 !
64	INTEGER, PARAMETER :: jpr_otz1 = 3 !
65	INTEGER, PARAMETER :: jpr_otx2 = 4 ! 3 atmosphere-ocean stress components on grid 2
66	INTEGER, PARAMETER :: jpr_oty2 = 5 !
67	INTEGER, PARAMETER :: jpr_otz2 = 6 !
68	INTEGER, PARAMETER :: jpr_itx1 = 7 ! 3 atmosphere-ice stress components on grid 1
69	INTEGER, PARAMETER :: jpr_ity1 = 8 !
70	INTEGER, PARAMETER :: jpr_itz1 = 9 !
71	INTEGER, PARAMETER :: jpr_itx2 = 10 ! 3 atmosphere-ice stress components on grid 2
72	INTEGER, PARAMETER :: jpr_ity2 = 11 !
73	INTEGER, PARAMETER :: jpr_itz2 = 12 !
74	INTEGER, PARAMETER :: jpr_qsroce = 13 ! Qsr above the ocean
75	INTEGER, PARAMETER :: jpr_qsrice = 14 ! Qsr above the ice
76	INTEGER, PARAMETER :: jpr_qsrmix = 15
77	INTEGER, PARAMETER :: jpr_qnsoce = 16 ! Qns above the ocean
78	INTEGER, PARAMETER :: jpr_qnsice = 17 ! Qns above the ice
79	INTEGER, PARAMETER :: jpr_qnsmix = 18
80	INTEGER, PARAMETER :: jpr_rain = 19 ! total liquid precipitation (rain)
81	INTEGER, PARAMETER :: jpr_snow = 20 ! solid precipitation over the ocean (snow)
82	INTEGER, PARAMETER :: jpr_tevp = 21 ! total evaporation
83	INTEGER, PARAMETER :: jpr_ievp = 22 ! solid evaporation (sublimation)
84	INTEGER, PARAMETER :: jpr_sbpr = 23 ! sublimation - liquid precipitation - solid precipitation
85	INTEGER, PARAMETER :: jpr_semp = 24 ! solid freshwater budget (sublimation - snow)
86	INTEGER, PARAMETER :: jpr_oemp = 25 ! ocean freshwater budget (evap - precip)
87	INTEGER, PARAMETER :: jpr_w10m = 26 ! 10m wind
88	INTEGER, PARAMETER :: jpr_dqnsdt = 27 ! d(Q non solar)/d(temperature)
89	INTEGER, PARAMETER :: jpr_rnf = 28 ! runoffs
90	INTEGER, PARAMETER :: jpr_cal = 29 ! calving
91	INTEGER, PARAMETER :: jpr_taum = 30 ! wind stress module
92	INTEGER, PARAMETER :: jpr_co2 = 31
93	INTEGER, PARAMETER :: jpr_topm = 32 ! topmeltn
94	INTEGER, PARAMETER :: jpr_botm = 33 ! botmeltn
95	INTEGER, PARAMETER :: jpr_sflx = 34 ! salt flux
96	INTEGER, PARAMETER :: jpr_toce = 35 ! ocean temperature
97	INTEGER, PARAMETER :: jpr_soce = 36 ! ocean salinity
98	INTEGER, PARAMETER :: jpr_ocx1 = 37 ! ocean current on grid 1
99	INTEGER, PARAMETER :: jpr_ocy1 = 38 !
100	INTEGER, PARAMETER :: jpr_ssh = 39 ! sea surface height
101	INTEGER, PARAMETER :: jpr_fice = 40 ! ice fraction
102	INTEGER, PARAMETER :: jpr_e3t1st = 41 ! first T level thickness
103	INTEGER, PARAMETER :: jpr_fraqsr = 42 ! fraction of solar net radiation absorbed in the first ocean level
104	INTEGER, PARAMETER :: jpr_mslp = 43 ! mean sea level pressure
105	INTEGER, PARAMETER :: jpr_hsig = 44 ! Hsig
106	INTEGER, PARAMETER :: jpr_phioc = 45 ! Wave=>ocean energy flux
107	INTEGER, PARAMETER :: jpr_sdrftx = 46 ! Stokes drift on grid 1
108	INTEGER, PARAMETER :: jpr_sdrfty = 47 ! Stokes drift on grid 2
109	INTEGER, PARAMETER :: jpr_wper = 48 ! Mean wave period
110	INTEGER, PARAMETER :: jpr_wnum = 49 ! Mean wavenumber
111	INTEGER, PARAMETER :: jpr_wstrf = 50 ! Stress fraction adsorbed by waves
112	INTEGER, PARAMETER :: jpr_wdrag = 51 ! Neutral surface drag coefficient
113	INTEGER, PARAMETER :: jpr_isf = 52
114	INTEGER, PARAMETER :: jpr_icb = 53
115	INTEGER, PARAMETER :: jpr_ts_ice = 54 ! Sea ice surface temp
116	INTEGER, PARAMETER :: jpr_rcv = 55
117
118	INTEGER, PARAMETER :: jprcv = 55 ! total number of fields received
119
120	INTEGER, PARAMETER :: jps_fice = 1 ! ice fraction sent to the atmosphere
121	INTEGER, PARAMETER :: jps_toce = 2 ! ocean temperature
122	INTEGER, PARAMETER :: jps_tice = 3 ! ice temperature
123	INTEGER, PARAMETER :: jps_tmix = 4 ! mixed temperature (ocean+ice)
124	INTEGER, PARAMETER :: jps_albice = 5 ! ice albedo
125	INTEGER, PARAMETER :: jps_albmix = 6 ! mixed albedo
126	INTEGER, PARAMETER :: jps_hice = 7 ! ice thickness
127	INTEGER, PARAMETER :: jps_hsnw = 8 ! snow thickness
128	INTEGER, PARAMETER :: jps_ocx1 = 9 ! ocean current on grid 1
129	INTEGER, PARAMETER :: jps_ocy1 = 10 !
130	INTEGER, PARAMETER :: jps_ocz1 = 11 !
131	INTEGER, PARAMETER :: jps_ivx1 = 12 ! ice current on grid 1
132	INTEGER, PARAMETER :: jps_ivy1 = 13 !
133	INTEGER, PARAMETER :: jps_ivz1 = 14 !
134	INTEGER, PARAMETER :: jps_co2 = 15
135	INTEGER, PARAMETER :: jps_soce = 16 ! ocean salinity
136	INTEGER, PARAMETER :: jps_ssh = 17 ! sea surface height
137	INTEGER, PARAMETER :: jps_qsroce = 18 ! Qsr above the ocean
138	INTEGER, PARAMETER :: jps_qnsoce = 19 ! Qns above the ocean
139	INTEGER, PARAMETER :: jps_oemp = 20 ! ocean freshwater budget (evap - precip)
140	INTEGER, PARAMETER :: jps_sflx = 21 ! salt flux
141	INTEGER, PARAMETER :: jps_otx1 = 22 ! 2 atmosphere-ocean stress components on grid 1
142	INTEGER, PARAMETER :: jps_oty1 = 23 !
143	INTEGER, PARAMETER :: jps_rnf = 24 ! runoffs
144	INTEGER, PARAMETER :: jps_taum = 25 ! wind stress module
145	INTEGER, PARAMETER :: jps_fice2 = 26 ! ice fraction sent to OPA (by SAS when doing SAS-OPA coupling)
146	INTEGER, PARAMETER :: jps_e3t1st = 27 ! first level depth (vvl)
147	INTEGER, PARAMETER :: jps_fraqsr = 28 ! fraction of solar net radiation absorbed in the first ocean level
148	INTEGER, PARAMETER :: jps_ficet = 29 ! total ice fraction
149	INTEGER, PARAMETER :: jps_ocxw = 30 ! currents on grid 1
150	INTEGER, PARAMETER :: jps_ocyw = 31 ! currents on grid 2
151	INTEGER, PARAMETER :: jps_wlev = 32 ! water level
152	INTEGER, PARAMETER :: jps_fice1 = 33 ! first-order ice concentration (for semi-implicit coupling of atmos-ice fluxes)
153	INTEGER, PARAMETER :: jps_a_p = 34 ! meltpond area
154	INTEGER, PARAMETER :: jps_ht_p = 35 ! meltpond thickness
155	INTEGER, PARAMETER :: jps_kice = 36 ! sea ice effective conductivity
156	INTEGER, PARAMETER :: jps_sstfrz = 37 ! sea surface freezing temperature
157	INTEGER, PARAMETER :: jps_ttilyr = 38 ! sea ice top layer temp
158	INTEGER, PARAMETER :: jpsnd = 38 ! total number of fields sent
159
160	! !! namelist namsbc_cpl
161	TYPE :: FLD_C !
162	CHARACTER(len = 32) :: cldes ! desciption of the coupling strategy
163	CHARACTER(len = 32) :: clcat ! multiple ice categories strategy
164	CHARACTER(len = 32) :: clvref ! reference of vector ('spherical' or 'cartesian')
165	CHARACTER(len = 32) :: clvor ! orientation of vector fields ('eastward-northward' or 'local grid')
166	CHARACTER(len = 32) :: clvgrd ! grids on which is located the vector fields
167	END TYPE FLD_C
168	! ! Send to the atmosphere
169	TYPE(FLD_C) :: sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2, sn_snd_thick1, sn_snd_cond, sn_snd_mpnd, sn_snd_sstfrz, sn_snd_ttilyr
170	! ! Received from the atmosphere
171	TYPE(FLD_C) :: sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf, sn_rcv_ts_ice
172	TYPE(FLD_C) :: sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp, sn_rcv_icb, sn_rcv_isf
173	! Send to waves
174	TYPE(FLD_C) :: sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev
175	! Received from waves
176	TYPE(FLD_C) :: sn_rcv_hsig,sn_rcv_phioc,sn_rcv_sdrfx,sn_rcv_sdrfy,sn_rcv_wper,sn_rcv_wnum,sn_rcv_wstrf,sn_rcv_wdrag
177	! ! Other namelist parameters
178	INTEGER :: nn_cplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data
179	LOGICAL :: ln_usecplmask ! use a coupling mask file to merge data received from several models
180	! -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel)
181	TYPE :: DYNARR
182	REAL(wp), POINTER, DIMENSION(:,:,:) :: z3
183	END TYPE DYNARR
184
185	TYPE( DYNARR ), SAVE, DIMENSION(jprcv) :: frcv ! all fields recieved from the atmosphere
186
187	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: alb_oce_mix ! ocean albedo sent to atmosphere (mix clear/overcast sky)
188
189	REAL(wp) :: rpref = 101000._wp ! reference atmospheric pressure[N/m2]
190	REAL(wp) :: r1_grau ! = 1.e0 / (grav * rau0)
191
192	INTEGER , ALLOCATABLE, SAVE, DIMENSION( :) :: nrcvinfo ! OASIS info argument
193
194	!! Substitution
195	# include "vectopt_loop_substitute.h90"
196	!!----------------------------------------------------------------------
197	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
198	!! $Id$
199	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
200	!!----------------------------------------------------------------------
201	CONTAINS
202
203	INTEGER FUNCTION sbc_cpl_alloc()
204	!!----------------------------------------------------------------------
205	!! * FUNCTION sbc_cpl_alloc *
206	!!----------------------------------------------------------------------
207	INTEGER :: ierr(4)
208	!!----------------------------------------------------------------------
209	ierr(:) = 0
210	!
211	ALLOCATE( alb_oce_mix(jpi,jpj), nrcvinfo(jprcv), STAT=ierr(1) )
212
213	#if ! defined key_lim3 && ! defined key_cice
214	ALLOCATE( a_i(jpi,jpj,1) , STAT=ierr(2) ) ! used in sbcice_if.F90 (done here as there is no sbc_ice_if_init)
215	#endif
216	ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) )
217	!
218	IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(4) )
219
220	sbc_cpl_alloc = MAXVAL( ierr )
221	IF( lk_mpp ) CALL mpp_sum ( sbc_cpl_alloc )
222	IF( sbc_cpl_alloc > 0 ) CALL ctl_warn('sbc_cpl_alloc: allocation of arrays failed')
223	!
224	END FUNCTION sbc_cpl_alloc
225
226
227	SUBROUTINE sbc_cpl_init( k_ice )
228	!!----------------------------------------------------------------------
229	!! * ROUTINE sbc_cpl_init *
230	!!
231	!! ** Purpose : Initialisation of send and received information from
232	!! the atmospheric component
233	!!
234	!! ** Method : * Read namsbc_cpl namelist
235	!! * define the receive interface
236	!! * define the send interface
237	!! * initialise the OASIS coupler
238	!!----------------------------------------------------------------------
239	INTEGER, INTENT(in) :: k_ice ! ice management in the sbc (=0/1/2/3)
240	!
241	INTEGER :: jn ! dummy loop index
242	INTEGER :: ios, inum ! Local integer
243	REAL(wp), POINTER, DIMENSION(:,:) :: zacs, zaos
244	!!
245	LOGICAL :: ln_iceshelf_init_atmos
246	NAMELIST/namsbc_cpl/ sn_snd_temp , sn_snd_alb , sn_snd_thick , sn_snd_crt , sn_snd_co2, &
247	& sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau , sn_rcv_dqnsdt, sn_rcv_qsr, &
248	& sn_snd_thick1, sn_snd_cond, sn_snd_mpnd, sn_snd_sstfrz , sn_rcv_ts_ice, sn_snd_ttilyr, &
249	& sn_snd_ifrac, sn_snd_crtw , sn_snd_wlev , sn_rcv_hsig , sn_rcv_phioc , &
250	& sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper , sn_rcv_wnum , sn_rcv_wstrf , &
251	& sn_rcv_wdrag, sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal , &
252	& sn_rcv_iceflx,sn_rcv_co2 , nn_cplmodel , ln_usecplmask, sn_rcv_mslp , &
253	& sn_rcv_icb , sn_rcv_isf, ln_iceshelf_init_atmos
254
255	!!---------------------------------------------------------------------
256	!
257	IF( nn_timing == 1 ) CALL timing_start('sbc_cpl_init')
258	!
259	CALL wrk_alloc( jpi,jpj, zacs, zaos )
260
261	! ================================ !
262	! Namelist informations !
263	! ================================ !
264	!
265	REWIND( numnam_ref ) ! Namelist namsbc_cpl in reference namelist : Variables for OASIS coupling
266	READ ( numnam_ref, namsbc_cpl, IOSTAT = ios, ERR = 901)
267	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_cpl in reference namelist', lwp )
268	!
269	REWIND( numnam_cfg ) ! Namelist namsbc_cpl in configuration namelist : Variables for OASIS coupling
270	READ ( numnam_cfg, namsbc_cpl, IOSTAT = ios, ERR = 902 )
271	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_cpl in configuration namelist', lwp )
272	IF(lwm) WRITE ( numond, namsbc_cpl )
273	!
274	IF(lwp) THEN ! control print
275	WRITE(numout,*)
276	WRITE(numout,*)'sbc_cpl_init : namsbc_cpl namelist '
277	WRITE(numout,*)'~~~~~~~~~~~~'
278	ENDIF
279
280	!!!!! Getting NEMO4-LIM working at the Met Office: Hardcode number of ice cats to 5 during the initialisation
281	jpl = 5
282	!!!!!
283
284	IF( lwp .AND. ln_cpl ) THEN ! control print
285	WRITE(numout,*)' received fields (mutiple ice categogies)'
286	WRITE(numout,*)' 10m wind module = ', TRIM(sn_rcv_w10m%cldes ), ' (', TRIM(sn_rcv_w10m%clcat ), ')'
287	WRITE(numout,*)' stress module = ', TRIM(sn_rcv_taumod%cldes), ' (', TRIM(sn_rcv_taumod%clcat), ')'
288	WRITE(numout,*)' surface stress = ', TRIM(sn_rcv_tau%cldes ), ' (', TRIM(sn_rcv_tau%clcat ), ')'
289	WRITE(numout,*)' - referential = ', sn_rcv_tau%clvref
290	WRITE(numout,*)' - orientation = ', sn_rcv_tau%clvor
291	WRITE(numout,*)' - mesh = ', sn_rcv_tau%clvgrd
292	WRITE(numout,*)' non-solar heat flux sensitivity = ', TRIM(sn_rcv_dqnsdt%cldes), ' (', TRIM(sn_rcv_dqnsdt%clcat), ')'
293	WRITE(numout,*)' solar heat flux = ', TRIM(sn_rcv_qsr%cldes ), ' (', TRIM(sn_rcv_qsr%clcat ), ')'
294	WRITE(numout,*)' non-solar heat flux = ', TRIM(sn_rcv_qns%cldes ), ' (', TRIM(sn_rcv_qns%clcat ), ')'
295	WRITE(numout,*)' freshwater budget = ', TRIM(sn_rcv_emp%cldes ), ' (', TRIM(sn_rcv_emp%clcat ), ')'
296	WRITE(numout,*)' runoffs = ', TRIM(sn_rcv_rnf%cldes ), ' (', TRIM(sn_rcv_rnf%clcat ), ')'
297	WRITE(numout,*)' calving = ', TRIM(sn_rcv_cal%cldes ), ' (', TRIM(sn_rcv_cal%clcat ), ')'
298	WRITE(numout,*)' iceberg = ', TRIM(sn_rcv_icb%cldes ), ' (', TRIM(sn_rcv_icb%clcat ), ')'
299	WRITE(numout,*)' ice shelf = ', TRIM(sn_rcv_isf%cldes ), ' (', TRIM(sn_rcv_isf%clcat ), ')'
300	WRITE(numout,*)' sea ice heat fluxes = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')'
301	WRITE(numout,*)' atm co2 = ', TRIM(sn_rcv_co2%cldes ), ' (', TRIM(sn_rcv_co2%clcat ), ')'
302	WRITE(numout,*)' significant wave heigth = ', TRIM(sn_rcv_hsig%cldes ), ' (', TRIM(sn_rcv_hsig%clcat ), ')'
303	WRITE(numout,*)' wave to oce energy flux = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')'
304	WRITE(numout,*)' Surface Stokes drift grid u = ', TRIM(sn_rcv_sdrfx%cldes ), ' (', TRIM(sn_rcv_sdrfx%clcat ), ')'
305	WRITE(numout,*)' Surface Stokes drift grid v = ', TRIM(sn_rcv_sdrfy%cldes ), ' (', TRIM(sn_rcv_sdrfy%clcat ), ')'
306	WRITE(numout,*)' Mean wave period = ', TRIM(sn_rcv_wper%cldes ), ' (', TRIM(sn_rcv_wper%clcat ), ')'
307	WRITE(numout,*)' Mean wave number = ', TRIM(sn_rcv_wnum%cldes ), ' (', TRIM(sn_rcv_wnum%clcat ), ')'
308	WRITE(numout,*)' Stress frac adsorbed by waves = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')'
309	WRITE(numout,*)' Neutral surf drag coefficient = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')'
310	WRITE(numout,*)' Sea ice surface skin temperature = ', TRIM(sn_rcv_ts_ice%cldes ), ' (', TRIM(sn_rcv_ts_ice%clcat ), ')'
311	WRITE(numout,*)' sent fields (multiple ice categories)'
312	WRITE(numout,*)' surface temperature = ', TRIM(sn_snd_temp%cldes ), ' (', TRIM(sn_snd_temp%clcat ), ')'
313	WRITE(numout,*)' top ice layer temperature = ', TRIM(sn_snd_ttilyr%cldes ), ' (', TRIM(sn_snd_ttilyr%clcat ), ')'
314	WRITE(numout,*)' albedo = ', TRIM(sn_snd_alb%cldes ), ' (', TRIM(sn_snd_alb%clcat ), ')'
315	WRITE(numout,*)' ice/snow thickness = ', TRIM(sn_snd_thick%cldes ), ' (', TRIM(sn_snd_thick%clcat ), ')'
316	WRITE(numout,*)' total ice fraction = ', TRIM(sn_snd_ifrac%cldes ), ' (', TRIM(sn_snd_ifrac%clcat ), ')'
317	WRITE(numout,*)' surface current = ', TRIM(sn_snd_crt%cldes ), ' (', TRIM(sn_snd_crt%clcat ), ')'
318	WRITE(numout,*)' - referential = ', sn_snd_crt%clvref
319	WRITE(numout,*)' - orientation = ', sn_snd_crt%clvor
320	WRITE(numout,*)' - mesh = ', sn_snd_crt%clvgrd
321	WRITE(numout,*)' oce co2 flux = ', TRIM(sn_snd_co2%cldes ), ' (', TRIM(sn_snd_co2%clcat ), ')'
322	WRITE(numout,*)' ice effective conductivity = ', TRIM(sn_snd_cond%cldes ), ' (', TRIM(sn_snd_cond%clcat ), ')'
323	WRITE(numout,*)' meltponds fraction and depth = ', TRIM(sn_snd_mpnd%cldes ), ' (', TRIM(sn_snd_mpnd%clcat ), ')'
324	WRITE(numout,*)' sea surface freezing temp = ', TRIM(sn_snd_sstfrz%cldes ), ' (', TRIM(sn_snd_sstfrz%clcat ), ')'
325	WRITE(numout,*)' water level = ', TRIM(sn_snd_wlev%cldes ), ' (', TRIM(sn_snd_wlev%clcat ), ')'
326	WRITE(numout,*)' mean sea level pressure = ', TRIM(sn_rcv_mslp%cldes ), ' (', TRIM(sn_rcv_mslp%clcat ), ')'
327	WRITE(numout,*)' surface current to waves = ', TRIM(sn_snd_crtw%cldes ), ' (', TRIM(sn_snd_crtw%clcat ), ')'
328	WRITE(numout,*)' - referential = ', sn_snd_crtw%clvref
329	WRITE(numout,*)' - orientation = ', sn_snd_crtw%clvor
330	WRITE(numout,*)' - mesh = ', sn_snd_crtw%clvgrd
331	WRITE(numout,*)' nn_cplmodel = ', nn_cplmodel
332	WRITE(numout,*)' ln_usecplmask = ', ln_usecplmask
333	ENDIF
334
335	! ! allocate sbccpl arrays
336	IF( sbc_cpl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' )
337
338	! ================================ !
339	! Define the receive interface !
340	! ================================ !
341	nrcvinfo(:) = OASIS_idle ! needed by nrcvinfo(jpr_otx1) if we do not receive ocean stress
342
343	! for each field: define the OASIS name (srcv(:)%clname)
344	! define receive or not from the namelist parameters (srcv(:)%laction)
345	! define the north fold type of lbc (srcv(:)%nsgn)
346
347	! default definitions of srcv
348	srcv(:)%laction = .FALSE. ; srcv(:)%clgrid = 'T' ; srcv(:)%nsgn = 1. ; srcv(:)%nct = 1
349
350	! ! ------------------------- !
351	! ! ice and ocean wind stress !
352	! ! ------------------------- !
353	! ! Name
354	srcv(jpr_otx1)%clname = 'O_OTaux1' ! 1st ocean component on grid ONE (T or U)
355	srcv(jpr_oty1)%clname = 'O_OTauy1' ! 2nd - - - -
356	srcv(jpr_otz1)%clname = 'O_OTauz1' ! 3rd - - - -
357	srcv(jpr_otx2)%clname = 'O_OTaux2' ! 1st ocean component on grid TWO (V)
358	srcv(jpr_oty2)%clname = 'O_OTauy2' ! 2nd - - - -
359	srcv(jpr_otz2)%clname = 'O_OTauz2' ! 3rd - - - -
360	!
361	srcv(jpr_itx1)%clname = 'O_ITaux1' ! 1st ice component on grid ONE (T, F, I or U)
362	srcv(jpr_ity1)%clname = 'O_ITauy1' ! 2nd - - - -
363	srcv(jpr_itz1)%clname = 'O_ITauz1' ! 3rd - - - -
364	srcv(jpr_itx2)%clname = 'O_ITaux2' ! 1st ice component on grid TWO (V)
365	srcv(jpr_ity2)%clname = 'O_ITauy2' ! 2nd - - - -
366	srcv(jpr_itz2)%clname = 'O_ITauz2' ! 3rd - - - -
367	!
368	! Vectors: change of sign at north fold ONLY if on the local grid
369	IF( TRIM( sn_rcv_tau%cldes ) == 'oce only' .OR. TRIM(sn_rcv_tau%cldes ) == 'oce and ice') THEN ! avoid working with the atmospheric fields if they are not coupled
370	IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' ) srcv(jpr_otx1:jpr_itz2)%nsgn = -1.
371
372	! ! Set grid and action
373	SELECT CASE( TRIM( sn_rcv_tau%clvgrd ) ) ! 'T', 'U,V', 'U,V,I', 'U,V,F', 'T,I', 'T,F', or 'T,U,V'
374	CASE( 'T' )
375	srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point
376	srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1
377	srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1
378	CASE( 'U,V' )
379	srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point
380	srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point
381	srcv(jpr_itx1:jpr_itz1)%clgrid = 'U' ! ice components given at U-point
382	srcv(jpr_itx2:jpr_itz2)%clgrid = 'V' ! and V-point
383	srcv(jpr_otx1:jpr_itz2)%laction = .TRUE. ! receive oce and ice components on both grid 1 & 2
384	CASE( 'U,V,T' )
385	srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point
386	srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point
387	srcv(jpr_itx1:jpr_itz1)%clgrid = 'T' ! ice components given at T-point
388	srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2
389	srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only
390	CASE( 'U,V,I' )
391	srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point
392	srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point
393	srcv(jpr_itx1:jpr_itz1)%clgrid = 'I' ! ice components given at I-point
394	srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2
395	srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only
396	CASE( 'U,V,F' )
397	srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point
398	srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point
399	srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point
400	srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2
401	srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only
402	CASE( 'T,I' )
403	srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point
404	srcv(jpr_itx1:jpr_itz1)%clgrid = 'I' ! ice components given at I-point
405	srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1
406	srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1
407	CASE( 'T,F' )
408	srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point
409	srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point
410	srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1
411	srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1
412	CASE( 'T,U,V' )
413	srcv(jpr_otx1:jpr_otz1)%clgrid = 'T' ! oce components given at T-point
414	srcv(jpr_itx1:jpr_itz1)%clgrid = 'U' ! ice components given at U-point
415	srcv(jpr_itx2:jpr_itz2)%clgrid = 'V' ! and V-point
416	srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 only
417	srcv(jpr_itx1:jpr_itz2)%laction = .TRUE. ! receive ice components on grid 1 & 2
418	CASE default
419	CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_tau%clvgrd' )
420	END SELECT
421	!
422	IF( TRIM( sn_rcv_tau%clvref ) == 'spherical' ) & ! spherical: 3rd component not received
423	& srcv( (/jpr_otz1, jpr_otz2, jpr_itz1, jpr_itz2/) )%laction = .FALSE.
424	!
425	IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' ) THEN ! already on local grid -> no need of the second grid
426	srcv(jpr_otx2:jpr_otz2)%laction = .FALSE.
427	srcv(jpr_itx2:jpr_itz2)%laction = .FALSE.
428	srcv(jpr_oty1)%clgrid = srcv(jpr_oty2)%clgrid ! not needed but cleaner...
429	srcv(jpr_ity1)%clgrid = srcv(jpr_ity2)%clgrid ! not needed but cleaner...
430	ENDIF
431	!
432	IF( TRIM( sn_rcv_tau%cldes ) /= 'oce and ice' ) THEN ! 'oce and ice' case ocean stress on ocean mesh used
433	srcv(jpr_itx1:jpr_itz2)%laction = .FALSE. ! ice components not received
434	srcv(jpr_itx1)%clgrid = 'U' ! ocean stress used after its transformation
435	srcv(jpr_ity1)%clgrid = 'V' ! i.e. it is always at U- & V-points for i- & j-comp. resp.
436	ENDIF
437	ENDIF
438
439	! ! ------------------------- !
440	! ! freshwater budget ! E-P
441	! ! ------------------------- !
442	! we suppose that atmosphere modele do not make the difference between precipiration (liquide or solid)
443	! over ice of free ocean within the same atmospheric cell.cd
444	srcv(jpr_rain)%clname = 'OTotRain' ! Rain = liquid precipitation
445	srcv(jpr_snow)%clname = 'OTotSnow' ! Snow = solid precipitation
446	srcv(jpr_tevp)%clname = 'OTotEvap' ! total evaporation (over oce + ice sublimation)
447	srcv(jpr_ievp)%clname = 'OIceEvap' ! evaporation over ice = sublimation
448	srcv(jpr_sbpr)%clname = 'OSubMPre' ! sublimation - liquid precipitation - solid precipitation
449	srcv(jpr_semp)%clname = 'OISubMSn' ! ice solid water budget = sublimation - solid precipitation
450	srcv(jpr_oemp)%clname = 'OOEvaMPr' ! ocean water budget = ocean Evap - ocean precip
451	SELECT CASE( TRIM( sn_rcv_emp%cldes ) )
452	CASE( 'none' ) ! nothing to do
453	CASE( 'oce only' ) ; srcv( jpr_oemp )%laction = .TRUE.
454	CASE( 'conservative' )
455	srcv( (/jpr_rain, jpr_snow, jpr_ievp, jpr_tevp/) )%laction = .TRUE.
456	IF ( k_ice <= 1 ) srcv(jpr_ievp)%laction = .FALSE.
457	CASE( 'oce and ice' ) ; srcv( (/jpr_ievp, jpr_sbpr, jpr_semp, jpr_oemp/) )%laction = .TRUE.
458	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_emp%cldes' )
459	END SELECT
460	!
461	! ! ------------------------- !
462	! ! Runoffs & Calving !
463	! ! ------------------------- !
464	srcv(jpr_rnf )%clname = 'O_Runoff'
465	IF( TRIM( sn_rcv_rnf%cldes ) == 'coupled' ) THEN
466	srcv(jpr_rnf)%laction = .TRUE.
467	l_rnfcpl = .TRUE. ! -> no need to read runoffs in sbcrnf
468	ln_rnf = nn_components /= jp_iam_sas ! -> force to go through sbcrnf if not sas
469	IF(lwp) WRITE(numout,*)
470	IF(lwp) WRITE(numout,*) ' runoffs received from oasis -> force ln_rnf = ', ln_rnf
471	ENDIF
472	!
473	srcv(jpr_cal)%clname = 'OCalving' ; IF( TRIM( sn_rcv_cal%cldes) == 'coupled' ) srcv(jpr_cal)%laction = .TRUE.
474	srcv(jpr_isf)%clname = 'OIcshelf' ; IF( TRIM( sn_rcv_isf%cldes) == 'coupled' ) srcv(jpr_isf)%laction = .TRUE.
475	srcv(jpr_icb)%clname = 'OIceberg' ; IF( TRIM( sn_rcv_icb%cldes) == 'coupled' ) srcv(jpr_icb)%laction = .TRUE.
476
477	IF( srcv(jpr_isf)%laction .AND. ln_isf ) THEN
478	l_isfcpl = .TRUE. ! -> no need to read isf in sbcisf
479	IF(lwp) WRITE(numout,*)
480	IF(lwp) WRITE(numout,*) ' iceshelf received from oasis '
481	ENDIF
482	!
483	! ! ------------------------- !
484	! ! non solar radiation ! Qns
485	! ! ------------------------- !
486	srcv(jpr_qnsoce)%clname = 'O_QnsOce'
487	srcv(jpr_qnsice)%clname = 'O_QnsIce'
488	srcv(jpr_qnsmix)%clname = 'O_QnsMix'
489	SELECT CASE( TRIM( sn_rcv_qns%cldes ) )
490	CASE( 'none' ) ! nothing to do
491	CASE( 'oce only' ) ; srcv( jpr_qnsoce )%laction = .TRUE.
492	CASE( 'conservative' ) ; srcv( (/jpr_qnsice, jpr_qnsmix/) )%laction = .TRUE.
493	CASE( 'oce and ice' ) ; srcv( (/jpr_qnsice, jpr_qnsoce/) )%laction = .TRUE.
494	CASE( 'mixed oce-ice' ) ; srcv( jpr_qnsmix )%laction = .TRUE.
495	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_qns%cldes' )
496	END SELECT
497	IF( TRIM( sn_rcv_qns%cldes ) == 'mixed oce-ice' .AND. jpl > 1 ) &
498	CALL ctl_stop( 'sbc_cpl_init: sn_rcv_qns%cldes not currently allowed to be mixed oce-ice for multi-category ice' )
499	! ! ------------------------- !
500	! ! solar radiation ! Qsr
501	! ! ------------------------- !
502	srcv(jpr_qsroce)%clname = 'O_QsrOce'
503	srcv(jpr_qsrice)%clname = 'O_QsrIce'
504	srcv(jpr_qsrmix)%clname = 'O_QsrMix'
505	SELECT CASE( TRIM( sn_rcv_qsr%cldes ) )
506	CASE( 'none' ) ! nothing to do
507	CASE( 'oce only' ) ; srcv( jpr_qsroce )%laction = .TRUE.
508	CASE( 'conservative' ) ; srcv( (/jpr_qsrice, jpr_qsrmix/) )%laction = .TRUE.
509	CASE( 'oce and ice' ) ; srcv( (/jpr_qsrice, jpr_qsroce/) )%laction = .TRUE.
510	CASE( 'mixed oce-ice' ) ; srcv( jpr_qsrmix )%laction = .TRUE.
511	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_qsr%cldes' )
512	END SELECT
513	IF( TRIM( sn_rcv_qsr%cldes ) == 'mixed oce-ice' .AND. jpl > 1 ) &
514	CALL ctl_stop( 'sbc_cpl_init: sn_rcv_qsr%cldes not currently allowed to be mixed oce-ice for multi-category ice' )
515	! ! ------------------------- !
516	! ! non solar sensitivity ! d(Qns)/d(T)
517	! ! ------------------------- !
518	srcv(jpr_dqnsdt)%clname = 'O_dQnsdT'
519	IF( TRIM( sn_rcv_dqnsdt%cldes ) == 'coupled' ) srcv(jpr_dqnsdt)%laction = .TRUE.
520	!
521	! non solar sensitivity mandatory for LIM ice model
522
523	!!!!! Getting NEMO4-LIM working at Met Office: Disable this check because we don't need dqnsdt for JULES-style coupling
524	!IF( TRIM( sn_rcv_dqnsdt%cldes ) == 'none' .AND. k_ice /= 0 .AND. k_ice /= 3 .AND. nn_components /= jp_iam_sas ) &
525	! CALL ctl_stop( 'sbc_cpl_init: sn_rcv_dqnsdt%cldes must be coupled in namsbc_cpl namelist' )
526	!!!!!
527
528	! non solar sensitivity mandatory for mixed oce-ice solar radiation coupling technique
529	IF( TRIM( sn_rcv_dqnsdt%cldes ) == 'none' .AND. TRIM( sn_rcv_qns%cldes ) == 'mixed oce-ice' ) &
530	CALL ctl_stop( 'sbc_cpl_init: namsbc_cpl namelist mismatch between sn_rcv_qns%cldes and sn_rcv_dqnsdt%cldes' )
531	! ! ------------------------- !
532	! ! 10m wind module !
533	! ! ------------------------- !
534	srcv(jpr_w10m)%clname = 'O_Wind10' ; IF( TRIM(sn_rcv_w10m%cldes ) == 'coupled' ) srcv(jpr_w10m)%laction = .TRUE.
535	!
536	! ! ------------------------- !
537	! ! wind stress module !
538	! ! ------------------------- !
539	srcv(jpr_taum)%clname = 'O_TauMod' ; IF( TRIM(sn_rcv_taumod%cldes) == 'coupled' ) srcv(jpr_taum)%laction = .TRUE.
540	lhftau = srcv(jpr_taum)%laction
541
542	! ! ------------------------- !
543	! ! Atmospheric CO2 !
544	! ! ------------------------- !
545	srcv(jpr_co2 )%clname = 'O_AtmCO2'
546	IF( TRIM(sn_rcv_co2%cldes ) == 'coupled' ) THEN
547	srcv(jpr_co2 )%laction = .TRUE.
548	l_co2cpl = .TRUE.
549	IF(lwp) WRITE(numout,*)
550	IF(lwp) WRITE(numout,*) ' Atmospheric pco2 received from oasis '
551	IF(lwp) WRITE(numout,*)
552	ENDIF
553
554	! ! ------------------------- !
555	! ! Mean Sea Level Pressure !
556	! ! ------------------------- !
557	srcv(jpr_mslp)%clname = 'O_MSLP' ; IF( TRIM(sn_rcv_mslp%cldes ) == 'coupled' ) srcv(jpr_mslp)%laction = .TRUE.
558
559	! ! ------------------------- !
560	! ! topmelt and botmelt !
561	! ! ------------------------- !
562	srcv(jpr_topm )%clname = 'OTopMlt'
563	srcv(jpr_botm )%clname = 'OBotMlt'
564	IF( TRIM(sn_rcv_iceflx%cldes) == 'coupled' ) THEN
565	IF ( TRIM( sn_rcv_iceflx%clcat ) == 'yes' ) THEN
566	srcv(jpr_topm:jpr_botm)%nct = jpl
567	ELSE
568	CALL ctl_stop( 'sbc_cpl_init: sn_rcv_iceflx%clcat should always be set to yes currently' )
569	ENDIF
570	srcv(jpr_topm:jpr_botm)%laction = .TRUE.
571	ENDIF
572	! ! ----------------------------- !
573
574	!!!!! To get NEMO4-LIM working at Met Office
575	srcv(jpr_ts_ice)%clname = 'OTsfIce'
576	IF ( TRIM( sn_rcv_ts_ice%cldes ) == 'ice' ) srcv(jpr_ts_ice)%laction = .TRUE.
577	IF ( TRIM( sn_rcv_ts_ice%clcat ) == 'yes' ) srcv(jpr_ts_ice)%nct = jpl
578	IF ( TRIM( sn_rcv_emp%clcat ) == 'yes' ) srcv(jpr_ievp)%nct = jpl
579	!!!!!
580
581	! ! ------------------------- !
582	! ! Wave breaking !
583	! ! ------------------------- !
584	srcv(jpr_hsig)%clname = 'O_Hsigwa' ! significant wave height
585	IF( TRIM(sn_rcv_hsig%cldes ) == 'coupled' ) THEN
586	srcv(jpr_hsig)%laction = .TRUE.
587	cpl_hsig = .TRUE.
588	ENDIF
589	srcv(jpr_phioc)%clname = 'O_PhiOce' ! wave to ocean energy
590	IF( TRIM(sn_rcv_phioc%cldes ) == 'coupled' ) THEN
591	srcv(jpr_phioc)%laction = .TRUE.
592	cpl_phioc = .TRUE.
593	ENDIF
594	srcv(jpr_sdrftx)%clname = 'O_Sdrfx' ! Stokes drift in the u direction
595	IF( TRIM(sn_rcv_sdrfx%cldes ) == 'coupled' ) THEN
596	srcv(jpr_sdrftx)%laction = .TRUE.
597	cpl_sdrftx = .TRUE.
598	ENDIF
599	srcv(jpr_sdrfty)%clname = 'O_Sdrfy' ! Stokes drift in the v direction
600	IF( TRIM(sn_rcv_sdrfy%cldes ) == 'coupled' ) THEN
601	srcv(jpr_sdrfty)%laction = .TRUE.
602	cpl_sdrfty = .TRUE.
603	ENDIF
604	srcv(jpr_wper)%clname = 'O_WPer' ! mean wave period
605	IF( TRIM(sn_rcv_wper%cldes ) == 'coupled' ) THEN
606	srcv(jpr_wper)%laction = .TRUE.
607	cpl_wper = .TRUE.
608	ENDIF
609	srcv(jpr_wnum)%clname = 'O_WNum' ! mean wave number
610	IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' ) THEN
611	srcv(jpr_wnum)%laction = .TRUE.
612	cpl_wnum = .TRUE.
613	ENDIF
614	srcv(jpr_wstrf)%clname = 'O_WStrf' ! stress fraction adsorbed by the wave
615	IF( TRIM(sn_rcv_wstrf%cldes ) == 'coupled' ) THEN
616	srcv(jpr_wstrf)%laction = .TRUE.
617	cpl_wstrf = .TRUE.
618	ENDIF
619	srcv(jpr_wdrag)%clname = 'O_WDrag' ! neutral surface drag coefficient
620	IF( TRIM(sn_rcv_wdrag%cldes ) == 'coupled' ) THEN
621	srcv(jpr_wdrag)%laction = .TRUE.
622	cpl_wdrag = .TRUE.
623	ENDIF
624	!
625	! ! ------------------------------- !
626	! ! OPA-SAS coupling - rcv by opa !
627	! ! ------------------------------- !
628	srcv(jpr_sflx)%clname = 'O_SFLX'
629	srcv(jpr_fice)%clname = 'RIceFrc'
630	!
631	IF( nn_components == jp_iam_opa ) THEN ! OPA coupled to SAS via OASIS: force received field by OPA (sent by SAS)
632	srcv(:)%laction = .FALSE. ! force default definition in case of opa <-> sas coupling
633	srcv(:)%clgrid = 'T' ! force default definition in case of opa <-> sas coupling
634	srcv(:)%nsgn = 1. ! force default definition in case of opa <-> sas coupling
635	srcv( (/jpr_qsroce, jpr_qnsoce, jpr_oemp, jpr_sflx, jpr_fice, jpr_otx1, jpr_oty1, jpr_taum/) )%laction = .TRUE.
636	srcv(jpr_otx1)%clgrid = 'U' ! oce components given at U-point
637	srcv(jpr_oty1)%clgrid = 'V' ! and V-point
638	! Vectors: change of sign at north fold ONLY if on the local grid
639	srcv( (/jpr_otx1,jpr_oty1/) )%nsgn = -1.
640	sn_rcv_tau%clvgrd = 'U,V'
641	sn_rcv_tau%clvor = 'local grid'
642	sn_rcv_tau%clvref = 'spherical'
643	sn_rcv_emp%cldes = 'oce only'
644	!
645	IF(lwp) THEN ! control print
646	WRITE(numout,*)
647	WRITE(numout,*)' Special conditions for SAS-OPA coupling '
648	WRITE(numout,*)' OPA component '
649	WRITE(numout,*)
650	WRITE(numout,*)' received fields from SAS component '
651	WRITE(numout,*)' ice cover '
652	WRITE(numout,*)' oce only EMP '
653	WRITE(numout,*)' salt flux '
654	WRITE(numout,*)' mixed oce-ice solar flux '
655	WRITE(numout,*)' mixed oce-ice non solar flux '
656	WRITE(numout,*)' wind stress U,V on local grid and sperical coordinates '
657	WRITE(numout,*)' wind stress module'
658	WRITE(numout,*)
659	ENDIF
660	ENDIF
661	! ! -------------------------------- !
662	! ! OPA-SAS coupling - rcv by sas !
663	! ! -------------------------------- !
664	srcv(jpr_toce )%clname = 'I_SSTSST'
665	srcv(jpr_soce )%clname = 'I_SSSal'
666	srcv(jpr_ocx1 )%clname = 'I_OCurx1'
667	srcv(jpr_ocy1 )%clname = 'I_OCury1'
668	srcv(jpr_ssh )%clname = 'I_SSHght'
669	srcv(jpr_e3t1st)%clname = 'I_E3T1st'
670	srcv(jpr_fraqsr)%clname = 'I_FraQsr'
671	!
672	IF( nn_components == jp_iam_sas ) THEN
673	IF( .NOT. ln_cpl ) srcv(:)%laction = .FALSE. ! force default definition in case of opa <-> sas coupling
674	IF( .NOT. ln_cpl ) srcv(:)%clgrid = 'T' ! force default definition in case of opa <-> sas coupling
675	IF( .NOT. ln_cpl ) srcv(:)%nsgn = 1. ! force default definition in case of opa <-> sas coupling
676	srcv( (/jpr_toce, jpr_soce, jpr_ssh, jpr_fraqsr, jpr_ocx1, jpr_ocy1/) )%laction = .TRUE.
677	srcv( jpr_e3t1st )%laction = .NOT.ln_linssh
678	srcv(jpr_ocx1)%clgrid = 'U' ! oce components given at U-point
679	srcv(jpr_ocy1)%clgrid = 'V' ! and V-point
680	! Vectors: change of sign at north fold ONLY if on the local grid
681	srcv(jpr_ocx1:jpr_ocy1)%nsgn = -1.
682	! Change first letter to couple with atmosphere if already coupled OPA
683	! this is nedeed as each variable name used in the namcouple must be unique:
684	! for example O_Runoff received by OPA from SAS and therefore O_Runoff received by SAS from the Atmosphere
685	DO jn = 1, jprcv
686	IF ( srcv(jn)%clname(1:1) == "O" ) srcv(jn)%clname = "S"//srcv(jn)%clname(2:LEN(srcv(jn)%clname))
687	END DO
688	!
689	IF(lwp) THEN ! control print
690	WRITE(numout,*)
691	WRITE(numout,*)' Special conditions for SAS-OPA coupling '
692	WRITE(numout,*)' SAS component '
693	WRITE(numout,*)
694	IF( .NOT. ln_cpl ) THEN
695	WRITE(numout,*)' received fields from OPA component '
696	ELSE
697	WRITE(numout,*)' Additional received fields from OPA component : '
698	ENDIF
699	WRITE(numout,*)' sea surface temperature (Celsius) '
700	WRITE(numout,*)' sea surface salinity '
701	WRITE(numout,*)' surface currents '
702	WRITE(numout,*)' sea surface height '
703	WRITE(numout,*)' thickness of first ocean T level '
704	WRITE(numout,*)' fraction of solar net radiation absorbed in the first ocean level'
705	WRITE(numout,*)
706	ENDIF
707	ENDIF
708
709	! =================================================== !
710	! Allocate all parts of frcv used for received fields !
711	! =================================================== !
712	DO jn = 1, jprcv
713	IF ( srcv(jn)%laction ) ALLOCATE( frcv(jn)%z3(jpi,jpj,srcv(jn)%nct) )
714	END DO
715	! Allocate taum part of frcv which is used even when not received as coupling field
716	IF ( .NOT. srcv(jpr_taum)%laction ) ALLOCATE( frcv(jpr_taum)%z3(jpi,jpj,srcv(jpr_taum)%nct) )
717	! Allocate w10m part of frcv which is used even when not received as coupling field
718	IF ( .NOT. srcv(jpr_w10m)%laction ) ALLOCATE( frcv(jpr_w10m)%z3(jpi,jpj,srcv(jpr_w10m)%nct) )
719	! Allocate jpr_otx1 part of frcv which is used even when not received as coupling field
720	IF ( .NOT. srcv(jpr_otx1)%laction ) ALLOCATE( frcv(jpr_otx1)%z3(jpi,jpj,srcv(jpr_otx1)%nct) )
721	IF ( .NOT. srcv(jpr_oty1)%laction ) ALLOCATE( frcv(jpr_oty1)%z3(jpi,jpj,srcv(jpr_oty1)%nct) )
722	! Allocate itx1 and ity1 as they are used in sbc_cpl_ice_tau even if srcv(jpr_itx1)%laction = .FALSE.
723	IF( k_ice /= 0 ) THEN
724	IF ( .NOT. srcv(jpr_itx1)%laction ) ALLOCATE( frcv(jpr_itx1)%z3(jpi,jpj,srcv(jpr_itx1)%nct) )
725	IF ( .NOT. srcv(jpr_ity1)%laction ) ALLOCATE( frcv(jpr_ity1)%z3(jpi,jpj,srcv(jpr_ity1)%nct) )
726	END IF
727
728	! ================================ !
729	! Define the send interface !
730	! ================================ !
731	! for each field: define the OASIS name (ssnd(:)%clname)
732	! define send or not from the namelist parameters (ssnd(:)%laction)
733	! define the north fold type of lbc (ssnd(:)%nsgn)
734
735	! default definitions of nsnd
736	ssnd(:)%laction = .FALSE. ; ssnd(:)%clgrid = 'T' ; ssnd(:)%nsgn = 1. ; ssnd(:)%nct = 1
737
738	! ! ------------------------- !
739	! ! Surface temperature !
740	! ! ------------------------- !
741	ssnd(jps_toce)%clname = 'O_SSTSST'
742	ssnd(jps_tice)%clname = 'O_TepIce'
743	ssnd(jps_ttilyr)%clname = 'O_TtiLyr'
744	ssnd(jps_tmix)%clname = 'O_TepMix'
745	SELECT CASE( TRIM( sn_snd_temp%cldes ) )
746	CASE( 'none' ) ! nothing to do
747	CASE( 'oce only' ) ; ssnd( jps_toce )%laction = .TRUE.
748	CASE( 'oce and ice' , 'weighted oce and ice' , 'oce and weighted ice')
749	ssnd( (/jps_toce, jps_tice/) )%laction = .TRUE.
750	IF ( TRIM( sn_snd_temp%clcat ) == 'yes' ) ssnd(jps_tice)%nct = jpl
751	CASE( 'mixed oce-ice' ) ; ssnd( jps_tmix )%laction = .TRUE.
752	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_temp%cldes' )
753	END SELECT
754
755	! ! ------------------------- !
756	! ! Albedo !
757	! ! ------------------------- !
758	ssnd(jps_albice)%clname = 'O_AlbIce'
759	ssnd(jps_albmix)%clname = 'O_AlbMix'
760	SELECT CASE( TRIM( sn_snd_alb%cldes ) )
761	CASE( 'none' ) ! nothing to do
762	CASE( 'ice' , 'weighted ice' ) ; ssnd(jps_albice)%laction = .TRUE.
763	CASE( 'mixed oce-ice' ) ; ssnd(jps_albmix)%laction = .TRUE.
764	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_alb%cldes' )
765	END SELECT
766	!
767	! Need to calculate oceanic albedo if
768	! 1. sending mixed oce-ice albedo or
769	! 2. receiving mixed oce-ice solar radiation
770	IF ( TRIM ( sn_snd_alb%cldes ) == 'mixed oce-ice' .OR. TRIM ( sn_rcv_qsr%cldes ) == 'mixed oce-ice' ) THEN
771	CALL oce_alb( zaos, zacs )
772	! Due to lack of information on nebulosity : mean clear/overcast sky
773	alb_oce_mix(:,:) = ( zacs(:,:) + zaos(:,:) ) * 0.5
774	ENDIF
775
776	! ! ------------------------- !
777	! ! Ice fraction & Thickness !
778	! ! ------------------------- !
779	ssnd(jps_fice)%clname = 'OIceFrc'
780	ssnd(jps_ficet)%clname = 'OIceFrcT'
781	ssnd(jps_hice)%clname = 'OIceTck'
782	ssnd(jps_a_p)%clname = 'OPndFrc'
783	ssnd(jps_ht_p)%clname = 'OPndTck'
784	ssnd(jps_hsnw)%clname = 'OSnwTck'
785	ssnd(jps_fice1)%clname = 'OIceFrd'
786	IF( k_ice /= 0 ) THEN
787	ssnd(jps_fice)%laction = .TRUE. ! if ice treated in the ocean (even in climato case)
788	ssnd(jps_fice1)%laction = .TRUE. ! First-order regridded ice concentration, to be used
789	! producing atmos-to-ice fluxes
790	! Currently no namelist entry to determine sending of multi-category ice fraction so use the thickness entry for now
791	IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_fice)%nct = jpl
792	IF ( TRIM( sn_snd_thick1%clcat ) == 'yes' ) ssnd(jps_fice1)%nct = jpl
793	ENDIF
794
795	IF (TRIM( sn_snd_ifrac%cldes ) == 'coupled') ssnd(jps_ficet)%laction = .TRUE.
796
797	SELECT CASE ( TRIM( sn_snd_thick%cldes ) )
798	CASE( 'none' ) ! nothing to do
799	CASE( 'ice and snow' )
800	ssnd(jps_hice:jps_hsnw)%laction = .TRUE.
801	IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) THEN
802	ssnd(jps_hice:jps_hsnw)%nct = jpl
803	ENDIF
804	CASE ( 'weighted ice and snow' )
805	ssnd(jps_hice:jps_hsnw)%laction = .TRUE.
806	IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_hice:jps_hsnw)%nct = jpl
807	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_thick%cldes' )
808	END SELECT
809
810	! ! ------------------------- !
811	! ! Ice Meltponds !
812	! ! ------------------------- !
813
814
815	!!!!! Getting NEMO4-LIM to work at Met Office
816	ssnd(jps_a_p)%clname = 'OPndFrc'
817	ssnd(jps_ht_p)%clname = 'OPndTck'
818	SELECT CASE ( TRIM( sn_snd_mpnd%cldes ) )
819	CASE ( 'none' )
820	ssnd(jps_a_p)%laction = .FALSE.
821	ssnd(jps_ht_p)%laction = .FALSE.
822	CASE ( 'ice only' )
823	ssnd(jps_a_p)%laction = .TRUE.
824	ssnd(jps_ht_p)%laction = .TRUE.
825	IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN
826	ssnd(jps_a_p)%nct = jpl
827	ssnd(jps_ht_p)%nct = jpl
828	ELSE
829	IF ( jpl > 1 ) THEN
830	CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_mpnd%cldes if not exchanging category fields' )
831	ENDIF
832	ENDIF
833	CASE ( 'weighted ice' )
834	ssnd(jps_a_p)%laction = .TRUE.
835	ssnd(jps_ht_p)%laction = .TRUE.
836	IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN
837	ssnd(jps_a_p)%nct = jpl
838	ssnd(jps_ht_p)%nct = jpl
839	ENDIF
840	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_mpnd%cldes; '//sn_snd_mpnd%cldes )
841	END SELECT
842	!!!!!
843
844	! ! ------------------------- !
845	! ! Surface current !
846	! ! ------------------------- !
847	! ocean currents ! ice velocities
848	ssnd(jps_ocx1)%clname = 'O_OCurx1' ; ssnd(jps_ivx1)%clname = 'O_IVelx1'
849	ssnd(jps_ocy1)%clname = 'O_OCury1' ; ssnd(jps_ivy1)%clname = 'O_IVely1'
850	ssnd(jps_ocz1)%clname = 'O_OCurz1' ; ssnd(jps_ivz1)%clname = 'O_IVelz1'
851	ssnd(jps_ocxw)%clname = 'O_OCurxw'
852	ssnd(jps_ocyw)%clname = 'O_OCuryw'
853	!
854	ssnd(jps_ocx1:jps_ivz1)%nsgn = -1. ! vectors: change of the sign at the north fold
855
856	IF( sn_snd_crt%clvgrd == 'U,V' ) THEN
857	ssnd(jps_ocx1)%clgrid = 'U' ; ssnd(jps_ocy1)%clgrid = 'V'
858	ELSE IF( sn_snd_crt%clvgrd /= 'T' ) THEN
859	CALL ctl_stop( 'sn_snd_crt%clvgrd must be equal to T' )
860	ssnd(jps_ocx1:jps_ivz1)%clgrid = 'T' ! all oce and ice components on the same unique grid
861	ENDIF
862	ssnd(jps_ocx1:jps_ivz1)%laction = .TRUE. ! default: all are send
863	IF( TRIM( sn_snd_crt%clvref ) == 'spherical' ) ssnd( (/jps_ocz1, jps_ivz1/) )%laction = .FALSE.
864	IF( TRIM( sn_snd_crt%clvor ) == 'eastward-northward' ) ssnd(jps_ocx1:jps_ivz1)%nsgn = 1.
865	SELECT CASE( TRIM( sn_snd_crt%cldes ) )
866	CASE( 'none' ) ; ssnd(jps_ocx1:jps_ivz1)%laction = .FALSE.
867	CASE( 'oce only' ) ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.
868	CASE( 'weighted oce and ice' ) ! nothing to do
869	CASE( 'mixed oce-ice' ) ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.
870	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crt%cldes' )
871	END SELECT
872
873	ssnd(jps_ocxw:jps_ocyw)%nsgn = -1. ! vectors: change of the sign at the north fold
874
875	IF( sn_snd_crtw%clvgrd == 'U,V' ) THEN
876	ssnd(jps_ocxw)%clgrid = 'U' ; ssnd(jps_ocyw)%clgrid = 'V'
877	ELSE IF( sn_snd_crtw%clvgrd /= 'T' ) THEN
878	CALL ctl_stop( 'sn_snd_crtw%clvgrd must be equal to T' )
879	ENDIF
880	IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) ssnd(jps_ocxw:jps_ocyw)%nsgn = 1.
881	SELECT CASE( TRIM( sn_snd_crtw%cldes ) )
882	CASE( 'none' ) ; ssnd(jps_ocxw:jps_ocyw)%laction = .FALSE.
883	CASE( 'oce only' ) ; ssnd(jps_ocxw:jps_ocyw)%laction = .TRUE.
884	CASE( 'weighted oce and ice' ) ! nothing to do
885	CASE( 'mixed oce-ice' ) ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.
886	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crtw%cldes' )
887	END SELECT
888
889	! ! ------------------------- !
890	! ! CO2 flux !
891	! ! ------------------------- !
892	ssnd(jps_co2)%clname = 'O_CO2FLX' ; IF( TRIM(sn_snd_co2%cldes) == 'coupled' ) ssnd(jps_co2 )%laction = .TRUE.
893	!
894
895	!!!!! Getting NEMO4-LIM to work at the Met Office
896	! ! ------------------------- !
897	! ! Sea surface freezing temp !
898	! ! ------------------------- !
899	ssnd(jps_sstfrz)%clname = 'O_SSTFrz' ; IF( TRIM(sn_snd_sstfrz%cldes) == 'coupled' ) ssnd(jps_sstfrz)%laction = .TRUE.
900	!!!!!
901
902	!
903	! ! ------------------------- !
904	! ! Ice conductivity !
905	! ! ------------------------- !
906	! Note that ultimately we will move to passing an ocean effective conductivity as well so there
907	! will be some changes to the parts of the code which currently relate only to ice conductivity
908
909	ssnd(jps_ttilyr )%clname = 'O_TtiLyr'
910	SELECT CASE ( TRIM( sn_snd_ttilyr%cldes ) )
911	CASE ( 'none' )
912	ssnd(jps_ttilyr)%laction = .FALSE.
913	CASE ( 'ice only' )
914	ssnd(jps_ttilyr)%laction = .TRUE.
915	IF ( TRIM( sn_snd_ttilyr%clcat ) == 'yes' ) THEN
916	ssnd(jps_ttilyr)%nct = jpl
917	ELSE
918	IF ( jpl > 1 ) THEN
919	CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_ttilyr%cldes if not exchanging category fields' )
920	ENDIF
921	ENDIF
922	CASE ( 'weighted ice' )
923	ssnd(jps_ttilyr)%laction = .TRUE.
924	IF ( TRIM( sn_snd_ttilyr%clcat ) == 'yes' ) ssnd(jps_ttilyr)%nct = jpl
925	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_ttilyr%cldes;'//sn_snd_ttilyr%cldes )
926	END SELECT
927
928	ssnd(jps_kice )%clname = 'OIceKn'
929	SELECT CASE ( TRIM( sn_snd_cond%cldes ) )
930	CASE ( 'none' )
931	ssnd(jps_kice)%laction = .FALSE.
932	CASE ( 'ice only' )
933	ssnd(jps_kice)%laction = .TRUE.
934	IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) THEN
935	ssnd(jps_kice)%nct = jpl
936	ELSE
937	IF ( jpl > 1 ) THEN
938	CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_cond%cldes if not exchanging category fields' )
939	ENDIF
940	ENDIF
941	CASE ( 'weighted ice' )
942	ssnd(jps_kice)%laction = .TRUE.
943	IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) ssnd(jps_kice)%nct = jpl
944	CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_cond%cldes;'//sn_snd_cond%cldes )
945	END SELECT
946	!
947
948	! ! ------------------------- !
949	! ! Sea surface height !
950	! ! ------------------------- !
951	ssnd(jps_wlev)%clname = 'O_Wlevel' ; IF( TRIM(sn_snd_wlev%cldes) == 'coupled' ) ssnd(jps_wlev)%laction = .TRUE.
952
953	! ! ------------------------------- !
954	! ! OPA-SAS coupling - snd by opa !
955	! ! ------------------------------- !
956	ssnd(jps_ssh )%clname = 'O_SSHght'
957	ssnd(jps_soce )%clname = 'O_SSSal'
958	ssnd(jps_e3t1st)%clname = 'O_E3T1st'
959	ssnd(jps_fraqsr)%clname = 'O_FraQsr'
960	!
961	IF( nn_components == jp_iam_opa ) THEN
962	ssnd(:)%laction = .FALSE. ! force default definition in case of opa <-> sas coupling
963	ssnd( (/jps_toce, jps_soce, jps_ssh, jps_fraqsr, jps_ocx1, jps_ocy1/) )%laction = .TRUE.
964	ssnd( jps_e3t1st )%laction = .NOT.ln_linssh
965	! vector definition: not used but cleaner...
966	ssnd(jps_ocx1)%clgrid = 'U' ! oce components given at U-point
967	ssnd(jps_ocy1)%clgrid = 'V' ! and V-point
968	sn_snd_crt%clvgrd = 'U,V'
969	sn_snd_crt%clvor = 'local grid'
970	sn_snd_crt%clvref = 'spherical'
971	!
972	IF(lwp) THEN ! control print
973	WRITE(numout,*)
974	WRITE(numout,*)' sent fields to SAS component '
975	WRITE(numout,*)' sea surface temperature (T before, Celsius) '
976	WRITE(numout,*)' sea surface salinity '
977	WRITE(numout,*)' surface currents U,V on local grid and spherical coordinates'
978	WRITE(numout,*)' sea surface height '
979	WRITE(numout,*)' thickness of first ocean T level '
980	WRITE(numout,*)' fraction of solar net radiation absorbed in the first ocean level'
981	WRITE(numout,*)
982	ENDIF
983	ENDIF
984	! ! ------------------------------- !
985	! ! OPA-SAS coupling - snd by sas !
986	! ! ------------------------------- !
987	ssnd(jps_sflx )%clname = 'I_SFLX'
988	ssnd(jps_fice2 )%clname = 'IIceFrc'
989	ssnd(jps_qsroce)%clname = 'I_QsrOce'
990	ssnd(jps_qnsoce)%clname = 'I_QnsOce'
991	ssnd(jps_oemp )%clname = 'IOEvaMPr'
992	ssnd(jps_otx1 )%clname = 'I_OTaux1'
993	ssnd(jps_oty1 )%clname = 'I_OTauy1'
994	ssnd(jps_rnf )%clname = 'I_Runoff'
995	ssnd(jps_taum )%clname = 'I_TauMod'
996	!
997	IF( nn_components == jp_iam_sas ) THEN
998	IF( .NOT. ln_cpl ) ssnd(:)%laction = .FALSE. ! force default definition in case of opa <-> sas coupling
999	ssnd( (/jps_qsroce, jps_qnsoce, jps_oemp, jps_fice2, jps_sflx, jps_otx1, jps_oty1, jps_taum/) )%laction = .TRUE.
1000	!
1001	! Change first letter to couple with atmosphere if already coupled with sea_ice
1002	! this is nedeed as each variable name used in the namcouple must be unique:
1003	! for example O_SSTSST sent by OPA to SAS and therefore S_SSTSST sent by SAS to the Atmosphere
1004	DO jn = 1, jpsnd
1005	IF ( ssnd(jn)%clname(1:1) == "O" ) ssnd(jn)%clname = "S"//ssnd(jn)%clname(2:LEN(ssnd(jn)%clname))
1006	END DO
1007	!
1008	IF(lwp) THEN ! control print
1009	WRITE(numout,*)
1010	IF( .NOT. ln_cpl ) THEN
1011	WRITE(numout,*)' sent fields to OPA component '
1012	ELSE
1013	WRITE(numout,*)' Additional sent fields to OPA component : '
1014	ENDIF
1015	WRITE(numout,*)' ice cover '
1016	WRITE(numout,*)' oce only EMP '
1017	WRITE(numout,*)' salt flux '
1018	WRITE(numout,*)' mixed oce-ice solar flux '
1019	WRITE(numout,*)' mixed oce-ice non solar flux '
1020	WRITE(numout,*)' wind stress U,V components'
1021	WRITE(numout,*)' wind stress module'
1022	ENDIF
1023	ENDIF
1024
1025	!
1026	! ================================ !
1027	! initialisation of the coupler !
1028	! ================================ !
1029
1030	CALL cpl_define(jprcv, jpsnd, nn_cplmodel)
1031
1032	IF (ln_usecplmask) THEN
1033	xcplmask(:,:,:) = 0.
1034	CALL iom_open( 'cplmask', inum )
1035	CALL iom_get( inum, jpdom_unknown, 'cplmask', xcplmask(1:nlci,1:nlcj,1:nn_cplmodel), &
1036	& kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,nn_cplmodel /) )
1037	CALL iom_close( inum )
1038	ELSE
1039	xcplmask(:,:,:) = 1.
1040	ENDIF
1041	xcplmask(:,:,0) = 1. - SUM( xcplmask(:,:,1:nn_cplmodel), dim = 3 )
1042	!
1043	ncpl_qsr_freq = cpl_freq( 'O_QsrOce' ) + cpl_freq( 'O_QsrMix' ) + cpl_freq( 'I_QsrOce' ) + cpl_freq( 'I_QsrMix' )
1044	IF( ln_dm2dc .AND. ln_cpl .AND. ncpl_qsr_freq /= 86400 ) &
1045	& CALL ctl_stop( 'sbc_cpl_init: diurnal cycle reconstruction (ln_dm2dc) needs daily couping for solar radiation' )
1046	IF( ln_dm2dc .AND. ln_cpl ) ncpl_qsr_freq = 86400 / ncpl_qsr_freq
1047
1048	CALL wrk_dealloc( jpi,jpj, zacs, zaos )
1049	!
1050	IF( nn_timing == 1 ) CALL timing_stop('sbc_cpl_init')
1051	!
1052	END SUBROUTINE sbc_cpl_init
1053
1054
1055	SUBROUTINE sbc_cpl_rcv( kt, k_fsbc, k_ice )
1056	!!----------------------------------------------------------------------
1057	!! * ROUTINE sbc_cpl_rcv *
1058	!!
1059	!! ** Purpose : provide the stress over the ocean and, if no sea-ice,
1060	!! provide the ocean heat and freshwater fluxes.
1061	!!
1062	!! ** Method : - Receive all the atmospheric fields (stored in frcv array). called at each time step.
1063	!! OASIS controls if there is something do receive or not. nrcvinfo contains the info
1064	!! to know if the field was really received or not
1065	!!
1066	!! --> If ocean stress was really received:
1067	!!
1068	!! - transform the received ocean stress vector from the received
1069	!! referential and grid into an atmosphere-ocean stress in
1070	!! the (i,j) ocean referencial and at the ocean velocity point.
1071	!! The received stress are :
1072	!! - defined by 3 components (if cartesian coordinate)
1073	!! or by 2 components (if spherical)
1074	!! - oriented along geographical coordinate (if eastward-northward)
1075	!! or along the local grid coordinate (if local grid)
1076	!! - given at U- and V-point, resp. if received on 2 grids
1077	!! or at T-point if received on 1 grid
1078	!! Therefore and if necessary, they are successively
1079	!! processed in order to obtain them
1080	!! first as 2 components on the sphere
1081	!! second as 2 components oriented along the local grid
1082	!! third as 2 components on the U,V grid
1083	!!
1084	!! -->
1085	!!
1086	!! - In 'ocean only' case, non solar and solar ocean heat fluxes
1087	!! and total ocean freshwater fluxes
1088	!!
1089	!! ** Method : receive all fields from the atmosphere and transform
1090	!! them into ocean surface boundary condition fields
1091	!!
1092	!! ** Action : update utau, vtau ocean stress at U,V grid
1093	!! taum wind stress module at T-point
1094	!! wndm wind speed module at T-point over free ocean or leads in presence of sea-ice
1095	!! qns non solar heat fluxes including emp heat content (ocean only case)
1096	!! and the latent heat flux of solid precip. melting
1097	!! qsr solar ocean heat fluxes (ocean only case)
1098	!! emp upward mass flux [evap. - precip. (- runoffs) (- calving)] (ocean only case)
1099	!!----------------------------------------------------------------------
1100	USE zdf_oce, ONLY : ln_zdfqiao
1101
1102	IMPLICIT NONE
1103
1104	INTEGER, INTENT(in) :: kt ! ocean model time step index
1105	INTEGER, INTENT(in) :: k_fsbc ! frequency of sbc (-> ice model) computation
1106	INTEGER, INTENT(in) :: k_ice ! ice management in the sbc (=0/1/2/3)
1107	!!
1108	LOGICAL :: llnewtx, llnewtau ! update wind stress components and module??
1109	INTEGER :: ji, jj, jn ! dummy loop indices
1110	INTEGER :: isec ! number of seconds since nit000 (assuming rdt did not change since nit000)
1111	REAL(wp) :: zcumulneg, zcumulpos ! temporary scalars
1112	REAL(wp) :: zcoef ! temporary scalar
1113	REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3
1114	REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient
1115	REAL(wp) :: zzx, zzy ! temporary variables
1116	REAL(wp), POINTER, DIMENSION(:,:) :: ztx, zty, zmsk, zemp, zqns, zqsr
1117	!!----------------------------------------------------------------------
1118	!
1119	IF( nn_timing == 1 ) CALL timing_start('sbc_cpl_rcv')
1120	!
1121	CALL wrk_alloc( jpi,jpj, ztx, zty, zmsk, zemp, zqns, zqsr )
1122	!
1123	IF( ln_mixcpl ) zmsk(:,:) = 1. - xcplmask(:,:,0)
1124	!
1125	! ! ======================================================= !
1126	! ! Receive all the atmos. fields (including ice information)
1127	! ! ======================================================= !
1128	isec = ( kt - nit000 ) * NINT( rdt ) ! date of exchanges
1129	DO jn = 1, jprcv ! received fields sent by the atmosphere
1130	IF( srcv(jn)%laction ) CALL cpl_rcv( jn, isec, frcv(jn)%z3, xcplmask(:,:,1:nn_cplmodel), nrcvinfo(jn) )
1131	END DO
1132
1133	! ! ========================= !
1134	IF( srcv(jpr_otx1)%laction ) THEN ! ocean stress components !
1135	! ! ========================= !
1136	! define frcv(jpr_otx1)%z3(:,:,1) and frcv(jpr_oty1)%z3(:,:,1): stress at U/V point along model grid
1137	! => need to be done only when we receive the field
1138	IF( nrcvinfo(jpr_otx1) == OASIS_Rcv ) THEN
1139	!
1140	IF( TRIM( sn_rcv_tau%clvref ) == 'cartesian' ) THEN ! 2 components on the sphere
1141	! ! (cartesian to spherical -> 3 to 2 components)
1142	!
1143	CALL geo2oce( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), frcv(jpr_otz1)%z3(:,:,1), &
1144	& srcv(jpr_otx1)%clgrid, ztx, zty )
1145	frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 1st grid
1146	frcv(jpr_oty1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 1st grid
1147	!
1148	IF( srcv(jpr_otx2)%laction ) THEN
1149	CALL geo2oce( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), frcv(jpr_otz2)%z3(:,:,1), &
1150	& srcv(jpr_otx2)%clgrid, ztx, zty )
1151	frcv(jpr_otx2)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 2nd grid
1152	frcv(jpr_oty2)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 2nd grid
1153	ENDIF
1154	!
1155	ENDIF
1156	!
1157	IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid
1158	! ! (geographical to local grid -> rotate the components)
1159	CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx )
1160	IF( srcv(jpr_otx2)%laction ) THEN
1161	CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty )
1162	ELSE
1163	CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty )
1164	ENDIF
1165	frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st component on the 1st grid
1166	frcv(jpr_oty1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd component on the 2nd grid
1167	ENDIF
1168	!
1169	IF( srcv(jpr_otx1)%clgrid == 'T' ) THEN
1170	DO jj = 2, jpjm1 ! T ==> (U,V)
1171	DO ji = fs_2, fs_jpim1 ! vector opt.
1172	frcv(jpr_otx1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_otx1)%z3(ji+1,jj ,1) + frcv(jpr_otx1)%z3(ji,jj,1) )
1173	frcv(jpr_oty1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_oty1)%z3(ji ,jj+1,1) + frcv(jpr_oty1)%z3(ji,jj,1) )
1174	END DO
1175	END DO
1176	CALL lbc_lnk( frcv(jpr_otx1)%z3(:,:,1), 'U', -1. ) ; CALL lbc_lnk( frcv(jpr_oty1)%z3(:,:,1), 'V', -1. )
1177	ENDIF
1178	llnewtx = .TRUE.
1179	ELSE
1180	llnewtx = .FALSE.
1181	ENDIF
1182	! ! ========================= !
1183	ELSE ! No dynamical coupling !
1184	! ! ========================= !
1185	frcv(jpr_otx1)%z3(:,:,1) = 0.e0 ! here simply set to zero
1186	frcv(jpr_oty1)%z3(:,:,1) = 0.e0 ! an external read in a file can be added instead
1187	llnewtx = .TRUE.
1188	!
1189	ENDIF
1190	! ! ========================= !
1191	! ! wind stress module ! (taum)
1192	! ! ========================= !
1193	!
1194	IF( .NOT. srcv(jpr_taum)%laction ) THEN ! compute wind stress module from its components if not received
1195	! => need to be done only when otx1 was changed
1196	IF( llnewtx ) THEN
1197	DO jj = 2, jpjm1
1198	DO ji = fs_2, fs_jpim1 ! vect. opt.
1199	zzx = frcv(jpr_otx1)%z3(ji-1,jj ,1) + frcv(jpr_otx1)%z3(ji,jj,1)
1200	zzy = frcv(jpr_oty1)%z3(ji ,jj-1,1) + frcv(jpr_oty1)%z3(ji,jj,1)
1201	frcv(jpr_taum)%z3(ji,jj,1) = 0.5 * SQRT( zzx * zzx + zzy * zzy )
1202	END DO
1203	END DO
1204	CALL lbc_lnk( frcv(jpr_taum)%z3(:,:,1), 'T', 1. )
1205	llnewtau = .TRUE.
1206	ELSE
1207	llnewtau = .FALSE.
1208	ENDIF
1209	ELSE
1210	llnewtau = nrcvinfo(jpr_taum) == OASIS_Rcv
1211	! Stress module can be negative when received (interpolation problem)
1212	IF( llnewtau ) THEN
1213	frcv(jpr_taum)%z3(:,:,1) = MAX( 0._wp, frcv(jpr_taum)%z3(:,:,1) )
1214	ENDIF
1215	ENDIF
1216	!
1217	! ! ========================= !
1218	! ! 10 m wind speed ! (wndm)
1219	! ! ========================= !
1220	!
1221	IF( .NOT. srcv(jpr_w10m)%laction ) THEN ! compute wind spreed from wind stress module if not received
1222	! => need to be done only when taumod was changed
1223	IF( llnewtau ) THEN
1224	zcoef = 1. / ( zrhoa * zcdrag )
1225	DO jj = 1, jpj
1226	DO ji = 1, jpi
1227	frcv(jpr_w10m)%z3(ji,jj,1) = SQRT( frcv(jpr_taum)%z3(ji,jj,1) * zcoef )
1228	END DO
1229	END DO
1230	ENDIF
1231	ENDIF
1232
1233	! u(v)tau and taum will be modified by ice model
1234	! -> need to be reset before each call of the ice/fsbc
1235	IF( MOD( kt-1, k_fsbc ) == 0 ) THEN
1236	!
1237	IF( ln_mixcpl ) THEN
1238	utau(:,:) = utau(:,:) * xcplmask(:,:,0) + frcv(jpr_otx1)%z3(:,:,1) * zmsk(:,:)
1239	vtau(:,:) = vtau(:,:) * xcplmask(:,:,0) + frcv(jpr_oty1)%z3(:,:,1) * zmsk(:,:)
1240	taum(:,:) = taum(:,:) * xcplmask(:,:,0) + frcv(jpr_taum)%z3(:,:,1) * zmsk(:,:)
1241	wndm(:,:) = wndm(:,:) * xcplmask(:,:,0) + frcv(jpr_w10m)%z3(:,:,1) * zmsk(:,:)
1242	ELSE
1243	utau(:,:) = frcv(jpr_otx1)%z3(:,:,1)
1244	vtau(:,:) = frcv(jpr_oty1)%z3(:,:,1)
1245	taum(:,:) = frcv(jpr_taum)%z3(:,:,1)
1246	wndm(:,:) = frcv(jpr_w10m)%z3(:,:,1)
1247	ENDIF
1248	CALL iom_put( "taum_oce", taum ) ! output wind stress module
1249	!
1250	ENDIF
1251
1252	! ! ================== !
1253	! ! atmosph. CO2 (ppm) !
1254	! ! ================== !
1255	IF( srcv(jpr_co2)%laction ) atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1)
1256	!
1257
1258	!!!!! Getting NEMO4-LIM to work at the Met Office
1259	! ! Sea ice surface skin temp:
1260	IF( srcv(jpr_ts_ice)%laction ) THEN
1261	DO jn = 1, jpl
1262	DO jj = 1, jpj
1263	DO ji = 1, jpi
1264	IF (frcv(jpr_ts_ice)%z3(ji,jj,jn) > 0.0) THEN
1265	tsfc_ice(ji,jj,jn) = 0.0
1266	ELSE IF (frcv(jpr_ts_ice)%z3(ji,jj,jn) < -60.0) THEN
1267	tsfc_ice(ji,jj,jn) = -60.0
1268	ELSE
1269	tsfc_ice(ji,jj,jn) = frcv(jpr_ts_ice)%z3(ji,jj,jn)
1270	ENDIF
1271	END DO
1272	END DO
1273	END DO
1274	ENDIF
1275	!!!!!
1276
1277
1278	! ! ========================= !
1279	! ! Mean Sea Level Pressure ! (taum)
1280	! ! ========================= !
1281	!
1282	IF( srcv(jpr_mslp)%laction ) THEN ! UKMO SHELF effect of atmospheric pressure on SSH
1283	IF( kt /= nit000 ) ssh_ibb(:,:) = ssh_ib(:,:) !* Swap of ssh_ib fields
1284
1285	r1_grau = 1.e0 / (grav * rau0) !* constant for optimization
1286	ssh_ib(:,:) = - ( frcv(jpr_mslp)%z3(:,:,1) - rpref ) * r1_grau ! equivalent ssh (inverse barometer)
1287	apr (:,:) = frcv(jpr_mslp)%z3(:,:,1) !atmospheric pressure
1288
1289	IF( kt == nit000 ) ssh_ibb(:,:) = ssh_ib(:,:) ! correct this later (read from restart if possible)
1290	END IF
1291	!
1292	IF( ln_sdw ) THEN ! Stokes Drift correction activated
1293	! ! ========================= !
1294	! ! Stokes drift u !
1295	! ! ========================= !
1296	IF( srcv(jpr_sdrftx)%laction ) ut0sd(:,:) = frcv(jpr_sdrftx)%z3(:,:,1)
1297	!
1298	! ! ========================= !
1299	! ! Stokes drift v !
1300	! ! ========================= !
1301	IF( srcv(jpr_sdrfty)%laction ) vt0sd(:,:) = frcv(jpr_sdrfty)%z3(:,:,1)
1302	!
1303	! ! ========================= !
1304	! ! Wave mean period !
1305	! ! ========================= !
1306	IF( srcv(jpr_wper)%laction ) wmp(:,:) = frcv(jpr_wper)%z3(:,:,1)
1307	!
1308	! ! ========================= !
1309	! ! Significant wave height !
1310	! ! ========================= !
1311	IF( srcv(jpr_hsig)%laction ) hsw(:,:) = frcv(jpr_hsig)%z3(:,:,1)
1312	!
1313	! ! ========================= !
1314	! ! Vertical mixing Qiao !
1315	! ! ========================= !
1316	IF( srcv(jpr_wnum)%laction .AND. ln_zdfqiao ) wnum(:,:) = frcv(jpr_wnum)%z3(:,:,1)
1317
1318	! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode
1319	IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. srcv(jpr_wper)%laction &
1320	.OR. srcv(jpr_hsig)%laction ) THEN
1321	CALL sbc_stokes()
1322	ENDIF
1323	ENDIF
1324	! ! ========================= !
1325	! ! Stress adsorbed by waves !
1326	! ! ========================= !
1327	IF( srcv(jpr_wstrf)%laction .AND. ln_tauoc ) tauoc_wave(:,:) = frcv(jpr_wstrf)%z3(:,:,1)
1328
1329	! ! ========================= !
1330	! ! Wave drag coefficient !
1331	! ! ========================= !
1332	IF( srcv(jpr_wdrag)%laction .AND. ln_cdgw ) cdn_wave(:,:) = frcv(jpr_wdrag)%z3(:,:,1)
1333
1334	! Fields received by SAS when OASIS coupling
1335	! (arrays no more filled at sbcssm stage)
1336	! ! ================== !
1337	! ! SSS !
1338	! ! ================== !
1339	IF( srcv(jpr_soce)%laction ) THEN ! received by sas in case of opa <-> sas coupling
1340	sss_m(:,:) = frcv(jpr_soce)%z3(:,:,1)
1341	CALL iom_put( 'sss_m', sss_m )
1342	ENDIF
1343	!
1344	! ! ================== !
1345	! ! SST !
1346	! ! ================== !
1347	IF( srcv(jpr_toce)%laction ) THEN ! received by sas in case of opa <-> sas coupling
1348	sst_m(:,:) = frcv(jpr_toce)%z3(:,:,1)
1349	IF( srcv(jpr_soce)%laction .AND. l_useCT ) THEN ! make sure that sst_m is the potential temperature
1350	sst_m(:,:) = eos_pt_from_ct( sst_m(:,:), sss_m(:,:) )
1351	ENDIF
1352	ENDIF
1353	! ! ================== !
1354	! ! SSH !
1355	! ! ================== !
1356	IF( srcv(jpr_ssh )%laction ) THEN ! received by sas in case of opa <-> sas coupling
1357	ssh_m(:,:) = frcv(jpr_ssh )%z3(:,:,1)
1358	CALL iom_put( 'ssh_m', ssh_m )
1359	ENDIF
1360	! ! ================== !
1361	! ! surface currents !
1362	! ! ================== !
1363	IF( srcv(jpr_ocx1)%laction ) THEN ! received by sas in case of opa <-> sas coupling
1364	ssu_m(:,:) = frcv(jpr_ocx1)%z3(:,:,1)
1365	ub (:,:,1) = ssu_m(:,:) ! will be used in icestp in the call of lim_sbc_tau
1366	un (:,:,1) = ssu_m(:,:) ! will be used in sbc_cpl_snd if atmosphere coupling
1367	CALL iom_put( 'ssu_m', ssu_m )
1368	ENDIF
1369	IF( srcv(jpr_ocy1)%laction ) THEN
1370	ssv_m(:,:) = frcv(jpr_ocy1)%z3(:,:,1)
1371	vb (:,:,1) = ssv_m(:,:) ! will be used in icestp in the call of lim_sbc_tau
1372	vn (:,:,1) = ssv_m(:,:) ! will be used in sbc_cpl_snd if atmosphere coupling
1373	CALL iom_put( 'ssv_m', ssv_m )
1374	ENDIF
1375	! ! ======================== !
1376	! ! first T level thickness !
1377	! ! ======================== !
1378	IF( srcv(jpr_e3t1st )%laction ) THEN ! received by sas in case of opa <-> sas coupling
1379	e3t_m(:,:) = frcv(jpr_e3t1st )%z3(:,:,1)
1380	CALL iom_put( 'e3t_m', e3t_m(:,:) )
1381	ENDIF
1382	! ! ================================ !
1383	! ! fraction of solar net radiation !
1384	! ! ================================ !
1385	IF( srcv(jpr_fraqsr)%laction ) THEN ! received by sas in case of opa <-> sas coupling
1386	frq_m(:,:) = frcv(jpr_fraqsr)%z3(:,:,1)
1387	CALL iom_put( 'frq_m', frq_m )
1388	ENDIF
1389
1390	! ! ========================= !
1391	IF( k_ice <= 1 .AND. MOD( kt-1, k_fsbc ) == 0 ) THEN ! heat & freshwater fluxes ! (Ocean only case)
1392	! ! ========================= !
1393	!
1394	! ! total freshwater fluxes over the ocean (emp)
1395	IF( srcv(jpr_oemp)%laction .OR. srcv(jpr_rain)%laction ) THEN
1396	SELECT CASE( TRIM( sn_rcv_emp%cldes ) ) ! evaporation - precipitation
1397	CASE( 'conservative' )
1398	zemp(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ( frcv(jpr_rain)%z3(:,:,1) + frcv(jpr_snow)%z3(:,:,1) )
1399	CASE( 'oce only', 'oce and ice' )
1400	zemp(:,:) = frcv(jpr_oemp)%z3(:,:,1)
1401	CASE default
1402	CALL ctl_stop( 'sbc_cpl_rcv: wrong definition of sn_rcv_emp%cldes' )
1403	END SELECT
1404	ELSE
1405	zemp(:,:) = 0._wp
1406	ENDIF
1407	!
1408	! ! runoffs and calving (added in emp)
1409	IF( srcv(jpr_rnf)%laction ) rnf(:,:) = frcv(jpr_rnf)%z3(:,:,1)
1410	IF( srcv(jpr_cal)%laction ) zemp(:,:) = zemp(:,:) - frcv(jpr_cal)%z3(:,:,1)
1411
1412	IF( srcv(jpr_icb)%laction ) THEN
1413	fwficb(:,:) = frcv(jpr_icb)%z3(:,:,1)
1414	rnf(:,:) = rnf(:,:) + fwficb(:,:) ! iceberg added to runfofs
1415	ENDIF
1416	IF( srcv(jpr_isf)%laction ) fwfisf(:,:) = - frcv(jpr_isf)%z3(:,:,1) ! fresh water flux from the isf (fwfisf <0 mean melting)
1417
1418	IF( ln_mixcpl ) THEN ; emp(:,:) = emp(:,:) * xcplmask(:,:,0) + zemp(:,:) * zmsk(:,:)
1419	ELSE ; emp(:,:) = zemp(:,:)
1420	ENDIF
1421	!
1422	! ! non solar heat flux over the ocean (qns)
1423	IF( srcv(jpr_qnsoce)%laction ) THEN ; zqns(:,:) = frcv(jpr_qnsoce)%z3(:,:,1)
1424	ELSE IF( srcv(jpr_qnsmix)%laction ) THEN ; zqns(:,:) = frcv(jpr_qnsmix)%z3(:,:,1)
1425	ELSE ; zqns(:,:) = 0._wp
1426	END IF
1427	! update qns over the free ocean with:
1428	IF( nn_components /= jp_iam_opa ) THEN
1429	zqns(:,:) = zqns(:,:) - zemp(:,:) * sst_m(:,:) * rcp ! remove heat content due to mass flux (assumed to be at SST)
1430	IF( srcv(jpr_snow )%laction ) THEN
1431	zqns(:,:) = zqns(:,:) - frcv(jpr_snow)%z3(:,:,1) * lfus ! energy for melting solid precipitation over the free ocean
1432	ENDIF
1433	ENDIF
1434	!
1435	IF( srcv(jpr_icb)%laction ) zqns(:,:) = zqns(:,:) - frcv(jpr_icb)%z3(:,:,1) * lfus ! remove heat content associated to iceberg melting
1436	!
1437	IF( ln_mixcpl ) THEN ; qns(:,:) = qns(:,:) * xcplmask(:,:,0) + zqns(:,:) * zmsk(:,:)
1438	ELSE ; qns(:,:) = zqns(:,:)
1439	ENDIF
1440
1441	! ! solar flux over the ocean (qsr)
1442	IF ( srcv(jpr_qsroce)%laction ) THEN ; zqsr(:,:) = frcv(jpr_qsroce)%z3(:,:,1)
1443	ELSE IF( srcv(jpr_qsrmix)%laction ) then ; zqsr(:,:) = frcv(jpr_qsrmix)%z3(:,:,1)
1444	ELSE ; zqsr(:,:) = 0._wp
1445	ENDIF
1446	IF( ln_dm2dc .AND. ln_cpl ) zqsr(:,:) = sbc_dcy( zqsr ) ! modify qsr to include the diurnal cycle
1447	IF( ln_mixcpl ) THEN ; qsr(:,:) = qsr(:,:) * xcplmask(:,:,0) + zqsr(:,:) * zmsk(:,:)
1448	ELSE ; qsr(:,:) = zqsr(:,:)
1449	ENDIF
1450	!
1451	! salt flux over the ocean (received by opa in case of opa <-> sas coupling)
1452	IF( srcv(jpr_sflx )%laction ) sfx(:,:) = frcv(jpr_sflx )%z3(:,:,1)
1453	! Ice cover (received by opa in case of opa <-> sas coupling)
1454	IF( srcv(jpr_fice )%laction ) fr_i(:,:) = frcv(jpr_fice )%z3(:,:,1)
1455	!
1456	ENDIF
1457	!
1458	CALL wrk_dealloc( jpi,jpj, ztx, zty, zmsk, zemp, zqns, zqsr )
1459	!
1460	IF( nn_timing == 1 ) CALL timing_stop('sbc_cpl_rcv')
1461	!
1462	END SUBROUTINE sbc_cpl_rcv
1463
1464
1465	SUBROUTINE sbc_cpl_ice_tau( p_taui, p_tauj )
1466	!!----------------------------------------------------------------------
1467	!! * ROUTINE sbc_cpl_ice_tau *
1468	!!
1469	!! ** Purpose : provide the stress over sea-ice in coupled mode
1470	!!
1471	!! ** Method : transform the received stress from the atmosphere into
1472	!! an atmosphere-ice stress in the (i,j) ocean referencial
1473	!! and at the velocity point of the sea-ice model (cp_ice_msh):
1474	!! 'C'-grid : i- (j-) components given at U- (V-) point
1475	!! 'I'-grid : B-grid lower-left corner: both components given at I-point
1476	!!
1477	!! The received stress are :
1478	!! - defined by 3 components (if cartesian coordinate)
1479	!! or by 2 components (if spherical)
1480	!! - oriented along geographical coordinate (if eastward-northward)
1481	!! or along the local grid coordinate (if local grid)
1482	!! - given at U- and V-point, resp. if received on 2 grids
1483	!! or at a same point (T or I) if received on 1 grid
1484	!! Therefore and if necessary, they are successively
1485	!! processed in order to obtain them
1486	!! first as 2 components on the sphere
1487	!! second as 2 components oriented along the local grid
1488	!! third as 2 components on the cp_ice_msh point
1489	!!
1490	!! Except in 'oce and ice' case, only one vector stress field
1491	!! is received. It has already been processed in sbc_cpl_rcv
1492	!! so that it is now defined as (i,j) components given at U-
1493	!! and V-points, respectively. Therefore, only the third
1494	!! transformation is done and only if the ice-grid is a 'I'-grid.
1495	!!
1496	!! ** Action : return ptau_i, ptau_j, the stress over the ice at cp_ice_msh point
1497	!!----------------------------------------------------------------------
1498	REAL(wp), INTENT(out), DIMENSION(:,:) :: p_taui ! i- & j-components of atmos-ice stress [N/m2]
1499	REAL(wp), INTENT(out), DIMENSION(:,:) :: p_tauj ! at I-point (B-grid) or U & V-point (C-grid)
1500	!!
1501	INTEGER :: ji, jj ! dummy loop indices
1502	INTEGER :: itx ! index of taux over ice
1503	REAL(wp), POINTER, DIMENSION(:,:) :: ztx, zty
1504	!!----------------------------------------------------------------------
1505	!
1506	IF( nn_timing == 1 ) CALL timing_start('sbc_cpl_ice_tau')
1507	!
1508	CALL wrk_alloc( jpi,jpj, ztx, zty )
1509
1510	IF( srcv(jpr_itx1)%laction ) THEN ; itx = jpr_itx1
1511	ELSE ; itx = jpr_otx1
1512	ENDIF
1513
1514	! do something only if we just received the stress from atmosphere
1515	IF( nrcvinfo(itx) == OASIS_Rcv ) THEN
1516	! ! ======================= !
1517	IF( srcv(jpr_itx1)%laction ) THEN ! ice stress received !
1518	! ! ======================= !
1519	!
1520	IF( TRIM( sn_rcv_tau%clvref ) == 'cartesian' ) THEN ! 2 components on the sphere
1521	! ! (cartesian to spherical -> 3 to 2 components)
1522	CALL geo2oce( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), frcv(jpr_itz1)%z3(:,:,1), &
1523	& srcv(jpr_itx1)%clgrid, ztx, zty )
1524	frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 1st grid
1525	frcv(jpr_ity1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 1st grid
1526	!
1527	IF( srcv(jpr_itx2)%laction ) THEN
1528	CALL geo2oce( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), frcv(jpr_itz2)%z3(:,:,1), &
1529	& srcv(jpr_itx2)%clgrid, ztx, zty )
1530	frcv(jpr_itx2)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 2nd grid
1531	frcv(jpr_ity2)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 2nd grid
1532	ENDIF
1533	!
1534	ENDIF
1535	!
1536	IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid
1537	! ! (geographical to local grid -> rotate the components)
1538	CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->i', ztx )
1539	IF( srcv(jpr_itx2)%laction ) THEN
1540	CALL rot_rep( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), srcv(jpr_itx2)%clgrid, 'en->j', zty )
1541	ELSE
1542	CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->j', zty )
1543	ENDIF
1544	frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st component on the 1st grid
1545	frcv(jpr_ity1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd component on the 1st grid
1546	ENDIF
1547	! ! ======================= !
1548	ELSE ! use ocean stress !
1549	! ! ======================= !
1550	frcv(jpr_itx1)%z3(:,:,1) = frcv(jpr_otx1)%z3(:,:,1)
1551	frcv(jpr_ity1)%z3(:,:,1) = frcv(jpr_oty1)%z3(:,:,1)
1552	!
1553	ENDIF
1554	! ! ======================= !
1555	! ! put on ice grid !
1556	! ! ======================= !
1557	!
1558	! j+1 j -----V---F
1559	! ice stress on ice velocity point (cp_ice_msh) ! \|
1560	! (C-grid ==>(U,V) or B-grid ==> I or F) j \| T U
1561	! \| \|
1562	! j j-1 -I-------\|
1563	! (for I) \| \|
1564	! i-1 i i
1565	! i i+1 (for I)
1566	SELECT CASE ( cp_ice_msh )
1567	!
1568	CASE( 'I' ) ! B-grid ==> I
1569	SELECT CASE ( srcv(jpr_itx1)%clgrid )
1570	CASE( 'U' )
1571	DO jj = 2, jpjm1 ! (U,V) ==> I
1572	DO ji = 2, jpim1 ! NO vector opt.
1573	p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji-1,jj ,1) + frcv(jpr_itx1)%z3(ji-1,jj-1,1) )
1574	p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji ,jj-1,1) + frcv(jpr_ity1)%z3(ji-1,jj-1,1) )
1575	END DO
1576	END DO
1577	CASE( 'F' )
1578	DO jj = 2, jpjm1 ! F ==> I
1579	DO ji = 2, jpim1 ! NO vector opt.
1580	p_taui(ji,jj) = frcv(jpr_itx1)%z3(ji-1,jj-1,1)
1581	p_tauj(ji,jj) = frcv(jpr_ity1)%z3(ji-1,jj-1,1)
1582	END DO
1583	END DO
1584	CASE( 'T' )
1585	DO jj = 2, jpjm1 ! T ==> I
1586	DO ji = 2, jpim1 ! NO vector opt.
1587	p_taui(ji,jj) = 0.25 * ( frcv(jpr_itx1)%z3(ji,jj ,1) + frcv(jpr_itx1)%z3(ji-1,jj ,1) &
1588	& + frcv(jpr_itx1)%z3(ji,jj-1,1) + frcv(jpr_itx1)%z3(ji-1,jj-1,1) )
1589	p_tauj(ji,jj) = 0.25 * ( frcv(jpr_ity1)%z3(ji,jj ,1) + frcv(jpr_ity1)%z3(ji-1,jj ,1) &
1590	& + frcv(jpr_oty1)%z3(ji,jj-1,1) + frcv(jpr_ity1)%z3(ji-1,jj-1,1) )
1591	END DO
1592	END DO
1593	CASE( 'I' )
1594	p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1) ! I ==> I
1595	p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1)
1596	END SELECT
1597	IF( srcv(jpr_itx1)%clgrid /= 'I' ) THEN
1598	CALL lbc_lnk( p_taui, 'I', -1. ) ; CALL lbc_lnk( p_tauj, 'I', -1. )
1599	ENDIF
1600	!
1601	CASE( 'F' ) ! B-grid ==> F
1602	SELECT CASE ( srcv(jpr_itx1)%clgrid )
1603	CASE( 'U' )
1604	DO jj = 2, jpjm1 ! (U,V) ==> F
1605	DO ji = fs_2, fs_jpim1 ! vector opt.
1606	p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji,jj,1) + frcv(jpr_itx1)%z3(ji ,jj+1,1) )
1607	p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji,jj,1) + frcv(jpr_ity1)%z3(ji+1,jj ,1) )
1608	END DO
1609	END DO
1610	CASE( 'I' )
1611	DO jj = 2, jpjm1 ! I ==> F
1612	DO ji = 2, jpim1 ! NO vector opt.
1613	p_taui(ji,jj) = frcv(jpr_itx1)%z3(ji+1,jj+1,1)
1614	p_tauj(ji,jj) = frcv(jpr_ity1)%z3(ji+1,jj+1,1)
1615	END DO
1616	END DO
1617	CASE( 'T' )
1618	DO jj = 2, jpjm1 ! T ==> F
1619	DO ji = 2, jpim1 ! NO vector opt.
1620	p_taui(ji,jj) = 0.25 * ( frcv(jpr_itx1)%z3(ji,jj ,1) + frcv(jpr_itx1)%z3(ji+1,jj ,1) &
1621	& + frcv(jpr_itx1)%z3(ji,jj+1,1) + frcv(jpr_itx1)%z3(ji+1,jj+1,1) )
1622	p_tauj(ji,jj) = 0.25 * ( frcv(jpr_ity1)%z3(ji,jj ,1) + frcv(jpr_ity1)%z3(ji+1,jj ,1) &
1623	& + frcv(jpr_ity1)%z3(ji,jj+1,1) + frcv(jpr_ity1)%z3(ji+1,jj+1,1) )
1624	END DO
1625	END DO
1626	CASE( 'F' )
1627	p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1) ! F ==> F
1628	p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1)
1629	END SELECT
1630	IF( srcv(jpr_itx1)%clgrid /= 'F' ) THEN
1631	CALL lbc_lnk( p_taui, 'F', -1. ) ; CALL lbc_lnk( p_tauj, 'F', -1. )
1632	ENDIF
1633	!
1634	CASE( 'C' ) ! C-grid ==> U,V
1635	SELECT CASE ( srcv(jpr_itx1)%clgrid )
1636	CASE( 'U' )
1637	p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1) ! (U,V) ==> (U,V)
1638	p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1)
1639	CASE( 'F' )
1640	DO jj = 2, jpjm1 ! F ==> (U,V)
1641	DO ji = fs_2, fs_jpim1 ! vector opt.
1642	p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji,jj,1) + frcv(jpr_itx1)%z3(ji ,jj-1,1) )
1643	p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(jj,jj,1) + frcv(jpr_ity1)%z3(ji-1,jj ,1) )
1644	END DO
1645	END DO
1646	CASE( 'T' )
1647	DO jj = 2, jpjm1 ! T ==> (U,V)
1648	DO ji = fs_2, fs_jpim1 ! vector opt.
1649	p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji+1,jj ,1) + frcv(jpr_itx1)%z3(ji,jj,1) )
1650	p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji ,jj+1,1) + frcv(jpr_ity1)%z3(ji,jj,1) )
1651	END DO
1652	END DO
1653	CASE( 'I' )
1654	DO jj = 2, jpjm1 ! I ==> (U,V)
1655	DO ji = 2, jpim1 ! NO vector opt.
1656	p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji+1,jj+1,1) + frcv(jpr_itx1)%z3(ji+1,jj ,1) )
1657	p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji+1,jj+1,1) + frcv(jpr_ity1)%z3(ji ,jj+1,1) )
1658	END DO
1659	END DO
1660	END SELECT
1661	IF( srcv(jpr_itx1)%clgrid /= 'U' ) THEN
1662	CALL lbc_lnk( p_taui, 'U', -1. ) ; CALL lbc_lnk( p_tauj, 'V', -1. )
1663	ENDIF
1664	END SELECT
1665
1666	ENDIF
1667	!
1668	CALL wrk_dealloc( jpi,jpj, ztx, zty )
1669	!
1670	IF( nn_timing == 1 ) CALL timing_stop('sbc_cpl_ice_tau')
1671	!
1672	END SUBROUTINE sbc_cpl_ice_tau
1673
1674
1675	SUBROUTINE sbc_cpl_ice_flx( picefr, palbi, psst, pist, phs, phi )
1676	!!----------------------------------------------------------------------
1677	!! * ROUTINE sbc_cpl_ice_flx *
1678	!!
1679	!! ** Purpose : provide the heat and freshwater fluxes of the ocean-ice system
1680	!!
1681	!! ** Method : transform the fields received from the atmosphere into
1682	!! surface heat and fresh water boundary condition for the
1683	!! ice-ocean system. The following fields are provided:
1684	!! * total non solar, solar and freshwater fluxes (qns_tot,
1685	!! qsr_tot and emp_tot) (total means weighted ice-ocean flux)
1686	!! NB: emp_tot include runoffs and calving.
1687	!! * fluxes over ice (qns_ice, qsr_ice, emp_ice) where
1688	!! emp_ice = sublimation - solid precipitation as liquid
1689	!! precipitation are re-routed directly to the ocean and
1690	!! calving directly enter the ocean (runoffs are read but included in trasbc.F90)
1691	!! * solid precipitation (sprecip), used to add to qns_tot
1692	!! the heat lost associated to melting solid precipitation
1693	!! over the ocean fraction.
1694	!! * heat content of rain, snow and evap can also be provided,
1695	!! otherwise heat flux associated with these mass flux are
1696	!! guessed (qemp_oce, qemp_ice)
1697	!!
1698	!! - the fluxes have been separated from the stress as
1699	!! (a) they are updated at each ice time step compare to
1700	!! an update at each coupled time step for the stress, and
1701	!! (b) the conservative computation of the fluxes over the
1702	!! sea-ice area requires the knowledge of the ice fraction
1703	!! after the ice advection and before the ice thermodynamics,
1704	!! so that the stress is updated before the ice dynamics
1705	!! while the fluxes are updated after it.
1706	!!
1707	!! ** Details
1708	!! qns_tot = (1-a) * qns_oce + a * qns_ice => provided
1709	!! + qemp_oce + qemp_ice => recalculated and added up to qns
1710	!!
1711	!! qsr_tot = (1-a) * qsr_oce + a * qsr_ice => provided
1712	!!
1713	!! emp_tot = emp_oce + emp_ice => calving is provided and added to emp_tot (and emp_oce).
1714	!! runoff (which includes rivers+icebergs) and iceshelf
1715	!! are provided but not included in emp here. Only runoff will
1716	!! be included in emp in other parts of NEMO code
1717	!! ** Action : update at each nf_ice time step:
1718	!! qns_tot, qsr_tot non-solar and solar total heat fluxes
1719	!! qns_ice, qsr_ice non-solar and solar heat fluxes over the ice
1720	!! emp_tot total evaporation - precipitation(liquid and solid) (-calving)
1721	!! emp_ice ice sublimation - solid precipitation over the ice
1722	!! dqns_ice d(non-solar heat flux)/d(Temperature) over the ice
1723	!! sprecip solid precipitation over the ocean
1724	!!----------------------------------------------------------------------
1725	REAL(wp), INTENT(in), DIMENSION(:,:) :: picefr ! ice fraction [0 to 1]
1726	! optional arguments, used only in 'mixed oce-ice' case
1727	REAL(wp), INTENT(in), DIMENSION(:,:,:), OPTIONAL :: palbi ! all skies ice albedo
1728	REAL(wp), INTENT(in), DIMENSION(:,: ), OPTIONAL :: psst ! sea surface temperature [Celsius]
1729	REAL(wp), INTENT(in), DIMENSION(:,:,:), OPTIONAL :: pist ! ice surface temperature [Kelvin]
1730	REAL(wp), INTENT(in), DIMENSION(:,:,:), OPTIONAL :: phs ! snow depth [m]
1731	REAL(wp), INTENT(in), DIMENSION(:,:,:), OPTIONAL :: phi ! ice thickness [m]
1732	!
1733	INTEGER :: ji,jj,jl ! dummy loop index
1734	REAL(wp), POINTER, DIMENSION(:,: ) :: zcptn, zcptrain, zcptsnw, ziceld, zmsk, zsnw
1735	REAL(wp), POINTER, DIMENSION(:,: ) :: zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip, zevap_oce, zdevap_ice
1736	REAL(wp), POINTER, DIMENSION(:,: ) :: zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice
1737	REAL(wp), POINTER, DIMENSION(:,:,:) :: zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zfrqsr_tr_i, zevap_ice
1738	!!----------------------------------------------------------------------
1739	!
1740	IF( nn_timing == 1 ) CALL timing_start('sbc_cpl_ice_flx')
1741	!
1742	CALL wrk_alloc( jpi,jpj, zcptn, zcptrain, zcptsnw, ziceld, zmsk, zsnw )
1743	CALL wrk_alloc( jpi,jpj, zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip, zevap_oce, zdevap_ice )
1744	CALL wrk_alloc( jpi,jpj, zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice )
1745	CALL wrk_alloc( jpi,jpj,jpl, zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zfrqsr_tr_i )
1746
1747	IF (sn_rcv_emp%clcat == 'yes') THEN
1748	CALL wrk_alloc( jpi, jpj, jpl, zevap_ice)
1749	ELSE
1750	CALL wrk_alloc( jpi, jpj, 1, zevap_ice)
1751	ENDIF
1752
1753	IF( ln_mixcpl ) zmsk(:,:) = 1. - xcplmask(:,:,0)
1754	ziceld(:,:) = 1. - picefr(:,:)
1755	zcptn(:,:) = rcp * sst_m(:,:)
1756	!
1757	! ! ========================= !
1758	! ! freshwater budget ! (emp_tot)
1759	! ! ========================= !
1760	!
1761	! ! solid Precipitation (sprecip)
1762	! ! liquid + solid Precipitation (tprecip)
1763	! ! total Evaporation - total Precipitation (emp_tot)
1764	! ! sublimation - solid precipitation (cell average) (emp_ice)
1765	SELECT CASE( TRIM( sn_rcv_emp%cldes ) )
1766	CASE( 'conservative' ) ! received fields: jpr_rain, jpr_snow, jpr_ievp, jpr_tevp
1767	zsprecip(:,:) = frcv(jpr_snow)%z3(:,:,1) ! May need to ensure positive here
1768	ztprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:) ! May need to ensure positive here
1769	zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:)
1770	zemp_ice(:,:) = ( frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) ) * picefr(:,:)
1771	CASE( 'oce and ice' ) ! received fields: jpr_sbpr, jpr_semp, jpr_oemp, jpr_ievp
1772	zemp_tot(:,:) = ziceld(:,:) * frcv(jpr_oemp)%z3(:,:,1) + picefr(:,:) * frcv(jpr_sbpr)%z3(:,:,1)
1773	zemp_ice(:,:) = frcv(jpr_semp)%z3(:,:,1) * picefr(:,:)
1774	zsprecip(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_semp)%z3(:,:,1)
1775	ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:)
1776	END SELECT
1777
1778	#if defined key_lim3
1779	! zsnw = snow fraction over ice after wind blowing (=picefr if no blowing)
1780	zsnw(:,:) = 0._wp ; CALL ice_thd_snwblow( ziceld, zsnw )
1781
1782	! --- evaporation minus precipitation corrected (because of wind blowing on snow) --- !
1783	zemp_ice(:,:) = zemp_ice(:,:) + zsprecip(:,:) * ( picefr(:,:) - zsnw(:,:) ) ! emp_ice = A * sublimation - zsnw * sprecip
1784	zemp_oce(:,:) = zemp_tot(:,:) - zemp_ice(:,:) ! emp_oce = emp_tot - emp_ice
1785
1786	! --- evaporation over ocean (used later for qemp) --- !
1787	zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:)
1788
1789	! --- evaporation over ice (kg/m2/s) --- !
1790	IF (sn_rcv_emp%clcat == 'yes') THEN
1791	DO jl=1,jpl
1792	zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl)
1793	ENDDO
1794	ELSE
1795	zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1)
1796	ENDIF
1797
1798	! since the sensitivity of evap to temperature (devap/dT) is not prescribed by the atmosphere, we set it to 0
1799	! therefore, sublimation is not redistributed over the ice categories when no subgrid scale fluxes are provided by atm.
1800	zdevap_ice(:,:) = 0._wp
1801
1802	! --- Continental fluxes --- !
1803	IF( srcv(jpr_rnf)%laction ) THEN ! runoffs (included in emp later on)
1804	rnf(:,:) = frcv(jpr_rnf)%z3(:,:,1)
1805	ENDIF
1806	IF( srcv(jpr_cal)%laction ) THEN ! calving (put in emp_tot and emp_oce)
1807	zemp_tot(:,:) = zemp_tot(:,:) - frcv(jpr_cal)%z3(:,:,1)
1808	zemp_oce(:,:) = zemp_oce(:,:) - frcv(jpr_cal)%z3(:,:,1)
1809	ENDIF
1810	IF( srcv(jpr_icb)%laction ) THEN ! iceberg added to runoffs
1811	fwficb(:,:) = frcv(jpr_icb)%z3(:,:,1)
1812	rnf(:,:) = rnf(:,:) + fwficb(:,:)
1813	ENDIF
1814	IF( srcv(jpr_isf)%laction ) THEN ! iceshelf (fwfisf <0 mean melting)
1815	fwfisf(:,:) = - frcv(jpr_isf)%z3(:,:,1)
1816	ENDIF
1817
1818	IF( ln_mixcpl ) THEN
1819	emp_tot(:,:) = emp_tot(:,:) * xcplmask(:,:,0) + zemp_tot(:,:) * zmsk(:,:)
1820	emp_ice(:,:) = emp_ice(:,:) * xcplmask(:,:,0) + zemp_ice(:,:) * zmsk(:,:)
1821	emp_oce(:,:) = emp_oce(:,:) * xcplmask(:,:,0) + zemp_oce(:,:) * zmsk(:,:)
1822	sprecip(:,:) = sprecip(:,:) * xcplmask(:,:,0) + zsprecip(:,:) * zmsk(:,:)
1823	tprecip(:,:) = tprecip(:,:) * xcplmask(:,:,0) + ztprecip(:,:) * zmsk(:,:)
1824	DO jl=1,jpl
1825	evap_ice (:,:,jl) = evap_ice (:,:,jl) * xcplmask(:,:,0) + zevap_ice (:,:,1) * zmsk(:,:)
1826	devap_ice(:,:,jl) = devap_ice(:,:,jl) * xcplmask(:,:,0) + zdevap_ice(:,:) * zmsk(:,:)
1827	ENDDO
1828	ELSE
1829	emp_tot(:,:) = zemp_tot(:,:)
1830	emp_ice(:,:) = zemp_ice(:,:)
1831	emp_oce(:,:) = zemp_oce(:,:)
1832	sprecip(:,:) = zsprecip(:,:)
1833	tprecip(:,:) = ztprecip(:,:)
1834	DO jl=1,jpl
1835	IF (sn_rcv_emp%clcat == 'yes') THEN
1836	evap_ice (:,:,jl) = zevap_ice (:,:,jl)
1837	ELSE
1838	evap_ice (:,:,jl) = zevap_ice (:,:,1)
1839	ENDIF
1840
1841	!!!!! Getting NEMO4-LIM to work at the Met Office: Semi-implicit coupling
1842	evap_ice (:,:,:) = evap_ice (:,:,:) * a_i(:,:,:)
1843	!!!!!
1844
1845	devap_ice(:,:,jl) = zdevap_ice(:,:)
1846	ENDDO
1847	ENDIF
1848
1849	#else
1850	zsnw(:,:) = picefr(:,:)
1851	! --- Continental fluxes --- !
1852	IF( srcv(jpr_rnf)%laction ) THEN ! runoffs (included in emp later on)
1853	rnf(:,:) = frcv(jpr_rnf)%z3(:,:,1)
1854	ENDIF
1855	IF( srcv(jpr_cal)%laction ) THEN ! calving (put in emp_tot)
1856	zemp_tot(:,:) = zemp_tot(:,:) - frcv(jpr_cal)%z3(:,:,1)
1857	ENDIF
1858	IF( srcv(jpr_icb)%laction ) THEN ! iceberg added to runoffs
1859	fwficb(:,:) = frcv(jpr_icb)%z3(:,:,1)
1860	rnf(:,:) = rnf(:,:) + fwficb(:,:)
1861	ENDIF
1862	IF( srcv(jpr_isf)%laction ) THEN ! iceshelf (fwfisf <0 mean melting)
1863	fwfisf(:,:) = - frcv(jpr_isf)%z3(:,:,1)
1864	ENDIF
1865
1866	IF( ln_mixcpl ) THEN
1867	emp_tot(:,:) = emp_tot(:,:) * xcplmask(:,:,0) + zemp_tot(:,:) * zmsk(:,:)
1868	emp_ice(:,:) = emp_ice(:,:) * xcplmask(:,:,0) + zemp_ice(:,:) * zmsk(:,:)
1869	sprecip(:,:) = sprecip(:,:) * xcplmask(:,:,0) + zsprecip(:,:) * zmsk(:,:)
1870	tprecip(:,:) = tprecip(:,:) * xcplmask(:,:,0) + ztprecip(:,:) * zmsk(:,:)
1871	ELSE
1872	emp_tot(:,:) = zemp_tot(:,:)
1873	emp_ice(:,:) = zemp_ice(:,:)
1874	sprecip(:,:) = zsprecip(:,:)
1875	tprecip(:,:) = ztprecip(:,:)
1876	ENDIF
1877
1878	#endif
1879	! outputs
1880	!! IF( srcv(jpr_rnf)%laction ) CALL iom_put( 'runoffs' , rnf(:,:) * tmask(:,:,1) ) ! runoff
1881	!! IF( srcv(jpr_isf)%laction ) CALL iom_put( 'iceshelf_cea', -fwfisf(:,:) * tmask(:,:,1) ) ! iceshelf
1882	IF( srcv(jpr_cal)%laction ) CALL iom_put( 'calving_cea' , frcv(jpr_cal)%z3(:,:,1) * tmask(:,:,1) ) ! calving
1883	IF( srcv(jpr_icb)%laction ) CALL iom_put( 'iceberg_cea' , frcv(jpr_icb)%z3(:,:,1) * tmask(:,:,1) ) ! icebergs
1884	IF( iom_use('snowpre') ) CALL iom_put( 'snowpre' , sprecip(:,:) ) ! Snow
1885	IF( iom_use('precip') ) CALL iom_put( 'precip' , tprecip(:,:) ) ! total precipitation
1886	IF( iom_use('rain') ) CALL iom_put( 'rain' , tprecip(:,:) - sprecip(:,:) ) ! liquid precipitation
1887	IF( iom_use('snow_ao_cea') ) CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) ) ) ! Snow over ice-free ocean (cell average)
1888	IF( iom_use('snow_ai_cea') ) CALL iom_put( 'snow_ai_cea' , sprecip(:,:) * zsnw(:,:) ) ! Snow over sea-ice (cell average)
1889	IF( iom_use('subl_ai_cea') ) CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) ) ! Sublimation over sea-ice (cell average)
1890	IF( iom_use('evap_ao_cea') ) CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1) &
1891	& - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average)
1892	! note: runoff output is done in sbcrnf (which includes icebergs too) and iceshelf output is done in sbcisf
1893	!
1894	! ! ========================= !
1895	SELECT CASE( TRIM( sn_rcv_qns%cldes ) ) ! non solar heat fluxes ! (qns)
1896	! ! ========================= !
1897	CASE( 'oce only' ) ! the required field is directly provided
1898	zqns_tot(:,:) = frcv(jpr_qnsoce)%z3(:,:,1)
1899	CASE( 'conservative' ) ! the required fields are directly provided
1900	zqns_tot(:,:) = frcv(jpr_qnsmix)%z3(:,:,1)
1901	IF ( TRIM(sn_rcv_qns%clcat) == 'yes' ) THEN
1902	zqns_ice(:,:,1:jpl) = frcv(jpr_qnsice)%z3(:,:,1:jpl)
1903	ELSE
1904	DO jl=1,jpl
1905	zqns_ice(:,:,jl) = frcv(jpr_qnsice)%z3(:,:,1) ! Set all category values equal
1906	ENDDO
1907	ENDIF
1908	CASE( 'oce and ice' ) ! the total flux is computed from ocean and ice fluxes
1909	zqns_tot(:,:) = ziceld(:,:) * frcv(jpr_qnsoce)%z3(:,:,1)
1910	IF ( TRIM(sn_rcv_qns%clcat) == 'yes' ) THEN
1911	DO jl=1,jpl
1912	zqns_tot(:,: ) = zqns_tot(:,:) + a_i(:,:,jl) * frcv(jpr_qnsice)%z3(:,:,jl)
1913	zqns_ice(:,:,jl) = frcv(jpr_qnsice)%z3(:,:,jl)
1914	ENDDO
1915	ELSE
1916	qns_tot(:,:) = qns_tot(:,:) + picefr(:,:) * frcv(jpr_qnsice)%z3(:,:,1)
1917	DO jl=1,jpl
1918	zqns_tot(:,: ) = zqns_tot(:,:) + picefr(:,:) * frcv(jpr_qnsice)%z3(:,:,1)
1919	zqns_ice(:,:,jl) = frcv(jpr_qnsice)%z3(:,:,1)
1920	ENDDO
1921	ENDIF
1922	CASE( 'mixed oce-ice' ) ! the ice flux is cumputed from the total flux, the SST and ice informations
1923	! NEED TO SORT OUT HOW THIS SHOULD WORK IN THE MULTI-CATEGORY CASE - CURRENTLY NOT ALLOWED WHEN INTERFACE INITIALISED
1924	zqns_tot(:,: ) = frcv(jpr_qnsmix)%z3(:,:,1)
1925	zqns_ice(:,:,1) = frcv(jpr_qnsmix)%z3(:,:,1) &
1926	& + frcv(jpr_dqnsdt)%z3(:,:,1) * ( pist(:,:,1) - ( (rt0 + psst(:,: ) ) * ziceld(:,:) &
1927	& + pist(:,:,1) * picefr(:,:) ) )
1928	END SELECT
1929	!
1930	! --- calving (removed from qns_tot) --- !
1931	IF( srcv(jpr_cal)%laction ) zqns_tot(:,:) = zqns_tot(:,:) - frcv(jpr_cal)%z3(:,:,1) * lfus ! remove latent heat of calving
1932	! we suppose it melts at 0deg, though it should be temp. of surrounding ocean
1933	! --- iceberg (removed from qns_tot) --- !
1934	IF( srcv(jpr_icb)%laction ) zqns_tot(:,:) = zqns_tot(:,:) - frcv(jpr_icb)%z3(:,:,1) * lfus ! remove latent heat of iceberg melting
1935
1936	#if defined key_lim3
1937	! --- non solar flux over ocean --- !
1938	! note: ziceld cannot be = 0 since we limit the ice concentration to amax
1939	zqns_oce = 0._wp
1940	WHERE( ziceld /= 0._wp ) zqns_oce(:,:) = ( zqns_tot(:,:) - SUM( a_i * zqns_ice, dim=3 ) ) / ziceld(:,:)
1941
1942	! Heat content per unit mass of snow (J/kg)
1943	WHERE( SUM( a_i, dim=3 ) > 1.e-10 ) ; zcptsnw(:,:) = cpic * SUM( (tn_ice - rt0) * a_i, dim=3 ) / SUM( a_i, dim=3 )
1944	ELSEWHERE ; zcptsnw(:,:) = zcptn(:,:)
1945	ENDWHERE
1946	! Heat content per unit mass of rain (J/kg)
1947	zcptrain(:,:) = rcp * ( SUM( (tn_ice(:,:,:) - rt0) * a_i(:,:,:), dim=3 ) + sst_m(:,:) * ziceld(:,:) )
1948
1949	! --- enthalpy of snow precip over ice in J/m3 (to be used in 1D-thermo) --- !
1950	zqprec_ice(:,:) = rhosn * ( zcptsnw(:,:) - lfus )
1951
1952	! --- heat content of evap over ice in W/m2 (to be used in 1D-thermo) --- !
1953	DO jl = 1, jpl
1954	zqevap_ice(:,:,jl) = 0._wp ! should be -evap * ( ( Tice - rt0 ) * cpic ) but atm. does not take it into account
1955	END DO
1956
1957	! --- heat flux associated with emp (W/m2) --- !
1958	zqemp_oce(:,:) = - zevap_oce(:,:) * zcptn (:,:) & ! evap
1959	& + ( ztprecip(:,:) - zsprecip(:,:) ) * zcptrain(:,:) & ! liquid precip
1960	& + zsprecip(:,:) * ( 1._wp - zsnw ) * ( zcptsnw (:,:) - lfus ) ! solid precip over ocean + snow melting
1961	zqemp_ice(:,:) = zsprecip(:,:) * zsnw * ( zcptsnw (:,:) - lfus ) ! solid precip over ice (qevap_ice=0 since atm. does not take it into account)
1962	!! zqemp_ice(:,:) = - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * zcptsnw (:,:) & ! ice evap
1963	!! & + zsprecip(:,:) * zsnw * zqprec_ice(:,:) * r1_rhosn ! solid precip over ice
1964
1965	! --- total non solar flux (including evap/precip) --- !
1966	zqns_tot(:,:) = zqns_tot(:,:) + zqemp_ice(:,:) + zqemp_oce(:,:)
1967
1968	! --- in case both coupled/forced are active, we must mix values --- !
1969	IF( ln_mixcpl ) THEN
1970	qns_tot(:,:) = qns_tot(:,:) * xcplmask(:,:,0) + zqns_tot(:,:)* zmsk(:,:)
1971	qns_oce(:,:) = qns_oce(:,:) * xcplmask(:,:,0) + zqns_oce(:,:)* zmsk(:,:)
1972	DO jl=1,jpl
1973	qns_ice (:,:,jl) = qns_ice (:,:,jl) * xcplmask(:,:,0) + zqns_ice (:,:,jl)* zmsk(:,:)
1974	qevap_ice(:,:,jl) = qevap_ice(:,:,jl) * xcplmask(:,:,0) + zqevap_ice(:,:,jl)* zmsk(:,:)
1975	ENDDO
1976	qprec_ice(:,:) = qprec_ice(:,:) * xcplmask(:,:,0) + zqprec_ice(:,:)* zmsk(:,:)
1977	qemp_oce (:,:) = qemp_oce(:,:) * xcplmask(:,:,0) + zqemp_oce(:,:)* zmsk(:,:)
1978	qemp_ice (:,:) = qemp_ice(:,:) * xcplmask(:,:,0) + zqemp_ice(:,:)* zmsk(:,:)
1979	ELSE
1980	qns_tot (:,: ) = zqns_tot (:,: )
1981	qns_oce (:,: ) = zqns_oce (:,: )
1982	qns_ice (:,:,:) = zqns_ice (:,:,:)
1983	qevap_ice(:,:,:) = zqevap_ice(:,:,:)
1984	qprec_ice(:,: ) = zqprec_ice(:,: )
1985	qemp_oce (:,: ) = zqemp_oce (:,: )
1986	qemp_ice (:,: ) = zqemp_ice (:,: )
1987	ENDIF
1988
1989	#else
1990	zcptsnw (:,:) = zcptn(:,:)
1991	zcptrain(:,:) = zcptn(:,:)
1992
1993	! clem: this formulation is certainly wrong... but better than it was...
1994	zqns_tot(:,:) = zqns_tot(:,:) & ! zqns_tot update over free ocean with:
1995	& - ( ziceld(:,:) * zsprecip(:,:) * lfus ) & ! remove the latent heat flux of solid precip. melting
1996	& - ( zemp_tot(:,:) & ! remove the heat content of mass flux (assumed to be at SST)
1997	& - zemp_ice(:,:) ) * zcptn(:,:)
1998
1999	IF( ln_mixcpl ) THEN
2000	qns_tot(:,:) = qns(:,:) * ziceld(:,:) + SUM( qns_ice(:,:,:) * a_i(:,:,:), dim=3 ) ! total flux from blk
2001	qns_tot(:,:) = qns_tot(:,:) * xcplmask(:,:,0) + zqns_tot(:,:)* zmsk(:,:)
2002	DO jl=1,jpl
2003	qns_ice(:,:,jl) = qns_ice(:,:,jl) * xcplmask(:,:,0) + zqns_ice(:,:,jl)* zmsk(:,:)
2004	ENDDO
2005	ELSE
2006	qns_tot(:,: ) = zqns_tot(:,: )
2007	qns_ice(:,:,:) = zqns_ice(:,:,:)
2008	ENDIF
2009
2010	#endif
2011	! outputs
2012	IF( srcv(jpr_cal)%laction ) CALL iom_put('hflx_cal_cea' , - frcv(jpr_cal)%z3(:,:,1) * lfus ) ! latent heat from calving
2013	IF( srcv(jpr_icb)%laction ) CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * lfus ) ! latent heat from icebergs melting
2014	IF( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea', sprecip(:,:) * ( zcptsnw(:,:) - Lfus ) ) ! heat flux from snow (cell average)
2015	IF( iom_use('hflx_rain_cea') ) CALL iom_put('hflx_rain_cea',( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) ) ! heat flux from rain (cell average)
2016	IF( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea',(frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) & ! heat flux from from evap (cell average)
2017	& ) * zcptn(:,:) * tmask(:,:,1) )
2018	IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea',sprecip(:,:) * (zcptsnw(:,:) - Lfus) * (1._wp - zsnw(:,:)) ) ! heat flux from snow (over ocean)
2019	IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea',sprecip(:,:) * (zcptsnw(:,:) - Lfus) * zsnw(:,:) ) ! heat flux from snow (over ice)
2020	! note: hflx for runoff and iceshelf are done in sbcrnf and sbcisf resp.
2021	!
2022	! ! ========================= !
2023	SELECT CASE( TRIM( sn_rcv_qsr%cldes ) ) ! solar heat fluxes ! (qsr)
2024	! ! ========================= !
2025	CASE( 'oce only' )
2026	zqsr_tot(:,: ) = MAX( 0._wp , frcv(jpr_qsroce)%z3(:,:,1) )
2027	CASE( 'conservative' )
2028	zqsr_tot(:,: ) = frcv(jpr_qsrmix)%z3(:,:,1)
2029	IF ( TRIM(sn_rcv_qsr%clcat) == 'yes' ) THEN
2030	zqsr_ice(:,:,1:jpl) = frcv(jpr_qsrice)%z3(:,:,1:jpl)
2031	ELSE
2032	! Set all category values equal for the moment
2033	DO jl=1,jpl
2034	zqsr_ice(:,:,jl) = frcv(jpr_qsrice)%z3(:,:,1)
2035	ENDDO
2036	ENDIF
2037	zqsr_tot(:,: ) = frcv(jpr_qsrmix)%z3(:,:,1)
2038	zqsr_ice(:,:,1) = frcv(jpr_qsrice)%z3(:,:,1)
2039	CASE( 'oce and ice' )
2040	zqsr_tot(:,: ) = ziceld(:,:) * frcv(jpr_qsroce)%z3(:,:,1)
2041	IF ( TRIM(sn_rcv_qsr%clcat) == 'yes' ) THEN
2042	DO jl=1,jpl
2043	zqsr_tot(:,: ) = zqsr_tot(:,:) + a_i(:,:,jl) * frcv(jpr_qsrice)%z3(:,:,jl)
2044	zqsr_ice(:,:,jl) = frcv(jpr_qsrice)%z3(:,:,jl)
2045	ENDDO
2046	ELSE
2047	qsr_tot(:,: ) = qsr_tot(:,:) + picefr(:,:) * frcv(jpr_qsrice)%z3(:,:,1)
2048	DO jl=1,jpl
2049	zqsr_tot(:,: ) = zqsr_tot(:,:) + picefr(:,:) * frcv(jpr_qsrice)%z3(:,:,1)
2050	zqsr_ice(:,:,jl) = frcv(jpr_qsrice)%z3(:,:,1)
2051	ENDDO
2052	ENDIF
2053	CASE( 'mixed oce-ice' )
2054	zqsr_tot(:,: ) = frcv(jpr_qsrmix)%z3(:,:,1)
2055	! NEED TO SORT OUT HOW THIS SHOULD WORK IN THE MULTI-CATEGORY CASE - CURRENTLY NOT ALLOWED WHEN INTERFACE INITIALISED
2056	! Create solar heat flux over ice using incoming solar heat flux and albedos
2057	! ( see OASIS3 user guide, 5th edition, p39 )
2058	zqsr_ice(:,:,1) = frcv(jpr_qsrmix)%z3(:,:,1) * ( 1.- palbi(:,:,1) ) &
2059	& / ( 1.- ( alb_oce_mix(:,: ) * ziceld(:,:) &
2060	& + palbi (:,:,1) * picefr(:,:) ) )
2061	END SELECT
2062	IF( ln_dm2dc .AND. ln_cpl ) THEN ! modify qsr to include the diurnal cycle
2063	zqsr_tot(:,: ) = sbc_dcy( zqsr_tot(:,: ) )
2064	DO jl=1,jpl
2065	zqsr_ice(:,:,jl) = sbc_dcy( zqsr_ice(:,:,jl) )
2066	ENDDO
2067	ENDIF
2068
2069	#if defined key_lim3
2070	! --- solar flux over ocean --- !
2071	! note: ziceld cannot be = 0 since we limit the ice concentration to amax
2072	zqsr_oce = 0._wp
2073	WHERE( ziceld /= 0._wp ) zqsr_oce(:,:) = ( zqsr_tot(:,:) - SUM( a_i * zqsr_ice, dim=3 ) ) / ziceld(:,:)
2074
2075	IF( ln_mixcpl ) THEN ; qsr_oce(:,:) = qsr_oce(:,:) * xcplmask(:,:,0) + zqsr_oce(:,:)* zmsk(:,:)
2076	ELSE ; qsr_oce(:,:) = zqsr_oce(:,:) ; ENDIF
2077	#endif
2078
2079	IF( ln_mixcpl ) THEN
2080	qsr_tot(:,:) = qsr(:,:) * ziceld(:,:) + SUM( qsr_ice(:,:,:) * a_i(:,:,:), dim=3 ) ! total flux from blk
2081	qsr_tot(:,:) = qsr_tot(:,:) * xcplmask(:,:,0) + zqsr_tot(:,:)* zmsk(:,:)
2082	DO jl=1,jpl
2083	qsr_ice(:,:,jl) = qsr_ice(:,:,jl) * xcplmask(:,:,0) + zqsr_ice(:,:,jl)* zmsk(:,:)
2084	ENDDO
2085	ELSE
2086	qsr_tot(:,: ) = zqsr_tot(:,: )
2087	qsr_ice(:,:,:) = zqsr_ice(:,:,:)
2088	ENDIF
2089
2090	! ! ========================= !
2091	SELECT CASE( TRIM( sn_rcv_dqnsdt%cldes ) ) ! d(qns)/dt !
2092	! ! ========================= !
2093	CASE ('coupled')
2094	IF ( TRIM(sn_rcv_dqnsdt%clcat) == 'yes' ) THEN
2095	zdqns_ice(:,:,1:jpl) = frcv(jpr_dqnsdt)%z3(:,:,1:jpl)
2096	ELSE
2097	! Set all category values equal for the moment
2098	DO jl=1,jpl
2099	zdqns_ice(:,:,jl) = frcv(jpr_dqnsdt)%z3(:,:,1)
2100	ENDDO
2101	ENDIF
2102	END SELECT
2103
2104	IF( ln_mixcpl ) THEN
2105	DO jl=1,jpl
2106	dqns_ice(:,:,jl) = dqns_ice(:,:,jl) * xcplmask(:,:,0) + zdqns_ice(:,:,jl) * zmsk(:,:)
2107	ENDDO
2108	ELSE
2109	dqns_ice(:,:,:) = zdqns_ice(:,:,:)
2110	ENDIF
2111
2112	! ! ========================= !
2113	SELECT CASE( TRIM( sn_rcv_iceflx%cldes ) ) ! topmelt and botmelt !
2114	! ! ========================= !
2115	CASE ('coupled')
2116	qml_ice(:,:,:)=frcv(jpr_topm)%z3(:,:,:) * a_i(:,:,:)
2117	qcn_ice(:,:,:)=frcv(jpr_botm)%z3(:,:,:) * a_i(:,:,:)
2118
2119	!!!!! Getting NEMO4-LIM to work at Met Office: Semi-implicit coupling
2120	qml_ice(:,:,:)=qml_ice(:,:,:) * a_i(:,:,:)
2121	qcn_ice(:,:,:)=qcn_ice(:,:,:) * a_i(:,:,:)
2122	END SELECT
2123
2124	! --- Transmitted shortwave radiation (W/m2) --- !
2125
2126	IF ( nice_jules == 0 ) THEN
2127
2128	zfrqsr_tr_i(:,:,:) = 0._wp ! surface transmission parameter
2129
2130	! former coding was
2131	! Coupled case: since cloud cover is not received from atmosphere
2132	! ===> used prescribed cloud fraction representative for polar oceans in summer (0.81)
2133	! fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice )
2134	! fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice )
2135
2136	! to retrieve that coding, we needed to access h_i & h_s from here
2137	! we could even retrieve cloud fraction from the coupler
2138
2139	DO jl = 1, jpl
2140	DO jj = 1 , jpj
2141	DO ji = 1, jpi
2142
2143	!--- surface transmission parameter (Grenfell Maykut 77) --- !
2144	zfrqsr_tr_i(ji,jj,jl) = 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice
2145
2146	! --- influence of snow and thin ice --- !
2147	IF ( phs(ji,jj,jl) >= 0.0_wp ) zfrqsr_tr_i(ji,jj,jl) = 0._wp ! snow fully opaque
2148	IF ( phi(ji,jj,jl) <= 0.1_wp ) zfrqsr_tr_i(ji,jj,jl) = 1._wp ! thin ice transmits all solar radiation
2149	END DO
2150	END DO
2151	END DO
2152
2153	qsr_ice_tr(:,:,:) = zfrqsr_tr_i(:,:,:) * qsr_ice(:,:,:) ! transmitted solar radiation
2154
2155	ENDIF
2156
2157	IF ( nice_jules == 2 ) THEN
2158
2159	! here we must receive the qsr_ice_tr array from the coupler
2160	! for now just assume zero
2161
2162	qsr_ice_tr(:,:,:) = 0.0_wp
2163
2164	ENDIF
2165
2166
2167
2168	CALL wrk_dealloc( jpi,jpj, zcptn, zcptrain, zcptsnw, ziceld, zmsk, zsnw )
2169	CALL wrk_dealloc( jpi,jpj, zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip, zevap_oce, zdevap_ice )
2170	CALL wrk_dealloc( jpi,jpj, zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice )
2171	CALL wrk_dealloc( jpi,jpj,jpl, zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zfrqsr_tr_i )
2172	IF (sn_rcv_emp%clcat == 'yes') THEN
2173	CALL wrk_dealloc( jpi,jpj,jpl,zevap_ice)
2174	ELSE
2175	CALL wrk_dealloc( jpi,jpj,1,zevap_ice)
2176	ENDIF
2177	!
2178	IF( nn_timing == 1 ) CALL timing_stop('sbc_cpl_ice_flx')
2179	!
2180	END SUBROUTINE sbc_cpl_ice_flx
2181
2182
2183	SUBROUTINE sbc_cpl_snd( kt )
2184	!!----------------------------------------------------------------------
2185	!! * ROUTINE sbc_cpl_snd *
2186	!!
2187	!! ** Purpose : provide the ocean-ice informations to the atmosphere
2188	!!
2189	!! ** Method : send to the atmosphere through a call to cpl_snd
2190	!! all the needed fields (as defined in sbc_cpl_init)
2191	!!----------------------------------------------------------------------
2192	INTEGER, INTENT(in) :: kt
2193	!
2194	INTEGER :: ji, jj, jl ! dummy loop indices
2195	INTEGER :: isec, info ! local integer
2196	REAL(wp) :: zumax, zvmax
2197	REAL(wp), POINTER, DIMENSION(:,:) :: zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1
2198	REAL(wp), POINTER, DIMENSION(:,:,:) :: ztmp3, ztmp4
2199	!!----------------------------------------------------------------------
2200	!
2201	IF( nn_timing == 1 ) CALL timing_start('sbc_cpl_snd')
2202	!
2203	CALL wrk_alloc( jpi,jpj, zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 )
2204	CALL wrk_alloc( jpi,jpj,jpl, ztmp3, ztmp4 )
2205
2206	isec = ( kt - nit000 ) * NINT( rdt ) ! date of exchanges
2207
2208	zfr_l(:,:) = 1.- fr_i(:,:)
2209	! ! ------------------------- !
2210	! ! Surface temperature ! in Kelvin
2211	! ! ------------------------- !
2212	IF( ssnd(jps_toce)%laction .OR. ssnd(jps_tice)%laction .OR. ssnd(jps_tmix)%laction ) THEN
2213
2214	IF ( nn_components == jp_iam_opa ) THEN
2215	ztmp1(:,:) = tsn(:,:,1,jp_tem) ! send temperature as it is (potential or conservative) -> use of l_useCT on the received part
2216	ELSE
2217	! we must send the surface potential temperature
2218	IF( l_useCT ) THEN ; ztmp1(:,:) = eos_pt_from_ct( tsn(:,:,1,jp_tem), tsn(:,:,1,jp_sal) )
2219	ELSE ; ztmp1(:,:) = tsn(:,:,1,jp_tem)
2220	ENDIF
2221	!
2222	SELECT CASE( sn_snd_temp%cldes)
2223	CASE( 'oce only' ) ; ztmp1(:,:) = ztmp1(:,:) + rt0
2224	CASE( 'oce and ice' ) ; ztmp1(:,:) = ztmp1(:,:) + rt0
2225	SELECT CASE( sn_snd_temp%clcat )
2226	CASE( 'yes' )
2227	ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl)
2228	CASE( 'no' )
2229	WHERE( SUM( a_i, dim=3 ) /= 0. )
2230	ztmp3(:,:,1) = SUM( tn_ice * a_i, dim=3 ) / SUM( a_i, dim=3 )
2231	ELSEWHERE
2232	ztmp3(:,:,1) = rt0
2233	END WHERE
2234	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
2235	END SELECT
2236	CASE( 'weighted oce and ice' ) ; ztmp1(:,:) = ( ztmp1(:,:) + rt0 ) * zfr_l(:,:)
2237	SELECT CASE( sn_snd_temp%clcat )
2238	CASE( 'yes' )
2239	ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
2240	CASE( 'no' )
2241	ztmp3(:,:,:) = 0.0
2242	DO jl=1,jpl
2243	ztmp3(:,:,1) = ztmp3(:,:,1) + tn_ice(:,:,jl) * a_i(:,:,jl)
2244	ENDDO
2245	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
2246	END SELECT
2247	CASE( 'oce and weighted ice') ; ztmp1(:,:) = tsn(:,:,1,jp_tem) + rt0
2248	SELECT CASE( sn_snd_temp%clcat )
2249	CASE( 'yes' )
2250	ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
2251	CASE( 'no' )
2252	ztmp3(:,:,:) = 0.0
2253	DO jl=1,jpl
2254	ztmp3(:,:,1) = ztmp3(:,:,1) + tn_ice(:,:,jl) * a_i(:,:,jl)
2255	ENDDO
2256	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
2257	END SELECT
2258	CASE( 'mixed oce-ice' )
2259	ztmp1(:,:) = ( ztmp1(:,:) + rt0 ) * zfr_l(:,:)
2260	DO jl=1,jpl
2261	ztmp1(:,:) = ztmp1(:,:) + tn_ice(:,:,jl) * a_i(:,:,jl)
2262	ENDDO
2263	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%cldes' )
2264	END SELECT
2265	ENDIF
2266	IF( ssnd(jps_toce)%laction ) CALL cpl_snd( jps_toce, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
2267	IF( ssnd(jps_tice)%laction ) CALL cpl_snd( jps_tice, isec, ztmp3, info )
2268	IF( ssnd(jps_tmix)%laction ) CALL cpl_snd( jps_tmix, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
2269	ENDIF
2270
2271	!!!!! Getting NEMO4-LIM working at Met Office
2272	! Top layer ice temperature
2273	IF( ssnd(jps_ttilyr)%laction) THEN
2274	SELECT CASE( sn_snd_ttilyr%cldes)
2275	CASE ('weighted ice')
2276	ttlyr_ice(:,:,:) = t1_ice(:,:,:) + rt0
2277	ztmp3(:,:,1:jpl) = ttlyr_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
2278	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_ttilyr%cldes' )
2279	END SELECT
2280	IF( ssnd(jps_ttilyr)%laction ) CALL cpl_snd( jps_ttilyr, isec, ztmp3, info )
2281	ENDIF
2282	!!!!!
2283
2284
2285	! ! ------------------------- !
2286	! ! Albedo !
2287	! ! ------------------------- !
2288	IF( ssnd(jps_albice)%laction ) THEN ! ice
2289	SELECT CASE( sn_snd_alb%cldes )
2290	CASE( 'ice' )
2291	SELECT CASE( sn_snd_alb%clcat )
2292	CASE( 'yes' )
2293	ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl)
2294	CASE( 'no' )
2295	WHERE( SUM( a_i, dim=3 ) /= 0. )
2296	ztmp1(:,:) = SUM( alb_ice (:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 ) / SUM( a_i(:,:,1:jpl), dim=3 )
2297	ELSEWHERE
2298	ztmp1(:,:) = alb_oce_mix(:,:)
2299	END WHERE
2300	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_alb%clcat' )
2301	END SELECT
2302	CASE( 'weighted ice' ) ;
2303	SELECT CASE( sn_snd_alb%clcat )
2304	CASE( 'yes' )
2305	ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
2306	CASE( 'no' )
2307	WHERE( fr_i (:,:) > 0. )
2308	ztmp1(:,:) = SUM ( alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 )
2309	ELSEWHERE
2310	ztmp1(:,:) = 0.
2311	END WHERE
2312	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_ice%clcat' )
2313	END SELECT
2314	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_alb%cldes' )
2315	END SELECT
2316
2317	SELECT CASE( sn_snd_alb%clcat )
2318	CASE( 'yes' )
2319	CALL cpl_snd( jps_albice, isec, ztmp3, info ) !-> MV this has never been checked in coupled mode
2320	CASE( 'no' )
2321	CALL cpl_snd( jps_albice, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
2322	END SELECT
2323	ENDIF
2324
2325	IF( ssnd(jps_albmix)%laction ) THEN ! mixed ice-ocean
2326	ztmp1(:,:) = alb_oce_mix(:,:) * zfr_l(:,:)
2327	DO jl=1,jpl
2328	ztmp1(:,:) = ztmp1(:,:) + alb_ice(:,:,jl) * a_i(:,:,jl)
2329	ENDDO
2330	CALL cpl_snd( jps_albmix, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
2331	ENDIF
2332	! ! ------------------------- !
2333	! ! Ice fraction & Thickness !
2334	! ! ------------------------- !
2335	! Send ice fraction field to atmosphere
2336	IF( ssnd(jps_fice)%laction ) THEN
2337	SELECT CASE( sn_snd_thick%clcat )
2338	CASE( 'yes' ) ; ztmp3(:,:,1:jpl) = a_i(:,:,1:jpl)
2339	CASE( 'no' ) ; ztmp3(:,:,1 ) = fr_i(:,: )
2340	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%clcat' )
2341	END SELECT
2342	IF( ssnd(jps_fice)%laction ) CALL cpl_snd( jps_fice, isec, ztmp3, info )
2343	ENDIF
2344
2345	IF( ssnd(jps_fice1)%laction ) THEN
2346	SELECT CASE( sn_snd_thick1%clcat )
2347	CASE( 'yes' ) ; ztmp3(:,:,1:jpl) = a_i(:,:,1:jpl)
2348	CASE( 'no' ) ; ztmp3(:,:,1 ) = fr_i(:,: )
2349	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick1%clcat' )
2350	END SELECT
2351	CALL cpl_snd( jps_fice1, isec, ztmp3, info )
2352	ENDIF
2353
2354	! Send ice fraction field to OPA (sent by SAS in SAS-OPA coupling)
2355	IF( ssnd(jps_fice2)%laction ) THEN
2356	ztmp3(:,:,1) = fr_i(:,:)
2357	IF( ssnd(jps_fice2)%laction ) CALL cpl_snd( jps_fice2, isec, ztmp3, info )
2358	ENDIF
2359
2360	! Send ice and snow thickness field
2361	IF( ssnd(jps_hice)%laction .OR. ssnd(jps_hsnw)%laction ) THEN
2362	SELECT CASE( sn_snd_thick%cldes)
2363	CASE( 'none' ) ! nothing to do
2364	CASE( 'weighted ice and snow' )
2365	SELECT CASE( sn_snd_thick%clcat )
2366	CASE( 'yes' )
2367	ztmp3(:,:,1:jpl) = h_i(:,:,1:jpl) * a_i(:,:,1:jpl)
2368	ztmp4(:,:,1:jpl) = h_s(:,:,1:jpl) * a_i(:,:,1:jpl)
2369	CASE( 'no' )
2370	ztmp3(:,:,:) = 0.0 ; ztmp4(:,:,:) = 0.0
2371	DO jl=1,jpl
2372	ztmp3(:,:,1) = ztmp3(:,:,1) + h_i(:,:,jl) * a_i(:,:,jl)
2373	ztmp4(:,:,1) = ztmp4(:,:,1) + h_s(:,:,jl) * a_i(:,:,jl)
2374	ENDDO
2375	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%clcat' )
2376	END SELECT
2377	CASE( 'ice and snow' )
2378	SELECT CASE( sn_snd_thick%clcat )
2379	CASE( 'yes' )
2380	ztmp3(:,:,1:jpl) = h_i(:,:,1:jpl)
2381	ztmp4(:,:,1:jpl) = h_s(:,:,1:jpl)
2382	CASE( 'no' )
2383	WHERE( SUM( a_i, dim=3 ) /= 0. )
2384	ztmp3(:,:,1) = SUM( h_i * a_i, dim=3 ) / SUM( a_i, dim=3 )
2385	ztmp4(:,:,1) = SUM( h_s * a_i, dim=3 ) / SUM( a_i, dim=3 )
2386	ELSEWHERE
2387	ztmp3(:,:,1) = 0.
2388	ztmp4(:,:,1) = 0.
2389	END WHERE
2390	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%clcat' )
2391	END SELECT
2392	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%cldes' )
2393	END SELECT
2394	IF( ssnd(jps_hice)%laction ) CALL cpl_snd( jps_hice, isec, ztmp3, info )
2395	IF( ssnd(jps_hsnw)%laction ) CALL cpl_snd( jps_hsnw, isec, ztmp4, info )
2396	ENDIF
2397
2398	!!!!! Getting NEMO4-LIM to work
2399	!!!!! Temporary code while we're not modelling meltponds
2400	a_p(:,:,1:jpl) = 0.0
2401	ht_p(:,:,1:jpl) = 0.0
2402
2403	!
2404	! Send meltpond fields
2405	IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN
2406	SELECT CASE( sn_snd_mpnd%cldes)
2407	CASE( 'weighted ice' )
2408	SELECT CASE( sn_snd_mpnd%clcat )
2409	CASE( 'yes' )
2410	ztmp3(:,:,1:jpl) = a_p(:,:,1:jpl) * a_i(:,:,1:jpl)
2411	ztmp4(:,:,1:jpl) = ht_p(:,:,1:jpl) * a_i(:,:,1:jpl)
2412	CASE( 'no' )
2413	ztmp3(:,:,:) = 0.0
2414	ztmp4(:,:,:) = 0.0
2415	DO jl=1,jpl
2416	ztmp3(:,:,1) = ztmp3(:,:,1) + a_p(:,:,jpl) * a_i(:,:,jpl)
2417	ztmp4(:,:,1) = ztmp4(:,:,1) + ht_p(:,:,jpl) * a_i(:,:,jpl)
2418	ENDDO
2419	CASE default ; CALL ctl_stop( 'sbc_cpl_mpd: wrong definition of sn_snd_mpnd%clcat' )
2420	END SELECT
2421	CASE( 'ice only' )
2422	ztmp3(:,:,1:jpl) = a_p(:,:,1:jpl)
2423	ztmp4(:,:,1:jpl) = ht_p(:,:,1:jpl)
2424	END SELECT
2425	IF( ssnd(jps_a_p)%laction ) CALL cpl_snd( jps_a_p, isec, ztmp3, info )
2426	IF( ssnd(jps_ht_p)%laction ) CALL cpl_snd( jps_ht_p, isec, ztmp4, info )
2427	!
2428	! Send ice effective conductivity
2429	SELECT CASE( sn_snd_cond%cldes)
2430	CASE( 'weighted ice' )
2431	SELECT CASE( sn_snd_cond%clcat )
2432	CASE( 'yes' )
2433	ztmp3(:,:,1:jpl) = cnd_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
2434	CASE( 'no' )
2435	ztmp3(:,:,:) = 0.0
2436	DO jl=1,jpl
2437	ztmp3(:,:,1) = ztmp3(:,:,1) + cnd_ice(:,:,jl) * a_i(:,:,jl)
2438	ENDDO
2439	CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_cond%clcat' )
2440	END SELECT
2441	CASE( 'ice only' )
2442	ztmp3(:,:,1:jpl) = cnd_ice(:,:,1:jpl)
2443	END SELECT
2444	IF( ssnd(jps_kice)%laction ) CALL cpl_snd( jps_kice, isec, ztmp3, info )
2445	ENDIF
2446	!
2447	!!!!!
2448
2449
2450	! ! ------------------------- !
2451	! ! CO2 flux from PISCES !
2452	! ! ------------------------- !
2453	IF( ssnd(jps_co2)%laction .AND. l_co2cpl ) CALL cpl_snd( jps_co2, isec, RESHAPE ( oce_co2, (/jpi,jpj,1/) ) , info )
2454	!
2455	! ! ------------------------- !
2456	IF( ssnd(jps_ocx1)%laction ) THEN ! Surface current !
2457	! ! ------------------------- !
2458	!
2459	! j+1 j -----V---F
2460	! surface velocity always sent from T point ! \|
2461	! j \| T U
2462	! \| \|
2463	! j j-1 -I-------\|
2464	! (for I) \| \|
2465	! i-1 i i
2466	! i i+1 (for I)
2467	IF( nn_components == jp_iam_opa ) THEN
2468	zotx1(:,:) = un(:,:,1)
2469	zoty1(:,:) = vn(:,:,1)
2470	ELSE
2471	SELECT CASE( TRIM( sn_snd_crt%cldes ) )
2472	CASE( 'oce only' ) ! C-grid ==> T
2473	DO jj = 2, jpjm1
2474	DO ji = fs_2, fs_jpim1 ! vector opt.
2475	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) )
2476	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) )
2477	END DO
2478	END DO
2479	CASE( 'weighted oce and ice' )
2480	SELECT CASE ( cp_ice_msh )
2481	CASE( 'C' ) ! Ocean and Ice on C-grid ==> T
2482	DO jj = 2, jpjm1
2483	DO ji = fs_2, fs_jpim1 ! vector opt.
2484	zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(ji,jj)
2485	zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + vn (ji ,jj-1,1) ) * zfr_l(ji,jj)
2486	zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj)
2487	zity1(ji,jj) = 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj)
2488	END DO
2489	END DO
2490	CASE( 'I' ) ! Ocean on C grid, Ice on I-point (B-grid) ==> T
2491	DO jj = 2, jpjm1
2492	DO ji = 2, jpim1 ! NO vector opt.
2493	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj)
2494	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj)
2495	zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1) &
2496	& + u_ice(ji+1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2497	zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1) &
2498	& + v_ice(ji+1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2499	END DO
2500	END DO
2501	CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T
2502	DO jj = 2, jpjm1
2503	DO ji = 2, jpim1 ! NO vector opt.
2504	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj)
2505	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj)
2506	zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) &
2507	& + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2508	zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) &
2509	& + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2510	END DO
2511	END DO
2512	END SELECT
2513	CALL lbc_lnk( zitx1, 'T', -1. ) ; CALL lbc_lnk( zity1, 'T', -1. )
2514	CASE( 'mixed oce-ice' )
2515	SELECT CASE ( cp_ice_msh )
2516	CASE( 'C' ) ! Ocean and Ice on C-grid ==> T
2517	DO jj = 2, jpjm1
2518	DO ji = fs_2, fs_jpim1 ! vector opt.
2519	zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(ji,jj) &
2520	& + 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj)
2521	zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + vn (ji ,jj-1,1) ) * zfr_l(ji,jj) &
2522	& + 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj)
2523	END DO
2524	END DO
2525	CASE( 'I' ) ! Ocean on C grid, Ice on I-point (B-grid) ==> T
2526	DO jj = 2, jpjm1
2527	DO ji = 2, jpim1 ! NO vector opt.
2528	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) &
2529	& + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1) &
2530	& + u_ice(ji+1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2531	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) &
2532	& + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1) &
2533	& + v_ice(ji+1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2534	END DO
2535	END DO
2536	CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T
2537	DO jj = 2, jpjm1
2538	DO ji = 2, jpim1 ! NO vector opt.
2539	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) &
2540	& + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) &
2541	& + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2542	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) &
2543	& + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) &
2544	& + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2545	END DO
2546	END DO
2547	END SELECT
2548	END SELECT
2549	CALL lbc_lnk( zotx1, ssnd(jps_ocx1)%clgrid, -1. ) ; CALL lbc_lnk( zoty1, ssnd(jps_ocy1)%clgrid, -1. )
2550	!
2551	ENDIF
2552	!
2553	!
2554	IF( TRIM( sn_snd_crt%clvor ) == 'eastward-northward' ) THEN ! Rotation of the components
2555	! ! Ocean component
2556	CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component
2557	CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component
2558	zotx1(:,:) = ztmp1(:,:) ! overwrite the components
2559	zoty1(:,:) = ztmp2(:,:)
2560	IF( ssnd(jps_ivx1)%laction ) THEN ! Ice component
2561	CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component
2562	CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component
2563	zitx1(:,:) = ztmp1(:,:) ! overwrite the components
2564	zity1(:,:) = ztmp2(:,:)
2565	ENDIF
2566	ENDIF
2567	!
2568	! spherical coordinates to cartesian -> 2 components to 3 components
2569	IF( TRIM( sn_snd_crt%clvref ) == 'cartesian' ) THEN
2570	ztmp1(:,:) = zotx1(:,:) ! ocean currents
2571	ztmp2(:,:) = zoty1(:,:)
2572	CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 )
2573	!
2574	IF( ssnd(jps_ivx1)%laction ) THEN ! ice velocities
2575	ztmp1(:,:) = zitx1(:,:)
2576	ztmp1(:,:) = zity1(:,:)
2577	CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 )
2578	ENDIF
2579	ENDIF
2580	!
2581	IF( ssnd(jps_ocx1)%laction ) CALL cpl_snd( jps_ocx1, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info ) ! ocean x current 1st grid
2582	IF( ssnd(jps_ocy1)%laction ) CALL cpl_snd( jps_ocy1, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info ) ! ocean y current 1st grid
2583	IF( ssnd(jps_ocz1)%laction ) CALL cpl_snd( jps_ocz1, isec, RESHAPE ( zotz1, (/jpi,jpj,1/) ), info ) ! ocean z current 1st grid
2584	!
2585	IF( ssnd(jps_ivx1)%laction ) CALL cpl_snd( jps_ivx1, isec, RESHAPE ( zitx1, (/jpi,jpj,1/) ), info ) ! ice x current 1st grid
2586	IF( ssnd(jps_ivy1)%laction ) CALL cpl_snd( jps_ivy1, isec, RESHAPE ( zity1, (/jpi,jpj,1/) ), info ) ! ice y current 1st grid
2587	IF( ssnd(jps_ivz1)%laction ) CALL cpl_snd( jps_ivz1, isec, RESHAPE ( zitz1, (/jpi,jpj,1/) ), info ) ! ice z current 1st grid
2588	!
2589	ENDIF
2590	!
2591	! ! ------------------------- !
2592	! ! Surface current to waves !
2593	! ! ------------------------- !
2594	IF( ssnd(jps_ocxw)%laction .OR. ssnd(jps_ocyw)%laction ) THEN
2595	!
2596	! j+1 j -----V---F
2597	! surface velocity always sent from T point ! \|
2598	! j \| T U
2599	! \| \|
2600	! j j-1 -I-------\|
2601	! (for I) \| \|
2602	! i-1 i i
2603	! i i+1 (for I)
2604	SELECT CASE( TRIM( sn_snd_crtw%cldes ) )
2605	CASE( 'oce only' ) ! C-grid ==> T
2606	DO jj = 2, jpjm1
2607	DO ji = fs_2, fs_jpim1 ! vector opt.
2608	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) )
2609	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji , jj-1,1) )
2610	END DO
2611	END DO
2612	CASE( 'weighted oce and ice' )
2613	SELECT CASE ( cp_ice_msh )
2614	CASE( 'C' ) ! Ocean and Ice on C-grid ==> T
2615	DO jj = 2, jpjm1
2616	DO ji = fs_2, fs_jpim1 ! vector opt.
2617	zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(ji,jj)
2618	zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + vn (ji ,jj-1,1) ) * zfr_l(ji,jj)
2619	zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj)
2620	zity1(ji,jj) = 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj)
2621	END DO
2622	END DO
2623	CASE( 'I' ) ! Ocean on C grid, Ice on I-point (B-grid) ==> T
2624	DO jj = 2, jpjm1
2625	DO ji = 2, jpim1 ! NO vector opt.
2626	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj)
2627	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj)
2628	zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1) &
2629	& + u_ice(ji+1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2630	zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1) &
2631	& + v_ice(ji+1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2632	END DO
2633	END DO
2634	CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T
2635	DO jj = 2, jpjm1
2636	DO ji = 2, jpim1 ! NO vector opt.
2637	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj)
2638	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj)
2639	zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) &
2640	& + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2641	zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) &
2642	& + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2643	END DO
2644	END DO
2645	END SELECT
2646	CALL lbc_lnk( zitx1, 'T', -1. ) ; CALL lbc_lnk( zity1, 'T', -1. )
2647	CASE( 'mixed oce-ice' )
2648	SELECT CASE ( cp_ice_msh )
2649	CASE( 'C' ) ! Ocean and Ice on C-grid ==> T
2650	DO jj = 2, jpjm1
2651	DO ji = fs_2, fs_jpim1 ! vector opt.
2652	zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(ji,jj) &
2653	& + 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj)
2654	zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + vn (ji ,jj-1,1) ) * zfr_l(ji,jj) &
2655	& + 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj)
2656	END DO
2657	END DO
2658	CASE( 'I' ) ! Ocean on C grid, Ice on I-point (B-grid) ==> T
2659	DO jj = 2, jpjm1
2660	DO ji = 2, jpim1 ! NO vector opt.
2661	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) &
2662	& + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1) &
2663	& + u_ice(ji+1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2664	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) &
2665	& + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1) &
2666	& + v_ice(ji+1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2667	END DO
2668	END DO
2669	CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T
2670	DO jj = 2, jpjm1
2671	DO ji = 2, jpim1 ! NO vector opt.
2672	zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) &
2673	& + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) &
2674	& + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj)
2675	zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) &
2676	& + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) &
2677	& + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj)
2678	END DO
2679	END DO
2680	END SELECT
2681	END SELECT
2682	CALL lbc_lnk( zotx1, ssnd(jps_ocxw)%clgrid, -1. ) ; CALL lbc_lnk( zoty1, ssnd(jps_ocyw)%clgrid, -1. )
2683	!
2684	!
2685	IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) THEN ! Rotation of the components
2686	! ! Ocean component
2687	CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->e', ztmp1 ) ! 1st component
2688	CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->n', ztmp2 ) ! 2nd component
2689	zotx1(:,:) = ztmp1(:,:) ! overwrite the components
2690	zoty1(:,:) = ztmp2(:,:)
2691	IF( ssnd(jps_ivx1)%laction ) THEN ! Ice component
2692	CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component
2693	CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component
2694	zitx1(:,:) = ztmp1(:,:) ! overwrite the components
2695	zity1(:,:) = ztmp2(:,:)
2696	ENDIF
2697	ENDIF
2698	!
2699	! ! spherical coordinates to cartesian -> 2 components to 3 components
2700	! IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN
2701	! ztmp1(:,:) = zotx1(:,:) ! ocean currents
2702	! ztmp2(:,:) = zoty1(:,:)
2703	! CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 )
2704	! !
2705	! IF( ssnd(jps_ivx1)%laction ) THEN ! ice velocities
2706	! ztmp1(:,:) = zitx1(:,:)
2707	! ztmp1(:,:) = zity1(:,:)
2708	! CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 )
2709	! ENDIF
2710	! ENDIF
2711	!
2712	IF( ssnd(jps_ocxw)%laction ) CALL cpl_snd( jps_ocxw, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info ) ! ocean x current 1st grid
2713	IF( ssnd(jps_ocyw)%laction ) CALL cpl_snd( jps_ocyw, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info ) ! ocean y current 1st grid
2714	!
2715	ENDIF
2716	!
2717	IF( ssnd(jps_ficet)%laction ) THEN
2718	CALL cpl_snd( jps_ficet, isec, RESHAPE ( fr_i, (/jpi,jpj,1/) ), info )
2719	END IF
2720	! ! ------------------------- !
2721	! ! Water levels to waves !
2722	! ! ------------------------- !
2723	IF( ssnd(jps_wlev)%laction ) THEN
2724	IF( ln_apr_dyn ) THEN
2725	IF( kt /= nit000 ) THEN
2726	ztmp1(:,:) = sshb(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) )
2727	ELSE
2728	ztmp1(:,:) = sshb(:,:)
2729	ENDIF
2730	ELSE
2731	ztmp1(:,:) = sshn(:,:)
2732	ENDIF
2733	CALL cpl_snd( jps_wlev , isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
2734	END IF
2735	!
2736	! Fields sent by OPA to SAS when doing OPA<->SAS coupling
2737	! ! SSH
2738	IF( ssnd(jps_ssh )%laction ) THEN
2739	! ! removed inverse barometer ssh when Patm
2740	! forcing is used (for sea-ice dynamics)
2741	IF( ln_apr_dyn ) THEN ; ztmp1(:,:) = sshb(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) )
2742	ELSE ; ztmp1(:,:) = sshn(:,:)
2743	ENDIF
2744	CALL cpl_snd( jps_ssh , isec, RESHAPE ( ztmp1 , (/jpi,jpj,1/) ), info )
2745
2746	ENDIF
2747	! ! SSS
2748	IF( ssnd(jps_soce )%laction ) THEN
2749	CALL cpl_snd( jps_soce , isec, RESHAPE ( tsn(:,:,1,jp_sal), (/jpi,jpj,1/) ), info )
2750	ENDIF
2751	! ! first T level thickness
2752	IF( ssnd(jps_e3t1st )%laction ) THEN
2753	CALL cpl_snd( jps_e3t1st, isec, RESHAPE ( e3t_n(:,:,1) , (/jpi,jpj,1/) ), info )
2754	ENDIF
2755	! ! Qsr fraction
2756	IF( ssnd(jps_fraqsr)%laction ) THEN
2757	CALL cpl_snd( jps_fraqsr, isec, RESHAPE ( fraqsr_1lev(:,:) , (/jpi,jpj,1/) ), info )
2758	ENDIF
2759	!
2760	! Fields sent by SAS to OPA when OASIS coupling
2761	! ! Solar heat flux
2762	IF( ssnd(jps_qsroce)%laction ) CALL cpl_snd( jps_qsroce, isec, RESHAPE ( qsr , (/jpi,jpj,1/) ), info )
2763	IF( ssnd(jps_qnsoce)%laction ) CALL cpl_snd( jps_qnsoce, isec, RESHAPE ( qns , (/jpi,jpj,1/) ), info )
2764	IF( ssnd(jps_oemp )%laction ) CALL cpl_snd( jps_oemp , isec, RESHAPE ( emp , (/jpi,jpj,1/) ), info )
2765	IF( ssnd(jps_sflx )%laction ) CALL cpl_snd( jps_sflx , isec, RESHAPE ( sfx , (/jpi,jpj,1/) ), info )
2766	IF( ssnd(jps_otx1 )%laction ) CALL cpl_snd( jps_otx1 , isec, RESHAPE ( utau, (/jpi,jpj,1/) ), info )
2767	IF( ssnd(jps_oty1 )%laction ) CALL cpl_snd( jps_oty1 , isec, RESHAPE ( vtau, (/jpi,jpj,1/) ), info )
2768	IF( ssnd(jps_rnf )%laction ) CALL cpl_snd( jps_rnf , isec, RESHAPE ( rnf , (/jpi,jpj,1/) ), info )
2769	IF( ssnd(jps_taum )%laction ) CALL cpl_snd( jps_taum , isec, RESHAPE ( taum, (/jpi,jpj,1/) ), info )
2770
2771	ztmp1(:,:) = sstfrz(:,:) + rt0
2772	IF( ssnd(jps_sstfrz)%laction ) CALL cpl_snd( jps_sstfrz, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info)
2773
2774	CALL wrk_dealloc( jpi,jpj, zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 )
2775	CALL wrk_dealloc( jpi,jpj,jpl, ztmp3, ztmp4 )
2776	!
2777	IF( nn_timing == 1 ) CALL timing_stop('sbc_cpl_snd')
2778	!
2779	END SUBROUTINE sbc_cpl_snd
2780
2781	!!======================================================================
2782	END MODULE sbccpl

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

Download in other formats: