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_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 6668

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

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

Download in other formats: