source: branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 11107

Last change on this file since 11107 was 11107, checked in by frrh, 18 months ago

Commit changes from Dan Copsey's sea ice heat coupling
flux fixes in branch:
branches/UKMO/dev_r5518_GO6_fix_zemp_ice_10681
revisions 11028:11088.

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