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_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

source: branches/UKMO/dev_r5518_GSI7_GSI8_landice_bitcomp_medusa/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 6687

Last change on this file since 6687 was 6687, checked in by frrh, 8 years ago

Remove redundant line. Add comments for outgoing MEDUSA fields.

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