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

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

source: branches/UKMO/r5936_hadgem3_cplseq/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 7134

Last change on this file since 7134 was 7134, checked in by jcastill, 8 years ago

Version as in UKMO/dev_r5107_hadgem3_cplseq@5646

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