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

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

source: NEMO/branches/UKMO/r8395_coupling_sequence/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 10763

Last change on this file since 10763 was 10763, checked in by jcastill, 5 years ago

Remove svn keywords properly

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