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

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

source: branches/2015/dev_r5218_CNRS17_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 5403

Last change on this file since 5403 was 5403, checked in by smasson, 9 years ago

dev_r5218_CNRS17_coupling: sette ok

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