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

Last change on this file was 72, checked in by smasson, 11 years ago

put bugfixes from NEMO v3.4.1

File size: 96.7 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#if defined key_oasis3 || defined key_oasis4
12   !!----------------------------------------------------------------------
13   !!   'key_oasis3' or 'key_oasis4'   Coupled Ocean/Atmosphere formulation
14   !!----------------------------------------------------------------------
15   !!   namsbc_cpl      : coupled formulation namlist
16   !!   sbc_cpl_init    : initialisation of the coupled exchanges
17   !!   sbc_cpl_rcv     : receive fields from the atmosphere over the ocean (ocean only)
18   !!                     receive stress from the atmosphere over the ocean (ocean-ice case)
19   !!   sbc_cpl_ice_tau : receive stress from the atmosphere over ice
20   !!   sbc_cpl_ice_flx : receive fluxes from the atmosphere over ice
21   !!   sbc_cpl_snd     : send     fields to the atmosphere
22   !!----------------------------------------------------------------------
23   USE dom_oce         ! ocean space and time domain
24   USE sbc_oce         ! Surface boundary condition: ocean fields
25   USE sbc_ice         ! Surface boundary condition: ice fields
26   USE sbcdcy          ! surface boundary condition: diurnal cycle
27   USE phycst          ! physical constants
28#if defined key_lim3
29   USE par_ice         ! ice parameters
30   USE ice             ! ice variables
31#endif
32#if defined key_lim2
33   USE par_ice_2       ! ice parameters
34   USE ice_2           ! ice variables
35#endif
36#if defined key_oasis3
37   USE cpl_oasis3      ! OASIS3 coupling
38#endif
39#if defined key_oasis4
40   USE cpl_oasis4      ! OASIS4 coupling
41#endif
42   USE geo2ocean       !
43   USE oce   , ONLY : tsn, un, vn
44   USE albedo          !
45   USE in_out_manager  ! I/O manager
46   USE iom             ! NetCDF library
47   USE lib_mpp         ! distribued memory computing library
48   USE wrk_nemo        ! work arrays
49   USE timing          ! Timing
50   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
51#if defined key_cpl_carbon_cycle
52   USE p4zflx, ONLY : oce_co2
53#endif
54   USE diaar5, ONLY :   lk_diaar5
55#if defined key_cice
56   USE ice_domain_size, only: ncat
57#endif
58   IMPLICIT NONE
59   PRIVATE
60!EM XIOS-OASIS-MCT compliance
61   PUBLIC   sbc_cpl_init       ! routine called by sbcmod.F90
62   PUBLIC   sbc_cpl_rcv        ! routine called by sbc_ice_lim(_2).F90
63   PUBLIC   sbc_cpl_snd        ! routine called by step.F90
64   PUBLIC   sbc_cpl_ice_tau    ! routine called by sbc_ice_lim(_2).F90
65   PUBLIC   sbc_cpl_ice_flx    ! routine called by sbc_ice_lim(_2).F90
66
67   INTEGER, PARAMETER ::   jpr_otx1   =  1            ! 3 atmosphere-ocean stress components on grid 1
68   INTEGER, PARAMETER ::   jpr_oty1   =  2            !
69   INTEGER, PARAMETER ::   jpr_otz1   =  3            !
70   INTEGER, PARAMETER ::   jpr_otx2   =  4            ! 3 atmosphere-ocean stress components on grid 2
71   INTEGER, PARAMETER ::   jpr_oty2   =  5            !
72   INTEGER, PARAMETER ::   jpr_otz2   =  6            !
73   INTEGER, PARAMETER ::   jpr_itx1   =  7            ! 3 atmosphere-ice   stress components on grid 1
74   INTEGER, PARAMETER ::   jpr_ity1   =  8            !
75   INTEGER, PARAMETER ::   jpr_itz1   =  9            !
76   INTEGER, PARAMETER ::   jpr_itx2   = 10            ! 3 atmosphere-ice   stress components on grid 2
77   INTEGER, PARAMETER ::   jpr_ity2   = 11            !
78   INTEGER, PARAMETER ::   jpr_itz2   = 12            !
79   INTEGER, PARAMETER ::   jpr_qsroce = 13            ! Qsr above the ocean
80   INTEGER, PARAMETER ::   jpr_qsrice = 14            ! Qsr above the ice
81   INTEGER, PARAMETER ::   jpr_qsrmix = 15 
82   INTEGER, PARAMETER ::   jpr_qnsoce = 16            ! Qns above the ocean
83   INTEGER, PARAMETER ::   jpr_qnsice = 17            ! Qns above the ice
84   INTEGER, PARAMETER ::   jpr_qnsmix = 18
85   INTEGER, PARAMETER ::   jpr_rain   = 19            ! total liquid precipitation (rain)
86   INTEGER, PARAMETER ::   jpr_snow   = 20            ! solid precipitation over the ocean (snow)
87   INTEGER, PARAMETER ::   jpr_tevp   = 21            ! total evaporation
88   INTEGER, PARAMETER ::   jpr_ievp   = 22            ! solid evaporation (sublimation)
89   INTEGER, PARAMETER ::   jpr_sbpr   = 23            ! sublimation - liquid precipitation - solid precipitation
90   INTEGER, PARAMETER ::   jpr_semp   = 24            ! solid freshwater budget (sublimation - snow)
91   INTEGER, PARAMETER ::   jpr_oemp   = 25            ! ocean freshwater budget (evap - precip)
92   INTEGER, PARAMETER ::   jpr_w10m   = 26            ! 10m wind
93   INTEGER, PARAMETER ::   jpr_dqnsdt = 27            ! d(Q non solar)/d(temperature)
94   INTEGER, PARAMETER ::   jpr_rnf    = 28            ! runoffs
95   INTEGER, PARAMETER ::   jpr_cal    = 29            ! calving
96   INTEGER, PARAMETER ::   jpr_taum   = 30            ! wind stress module
97   INTEGER, PARAMETER ::   jpr_co2    = 31
98   INTEGER, PARAMETER ::   jpr_topm   = 32            ! topmeltn
99   INTEGER, PARAMETER ::   jpr_botm   = 33            ! botmeltn
100   INTEGER, PARAMETER ::   jprcv      = 33            ! total number of fields received
101
102   INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction
103   INTEGER, PARAMETER ::   jps_toce   =  2            ! ocean temperature
104   INTEGER, PARAMETER ::   jps_tice   =  3            ! ice   temperature
105   INTEGER, PARAMETER ::   jps_tmix   =  4            ! mixed temperature (ocean+ice)
106   INTEGER, PARAMETER ::   jps_albice =  5            ! ice   albedo
107   INTEGER, PARAMETER ::   jps_albmix =  6            ! mixed albedo
108   INTEGER, PARAMETER ::   jps_hice   =  7            ! ice  thickness
109   INTEGER, PARAMETER ::   jps_hsnw   =  8            ! snow thickness
110   INTEGER, PARAMETER ::   jps_ocx1   =  9            ! ocean current on grid 1
111   INTEGER, PARAMETER ::   jps_ocy1   = 10            !
112   INTEGER, PARAMETER ::   jps_ocz1   = 11            !
113   INTEGER, PARAMETER ::   jps_ivx1   = 12            ! ice   current on grid 1
114   INTEGER, PARAMETER ::   jps_ivy1   = 13            !
115   INTEGER, PARAMETER ::   jps_ivz1   = 14            !
116   INTEGER, PARAMETER ::   jps_co2    = 15
117   INTEGER, PARAMETER ::   jpsnd      = 15            ! total number of fields sended
118
119   !                                                         !!** namelist namsbc_cpl **
120   TYPE ::   FLD_C
121      CHARACTER(len = 32) ::   cldes                  ! desciption of the coupling strategy
122      CHARACTER(len = 32) ::   clcat                  ! multiple ice categories strategy
123      CHARACTER(len = 32) ::   clvref                 ! reference of vector ('spherical' or 'cartesian')
124      CHARACTER(len = 32) ::   clvor                  ! orientation of vector fields ('eastward-northward' or 'local grid')
125      CHARACTER(len = 32) ::   clvgrd                 ! grids on which is located the vector fields
126   END TYPE FLD_C
127   ! Send to the atmosphere                           !
128   TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2                       
129   ! Received from the atmosphere                     !
130   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
131   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2                       
132
133   TYPE ::   DYNARR     
134      REAL(wp), POINTER, DIMENSION(:,:,:)    ::   z3   
135   END TYPE DYNARR
136
137   TYPE( DYNARR ), SAVE, DIMENSION(jprcv) ::   frcv                      ! all fields recieved from the atmosphere
138
139   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   albedo_oce_mix     ! ocean albedo sent to atmosphere (mix clear/overcast sky)
140
141   INTEGER , ALLOCATABLE, SAVE, DIMENSION(    :) ::   nrcvinfo           ! OASIS info argument
142
143#if ! defined key_lim2   &&   ! defined key_lim3
144   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   u_ice, v_ice,fr1_i0,fr2_i0          ! jpi, jpj
145   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tn_ice, alb_ice, qns_ice, dqns_ice  ! (jpi,jpj,jpl)
146#endif
147
148#if defined key_cice
149   INTEGER, PARAMETER ::   jpl = ncat
150#elif ! defined key_lim2   &&   ! defined key_lim3
151   INTEGER, PARAMETER ::   jpl = 1 
152   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   emp_ice
153   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qsr_ice
154#endif
155
156#if ! defined key_lim3   &&  ! defined key_cice
157   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_i
158#endif
159
160#if ! defined key_lim3
161   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  ht_i, ht_s
162#endif
163
164#if ! defined key_cice
165   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  topmelt, botmelt
166#endif
167
168   !! Substitution
169#  include "vectopt_loop_substitute.h90"
170   !!----------------------------------------------------------------------
171   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
172   !! $Id: sbccpl.F90 3294 2012-01-28 16:44:18Z rblod $
173   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
174   !!----------------------------------------------------------------------
175
176CONTAINS
177 
178   INTEGER FUNCTION sbc_cpl_alloc()
179      !!----------------------------------------------------------------------
180      !!             ***  FUNCTION sbc_cpl_alloc  ***
181      !!----------------------------------------------------------------------
182      INTEGER :: ierr(4),jn
183      !!----------------------------------------------------------------------
184      ierr(:) = 0
185      !
186      ALLOCATE( albedo_oce_mix(jpi,jpj), nrcvinfo(jprcv),  STAT=ierr(1) )
187      !
188#if ! defined key_lim2 && ! defined key_lim3
189      ! quick patch to be able to run the coupled model without sea-ice...
190      ALLOCATE( u_ice(jpi,jpj) , fr1_i0(jpi,jpj) , tn_ice (jpi,jpj,1) ,     &
191                v_ice(jpi,jpj) , fr2_i0(jpi,jpj) , alb_ice(jpi,jpj,1),      &
192                emp_ice(jpi,jpj) , qns_ice(jpi,jpj,1) , dqns_ice(jpi,jpj,1) , STAT=ierr(2) )
193#endif
194
195#if ! defined key_lim3 && ! defined key_cice
196      ALLOCATE( a_i(jpi,jpj,jpl) , STAT=ierr(3) )
197#endif
198
199#if defined key_cice || defined key_lim2
200      ALLOCATE( ht_i(jpi,jpj,jpl) , ht_s(jpi,jpj,jpl) , STAT=ierr(4) )
201#endif
202      sbc_cpl_alloc = MAXVAL( ierr )
203      IF( lk_mpp            )   CALL mpp_sum ( sbc_cpl_alloc )
204      IF( sbc_cpl_alloc > 0 )   CALL ctl_warn('sbc_cpl_alloc: allocation of arrays failed')
205      !
206   END FUNCTION sbc_cpl_alloc
207
208
209   SUBROUTINE sbc_cpl_init( k_ice )     
210      !!----------------------------------------------------------------------
211      !!             ***  ROUTINE sbc_cpl_init  ***
212      !!
213      !! ** Purpose :   Initialisation of send and recieved information from
214      !!                the atmospheric component
215      !!
216      !! ** Method  : * Read namsbc_cpl namelist
217      !!              * define the receive interface
218      !!              * define the send    interface
219      !!              * initialise the OASIS coupler
220      !!----------------------------------------------------------------------
221      INTEGER, INTENT(in) ::   k_ice    ! ice management in the sbc (=0/1/2/3)
222      !!
223      INTEGER ::   jn   ! dummy loop index
224      INTEGER ::   ierror 
225      REAL(wp), POINTER, DIMENSION(:,:) ::   zacs, zaos
226      !!
227      NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2,   &
228         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr,   &
229         &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal   , sn_rcv_iceflx  , sn_rcv_co2
230      !!---------------------------------------------------------------------
231      !
232      IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_init')
233      !
234      CALL wrk_alloc( jpi,jpj, zacs, zaos )
235
236      ALLOCATE(srcv(nmaxfld), stat = ierror)
237      IF( ierror > 0 ) CALL ctl_stop('Pb of alloction of srcv(nmaxfld) in sbc_cpl_init')
238      ALLOCATE(ssnd(nmaxfld), stat = ierror)
239      IF( ierror > 0 ) CALL ctl_stop('Pb of alloction of ssnd(nmaxfld) in sbc_cpl_init')
240
241      ! ================================ !
242      !      Namelist informations       !
243      ! ================================ !
244
245      ! default definitions
246      !                    !     description       !  multiple  !    vector   !      vector          ! vector !
247      !                    !                       ! categories !  reference  !    orientation       ! grids  !
248      ! send
249      sn_snd_temp   = FLD_C( 'weighted oce and ice',    'no'    ,     ''      ,         ''           ,   ''   ) 
250      sn_snd_alb    = FLD_C( 'weighted ice'        ,    'no'    ,     ''      ,         ''           ,   ''   ) 
251      sn_snd_thick  = FLD_C( 'none'                ,    'no'    ,     ''      ,         ''           ,   ''   ) 
252      sn_snd_crt    = FLD_C( 'none'                ,    'no'    , 'spherical' , 'eastward-northward' ,  'T'   )     
253      sn_snd_co2    = FLD_C( 'none'                ,    'no'    ,     ''      ,         ''           ,   ''   )     
254      ! receive
255      sn_rcv_w10m   = FLD_C( 'none'                ,    'no'    ,     ''      ,         ''          ,   ''    )
256      sn_rcv_taumod = FLD_C( 'coupled'             ,    'no'    ,     ''      ,         ''          ,   ''    )
257      sn_rcv_tau    = FLD_C( 'oce only'            ,    'no'    , 'cartesian' , 'eastward-northward',  'U,V'  ) 
258      sn_rcv_dqnsdt = FLD_C( 'coupled'             ,    'no'    ,     ''      ,         ''          ,   ''    )
259      sn_rcv_qsr    = FLD_C( 'oce and ice'         ,    'no'    ,     ''      ,         ''          ,   ''    )
260      sn_rcv_qns    = FLD_C( 'oce and ice'         ,    'no'    ,     ''      ,         ''          ,   ''    )
261      sn_rcv_emp    = FLD_C( 'conservative'        ,    'no'    ,     ''      ,         ''          ,   ''    )
262      sn_rcv_rnf    = FLD_C( 'coupled'             ,    'no'    ,     ''      ,         ''          ,   ''    )
263      sn_rcv_cal    = FLD_C( 'coupled'             ,    'no'    ,     ''      ,         ''          ,   ''    )
264      sn_rcv_iceflx = FLD_C( 'none'                ,    'no'    ,     ''      ,         ''          ,   ''    )
265      sn_rcv_co2    = FLD_C( 'none'                ,    'no'    ,     ''      ,         ''          ,   ''    )
266
267      REWIND( numnam )                    ! ... read namlist namsbc_cpl
268      READ  ( numnam, namsbc_cpl )
269
270      IF(lwp) THEN                        ! control print
271         WRITE(numout,*)
272         WRITE(numout,*)'sbc_cpl_init : namsbc_cpl namelist '
273         WRITE(numout,*)'~~~~~~~~~~~~'
274         WRITE(numout,*)'  received fields (mutiple ice categogies)'
275         WRITE(numout,*)'      10m wind module                 = ', TRIM(sn_rcv_w10m%cldes  ), ' (', TRIM(sn_rcv_w10m%clcat  ), ')'
276         WRITE(numout,*)'      stress module                   = ', TRIM(sn_rcv_taumod%cldes), ' (', TRIM(sn_rcv_taumod%clcat), ')'
277         WRITE(numout,*)'      surface stress                  = ', TRIM(sn_rcv_tau%cldes   ), ' (', TRIM(sn_rcv_tau%clcat   ), ')'
278         WRITE(numout,*)'                     - referential    = ', sn_rcv_tau%clvref
279         WRITE(numout,*)'                     - orientation    = ', sn_rcv_tau%clvor
280         WRITE(numout,*)'                     - mesh           = ', sn_rcv_tau%clvgrd
281         WRITE(numout,*)'      non-solar heat flux sensitivity = ', TRIM(sn_rcv_dqnsdt%cldes), ' (', TRIM(sn_rcv_dqnsdt%clcat), ')'
282         WRITE(numout,*)'      solar heat flux                 = ', TRIM(sn_rcv_qsr%cldes   ), ' (', TRIM(sn_rcv_qsr%clcat   ), ')'
283         WRITE(numout,*)'      non-solar heat flux             = ', TRIM(sn_rcv_qns%cldes   ), ' (', TRIM(sn_rcv_qns%clcat   ), ')'
284         WRITE(numout,*)'      freshwater budget               = ', TRIM(sn_rcv_emp%cldes   ), ' (', TRIM(sn_rcv_emp%clcat   ), ')'
285         WRITE(numout,*)'      runoffs                         = ', TRIM(sn_rcv_rnf%cldes   ), ' (', TRIM(sn_rcv_rnf%clcat   ), ')'
286         WRITE(numout,*)'      calving                         = ', TRIM(sn_rcv_cal%cldes   ), ' (', TRIM(sn_rcv_cal%clcat   ), ')'
287         WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')'
288         WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')'
289         WRITE(numout,*)'  sent fields (multiple ice categories)'
290         WRITE(numout,*)'      surface temperature             = ', TRIM(sn_snd_temp%cldes  ), ' (', TRIM(sn_snd_temp%clcat  ), ')'
291         WRITE(numout,*)'      albedo                          = ', TRIM(sn_snd_alb%cldes   ), ' (', TRIM(sn_snd_alb%clcat   ), ')'
292         WRITE(numout,*)'      ice/snow thickness              = ', TRIM(sn_snd_thick%cldes ), ' (', TRIM(sn_snd_thick%clcat ), ')'
293         WRITE(numout,*)'      surface current                 = ', TRIM(sn_snd_crt%cldes   ), ' (', TRIM(sn_snd_crt%clcat   ), ')'
294         WRITE(numout,*)'                      - referential   = ', sn_snd_crt%clvref 
295         WRITE(numout,*)'                      - orientation   = ', sn_snd_crt%clvor
296         WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd
297         WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')'
298      ENDIF
299
300      !                                   ! allocate sbccpl arrays
301      IF( sbc_cpl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' )
302     
303      ! ================================ !
304      !   Define the receive interface   !
305      ! ================================ !
306      nrcvinfo(:) = OASIS_idle   ! needed by nrcvinfo(jpr_otx1) if we do not receive ocean stress
307
308      ! for each field: define the OASIS name                              (srcv(:)%clname)
309      !                 define receive or not from the namelist parameters (srcv(:)%laction)
310      !                 define the north fold type of lbc                  (srcv(:)%nsgn)
311
312      ! default definitions of srcv
313      srcv(:)%laction = .FALSE.   ;   srcv(:)%clgrid = 'T'   ;   srcv(:)%nsgn = 1.   ;   srcv(:)%nct = 1
314
315      !                                                      ! ------------------------- !
316      !                                                      ! ice and ocean wind stress !   
317      !                                                      ! ------------------------- !
318      !                                                           ! Name
319      srcv(jpr_otx1)%clname = 'O_OTaux1'      ! 1st ocean component on grid ONE (T or U)
320      srcv(jpr_oty1)%clname = 'O_OTauy1'      ! 2nd   -      -         -     -
321      srcv(jpr_otz1)%clname = 'O_OTauz1'      ! 3rd   -      -         -     -
322      srcv(jpr_otx2)%clname = 'O_OTaux2'      ! 1st ocean component on grid TWO (V)
323      srcv(jpr_oty2)%clname = 'O_OTauy2'      ! 2nd   -      -         -     -
324      srcv(jpr_otz2)%clname = 'O_OTauz2'      ! 3rd   -      -         -     -
325      !
326      srcv(jpr_itx1)%clname = 'O_ITaux1'      ! 1st  ice  component on grid ONE (T, F, I or U)
327      srcv(jpr_ity1)%clname = 'O_ITauy1'      ! 2nd   -      -         -     -
328      srcv(jpr_itz1)%clname = 'O_ITauz1'      ! 3rd   -      -         -     -
329      srcv(jpr_itx2)%clname = 'O_ITaux2'      ! 1st  ice  component on grid TWO (V)
330      srcv(jpr_ity2)%clname = 'O_ITauy2'      ! 2nd   -      -         -     -
331      srcv(jpr_itz2)%clname = 'O_ITauz2'      ! 3rd   -      -         -     -
332      !
333      ! Vectors: change of sign at north fold ONLY if on the local grid
334      IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' )   srcv(jpr_otx1:jpr_itz2)%nsgn = -1.
335     
336      !                                                           ! Set grid and action
337      SELECT CASE( TRIM( sn_rcv_tau%clvgrd ) )      !  'T', 'U,V', 'U,V,I', 'U,V,F', 'T,I', 'T,F', or 'T,U,V'
338      CASE( 'T' ) 
339         srcv(jpr_otx1:jpr_itz2)%clgrid  = 'T'        ! oce and ice components given at T-point
340         srcv(jpr_otx1:jpr_otz1)%laction = .TRUE.     ! receive oce components on grid 1
341         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1
342      CASE( 'U,V' ) 
343         srcv(jpr_otx1:jpr_otz1)%clgrid  = 'U'        ! oce components given at U-point
344         srcv(jpr_otx2:jpr_otz2)%clgrid  = 'V'        !           and           V-point
345         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'U'        ! ice components given at U-point
346         srcv(jpr_itx2:jpr_itz2)%clgrid  = 'V'        !           and           V-point
347         srcv(jpr_otx1:jpr_itz2)%laction = .TRUE.     ! receive oce and ice components on both grid 1 & 2
348      CASE( 'U,V,T' )
349         srcv(jpr_otx1:jpr_otz1)%clgrid  = 'U'        ! oce components given at U-point
350         srcv(jpr_otx2:jpr_otz2)%clgrid  = 'V'        !           and           V-point
351         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'T'        ! ice components given at T-point
352         srcv(jpr_otx1:jpr_otz2)%laction = .TRUE.     ! receive oce components on grid 1 & 2
353         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1 only
354      CASE( 'U,V,I' )
355         srcv(jpr_otx1:jpr_otz1)%clgrid  = 'U'        ! oce components given at U-point
356         srcv(jpr_otx2:jpr_otz2)%clgrid  = 'V'        !           and           V-point
357         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'I'        ! ice components given at I-point
358         srcv(jpr_otx1:jpr_otz2)%laction = .TRUE.     ! receive oce components on grid 1 & 2
359         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1 only
360      CASE( 'U,V,F' )
361         srcv(jpr_otx1:jpr_otz1)%clgrid  = 'U'        ! oce components given at U-point
362         srcv(jpr_otx2:jpr_otz2)%clgrid  = 'V'        !           and           V-point
363         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'F'        ! ice components given at F-point
364         srcv(jpr_otx1:jpr_otz2)%laction = .TRUE.     ! receive oce components on grid 1 & 2
365         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1 only
366      CASE( 'T,I' ) 
367         srcv(jpr_otx1:jpr_itz2)%clgrid  = 'T'        ! oce and ice components given at T-point
368         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'I'        ! ice components given at I-point
369         srcv(jpr_otx1:jpr_otz1)%laction = .TRUE.     ! receive oce components on grid 1
370         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1
371      CASE( 'T,F' ) 
372         srcv(jpr_otx1:jpr_itz2)%clgrid  = 'T'        ! oce and ice components given at T-point
373         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'F'        ! ice components given at F-point
374         srcv(jpr_otx1:jpr_otz1)%laction = .TRUE.     ! receive oce components on grid 1
375         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1
376      CASE( 'T,U,V' )
377         srcv(jpr_otx1:jpr_otz1)%clgrid  = 'T'        ! oce components given at T-point
378         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'U'        ! ice components given at U-point
379         srcv(jpr_itx2:jpr_itz2)%clgrid  = 'V'        !           and           V-point
380         srcv(jpr_otx1:jpr_otz1)%laction = .TRUE.     ! receive oce components on grid 1 only
381         srcv(jpr_itx1:jpr_itz2)%laction = .TRUE.     ! receive ice components on grid 1 & 2
382      CASE default   
383         CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_tau%clvgrd' )
384      END SELECT
385      !
386      IF( TRIM( sn_rcv_tau%clvref ) == 'spherical'  )   &         ! spherical: 3rd component not received
387         &     srcv( (/jpr_otz1, jpr_otz2, jpr_itz1, jpr_itz2/) )%laction = .FALSE. 
388      !
389      IF( TRIM( sn_rcv_tau%clvor  ) == 'local grid' ) THEN        ! already on local grid -> no need of the second grid
390            srcv(jpr_otx2:jpr_otz2)%laction = .FALSE. 
391            srcv(jpr_itx2:jpr_itz2)%laction = .FALSE. 
392            srcv(jpr_oty1)%clgrid = srcv(jpr_oty2)%clgrid   ! not needed but cleaner...
393            srcv(jpr_ity1)%clgrid = srcv(jpr_ity2)%clgrid   ! not needed but cleaner...
394      ENDIF
395      !
396      IF( TRIM( sn_rcv_tau%cldes ) /= 'oce and ice' ) THEN        ! 'oce and ice' case ocean stress on ocean mesh used
397         srcv(jpr_itx1:jpr_itz2)%laction = .FALSE.    ! ice components not received
398         srcv(jpr_itx1)%clgrid = 'U'                  ! ocean stress used after its transformation
399         srcv(jpr_ity1)%clgrid = 'V'                  ! i.e. it is always at U- & V-points for i- & j-comp. resp.
400      ENDIF
401       
402      !                                                      ! ------------------------- !
403      !                                                      !    freshwater budget      !   E-P
404      !                                                      ! ------------------------- !
405      ! we suppose that atmosphere modele do not make the difference between precipiration (liquide or solid)
406      ! over ice of free ocean within the same atmospheric cell.cd
407      srcv(jpr_rain)%clname = 'OTotRain'      ! Rain = liquid precipitation
408      srcv(jpr_snow)%clname = 'OTotSnow'      ! Snow = solid precipitation
409      srcv(jpr_tevp)%clname = 'OTotEvap'      ! total evaporation (over oce + ice sublimation)
410      srcv(jpr_ievp)%clname = 'OIceEvap'      ! evaporation over ice = sublimation
411      srcv(jpr_sbpr)%clname = 'OSubMPre'      ! sublimation - liquid precipitation - solid precipitation
412      srcv(jpr_semp)%clname = 'OISubMSn'      ! ice solid water budget = sublimation - solid precipitation
413      srcv(jpr_oemp)%clname = 'OOEvaMPr'      ! ocean water budget = ocean Evap - ocean precip
414      SELECT CASE( TRIM( sn_rcv_emp%cldes ) )
415      CASE( 'oce only'      )   ;   srcv(                                 jpr_oemp   )%laction = .TRUE. 
416      CASE( 'conservative'  )   ;   srcv( (/jpr_rain, jpr_snow, jpr_ievp, jpr_tevp/) )%laction = .TRUE.
417      CASE( 'oce and ice'   )   ;   srcv( (/jpr_ievp, jpr_sbpr, jpr_semp, jpr_oemp/) )%laction = .TRUE.
418      CASE default              ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_emp%cldes' )
419      END SELECT
420
421      !                                                      ! ------------------------- !
422      !                                                      !     Runoffs & Calving     !   
423      !                                                      ! ------------------------- !
424      srcv(jpr_rnf   )%clname = 'O_Runoff'   ;   IF( TRIM( sn_rcv_rnf%cldes ) == 'coupled' )   srcv(jpr_rnf)%laction = .TRUE.
425! This isn't right - really just want ln_rnf_emp changed
426!                                                 IF( TRIM( sn_rcv_rnf%cldes ) == 'climato' )   THEN   ;   ln_rnf = .TRUE.
427!                                                 ELSE                                                 ;   ln_rnf = .FALSE.
428!                                                 ENDIF
429      srcv(jpr_cal   )%clname = 'OCalving'   ;   IF( TRIM( sn_rcv_cal%cldes ) == 'coupled' )   srcv(jpr_cal)%laction = .TRUE.
430
431      !                                                      ! ------------------------- !
432      !                                                      !    non solar radiation    !   Qns
433      !                                                      ! ------------------------- !
434      srcv(jpr_qnsoce)%clname = 'O_QnsOce'
435      srcv(jpr_qnsice)%clname = 'O_QnsIce'
436      srcv(jpr_qnsmix)%clname = 'O_QnsMix'
437      SELECT CASE( TRIM( sn_rcv_qns%cldes ) )
438      CASE( 'oce only'      )   ;   srcv(               jpr_qnsoce   )%laction = .TRUE.
439      CASE( 'conservative'  )   ;   srcv( (/jpr_qnsice, jpr_qnsmix/) )%laction = .TRUE.
440      CASE( 'oce and ice'   )   ;   srcv( (/jpr_qnsice, jpr_qnsoce/) )%laction = .TRUE.
441      CASE( 'mixed oce-ice' )   ;   srcv(               jpr_qnsmix   )%laction = .TRUE. 
442      CASE default              ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_qns%cldes' )
443      END SELECT
444      IF( TRIM( sn_rcv_qns%cldes ) == 'mixed oce-ice' .AND. jpl > 1 ) &
445         CALL ctl_stop( 'sbc_cpl_init: sn_rcv_qns%cldes not currently allowed to be mixed oce-ice for multi-category ice' )
446      !                                                      ! ------------------------- !
447      !                                                      !    solar radiation        !   Qsr
448      !                                                      ! ------------------------- !
449      srcv(jpr_qsroce)%clname = 'O_QsrOce'
450      srcv(jpr_qsrice)%clname = 'O_QsrIce'
451      srcv(jpr_qsrmix)%clname = 'O_QsrMix'
452      SELECT CASE( TRIM( sn_rcv_qsr%cldes ) )
453      CASE( 'oce only'      )   ;   srcv(               jpr_qsroce   )%laction = .TRUE.
454      CASE( 'conservative'  )   ;   srcv( (/jpr_qsrice, jpr_qsrmix/) )%laction = .TRUE.
455      CASE( 'oce and ice'   )   ;   srcv( (/jpr_qsrice, jpr_qsroce/) )%laction = .TRUE.
456      CASE( 'mixed oce-ice' )   ;   srcv(               jpr_qsrmix   )%laction = .TRUE. 
457      CASE default              ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_qsr%cldes' )
458      END SELECT
459      IF( TRIM( sn_rcv_qsr%cldes ) == 'mixed oce-ice' .AND. jpl > 1 ) &
460         CALL ctl_stop( 'sbc_cpl_init: sn_rcv_qsr%cldes not currently allowed to be mixed oce-ice for multi-category ice' )
461      !                                                      ! ------------------------- !
462      !                                                      !   non solar sensitivity   !   d(Qns)/d(T)
463      !                                                      ! ------------------------- !
464      srcv(jpr_dqnsdt)%clname = 'O_dQnsdT'   
465      IF( TRIM( sn_rcv_dqnsdt%cldes ) == 'coupled' )   srcv(jpr_dqnsdt)%laction = .TRUE.
466      !
467      ! non solar sensitivity mandatory for LIM ice model
468      IF( TRIM( sn_rcv_dqnsdt%cldes ) == 'none' .AND. k_ice /= 0 .AND. k_ice /= 4) &
469         CALL ctl_stop( 'sbc_cpl_init: sn_rcv_dqnsdt%cldes must be coupled in namsbc_cpl namelist' )
470      ! non solar sensitivity mandatory for mixed oce-ice solar radiation coupling technique
471      IF( TRIM( sn_rcv_dqnsdt%cldes ) == 'none' .AND. TRIM( sn_rcv_qns%cldes ) == 'mixed oce-ice' ) &
472         CALL ctl_stop( 'sbc_cpl_init: namsbc_cpl namelist mismatch between sn_rcv_qns%cldes and sn_rcv_dqnsdt%cldes' )
473      !                                                      ! ------------------------- !
474      !                                                      !    Ice Qsr penetration    !   
475      !                                                      ! ------------------------- !
476      ! fraction of net shortwave radiation which is not absorbed in the thin surface layer
477      ! and penetrates inside the ice cover ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 )
478      ! Coupled case: since cloud cover is not received from atmosphere
479      !               ===> defined as constant value -> definition done in sbc_cpl_init
480      fr1_i0(:,:) = 0.18
481      fr2_i0(:,:) = 0.82
482      !                                                      ! ------------------------- !
483      !                                                      !      10m wind module      !   
484      !                                                      ! ------------------------- !
485      srcv(jpr_w10m)%clname = 'O_Wind10'   ;   IF( TRIM(sn_rcv_w10m%cldes  ) == 'coupled' )   srcv(jpr_w10m)%laction = .TRUE. 
486      !
487      !                                                      ! ------------------------- !
488      !                                                      !   wind stress module      !   
489      !                                                      ! ------------------------- !
490      srcv(jpr_taum)%clname = 'O_TauMod'   ;   IF( TRIM(sn_rcv_taumod%cldes) == 'coupled' )   srcv(jpr_taum)%laction = .TRUE.
491      lhftau = srcv(jpr_taum)%laction
492
493      !                                                      ! ------------------------- !
494      !                                                      !      Atmospheric CO2      !
495      !                                                      ! ------------------------- !
496      srcv(jpr_co2 )%clname = 'O_AtmCO2'   ;   IF( TRIM(sn_rcv_co2%cldes   ) == 'coupled' )    srcv(jpr_co2 )%laction = .TRUE.
497      !                                                      ! ------------------------- !
498      !                                                      !   topmelt and botmelt     !   
499      !                                                      ! ------------------------- !
500      srcv(jpr_topm )%clname = 'OTopMlt'
501      srcv(jpr_botm )%clname = 'OBotMlt'
502      IF( TRIM(sn_rcv_iceflx%cldes) == 'coupled' ) THEN
503         IF ( TRIM( sn_rcv_iceflx%clcat ) == 'yes' ) THEN
504            srcv(jpr_topm:jpr_botm)%nct = jpl
505         ELSE
506            CALL ctl_stop( 'sbc_cpl_init: sn_rcv_iceflx%clcat should always be set to yes currently' )
507         ENDIF
508         srcv(jpr_topm:jpr_botm)%laction = .TRUE.
509      ENDIF
510
511      ! Allocate all parts of frcv used for received fields
512      DO jn = 1, jprcv
513         IF ( srcv(jn)%laction ) ALLOCATE( frcv(jn)%z3(jpi,jpj,srcv(jn)%nct) )
514      END DO
515      ! Allocate taum part of frcv which is used even when not received as coupling field
516      IF ( .NOT. srcv(jpr_taum)%laction ) ALLOCATE( frcv(jpr_taum)%z3(jpi,jpj,srcv(jn)%nct) )
517
518      ! ================================ !
519      !     Define the send interface    !
520      ! ================================ !
521      ! for each field: define the OASIS name                           (ssnd(:)%clname)
522      !                 define send or not from the namelist parameters (ssnd(:)%laction)
523      !                 define the north fold type of lbc               (ssnd(:)%nsgn)
524     
525      ! default definitions of nsnd
526      ssnd(:)%laction = .FALSE.   ;   ssnd(:)%clgrid = 'T'   ;   ssnd(:)%nsgn = 1.  ; ssnd(:)%nct = 1
527         
528      !                                                      ! ------------------------- !
529      !                                                      !    Surface temperature    !
530      !                                                      ! ------------------------- !
531      ssnd(jps_toce)%clname = 'O_SSTSST'
532      ssnd(jps_tice)%clname = 'O_TepIce'
533      ssnd(jps_tmix)%clname = 'O_TepMix'
534      SELECT CASE( TRIM( sn_snd_temp%cldes ) )
535      CASE( 'none'          )       ! nothing to do
536      CASE( 'oce only'             )   ;   ssnd(   jps_toce             )%laction = .TRUE.
537      CASE( 'weighted oce and ice' )
538         ssnd( (/jps_toce, jps_tice/) )%laction = .TRUE.
539         IF ( TRIM( sn_snd_temp%clcat ) == 'yes' )  ssnd(jps_tice)%nct = jpl
540      CASE( 'mixed oce-ice'        )   ;   ssnd(   jps_tmix             )%laction = .TRUE.
541      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_temp%cldes' )
542      END SELECT
543     
544      !                                                      ! ------------------------- !
545      !                                                      !          Albedo           !
546      !                                                      ! ------------------------- !
547      ssnd(jps_albice)%clname = 'O_AlbIce' 
548      ssnd(jps_albmix)%clname = 'O_AlbMix'
549      SELECT CASE( TRIM( sn_snd_alb%cldes ) )
550      CASE( 'none'          )       ! nothing to do
551      CASE( 'weighted ice'  )   ;   ssnd(jps_albice)%laction = .TRUE.
552      CASE( 'mixed oce-ice' )   ;   ssnd(jps_albmix)%laction = .TRUE.
553      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_alb%cldes' )
554      END SELECT
555      !
556      ! Need to calculate oceanic albedo if
557      !     1. sending mixed oce-ice albedo or
558      !     2. receiving mixed oce-ice solar radiation
559      IF ( TRIM ( sn_snd_alb%cldes ) == 'mixed oce-ice' .OR. TRIM ( sn_rcv_qsr%cldes ) == 'mixed oce-ice' ) THEN
560         CALL albedo_oce( zaos, zacs )
561         ! Due to lack of information on nebulosity : mean clear/overcast sky
562         albedo_oce_mix(:,:) = ( zacs(:,:) + zaos(:,:) ) * 0.5
563      ENDIF
564
565      !                                                      ! ------------------------- !
566      !                                                      !  Ice fraction & Thickness !
567      !                                                      ! ------------------------- !
568      ssnd(jps_fice)%clname = 'OIceFrc'
569      ssnd(jps_hice)%clname = 'OIceTck'
570      ssnd(jps_hsnw)%clname = 'OSnwTck'
571      IF( k_ice /= 0 ) THEN
572         ssnd(jps_fice)%laction = .TRUE.                  ! if ice treated in the ocean (even in climato case)
573! Currently no namelist entry to determine sending of multi-category ice fraction so use the thickness entry for now
574         IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_fice)%nct = jpl
575      ENDIF
576
577      SELECT CASE ( TRIM( sn_snd_thick%cldes ) )
578      CASE( 'none'         )       ! nothing to do
579      CASE( 'ice and snow' ) 
580         ssnd(jps_hice:jps_hsnw)%laction = .TRUE.
581         IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) THEN
582            ssnd(jps_hice:jps_hsnw)%nct = jpl
583         ELSE
584            IF ( jpl > 1 ) THEN
585CALL ctl_stop( 'sbc_cpl_init: use weighted ice and snow option for sn_snd_thick%cldes if not exchanging category fields' )
586            ENDIF
587         ENDIF
588      CASE( 'weighted ice and snow' ) 
589         ssnd(jps_hice:jps_hsnw)%laction = .TRUE.
590         IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_hice:jps_hsnw)%nct = jpl
591      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_thick%cldes' )
592      END SELECT
593
594      !                                                      ! ------------------------- !
595      !                                                      !      Surface current      !
596      !                                                      ! ------------------------- !
597      !        ocean currents              !            ice velocities
598      ssnd(jps_ocx1)%clname = 'O_OCurx1'   ;   ssnd(jps_ivx1)%clname = 'O_IVelx1'
599      ssnd(jps_ocy1)%clname = 'O_OCury1'   ;   ssnd(jps_ivy1)%clname = 'O_IVely1'
600      ssnd(jps_ocz1)%clname = 'O_OCurz1'   ;   ssnd(jps_ivz1)%clname = 'O_IVelz1'
601      !
602      ssnd(jps_ocx1:jps_ivz1)%nsgn = -1.   ! vectors: change of the sign at the north fold
603
604      IF( sn_snd_crt%clvgrd == 'U,V' ) THEN
605         ssnd(jps_ocx1)%clgrid = 'U' ; ssnd(jps_ocy1)%clgrid = 'V'
606      ELSE IF( sn_snd_crt%clvgrd /= 'T' ) THEN 
607         CALL ctl_stop( 'sn_snd_crt%clvgrd must be equal to T' )
608         ssnd(jps_ocx1:jps_ivz1)%clgrid  = 'T'      ! all oce and ice components on the same unique grid
609      ENDIF
610      ssnd(jps_ocx1:jps_ivz1)%laction = .TRUE.   ! default: all are send
611      IF( TRIM( sn_snd_crt%clvref ) == 'spherical' )   ssnd( (/jps_ocz1, jps_ivz1/) )%laction = .FALSE. 
612      IF( TRIM( sn_snd_crt%clvor ) == 'eastward-northward' ) ssnd(jps_ocx1:jps_ivz1)%nsgn = 1.
613      SELECT CASE( TRIM( sn_snd_crt%cldes ) )
614      CASE( 'none'                 )   ;   ssnd(jps_ocx1:jps_ivz1)%laction = .FALSE.
615      CASE( 'oce only'             )   ;   ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.
616      CASE( 'weighted oce and ice' )   !   nothing to do
617      CASE( 'mixed oce-ice'        )   ;   ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.
618      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crt%cldes' )
619      END SELECT
620
621      !                                                      ! ------------------------- !
622      !                                                      !          CO2 flux         !
623      !                                                      ! ------------------------- !
624      ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE.
625      !
626      ! ================================ !
627      !   initialisation of the coupler  !
628      ! ================================ !
629
630      CALL cpl_prism_define(jprcv, jpsnd)           
631      !
632      IF( ln_dm2dc .AND. ( cpl_prism_freq( jpr_qsroce ) + cpl_prism_freq( jpr_qsrmix ) /= 86400 ) )   &
633         &   CALL ctl_stop( 'sbc_cpl_init: diurnal cycle reconstruction (ln_dm2dc) needs daily couping for solar radiation' )
634
635      CALL wrk_dealloc( jpi,jpj, zacs, zaos )
636      !
637      IF( nn_timing == 1 )  CALL timing_stop('sbc_cpl_init')
638      !
639   END SUBROUTINE sbc_cpl_init
640
641
642   SUBROUTINE sbc_cpl_rcv( kt, k_fsbc, k_ice )     
643      !!----------------------------------------------------------------------
644      !!             ***  ROUTINE sbc_cpl_rcv  ***
645      !!
646      !! ** Purpose :   provide the stress over the ocean and, if no sea-ice,
647      !!                provide the ocean heat and freshwater fluxes.
648      !!
649      !! ** Method  : - Receive all the atmospheric fields (stored in frcv array). called at each time step.
650      !!                OASIS controls if there is something do receive or not. nrcvinfo contains the info
651      !!                to know if the field was really received or not
652      !!
653      !!              --> If ocean stress was really received:
654      !!
655      !!                  - transform the received ocean stress vector from the received
656      !!                 referential and grid into an atmosphere-ocean stress in
657      !!                 the (i,j) ocean referencial and at the ocean velocity point.
658      !!                    The received stress are :
659      !!                     - defined by 3 components (if cartesian coordinate)
660      !!                            or by 2 components (if spherical)
661      !!                     - oriented along geographical   coordinate (if eastward-northward)
662      !!                            or  along the local grid coordinate (if local grid)
663      !!                     - given at U- and V-point, resp.   if received on 2 grids
664      !!                            or at T-point               if received on 1 grid
665      !!                    Therefore and if necessary, they are successively
666      !!                  processed in order to obtain them
667      !!                     first  as  2 components on the sphere
668      !!                     second as  2 components oriented along the local grid
669      !!                     third  as  2 components on the U,V grid
670      !!
671      !!              -->
672      !!
673      !!              - In 'ocean only' case, non solar and solar ocean heat fluxes
674      !!             and total ocean freshwater fluxes 
675      !!
676      !! ** Method  :   receive all fields from the atmosphere and transform
677      !!              them into ocean surface boundary condition fields
678      !!
679      !! ** Action  :   update  utau, vtau   ocean stress at U,V grid
680      !!                        taum, wndm   wind stres and wind speed module at T-point
681      !!                        qns , qsr    non solar and solar ocean heat fluxes   ('ocean only case)
682      !!                        emp = emps   evap. - precip. (- runoffs) (- calving) ('ocean only case)
683      !!----------------------------------------------------------------------
684      INTEGER, INTENT(in) ::   kt       ! ocean model time step index
685      INTEGER, INTENT(in) ::   k_fsbc   ! frequency of sbc (-> ice model) computation
686      INTEGER, INTENT(in) ::   k_ice    ! ice management in the sbc (=0/1/2/3)
687      !!
688      LOGICAL ::    llnewtx, llnewtau      ! update wind stress components and module??
689      INTEGER  ::   ji, jj, jn             ! dummy loop indices
690      INTEGER  ::   isec                   ! number of seconds since nit000 (assuming rdttra did not change since nit000)
691      REAL(wp) ::   zcumulneg, zcumulpos   ! temporary scalars     
692      REAL(wp) ::   zcoef                  ! temporary scalar
693      REAL(wp) ::   zrhoa  = 1.22          ! Air density kg/m3
694      REAL(wp) ::   zcdrag = 1.5e-3        ! drag coefficient
695      REAL(wp) ::   zzx, zzy               ! temporary variables
696      REAL(wp), POINTER, DIMENSION(:,:) ::   ztx, zty 
697      !!----------------------------------------------------------------------
698      !
699      IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_rcv')
700      !
701      CALL wrk_alloc( jpi,jpj, ztx, zty )
702
703!EM XIOS/OASIS-MCT compliance
704!EM      IF( kt == nit000 )   CALL sbc_cpl_init( k_ice )          ! initialisation
705
706      !                                                 ! Receive all the atmos. fields (including ice information)
707      isec = ( kt - nit000 ) * NINT( rdttra(1) )             ! date of exchanges
708      DO jn = 1, jprcv                                       ! received fields sent by the atmosphere
709         IF( srcv(jn)%laction )   CALL cpl_prism_rcv( jn, isec, frcv(jn)%z3, nrcvinfo(jn) )
710      END DO
711
712      !                                                      ! ========================= !
713      IF( srcv(jpr_otx1)%laction ) THEN                      !  ocean stress components  !
714         !                                                   ! ========================= !
715         ! define frcv(jpr_otx1)%z3(:,:,1) and frcv(jpr_oty1)%z3(:,:,1): stress at U/V point along model grid
716         ! => need to be done only when we receive the field
717         IF(  nrcvinfo(jpr_otx1) == OASIS_Rcv ) THEN
718            !
719            IF( TRIM( sn_rcv_tau%clvref ) == 'cartesian' ) THEN            ! 2 components on the sphere
720               !                                                       ! (cartesian to spherical -> 3 to 2 components)
721               !
722               CALL geo2oce( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), frcv(jpr_otz1)%z3(:,:,1),   &
723                  &          srcv(jpr_otx1)%clgrid, ztx, zty )
724               frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:)   ! overwrite 1st comp. on the 1st grid
725               frcv(jpr_oty1)%z3(:,:,1) = zty(:,:)   ! overwrite 2nd comp. on the 1st grid
726               !
727               IF( srcv(jpr_otx2)%laction ) THEN
728                  CALL geo2oce( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), frcv(jpr_otz2)%z3(:,:,1),   &
729                     &          srcv(jpr_otx2)%clgrid, ztx, zty )
730                  frcv(jpr_otx2)%z3(:,:,1) = ztx(:,:)   ! overwrite 1st comp. on the 2nd grid
731                  frcv(jpr_oty2)%z3(:,:,1) = zty(:,:)   ! overwrite 2nd comp. on the 2nd grid
732               ENDIF
733               !
734            ENDIF
735            !
736            IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN   ! 2 components oriented along the local grid
737               !                                                       ! (geographical to local grid -> rotate the components)
738               CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx )   
739               IF( srcv(jpr_otx2)%laction ) THEN
740                  CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty )   
741               ELSE     
742                  CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty ) 
743               ENDIF
744               frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid
745               frcv(jpr_oty1)%z3(:,:,1) = zty(:,:)      ! overwrite 2nd component on the 2nd grid
746            ENDIF
747            !                             
748            IF( srcv(jpr_otx1)%clgrid == 'T' ) THEN
749               DO jj = 2, jpjm1                                          ! T ==> (U,V)
750                  DO ji = fs_2, fs_jpim1   ! vector opt.
751                     frcv(jpr_otx1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_otx1)%z3(ji+1,jj  ,1) + frcv(jpr_otx1)%z3(ji,jj,1) )
752                     frcv(jpr_oty1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_oty1)%z3(ji  ,jj+1,1) + frcv(jpr_oty1)%z3(ji,jj,1) )
753                  END DO
754               END DO
755               CALL lbc_lnk( frcv(jpr_otx1)%z3(:,:,1), 'U',  -1. )   ;   CALL lbc_lnk( frcv(jpr_oty1)%z3(:,:,1), 'V',  -1. )
756            ENDIF
757            llnewtx = .TRUE.
758         ELSE
759            llnewtx = .FALSE.
760         ENDIF
761         !                                                   ! ========================= !
762      ELSE                                                   !   No dynamical coupling   !
763         !                                                   ! ========================= !
764         frcv(jpr_otx1)%z3(:,:,1) = 0.e0                               ! here simply set to zero
765         frcv(jpr_oty1)%z3(:,:,1) = 0.e0                               ! an external read in a file can be added instead
766         llnewtx = .TRUE.
767         !
768      ENDIF
769     
770      !                                                      ! ========================= !
771      !                                                      !    wind stress module     !   (taum)
772      !                                                      ! ========================= !
773      !
774      IF( .NOT. srcv(jpr_taum)%laction ) THEN                    ! compute wind stress module from its components if not received
775         ! => need to be done only when otx1 was changed
776         IF( llnewtx ) THEN
777!CDIR NOVERRCHK
778            DO jj = 2, jpjm1
779!CDIR NOVERRCHK
780               DO ji = fs_2, fs_jpim1   ! vect. opt.
781                  zzx = frcv(jpr_otx1)%z3(ji-1,jj  ,1) + frcv(jpr_otx1)%z3(ji,jj,1)
782                  zzy = frcv(jpr_oty1)%z3(ji  ,jj-1,1) + frcv(jpr_oty1)%z3(ji,jj,1)
783                  frcv(jpr_taum)%z3(ji,jj,1) = 0.5 * SQRT( zzx * zzx + zzy * zzy )
784               END DO
785            END DO
786            CALL lbc_lnk( frcv(jpr_taum)%z3(:,:,1), 'T', 1. )
787            llnewtau = .TRUE.
788         ELSE
789            llnewtau = .FALSE.
790         ENDIF
791      ELSE
792         llnewtau = nrcvinfo(jpr_taum) == OASIS_Rcv
793         ! Stress module can be negative when received (interpolation problem)
794         IF( llnewtau ) THEN
795            frcv(jpr_taum)%z3(:,:,1) = MAX( 0.0e0, frcv(jpr_taum)%z3(:,:,1) )
796         ENDIF
797      ENDIF
798     
799      !                                                      ! ========================= !
800      !                                                      !      10 m wind speed      !   (wndm)
801      !                                                      ! ========================= !
802      !
803      IF( .NOT. srcv(jpr_w10m)%laction ) THEN                    ! compute wind spreed from wind stress module if not received 
804         ! => need to be done only when taumod was changed
805         IF( llnewtau ) THEN
806            zcoef = 1. / ( zrhoa * zcdrag ) 
807!CDIR NOVERRCHK
808            DO jj = 1, jpj
809!CDIR NOVERRCHK
810               DO ji = 1, jpi 
811                  wndm(ji,jj) = SQRT( frcv(jpr_taum)%z3(ji,jj,1) * zcoef )
812               END DO
813            END DO
814         ENDIF
815      ELSE
816         IF ( nrcvinfo(jpr_w10m) == OASIS_Rcv ) wndm(:,:) = frcv(jpr_w10m)%z3(:,:,1)
817      ENDIF
818
819      ! u(v)tau and taum will be modified by ice model
820      ! -> need to be reset before each call of the ice/fsbc     
821      IF( MOD( kt-1, k_fsbc ) == 0 ) THEN
822         !
823         utau(:,:) = frcv(jpr_otx1)%z3(:,:,1)
824         vtau(:,:) = frcv(jpr_oty1)%z3(:,:,1)
825         taum(:,:) = frcv(jpr_taum)%z3(:,:,1)
826         CALL iom_put( "taum_oce", taum )   ! output wind stress module
827         
828      ENDIF
829
830#if defined key_cpl_carbon_cycle
831      !                                                              ! atmosph. CO2 (ppm)
832      IF( srcv(jpr_co2)%laction )   atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1)
833#endif
834
835      !                                                      ! ========================= !
836      IF( k_ice <= 1 ) THEN                                  !  heat & freshwater fluxes ! (Ocean only case)
837         !                                                   ! ========================= !
838         !
839         !                                                       ! non solar heat flux over the ocean (qns)
840         IF( srcv(jpr_qnsoce)%laction )   qns(:,:) = frcv(jpr_qnsoce)%z3(:,:,1)
841         IF( srcv(jpr_qnsmix)%laction )   qns(:,:) = frcv(jpr_qnsmix)%z3(:,:,1)
842         ! add the latent heat of solid precip. melting
843         IF( srcv(jpr_snow  )%laction )   qns(:,:) = qns(:,:) - frcv(jpr_snow)%z3(:,:,1) * lfus             
844
845         !                                                       ! solar flux over the ocean          (qsr)
846         IF( srcv(jpr_qsroce)%laction )   qsr(:,:) = frcv(jpr_qsroce)%z3(:,:,1)
847         IF( srcv(jpr_qsrmix)%laction )   qsr(:,:) = frcv(jpr_qsrmix)%z3(:,:,1)
848         IF( ln_dm2dc )   qsr(:,:) = sbc_dcy( qsr )                           ! modify qsr to include the diurnal cycle
849         !
850         !                                                       ! total freshwater fluxes over the ocean (emp, emps)
851         SELECT CASE( TRIM( sn_rcv_emp%cldes ) )                                    ! evaporation - precipitation
852         CASE( 'conservative' )
853            emp(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ( frcv(jpr_rain)%z3(:,:,1) + frcv(jpr_snow)%z3(:,:,1) )
854         CASE( 'oce only', 'oce and ice' )
855            emp(:,:) = frcv(jpr_oemp)%z3(:,:,1)
856         CASE default
857            CALL ctl_stop( 'sbc_cpl_rcv: wrong definition of sn_rcv_emp%cldes' )
858         END SELECT
859         !
860         !                                                        ! runoffs and calving (added in emp)
861         IF( srcv(jpr_rnf)%laction )   emp(:,:) = emp(:,:) - frcv(jpr_rnf)%z3(:,:,1)
862         IF( srcv(jpr_cal)%laction )   emp(:,:) = emp(:,:) - frcv(jpr_cal)%z3(:,:,1)
863         !
864!!gm :  this seems to be internal cooking, not sure to need that in a generic interface
865!!gm                                       at least should be optional...
866!!         IF( TRIM( sn_rcv_rnf%cldes ) == 'coupled' ) THEN     ! add to the total freshwater budget
867!!            ! remove negative runoff
868!!            zcumulpos = SUM( MAX( frcv(jpr_rnf)%z3(:,:,1), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) )
869!!            zcumulneg = SUM( MIN( frcv(jpr_rnf)%z3(:,:,1), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) )
870!!            IF( lk_mpp )   CALL mpp_sum( zcumulpos )   ! sum over the global domain
871!!            IF( lk_mpp )   CALL mpp_sum( zcumulneg )
872!!            IF( zcumulpos /= 0. ) THEN                 ! distribute negative runoff on positive runoff grid points
873!!               zcumulneg = 1.e0 + zcumulneg / zcumulpos
874!!               frcv(jpr_rnf)%z3(:,:,1) = MAX( frcv(jpr_rnf)%z3(:,:,1), 0.e0 ) * zcumulneg
875!!            ENDIF     
876!!            ! add runoff to e-p
877!!            emp(:,:) = emp(:,:) - frcv(jpr_rnf)%z3(:,:,1)
878!!         ENDIF
879!!gm  end of internal cooking
880         !
881         emps(:,:) = emp(:,:)                                        ! concentration/dilution = emp
882 
883      ENDIF
884      !
885      CALL wrk_dealloc( jpi,jpj, ztx, zty )
886      !
887      IF( nn_timing == 1 )  CALL timing_stop('sbc_cpl_rcv')
888      !
889   END SUBROUTINE sbc_cpl_rcv
890   
891
892   SUBROUTINE sbc_cpl_ice_tau( p_taui, p_tauj )     
893      !!----------------------------------------------------------------------
894      !!             ***  ROUTINE sbc_cpl_ice_tau  ***
895      !!
896      !! ** Purpose :   provide the stress over sea-ice in coupled mode
897      !!
898      !! ** Method  :   transform the received stress from the atmosphere into
899      !!             an atmosphere-ice stress in the (i,j) ocean referencial
900      !!             and at the velocity point of the sea-ice model (cp_ice_msh):
901      !!                'C'-grid : i- (j-) components given at U- (V-) point
902      !!                'I'-grid : B-grid lower-left corner: both components given at I-point
903      !!
904      !!                The received stress are :
905      !!                 - defined by 3 components (if cartesian coordinate)
906      !!                        or by 2 components (if spherical)
907      !!                 - oriented along geographical   coordinate (if eastward-northward)
908      !!                        or  along the local grid coordinate (if local grid)
909      !!                 - given at U- and V-point, resp.   if received on 2 grids
910      !!                        or at a same point (T or I) if received on 1 grid
911      !!                Therefore and if necessary, they are successively
912      !!             processed in order to obtain them
913      !!                 first  as  2 components on the sphere
914      !!                 second as  2 components oriented along the local grid
915      !!                 third  as  2 components on the cp_ice_msh point
916      !!
917      !!                In 'oce and ice' case, only one vector stress field
918      !!             is received. It has already been processed in sbc_cpl_rcv
919      !!             so that it is now defined as (i,j) components given at U-
920      !!             and V-points, respectively. Therefore, here only the third
921      !!             transformation is done and only if the ice-grid is a 'I'-grid.
922      !!
923      !! ** Action  :   return ptau_i, ptau_j, the stress over the ice at cp_ice_msh point
924      !!----------------------------------------------------------------------
925      REAL(wp), INTENT(out), DIMENSION(:,:) ::   p_taui   ! i- & j-components of atmos-ice stress [N/m2]
926      REAL(wp), INTENT(out), DIMENSION(:,:) ::   p_tauj   ! at I-point (B-grid) or U & V-point (C-grid)
927      !!
928      INTEGER ::   ji, jj                          ! dummy loop indices
929      INTEGER ::   itx                             ! index of taux over ice
930      REAL(wp), POINTER, DIMENSION(:,:) ::   ztx, zty 
931      !!----------------------------------------------------------------------
932      !
933      IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_ice_tau')
934      !
935      CALL wrk_alloc( jpi,jpj, ztx, zty )
936
937      IF( srcv(jpr_itx1)%laction ) THEN   ;   itx =  jpr_itx1   
938      ELSE                                ;   itx =  jpr_otx1
939      ENDIF
940
941      ! do something only if we just received the stress from atmosphere
942      IF(  nrcvinfo(itx) == OASIS_Rcv ) THEN
943
944         !                                                      ! ======================= !
945         IF( srcv(jpr_itx1)%laction ) THEN                      !   ice stress received   !
946            !                                                   ! ======================= !
947           
948            IF( TRIM( sn_rcv_tau%clvref ) == 'cartesian' ) THEN            ! 2 components on the sphere
949               !                                                       ! (cartesian to spherical -> 3 to 2 components)
950               CALL geo2oce(  frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), frcv(jpr_itz1)%z3(:,:,1),   &
951                  &          srcv(jpr_itx1)%clgrid, ztx, zty )
952               frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:)   ! overwrite 1st comp. on the 1st grid
953               frcv(jpr_ity1)%z3(:,:,1) = zty(:,:)   ! overwrite 2nd comp. on the 1st grid
954               !
955               IF( srcv(jpr_itx2)%laction ) THEN
956                  CALL geo2oce( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), frcv(jpr_itz2)%z3(:,:,1),   &
957                     &          srcv(jpr_itx2)%clgrid, ztx, zty )
958                  frcv(jpr_itx2)%z3(:,:,1) = ztx(:,:)   ! overwrite 1st comp. on the 2nd grid
959                  frcv(jpr_ity2)%z3(:,:,1) = zty(:,:)   ! overwrite 2nd comp. on the 2nd grid
960               ENDIF
961               !
962            ENDIF
963            !
964            IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN   ! 2 components oriented along the local grid
965               !                                                       ! (geographical to local grid -> rotate the components)
966               CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->i', ztx )   
967               IF( srcv(jpr_itx2)%laction ) THEN
968                  CALL rot_rep( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), srcv(jpr_itx2)%clgrid, 'en->j', zty )   
969               ELSE
970                  CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->j', zty ) 
971               ENDIF
972               frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid
973               frcv(jpr_ity1)%z3(:,:,1) = zty(:,:)      ! overwrite 2nd component on the 1st grid
974            ENDIF
975            !                                                   ! ======================= !
976         ELSE                                                   !     use ocean stress    !
977            !                                                   ! ======================= !
978            frcv(jpr_itx1)%z3(:,:,1) = frcv(jpr_otx1)%z3(:,:,1)
979            frcv(jpr_ity1)%z3(:,:,1) = frcv(jpr_oty1)%z3(:,:,1)
980            !
981         ENDIF
982
983         !                                                      ! ======================= !
984         !                                                      !     put on ice grid     !
985         !                                                      ! ======================= !
986         !   
987         !                                                  j+1   j     -----V---F
988         ! ice stress on ice velocity point (cp_ice_msh)                 !       |
989         ! (C-grid ==>(U,V) or B-grid ==> I or F)                 j      |   T   U
990         !                                                               |       |
991         !                                                   j    j-1   -I-------|
992         !                                               (for I)         |       |
993         !                                                              i-1  i   i
994         !                                                               i      i+1 (for I)
995         SELECT CASE ( cp_ice_msh )
996            !
997         CASE( 'I' )                                         ! B-grid ==> I
998            SELECT CASE ( srcv(jpr_itx1)%clgrid )
999            CASE( 'U' )
1000               DO jj = 2, jpjm1                                   ! (U,V) ==> I
1001                  DO ji = 2, jpim1   ! NO vector opt.
1002                     p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji-1,jj  ,1) + frcv(jpr_itx1)%z3(ji-1,jj-1,1) )
1003                     p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji  ,jj-1,1) + frcv(jpr_ity1)%z3(ji-1,jj-1,1) )
1004                  END DO
1005               END DO
1006            CASE( 'F' )
1007               DO jj = 2, jpjm1                                   ! F ==> I
1008                  DO ji = 2, jpim1   ! NO vector opt.
1009                     p_taui(ji,jj) = frcv(jpr_itx1)%z3(ji-1,jj-1,1)
1010                     p_tauj(ji,jj) = frcv(jpr_ity1)%z3(ji-1,jj-1,1)
1011                  END DO
1012               END DO
1013            CASE( 'T' )
1014               DO jj = 2, jpjm1                                   ! T ==> I
1015                  DO ji = 2, jpim1   ! NO vector opt.
1016                     p_taui(ji,jj) = 0.25 * ( frcv(jpr_itx1)%z3(ji,jj  ,1) + frcv(jpr_itx1)%z3(ji-1,jj  ,1)   &
1017                        &                   + frcv(jpr_itx1)%z3(ji,jj-1,1) + frcv(jpr_itx1)%z3(ji-1,jj-1,1) ) 
1018                     p_tauj(ji,jj) = 0.25 * ( frcv(jpr_ity1)%z3(ji,jj  ,1) + frcv(jpr_ity1)%z3(ji-1,jj  ,1)   &
1019                        &                   + frcv(jpr_oty1)%z3(ji,jj-1,1) + frcv(jpr_ity1)%z3(ji-1,jj-1,1) )
1020                  END DO
1021               END DO
1022            CASE( 'I' )
1023               p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1)                   ! I ==> I
1024               p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1)
1025            END SELECT
1026            IF( srcv(jpr_itx1)%clgrid /= 'I' ) THEN
1027               CALL lbc_lnk( p_taui, 'I',  -1. )   ;   CALL lbc_lnk( p_tauj, 'I',  -1. )
1028            ENDIF
1029            !
1030         CASE( 'F' )                                         ! B-grid ==> F
1031            SELECT CASE ( srcv(jpr_itx1)%clgrid )
1032            CASE( 'U' )
1033               DO jj = 2, jpjm1                                   ! (U,V) ==> F
1034                  DO ji = fs_2, fs_jpim1   ! vector opt.
1035                     p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji,jj,1) + frcv(jpr_itx1)%z3(ji  ,jj+1,1) )
1036                     p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji,jj,1) + frcv(jpr_ity1)%z3(ji+1,jj  ,1) )
1037                  END DO
1038               END DO
1039            CASE( 'I' )
1040               DO jj = 2, jpjm1                                   ! I ==> F
1041                  DO ji = 2, jpim1   ! NO vector opt.
1042                     p_taui(ji,jj) = frcv(jpr_itx1)%z3(ji+1,jj+1,1)
1043                     p_tauj(ji,jj) = frcv(jpr_ity1)%z3(ji+1,jj+1,1)
1044                  END DO
1045               END DO
1046            CASE( 'T' )
1047               DO jj = 2, jpjm1                                   ! T ==> F
1048                  DO ji = 2, jpim1   ! NO vector opt.
1049                     p_taui(ji,jj) = 0.25 * ( frcv(jpr_itx1)%z3(ji,jj  ,1) + frcv(jpr_itx1)%z3(ji+1,jj  ,1)   &
1050                        &                   + frcv(jpr_itx1)%z3(ji,jj+1,1) + frcv(jpr_itx1)%z3(ji+1,jj+1,1) ) 
1051                     p_tauj(ji,jj) = 0.25 * ( frcv(jpr_ity1)%z3(ji,jj  ,1) + frcv(jpr_ity1)%z3(ji+1,jj  ,1)   &
1052                        &                   + frcv(jpr_ity1)%z3(ji,jj+1,1) + frcv(jpr_ity1)%z3(ji+1,jj+1,1) )
1053                  END DO
1054               END DO
1055            CASE( 'F' )
1056               p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1)                   ! F ==> F
1057               p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1)
1058            END SELECT
1059            IF( srcv(jpr_itx1)%clgrid /= 'F' ) THEN
1060               CALL lbc_lnk( p_taui, 'F',  -1. )   ;   CALL lbc_lnk( p_tauj, 'F',  -1. )
1061            ENDIF
1062            !
1063         CASE( 'C' )                                         ! C-grid ==> U,V
1064            SELECT CASE ( srcv(jpr_itx1)%clgrid )
1065            CASE( 'U' )
1066               p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1)                   ! (U,V) ==> (U,V)
1067               p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1)
1068            CASE( 'F' )
1069               DO jj = 2, jpjm1                                   ! F ==> (U,V)
1070                  DO ji = fs_2, fs_jpim1   ! vector opt.
1071                     p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji,jj,1) + frcv(jpr_itx1)%z3(ji  ,jj-1,1) )
1072                     p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(jj,jj,1) + frcv(jpr_ity1)%z3(ji-1,jj  ,1) )
1073                  END DO
1074               END DO
1075            CASE( 'T' )
1076               DO jj = 2, jpjm1                                   ! T ==> (U,V)
1077                  DO ji = fs_2, fs_jpim1   ! vector opt.
1078                     p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji+1,jj  ,1) + frcv(jpr_itx1)%z3(ji,jj,1) )
1079                     p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji  ,jj+1,1) + frcv(jpr_ity1)%z3(ji,jj,1) )
1080                  END DO
1081               END DO
1082            CASE( 'I' )
1083               DO jj = 2, jpjm1                                   ! I ==> (U,V)
1084                  DO ji = 2, jpim1   ! NO vector opt.
1085                     p_taui(ji,jj) = 0.5 * ( frcv(jpr_itx1)%z3(ji+1,jj+1,1) + frcv(jpr_itx1)%z3(ji+1,jj  ,1) )
1086                     p_tauj(ji,jj) = 0.5 * ( frcv(jpr_ity1)%z3(ji+1,jj+1,1) + frcv(jpr_ity1)%z3(ji  ,jj+1,1) )
1087                  END DO
1088               END DO
1089            END SELECT
1090            IF( srcv(jpr_itx1)%clgrid /= 'U' ) THEN
1091               CALL lbc_lnk( p_taui, 'U',  -1. )   ;   CALL lbc_lnk( p_tauj, 'V',  -1. )
1092            ENDIF
1093         END SELECT
1094
1095      ENDIF
1096      !   
1097      CALL wrk_dealloc( jpi,jpj, ztx, zty )
1098      !
1099      IF( nn_timing == 1 )  CALL timing_stop('sbc_cpl_ice_tau')
1100      !
1101   END SUBROUTINE sbc_cpl_ice_tau
1102   
1103
1104   SUBROUTINE sbc_cpl_ice_flx( p_frld  , palbi   , psst    , pist    )
1105      !!----------------------------------------------------------------------
1106      !!             ***  ROUTINE sbc_cpl_ice_flx  ***
1107      !!
1108      !! ** Purpose :   provide the heat and freshwater fluxes of the
1109      !!              ocean-ice system.
1110      !!
1111      !! ** Method  :   transform the fields received from the atmosphere into
1112      !!             surface heat and fresh water boundary condition for the
1113      !!             ice-ocean system. The following fields are provided:
1114      !!              * total non solar, solar and freshwater fluxes (qns_tot,
1115      !!             qsr_tot and emp_tot) (total means weighted ice-ocean flux)
1116      !!             NB: emp_tot include runoffs and calving.
1117      !!              * fluxes over ice (qns_ice, qsr_ice, emp_ice) where
1118      !!             emp_ice = sublimation - solid precipitation as liquid
1119      !!             precipitation are re-routed directly to the ocean and
1120      !!             runoffs and calving directly enter the ocean.
1121      !!              * solid precipitation (sprecip), used to add to qns_tot
1122      !!             the heat lost associated to melting solid precipitation
1123      !!             over the ocean fraction.
1124      !!       ===>> CAUTION here this changes the net heat flux received from
1125      !!             the atmosphere
1126      !!
1127      !!                  - the fluxes have been separated from the stress as
1128      !!                 (a) they are updated at each ice time step compare to
1129      !!                 an update at each coupled time step for the stress, and
1130      !!                 (b) the conservative computation of the fluxes over the
1131      !!                 sea-ice area requires the knowledge of the ice fraction
1132      !!                 after the ice advection and before the ice thermodynamics,
1133      !!                 so that the stress is updated before the ice dynamics
1134      !!                 while the fluxes are updated after it.
1135      !!
1136      !! ** Action  :   update at each nf_ice time step:
1137      !!                   qns_tot, qsr_tot  non-solar and solar total heat fluxes
1138      !!                   qns_ice, qsr_ice  non-solar and solar heat fluxes over the ice
1139      !!                   emp_tot            total evaporation - precipitation(liquid and solid) (-runoff)(-calving)
1140      !!                   emp_ice            ice sublimation - solid precipitation over the ice
1141      !!                   dqns_ice           d(non-solar heat flux)/d(Temperature) over the ice
1142      !!                   sprecip             solid precipitation over the ocean 
1143      !!----------------------------------------------------------------------
1144      REAL(wp), INTENT(in   ), DIMENSION(:,:)   ::   p_frld     ! lead fraction                [0 to 1]
1145      ! optional arguments, used only in 'mixed oce-ice' case
1146      REAL(wp), INTENT(in   ), DIMENSION(:,:,:), OPTIONAL ::   palbi   ! ice albedo
1147      REAL(wp), INTENT(in   ), DIMENSION(:,:  ), OPTIONAL ::   psst    ! sea surface temperature     [Celcius]
1148      REAL(wp), INTENT(in   ), DIMENSION(:,:,:), OPTIONAL ::   pist    ! ice surface temperature     [Kelvin]
1149      !
1150      INTEGER ::   jl   ! dummy loop index
1151      REAL(wp), POINTER, DIMENSION(:,:) ::   zcptn, ztmp, zicefr
1152      !!----------------------------------------------------------------------
1153      !
1154      IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_ice_flx')
1155      !
1156      CALL wrk_alloc( jpi,jpj, zcptn, ztmp, zicefr )
1157
1158      zicefr(:,:) = 1.- p_frld(:,:)
1159      IF( lk_diaar5 )   zcptn(:,:) = rcp * tsn(:,:,1,jp_tem)
1160      !
1161      !                                                      ! ========================= !
1162      !                                                      !    freshwater budget      !   (emp)
1163      !                                                      ! ========================= !
1164      !
1165      !                                                           ! total Precipitations - total Evaporation (emp_tot)
1166      !                                                           ! solid precipitation  - sublimation       (emp_ice)
1167      !                                                           ! solid Precipitation                      (sprecip)
1168      SELECT CASE( TRIM( sn_rcv_emp%cldes ) )
1169      CASE( 'conservative'  )   ! received fields: jpr_rain, jpr_snow, jpr_ievp, jpr_tevp
1170         sprecip(:,:) = frcv(jpr_snow)%z3(:,:,1)                 ! May need to ensure positive here
1171         tprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + sprecip (:,:) ! May need to ensure positive here
1172         emp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - tprecip(:,:)
1173         emp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1)
1174                           CALL iom_put( 'rain'         , frcv(jpr_rain)%z3(:,:,1)              )   ! liquid precipitation
1175         IF( lk_diaar5 )   CALL iom_put( 'hflx_rain_cea', frcv(jpr_rain)%z3(:,:,1) * zcptn(:,:) )   ! heat flux from liq. precip.
1176         ztmp(:,:) = frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * zicefr(:,:)
1177                           CALL iom_put( 'evap_ao_cea'  , ztmp                            )   ! ice-free oce evap (cell average)
1178         IF( lk_diaar5 )   CALL iom_put( 'hflx_evap_cea', ztmp(:,:         ) * zcptn(:,:) )   ! heat flux from from evap (cell ave)
1179      CASE( 'oce and ice'   )   ! received fields: jpr_sbpr, jpr_semp, jpr_oemp, jpr_ievp
1180         emp_tot(:,:) = p_frld(:,:) * frcv(jpr_oemp)%z3(:,:,1) + zicefr(:,:) * frcv(jpr_sbpr)%z3(:,:,1)
1181         emp_ice(:,:) = frcv(jpr_semp)%z3(:,:,1)
1182         sprecip(:,:) = - frcv(jpr_semp)%z3(:,:,1) + frcv(jpr_ievp)%z3(:,:,1)
1183      END SELECT
1184
1185      CALL iom_put( 'snowpre'    , sprecip                                )   ! Snow
1186      CALL iom_put( 'snow_ao_cea', sprecip(:,:         ) * p_frld(:,:)    )   ! Snow        over ice-free ocean  (cell average)
1187      CALL iom_put( 'snow_ai_cea', sprecip(:,:         ) * zicefr(:,:)    )   ! Snow        over sea-ice         (cell average)
1188      CALL iom_put( 'subl_ai_cea', frcv(jpr_ievp)%z3(:,:,1) * zicefr(:,:) )   ! Sublimation over sea-ice         (cell average)
1189      !   
1190      !                                                           ! runoffs and calving (put in emp_tot)
1191      IF( srcv(jpr_rnf)%laction ) THEN
1192         emp_tot(:,:) = emp_tot(:,:) - frcv(jpr_rnf)%z3(:,:,1)
1193                           CALL iom_put( 'runoffs'      , frcv(jpr_rnf)%z3(:,:,1)              )   ! rivers
1194         IF( lk_diaar5 )   CALL iom_put( 'hflx_rnf_cea' , frcv(jpr_rnf)%z3(:,:,1) * zcptn(:,:) )   ! heat flux from rivers
1195      ENDIF
1196      IF( srcv(jpr_cal)%laction ) THEN
1197         emp_tot(:,:) = emp_tot(:,:) - frcv(jpr_cal)%z3(:,:,1)
1198         CALL iom_put( 'calving', frcv(jpr_cal)%z3(:,:,1) )
1199      ENDIF
1200      !
1201!!gm :  this seems to be internal cooking, not sure to need that in a generic interface
1202!!gm                                       at least should be optional...
1203!!       ! remove negative runoff                            ! sum over the global domain
1204!!       zcumulpos = SUM( MAX( frcv(jpr_rnf)%z3(:,:,1), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) )
1205!!       zcumulneg = SUM( MIN( frcv(jpr_rnf)%z3(:,:,1), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) )
1206!!       IF( lk_mpp )   CALL mpp_sum( zcumulpos )
1207!!       IF( lk_mpp )   CALL mpp_sum( zcumulneg )
1208!!       IF( zcumulpos /= 0. ) THEN                          ! distribute negative runoff on positive runoff grid points
1209!!          zcumulneg = 1.e0 + zcumulneg / zcumulpos
1210!!          frcv(jpr_rnf)%z3(:,:,1) = MAX( frcv(jpr_rnf)%z3(:,:,1), 0.e0 ) * zcumulneg
1211!!       ENDIF     
1212!!       emp_tot(:,:) = emp_tot(:,:) - frcv(jpr_rnf)%z3(:,:,1)   ! add runoff to e-p
1213!!
1214!!gm  end of internal cooking
1215
1216      !                                                      ! ========================= !
1217      SELECT CASE( TRIM( sn_rcv_qns%cldes ) )                !   non solar heat fluxes   !   (qns)
1218      !                                                      ! ========================= !
1219      CASE( 'oce only' )                                     ! the required field is directly provided
1220         qns_tot(:,:  ) = frcv(jpr_qnsoce)%z3(:,:,1)
1221      CASE( 'conservative' )                                      ! the required fields are directly provided
1222         qns_tot(:,:  ) = frcv(jpr_qnsmix)%z3(:,:,1)
1223         IF ( TRIM(sn_rcv_qns%clcat) == 'yes' ) THEN
1224            qns_ice(:,:,1:jpl) = frcv(jpr_qnsice)%z3(:,:,1:jpl)
1225         ELSE
1226            ! Set all category values equal for the moment
1227            DO jl=1,jpl
1228               qns_ice(:,:,jl) = frcv(jpr_qnsice)%z3(:,:,1)
1229            ENDDO
1230         ENDIF
1231      CASE( 'oce and ice' )       ! the total flux is computed from ocean and ice fluxes
1232         qns_tot(:,:  ) =  p_frld(:,:) * frcv(jpr_qnsoce)%z3(:,:,1)
1233         IF ( TRIM(sn_rcv_qns%clcat) == 'yes' ) THEN
1234            DO jl=1,jpl
1235               qns_tot(:,:   ) = qns_tot(:,:) + a_i(:,:,jl) * frcv(jpr_qnsice)%z3(:,:,jl)   
1236               qns_ice(:,:,jl) = frcv(jpr_qnsice)%z3(:,:,jl)
1237            ENDDO
1238         ELSE
1239            DO jl=1,jpl
1240               qns_tot(:,:   ) = qns_tot(:,:) + zicefr(:,:) * frcv(jpr_qnsice)%z3(:,:,1)
1241               qns_ice(:,:,jl) = frcv(jpr_qnsice)%z3(:,:,1)
1242            ENDDO
1243         ENDIF
1244      CASE( 'mixed oce-ice' )     ! the ice flux is cumputed from the total flux, the SST and ice informations
1245! ** NEED TO SORT OUT HOW THIS SHOULD WORK IN THE MULTI-CATEGORY CASE - CURRENTLY NOT ALLOWED WHEN INTERFACE INITIALISED **
1246         qns_tot(:,:  ) = frcv(jpr_qnsmix)%z3(:,:,1)
1247         qns_ice(:,:,1) = frcv(jpr_qnsmix)%z3(:,:,1)    &
1248            &            + frcv(jpr_dqnsdt)%z3(:,:,1) * ( pist(:,:,1) - ( (rt0 + psst(:,:  ) ) * p_frld(:,:)   &
1249            &                                                   +          pist(:,:,1)   * zicefr(:,:) ) )
1250      END SELECT
1251      ztmp(:,:) = p_frld(:,:) * sprecip(:,:) * lfus               ! add the latent heat of solid precip. melting
1252      qns_tot(:,:) = qns_tot(:,:) - ztmp(:,:)                     ! over free ocean
1253      IF( lk_diaar5 )   CALL iom_put( 'hflx_snow_cea', ztmp + sprecip(:,:) * zcptn(:,:) )   ! heat flux from snow (cell average)
1254!!gm
1255!!    currently it is taken into account in leads budget but not in the qns_tot, and thus not in
1256!!    the flux that enter the ocean....
1257!!    moreover 1 - it is not diagnose anywhere....
1258!!             2 - it is unclear for me whether this heat lost is taken into account in the atmosphere or not...
1259!!
1260!! similar job should be done for snow and precipitation temperature
1261      !                                     
1262      IF( srcv(jpr_cal)%laction ) THEN                            ! Iceberg melting
1263         ztmp(:,:) = frcv(jpr_cal)%z3(:,:,1) * lfus               ! add the latent heat of iceberg melting
1264         qns_tot(:,:) = qns_tot(:,:) - ztmp(:,:)
1265         IF( lk_diaar5 )   CALL iom_put( 'hflx_cal_cea', ztmp + frcv(jpr_cal)%z3(:,:,1) * zcptn(:,:) )   ! heat flux from calving
1266      ENDIF
1267
1268      !                                                      ! ========================= !
1269      SELECT CASE( TRIM( sn_rcv_qsr%cldes ) )                !      solar heat fluxes    !   (qsr)
1270      !                                                      ! ========================= !
1271      CASE( 'oce only' )
1272         qsr_tot(:,:  ) = MAX(0.0,frcv(jpr_qsroce)%z3(:,:,1))
1273      CASE( 'conservative' )
1274         qsr_tot(:,:  ) = frcv(jpr_qsrmix)%z3(:,:,1)
1275         IF ( TRIM(sn_rcv_qsr%clcat) == 'yes' ) THEN
1276            qsr_ice(:,:,1:jpl) = frcv(jpr_qsrice)%z3(:,:,1:jpl)
1277         ELSE
1278            ! Set all category values equal for the moment
1279            DO jl=1,jpl
1280               qsr_ice(:,:,jl) = frcv(jpr_qsrice)%z3(:,:,1)
1281            ENDDO
1282         ENDIF
1283         qsr_tot(:,:  ) = frcv(jpr_qsrmix)%z3(:,:,1)
1284         qsr_ice(:,:,1) = frcv(jpr_qsrice)%z3(:,:,1)
1285      CASE( 'oce and ice' )
1286         qsr_tot(:,:  ) =  p_frld(:,:) * frcv(jpr_qsroce)%z3(:,:,1)
1287         IF ( TRIM(sn_rcv_qsr%clcat) == 'yes' ) THEN
1288            DO jl=1,jpl
1289               qsr_tot(:,:   ) = qsr_tot(:,:) + a_i(:,:,jl) * frcv(jpr_qsrice)%z3(:,:,jl)   
1290               qsr_ice(:,:,jl) = frcv(jpr_qsrice)%z3(:,:,jl)
1291            ENDDO
1292         ELSE
1293            DO jl=1,jpl
1294               qsr_tot(:,:   ) = qsr_tot(:,:) + zicefr(:,:) * frcv(jpr_qsrice)%z3(:,:,1)
1295               qsr_ice(:,:,jl) = frcv(jpr_qsrice)%z3(:,:,1)
1296            ENDDO
1297         ENDIF
1298      CASE( 'mixed oce-ice' )
1299         qsr_tot(:,:  ) = frcv(jpr_qsrmix)%z3(:,:,1)
1300! ** NEED TO SORT OUT HOW THIS SHOULD WORK IN THE MULTI-CATEGORY CASE - CURRENTLY NOT ALLOWED WHEN INTERFACE INITIALISED **
1301!       Create solar heat flux over ice using incoming solar heat flux and albedos
1302!       ( see OASIS3 user guide, 5th edition, p39 )
1303         qsr_ice(:,:,1) = frcv(jpr_qsrmix)%z3(:,:,1) * ( 1.- palbi(:,:,1) )   &
1304            &            / (  1.- ( albedo_oce_mix(:,:  ) * p_frld(:,:)       &
1305            &                     + palbi         (:,:,1) * zicefr(:,:) ) )
1306      END SELECT
1307      IF( ln_dm2dc ) THEN   ! modify qsr to include the diurnal cycle
1308         qsr_tot(:,:  ) = sbc_dcy( qsr_tot(:,:  ) )
1309         DO jl=1,jpl
1310            qsr_ice(:,:,jl) = sbc_dcy( qsr_ice(:,:,jl) )
1311         ENDDO
1312      ENDIF
1313
1314      SELECT CASE( TRIM( sn_rcv_dqnsdt%cldes ) )
1315      CASE ('coupled')
1316         IF ( TRIM(sn_rcv_dqnsdt%clcat) == 'yes' ) THEN
1317            dqns_ice(:,:,1:jpl) = frcv(jpr_dqnsdt)%z3(:,:,1:jpl)
1318         ELSE
1319            ! Set all category values equal for the moment
1320            DO jl=1,jpl
1321               dqns_ice(:,:,jl) = frcv(jpr_dqnsdt)%z3(:,:,1)
1322            ENDDO
1323         ENDIF
1324      END SELECT
1325
1326      SELECT CASE( TRIM( sn_rcv_iceflx%cldes ) )
1327      CASE ('coupled')
1328         topmelt(:,:,:)=frcv(jpr_topm)%z3(:,:,:)
1329         botmelt(:,:,:)=frcv(jpr_botm)%z3(:,:,:)
1330      END SELECT
1331
1332      CALL wrk_dealloc( jpi,jpj, zcptn, ztmp, zicefr )
1333      !
1334      IF( nn_timing == 1 )  CALL timing_stop('sbc_cpl_ice_flx')
1335      !
1336   END SUBROUTINE sbc_cpl_ice_flx
1337   
1338   
1339   SUBROUTINE sbc_cpl_snd( kt )
1340      !!----------------------------------------------------------------------
1341      !!             ***  ROUTINE sbc_cpl_snd  ***
1342      !!
1343      !! ** Purpose :   provide the ocean-ice informations to the atmosphere
1344      !!
1345      !! ** Method  :   send to the atmosphere through a call to cpl_prism_snd
1346      !!              all the needed fields (as defined in sbc_cpl_init)
1347      !!----------------------------------------------------------------------
1348      INTEGER, INTENT(in) ::   kt
1349      !
1350      INTEGER ::   ji, jj, jl   ! dummy loop indices
1351      INTEGER ::   isec, info   ! local integer
1352      REAL(wp), POINTER, DIMENSION(:,:)   ::   zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1
1353      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztmp3, ztmp4   
1354      !!----------------------------------------------------------------------
1355      !
1356      IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_snd')
1357      !
1358      CALL wrk_alloc( jpi,jpj, zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 )
1359      CALL wrk_alloc( jpi,jpj,jpl, ztmp3, ztmp4 )
1360
1361      isec = ( kt - nit000 ) * NINT(rdttra(1))        ! date of exchanges
1362
1363      zfr_l(:,:) = 1.- fr_i(:,:)
1364
1365      !                                                      ! ------------------------- !
1366      !                                                      !    Surface temperature    !   in Kelvin
1367      !                                                      ! ------------------------- !
1368      IF( ssnd(jps_toce)%laction .OR. ssnd(jps_tice)%laction .OR. ssnd(jps_tmix)%laction ) THEN
1369         SELECT CASE( sn_snd_temp%cldes)
1370         CASE( 'oce only'             )   ;   ztmp1(:,:) =   tsn(:,:,1,jp_tem) + rt0
1371         CASE( 'weighted oce and ice' )   ;   ztmp1(:,:) = ( tsn(:,:,1,jp_tem) + rt0 ) * zfr_l(:,:)   
1372            SELECT CASE( sn_snd_temp%clcat )
1373            CASE( 'yes' )   
1374               ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
1375            CASE( 'no' )
1376               ztmp3(:,:,:) = 0.0
1377               DO jl=1,jpl
1378                  ztmp3(:,:,1) = ztmp3(:,:,1) + tn_ice(:,:,jl) * a_i(:,:,jl)
1379               ENDDO
1380            CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
1381            END SELECT
1382         CASE( 'mixed oce-ice'        )   
1383            ztmp1(:,:) = ( tsn(:,:,1,1) + rt0 ) * zfr_l(:,:) 
1384            DO jl=1,jpl
1385               ztmp1(:,:) = ztmp1(:,:) + tn_ice(:,:,jl) * a_i(:,:,jl)
1386            ENDDO
1387         CASE default                     ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%cldes' )
1388         END SELECT
1389         IF( ssnd(jps_toce)%laction )   CALL cpl_prism_snd( jps_toce, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
1390         IF( ssnd(jps_tice)%laction )   CALL cpl_prism_snd( jps_tice, isec, ztmp3, info )
1391         IF( ssnd(jps_tmix)%laction )   CALL cpl_prism_snd( jps_tmix, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
1392      ENDIF
1393      !
1394      !                                                      ! ------------------------- !
1395      !                                                      !           Albedo          !
1396      !                                                      ! ------------------------- !
1397      IF( ssnd(jps_albice)%laction ) THEN                         ! ice
1398         ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
1399         CALL cpl_prism_snd( jps_albice, isec, ztmp3, info )
1400      ENDIF
1401      IF( ssnd(jps_albmix)%laction ) THEN                         ! mixed ice-ocean
1402         ztmp1(:,:) = albedo_oce_mix(:,:) * zfr_l(:,:)
1403         DO jl=1,jpl
1404            ztmp1(:,:) = ztmp1(:,:) + alb_ice(:,:,jl) * a_i(:,:,jl)
1405         ENDDO
1406         CALL cpl_prism_snd( jps_albmix, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
1407      ENDIF
1408      !                                                      ! ------------------------- !
1409      !                                                      !  Ice fraction & Thickness !
1410      !                                                      ! ------------------------- !
1411      ! Send ice fraction field
1412      IF( ssnd(jps_fice)%laction ) THEN
1413         SELECT CASE( sn_snd_thick%clcat )
1414         CASE( 'yes' )   ;   ztmp3(:,:,1:jpl) =  a_i(:,:,1:jpl)
1415         CASE( 'no'  )   ;   ztmp3(:,:,1    ) = fr_i(:,:      )
1416         CASE default    ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%clcat' )
1417         END SELECT
1418         CALL cpl_prism_snd( jps_fice, isec, ztmp3, info )
1419      ENDIF
1420
1421      ! Send ice and snow thickness field
1422      IF( ssnd(jps_hice)%laction .OR. ssnd(jps_hsnw)%laction ) THEN
1423         SELECT CASE( sn_snd_thick%cldes)
1424         CASE( 'none'                  )       ! nothing to do
1425         CASE( 'weighted ice and snow' )   
1426            SELECT CASE( sn_snd_thick%clcat )
1427            CASE( 'yes' )   
1428               ztmp3(:,:,1:jpl) =  ht_i(:,:,1:jpl) * a_i(:,:,1:jpl)
1429               ztmp4(:,:,1:jpl) =  ht_s(:,:,1:jpl) * a_i(:,:,1:jpl)
1430            CASE( 'no' )
1431               ztmp3(:,:,:) = 0.0   ;  ztmp4(:,:,:) = 0.0
1432               DO jl=1,jpl
1433                  ztmp3(:,:,1) = ztmp3(:,:,1) + ht_i(:,:,jl) * a_i(:,:,jl)
1434                  ztmp4(:,:,1) = ztmp4(:,:,1) + ht_s(:,:,jl) * a_i(:,:,jl)
1435               ENDDO
1436            CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%clcat' )
1437            END SELECT
1438         CASE( 'ice and snow'         )   
1439            ztmp3(:,:,1:jpl) = ht_i(:,:,1:jpl)
1440            ztmp4(:,:,1:jpl) = ht_s(:,:,1:jpl)
1441         CASE default                     ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%cldes' )
1442         END SELECT
1443         IF( ssnd(jps_hice)%laction )   CALL cpl_prism_snd( jps_hice, isec, ztmp3, info )
1444         IF( ssnd(jps_hsnw)%laction )   CALL cpl_prism_snd( jps_hsnw, isec, ztmp4, info )
1445      ENDIF
1446      !
1447#if defined key_cpl_carbon_cycle
1448      !                                                      ! ------------------------- !
1449      !                                                      !  CO2 flux from PISCES     !
1450      !                                                      ! ------------------------- !
1451      IF( ssnd(jps_co2)%laction )   CALL cpl_prism_snd( jps_co2, isec, RESHAPE ( oce_co2, (/jpi,jpj,1/) ) , info )
1452      !
1453#endif
1454      !                                                      ! ------------------------- !
1455      IF( ssnd(jps_ocx1)%laction ) THEN                      !      Surface current      !
1456         !                                                   ! ------------------------- !
1457         !   
1458         !                                                  j+1   j     -----V---F
1459         ! surface velocity always sent from T point                     !       |
1460         !                                                        j      |   T   U
1461         !                                                               |       |
1462         !                                                   j    j-1   -I-------|
1463         !                                               (for I)         |       |
1464         !                                                              i-1  i   i
1465         !                                                               i      i+1 (for I)
1466         SELECT CASE( TRIM( sn_snd_crt%cldes ) )
1467         CASE( 'oce only'             )      ! C-grid ==> T
1468            DO jj = 2, jpjm1
1469               DO ji = fs_2, fs_jpim1   ! vector opt.
1470                  zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj  ,1) )
1471                  zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji  ,jj-1,1) ) 
1472               END DO
1473            END DO
1474         CASE( 'weighted oce and ice' )   
1475            SELECT CASE ( cp_ice_msh )
1476            CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T
1477               DO jj = 2, jpjm1
1478                  DO ji = fs_2, fs_jpim1   ! vector opt.
1479                     zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un   (ji-1,jj  ,1) ) * zfr_l(ji,jj) 
1480                     zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn   (ji  ,jj-1,1) ) * zfr_l(ji,jj)
1481                     zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)
1482                     zity1(ji,jj) = 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)
1483                  END DO
1484               END DO
1485            CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T
1486               DO jj = 2, jpjm1
1487                  DO ji = 2, jpim1   ! NO vector opt.
1488                     zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj) 
1489                     zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj) 
1490                     zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &
1491                        &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)
1492                     zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &
1493                        &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)
1494                  END DO
1495               END DO
1496            CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T
1497               DO jj = 2, jpjm1
1498                  DO ji = 2, jpim1   ! NO vector opt.
1499                     zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj) 
1500                     zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj) 
1501                     zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &
1502                        &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)
1503                     zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &
1504                        &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)
1505                  END DO
1506               END DO
1507            END SELECT
1508            CALL lbc_lnk( zitx1, 'T', -1. )   ;   CALL lbc_lnk( zity1, 'T', -1. )
1509         CASE( 'mixed oce-ice'        )
1510            SELECT CASE ( cp_ice_msh )
1511            CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T
1512               DO jj = 2, jpjm1
1513                  DO ji = fs_2, fs_jpim1   ! vector opt.
1514                     zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un   (ji-1,jj  ,1) ) * zfr_l(ji,jj)   &
1515                        &         + 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)
1516                     zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn   (ji  ,jj-1,1) ) * zfr_l(ji,jj)   &
1517                        &         + 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)
1518                  END DO
1519               END DO
1520            CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T
1521               DO jj = 2, jpjm1
1522                  DO ji = 2, jpim1   ! NO vector opt.
1523                     zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &   
1524                        &         + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &
1525                        &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)
1526                     zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   & 
1527                        &         + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &
1528                        &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)
1529                  END DO
1530               END DO
1531            CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T
1532               DO jj = 2, jpjm1
1533                  DO ji = 2, jpim1   ! NO vector opt.
1534                     zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &   
1535                        &         + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &
1536                        &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)
1537                     zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   & 
1538                        &         + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &
1539                        &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)
1540                  END DO
1541               END DO
1542            END SELECT
1543         END SELECT
1544         CALL lbc_lnk( zotx1, ssnd(jps_ocx1)%clgrid, -1. )   ;   CALL lbc_lnk( zoty1, ssnd(jps_ocy1)%clgrid, -1. )
1545         !
1546         !
1547         IF( TRIM( sn_snd_crt%clvor ) == 'eastward-northward' ) THEN             ! Rotation of the components
1548            !                                                                     ! Ocean component
1549            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 )       ! 1st component
1550            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 )       ! 2nd component
1551            zotx1(:,:) = ztmp1(:,:)                                                   ! overwrite the components
1552            zoty1(:,:) = ztmp2(:,:)
1553            IF( ssnd(jps_ivx1)%laction ) THEN                                     ! Ice component
1554               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 )    ! 1st component
1555               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 )    ! 2nd component
1556               zitx1(:,:) = ztmp1(:,:)                                                ! overwrite the components
1557               zity1(:,:) = ztmp2(:,:)
1558            ENDIF
1559         ENDIF
1560         !
1561         ! spherical coordinates to cartesian -> 2 components to 3 components
1562         IF( TRIM( sn_snd_crt%clvref ) == 'cartesian' ) THEN
1563            ztmp1(:,:) = zotx1(:,:)                     ! ocean currents
1564            ztmp2(:,:) = zoty1(:,:)
1565            CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 )
1566            !
1567            IF( ssnd(jps_ivx1)%laction ) THEN           ! ice velocities
1568               ztmp1(:,:) = zitx1(:,:)
1569               ztmp1(:,:) = zity1(:,:)
1570               CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 )
1571            ENDIF
1572         ENDIF
1573         !
1574         IF( ssnd(jps_ocx1)%laction )   CALL cpl_prism_snd( jps_ocx1, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info )   ! ocean x current 1st grid
1575         IF( ssnd(jps_ocy1)%laction )   CALL cpl_prism_snd( jps_ocy1, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info )   ! ocean y current 1st grid
1576         IF( ssnd(jps_ocz1)%laction )   CALL cpl_prism_snd( jps_ocz1, isec, RESHAPE ( zotz1, (/jpi,jpj,1/) ), info )   ! ocean z current 1st grid
1577         !
1578         IF( ssnd(jps_ivx1)%laction )   CALL cpl_prism_snd( jps_ivx1, isec, RESHAPE ( zitx1, (/jpi,jpj,1/) ), info )   ! ice   x current 1st grid
1579         IF( ssnd(jps_ivy1)%laction )   CALL cpl_prism_snd( jps_ivy1, isec, RESHAPE ( zity1, (/jpi,jpj,1/) ), info )   ! ice   y current 1st grid
1580         IF( ssnd(jps_ivz1)%laction )   CALL cpl_prism_snd( jps_ivz1, isec, RESHAPE ( zitz1, (/jpi,jpj,1/) ), info )   ! ice   z current 1st grid
1581         !
1582      ENDIF
1583      !
1584      CALL wrk_dealloc( jpi,jpj, zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 )
1585      CALL wrk_dealloc( jpi,jpj,jpl, ztmp3, ztmp4 )
1586      !
1587      IF( nn_timing == 1 )  CALL timing_stop('sbc_cpl_snd')
1588      !
1589   END SUBROUTINE sbc_cpl_snd
1590   
1591#else
1592   !!----------------------------------------------------------------------
1593   !!   Dummy module                                            NO coupling
1594   !!----------------------------------------------------------------------
1595   USE par_kind        ! kind definition
1596CONTAINS
1597   SUBROUTINE sbc_cpl_snd( kt )
1598      WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?', kt
1599   END SUBROUTINE sbc_cpl_snd
1600   !
1601   SUBROUTINE sbc_cpl_rcv( kt, k_fsbc, k_ice )     
1602      WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?', kt, k_fsbc, k_ice
1603   END SUBROUTINE sbc_cpl_rcv
1604   !
1605   SUBROUTINE sbc_cpl_ice_tau( p_taui, p_tauj )     
1606      REAL(wp), INTENT(out), DIMENSION(:,:) ::   p_taui   ! i- & j-components of atmos-ice stress [N/m2]
1607      REAL(wp), INTENT(out), DIMENSION(:,:) ::   p_tauj   ! at I-point (B-grid) or U & V-point (C-grid)
1608      p_taui(:,:) = 0.   ;   p_tauj(:,:) = 0. ! stupid definition to avoid warning message when compiling...
1609      WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?'
1610   END SUBROUTINE sbc_cpl_ice_tau
1611   !
1612   SUBROUTINE sbc_cpl_ice_flx( p_frld , palbi   , psst    , pist  )
1613      REAL(wp), INTENT(in   ), DIMENSION(:,:  ) ::   p_frld     ! lead fraction                [0 to 1]
1614      REAL(wp), INTENT(in   ), DIMENSION(:,:,:), OPTIONAL ::   palbi   ! ice albedo
1615      REAL(wp), INTENT(in   ), DIMENSION(:,:  ), OPTIONAL ::   psst    ! sea surface temperature      [Celcius]
1616      REAL(wp), INTENT(in   ), DIMENSION(:,:,:), OPTIONAL ::   pist    ! ice surface temperature      [Kelvin]
1617      WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?', p_frld(1,1), palbi(1,1,1), psst(1,1), pist(1,1,1) 
1618   END SUBROUTINE sbc_cpl_ice_flx
1619   
1620#endif
1621
1622   !!======================================================================
1623END MODULE sbccpl
Note: See TracBrowser for help on using the repository browser.