New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
sbccpl.F90 in NEMO/branches/UKMO/r12083_cpl-pressure/src/OCE/SBC – NEMO

source: NEMO/branches/UKMO/r12083_cpl-pressure/src/OCE/SBC/sbccpl.F90 @ 12647

Last change on this file since 12647 was 12461, checked in by jcastill, 4 years ago

Changes as the original branch updated to vn4.1

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