source: branches/UKMO/AMM15_v3_6_STABLE_package_collate_utils305/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 12004

Last change on this file since 12004 was 12004, checked in by dford, 10 months ago

Correct merge errors from previous commit.

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