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

source: trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 7761

Last change on this file since 7761 was 7761, checked in by clem, 7 years ago

make AGRIF and LIM3 fully compatible

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