source: NEMO/branches/UKMO/r8395_cpl-pressure/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 10803

Last change on this file since 10803 was 10803, checked in by jcastill, 19 months ago

Complete set of changes to use pressure read from coupling when available rather than from file

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