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

Context Navigation

source: branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 7809

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

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

Download in other formats: