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

source: trunk/NEMO/OPA_SRC/SBC/sbccpl.F90 @ 1693

Last change on this file since 1693 was 1693, checked in by smasson, 15 years ago

slight bug when sending surface current in coupled mode, see ticket:575

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