source: trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 5407

Last change on this file since 5407 was 5407, checked in by smasson, 6 years ago

merge dev_r5218_CNRS17_coupling into the trunk

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