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

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

source: branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 8917

Last change on this file since 8917 was 8917, checked in by jcastill, 6 years ago

Some fixes to take into account that the received wave fields are on the T grid and winds on the U,V grids, so that variables are calculated in the right grid

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