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

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

source: branches/UKMO/dev_r8126_LIM3_couple/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 8878

Last change on this file since 8878 was 8878, checked in by frrh, 6 years ago

Merge in http://fcm3/projects/NEMO.xm/log/branches/UKMO/dev_r8183_GC_couple_pkg
revisions 8731:8734 inclusive

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