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 @ 5357

Last change on this file since 5357 was 5357, checked in by clem, 9 years ago

LIM3: change the interface between the ice and atm for both coupled and forced modes. Some work still needs to be done to deal with sublimation in coupled mode.

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