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

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

dev_r5218_CNRS17_coupling: update for fraqsr

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