source: NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl/src/OCE/SBC/sbccpl.F90 @ 12948

Last change on this file since 12948 was 12948, checked in by cetlod, 5 months ago

Minor correction to take into account the mixed oce-ice coupling type for ice and ocean wind stress

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