source: NEMO/branches/2020/dev_r12702_ASINTER-02_emanuelaclementi_Waves/src/OCE/SBC/sbccpl.F90 @ 12991

Last change on this file since 12991 was 12991, checked in by emanuelaclementi, 5 months ago

Included wave-current processes following Couvelard et al 2019 and Gurvan code -ticket #2155

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