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/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

source: branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 2633

Last change on this file since 2633 was 2633, checked in by trackstand2, 13 years ago

Renamed wrk_use => wrk_in_use and wrk_release => wrk_not_released

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