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

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 4428

Last change on this file since 4428 was 3211, checked in by spickles2, 12 years ago

Stephen Pickles, 11 Dec 2011

Commit to bring the rest of the DCSE NEMO development branch
in line with the latest development version. This includes
array index re-ordering of all OPA_SRC/.

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