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

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

source: branches/2014/dev_4728_CNRS04_coupled_interface/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 4857

Last change on this file since 4857 was 4857, checked in by smasson, 9 years ago

dev_4728_CNRS04_coupled_interface: cleaning of the coupling interface for OASIS3-MCT

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