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

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

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

Last change on this file since 5382 was 5382, checked in by smasson, 10 years ago

dev_r5218_CNRS17_coupling: merge with trunk rev 5381

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