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

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

source: branches/UKMO/dev_r6912_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 9132

Last change on this file since 9132 was 9132, checked in by andmirek, 6 years ago

#1868 changes enabling coupling

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