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.
sbcice_cice.F90 in branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

source: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 4616

Last change on this file since 4616 was 4616, checked in by gm, 10 years ago

#1260 : see the associated wiki page for explanation

File size: 40.1 KB
Line 
1MODULE sbcice_cice
2   !!======================================================================
3   !!                       ***  MODULE  sbcice_cice  ***
4   !! To couple with sea ice model CICE (LANL)
5   !!=====================================================================
6#if defined key_cice
7   !!----------------------------------------------------------------------
8   !!   'key_cice' :                                     CICE sea-ice model
9   !!----------------------------------------------------------------------
10   !!   sbc_ice_cice  : sea-ice model time-stepping and update ocean sbc over ice-covered area
11   !!   
12   !!   
13   !!----------------------------------------------------------------------
14   USE oce             ! ocean dynamics and tracers
15   USE dom_oce         ! ocean space and time domain
16   USE domvvl
17   USE phycst, only : rcp, rau0, r1_rau0, rhosn, rhoic
18   USE in_out_manager  ! I/O manager
19   USE lib_mpp         ! distributed memory computing library
20   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
21   USE wrk_nemo        ! work arrays
22   USE timing          ! Timing
23   USE daymod          ! calendar
24   USE fldread         ! read input fields
25
26   USE sbc_oce         ! Surface boundary condition: ocean fields
27   USE sbc_ice         ! Surface boundary condition: ice   fields
28   USE sbcblk_core     ! Surface boundary condition: CORE bulk
29   USE sbccpl
30
31   USE ice_kinds_mod
32   USE ice_blocks
33   USE ice_domain
34   USE ice_domain_size
35   USE ice_boundary
36   USE ice_constants
37   USE ice_gather_scatter
38   USE ice_calendar, only: dt
39   USE ice_state, only: aice,aicen,uvel,vvel,vsno,vsnon,vice,vicen
40   USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  &
41                sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_gbm,     &
42                fresh_gbm,fhocn_gbm,fswthru_gbm,frzmlt,          &
43                flatn_f,fsurfn_f,fcondtopn_f,                    &
44                uatm,vatm,wind,fsw,flw,Tair,potT,Qa,rhoa,zlvl,   &
45                swvdr,swvdf,swidr,swidf
46   USE ice_forcing, only: frcvdr,frcvdf,frcidr,frcidf
47   USE ice_atmo, only: calc_strair
48   USE ice_therm_vertical, only: calc_Tsfc
49
50   USE CICE_InitMod
51   USE CICE_RunMod
52   USE CICE_FinalMod
53
54   IMPLICIT NONE
55   PRIVATE
56
57   PUBLIC cice_sbc_init   ! routine called by sbc_init
58   PUBLIC cice_sbc_final  ! routine called by sbc_final
59   PUBLIC sbc_ice_cice    ! routine called by sbc
60
61   INTEGER , PARAMETER ::   ji_off = INT ( (jpiglo - nx_global) / 2 )
62   INTEGER , PARAMETER ::   jj_off = INT ( (jpjglo - ny_global) / 2 )
63
64   INTEGER , PARAMETER ::   jpfld   = 13   ! maximum number of files to read
65   INTEGER , PARAMETER ::   jp_snow = 1    ! index of snow file
66   INTEGER , PARAMETER ::   jp_rain = 2    ! index of rain file
67   INTEGER , PARAMETER ::   jp_sblm = 3    ! index of sublimation file
68   INTEGER , PARAMETER ::   jp_top1 = 4    ! index of category 1 topmelt file
69   INTEGER , PARAMETER ::   jp_top2 = 5    ! index of category 2 topmelt file
70   INTEGER , PARAMETER ::   jp_top3 = 6    ! index of category 3 topmelt file
71   INTEGER , PARAMETER ::   jp_top4 = 7    ! index of category 4 topmelt file
72   INTEGER , PARAMETER ::   jp_top5 = 8    ! index of category 5 topmelt file
73   INTEGER , PARAMETER ::   jp_bot1 = 9    ! index of category 1 botmelt file
74   INTEGER , PARAMETER ::   jp_bot2 = 10   ! index of category 2 botmelt file
75   INTEGER , PARAMETER ::   jp_bot3 = 11   ! index of category 3 botmelt file
76   INTEGER , PARAMETER ::   jp_bot4 = 12   ! index of category 4 botmelt file
77   INTEGER , PARAMETER ::   jp_bot5 = 13   ! index of category 5 botmelt file
78   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf    ! structure of input fields (file informations, fields read)
79
80   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:), PRIVATE ::   png     ! local array used in sbc_cice_ice
81
82   !! * Substitutions
83#  include "domzgr_substitute.h90"
84   !!----------------------------------------------------------------------
85CONTAINS
86
87   INTEGER FUNCTION sbc_ice_cice_alloc()
88      !!----------------------------------------------------------------------
89      !!                ***  FUNCTION sbc_ice_cice_alloc  ***
90      !!----------------------------------------------------------------------
91      ALLOCATE( png(jpi,jpj,jpnij), STAT=sbc_ice_cice_alloc )
92      IF( lk_mpp                 )   CALL mpp_sum ( sbc_ice_cice_alloc )
93      IF( sbc_ice_cice_alloc > 0 )   CALL ctl_warn('sbc_ice_cice_alloc: allocation of arrays failed.')
94   END FUNCTION sbc_ice_cice_alloc
95
96   SUBROUTINE sbc_ice_cice( kt, nsbc )
97      !!---------------------------------------------------------------------
98      !!                  ***  ROUTINE sbc_ice_cice  ***
99      !!                   
100      !! ** Purpose :   update the ocean surface boundary condition via the
101      !!                CICE Sea Ice Model time stepping
102      !!
103      !! ** Method  : - Get any extra forcing fields for CICE 
104      !!              - Prepare forcing fields
105      !!              - CICE model time stepping
106      !!              - call the routine that computes mass and
107      !!                heat fluxes at the ice/ocean interface
108      !!
109      !! ** Action  : - time evolution of the CICE sea-ice model
110      !!              - update all sbc variables below sea-ice:
111      !!                utau, vtau, qns , qsr, emp , sfx
112      !!---------------------------------------------------------------------
113      INTEGER, INTENT(in) ::   kt      ! ocean time step
114      INTEGER, INTENT(in) ::   nsbc    ! surface forcing type
115      !!----------------------------------------------------------------------
116      !
117      IF( nn_timing == 1 )  CALL timing_start('sbc_ice_cice')
118      !
119      !                                        !----------------------!
120      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !  Ice time-step only  !
121         !                                     !----------------------!
122
123         ! Make sure any fluxes required for CICE are set
124         IF ( nsbc == 2 )  THEN
125            CALL cice_sbc_force(kt)
126         ELSE IF ( nsbc == 5 ) THEN
127            CALL sbc_cpl_ice_flx( 1.0-fr_i  )
128         ENDIF
129
130         CALL cice_sbc_in ( kt, nsbc )
131         CALL CICE_Run
132         CALL cice_sbc_out ( kt, nsbc )
133
134         IF ( nsbc == 5 )  CALL cice_sbc_hadgam(kt+1)
135
136      ENDIF                                          ! End sea-ice time step only
137      !
138      IF( nn_timing == 1 )  CALL timing_stop('sbc_ice_cice')
139
140   END SUBROUTINE sbc_ice_cice
141
142   SUBROUTINE cice_sbc_init (nsbc)
143      !!---------------------------------------------------------------------
144      !!                    ***  ROUTINE cice_sbc_init  ***
145      !! ** Purpose: Initialise ice related fields for NEMO and coupling
146      !!
147      INTEGER, INTENT( in  ) ::   nsbc                ! surface forcing type
148      REAL(wp), DIMENSION(:,:), POINTER :: ztmp1, ztmp2
149      REAL(wp) ::   zcoefu, zcoefv, zcoeff            ! local scalar
150      INTEGER  ::   ji, jj, jl                        ! dummy loop indices
151      !!---------------------------------------------------------------------
152
153      IF( nn_timing == 1 )  CALL timing_start('cice_sbc_init')
154      !
155      CALL wrk_alloc( jpi,jpj, ztmp1, ztmp2 )
156      !
157      IF(lwp) WRITE(numout,*)'cice_sbc_init'
158
159! Initialize CICE
160      CALL CICE_Initialize
161
162! Do some CICE consistency checks
163      IF ( (nsbc == 2) .OR. (nsbc == 5) ) THEN
164         IF ( calc_strair .OR. calc_Tsfc ) THEN
165            CALL ctl_stop( 'STOP', 'cice_sbc_init : Forcing option requires calc_strair=F and calc_Tsfc=F in ice_in' )
166         ENDIF
167      ELSEIF (nsbc == 4) THEN
168         IF ( .NOT. (calc_strair .AND. calc_Tsfc) ) THEN
169            CALL ctl_stop( 'STOP', 'cice_sbc_init : Forcing option requires calc_strair=T and calc_Tsfc=T in ice_in' )
170         ENDIF
171      ENDIF
172
173
174! allocate sbc_ice and sbc_cice arrays
175      IF( sbc_ice_alloc()      /= 0 )   CALL ctl_stop( 'STOP', 'sbc_ice_cice_alloc : unable to allocate arrays' )
176      IF( sbc_ice_cice_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_ice_cice_alloc : unable to allocate cice arrays' )
177
178! Ensure ocean temperatures are nowhere below freezing if not a NEMO restart
179      IF( .NOT. ln_rstart ) THEN
180         tsn(:,:,:,jp_tem) = MAX (tsn(:,:,:,jp_tem),Tocnfrz)
181         tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem)
182      ENDIF
183
184      fr_iu(:,:)=0.0
185      fr_iv(:,:)=0.0
186
187      CALL cice2nemo(aice,fr_i, 'T', 1. )
188      IF ( (nsbc == 2).OR.(nsbc == 5) ) THEN
189         DO jl=1,ncat
190            CALL cice2nemo(aicen(:,:,jl,:),a_i(:,:,jl), 'T', 1. )
191         ENDDO
192      ENDIF
193
194! T point to U point
195! T point to V point
196      DO jj=1,jpjm1
197         DO ji=1,jpim1
198            fr_iu(ji,jj)=0.5*(fr_i(ji,jj)+fr_i(ji+1,jj))*umask(ji,jj,1)
199            fr_iv(ji,jj)=0.5*(fr_i(ji,jj)+fr_i(ji,jj+1))*vmask(ji,jj,1)
200         ENDDO
201      ENDDO
202
203      CALL lbc_lnk ( fr_iu , 'U', 1. )
204      CALL lbc_lnk ( fr_iv , 'V', 1. )
205
206      !                                      ! embedded sea ice
207      IF( nn_ice_embd /= 0 ) THEN            ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass
208         CALL cice2nemo(vsno(:,:,:),ztmp1,'T', 1. )
209         CALL cice2nemo(vice(:,:,:),ztmp2,'T', 1. )
210         snwice_mass  (:,:) = ( rhosn * ztmp1(:,:) + rhoic * ztmp2(:,:)  )
211         snwice_mass_b(:,:) = snwice_mass(:,:)
212      ELSE
213         snwice_mass  (:,:) = 0.0_wp         ! no mass exchanges
214         snwice_mass_b(:,:) = 0.0_wp         ! no mass exchanges
215      ENDIF
216      IF( nn_ice_embd == 2 .AND.          &  ! full embedment (case 2) & no restart :
217         &   .NOT.ln_rstart ) THEN           ! deplete the initial ssh belew sea-ice area
218         sshn(:,:) = sshn(:,:) - snwice_mass(:,:) * r1_rau0
219         sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0
220         !
221         ! Note: Changed the initial values of sshb and sshn=>  need to recompute ssh[u,v,f]_[b,n]
222         !       which were previously set in domvvl
223         IF ( lk_vvl ) THEN            ! Is this necessary? embd 2 should be restricted to vvl only???
224            DO jj = 1, jpjm1
225               DO ji = 1, jpim1                    ! caution: use of Vector Opt. not possible
226                  zcoefu = 0.5  * umask(ji,jj,1) * r1_e1e2u(ji,jj)
227                  zcoefv = 0.5  * vmask(ji,jj,1) * r1_e1e2v(ji,jj)
228                  zcoeff = 0.25 * umask(ji,jj,1) * umask(ji,jj+1,1)
229                  sshu_b(ji,jj) = zcoefu * ( e1e2t(ji,jj) * sshb(ji,jj) + e1e2t(ji+1,jj  ) * sshb(ji+1,jj  ) )
230                  sshv_b(ji,jj) = zcoefv * ( e1e2t(ji,jj) * sshb(ji,jj) + e1e2t(ji  ,jj+1) * sshb(ji  ,jj+1) )
231                  sshu_n(ji,jj) = zcoefu * ( e1e2t(ji,jj) * sshn(ji,jj) + e1e2t(ji+1,jj  ) * sshn(ji+1,jj  ) )
232                  sshv_n(ji,jj) = zcoefv * ( e1e2t(ji,jj) * sshn(ji,jj) + e1e2t(ji  ,jj+1) * sshn(ji  ,jj+1) )
233               END DO
234            END DO
235            CALL lbc_lnk( sshu_b, 'U', 1. )   ;   CALL lbc_lnk( sshu_n, 'U', 1. )
236            CALL lbc_lnk( sshv_b, 'V', 1. )   ;   CALL lbc_lnk( sshv_n, 'V', 1. )
237            DO jj = 1, jpjm1
238               DO ji = 1, jpim1      ! NO Vector Opt.
239                  sshf_n(ji,jj) = 0.5  * umask(ji,jj,1) * umask(ji,jj+1,1)      &
240                       &               * (  e1e2u(ji,jj  ) * sshu_n(ji,jj  )    &
241                       &                  + e1e2u(ji,jj+1) * sshu_n(ji,jj+1)  ) * r1_e1e2f(ji,jj)
242              END DO
243            END DO
244            CALL lbc_lnk( sshf_n, 'F', 1. )
245          ENDIF
246      ENDIF
247 
248      CALL wrk_dealloc( jpi,jpj, ztmp1, ztmp2 )
249      !
250      IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_init')
251      !
252   END SUBROUTINE cice_sbc_init
253
254   
255   SUBROUTINE cice_sbc_in (kt, nsbc)
256      !!---------------------------------------------------------------------
257      !!                    ***  ROUTINE cice_sbc_in  ***
258      !! ** Purpose: Set coupling fields and pass to CICE
259      !!---------------------------------------------------------------------
260      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
261      INTEGER, INTENT(in   ) ::   nsbc ! surface forcing type
262      !
263      INTEGER  ::   ji, jj, jl                   ! dummy loop indices     
264      REAL(wp), DIMENSION(:,:), POINTER :: ztmp, zpice
265      REAL(wp), DIMENSION(:,:,:), POINTER :: ztmpn
266      REAL(wp) ::   zintb, zintn  ! dummy argument
267      !!---------------------------------------------------------------------
268
269      IF( nn_timing == 1 )  CALL timing_start('cice_sbc_in')
270      !
271      CALL wrk_alloc( jpi,jpj, ztmp, zpice )
272      CALL wrk_alloc( jpi,jpj,ncat, ztmpn )
273
274      IF( kt == nit000 )  THEN
275         IF(lwp) WRITE(numout,*)'cice_sbc_in'
276      ENDIF
277
278      ztmp(:,:)=0.0
279
280! Aggregate ice concentration already set in cice_sbc_out (or cice_sbc_init on
281! the first time-step)
282
283! forced and coupled case
284
285      IF ( (nsbc == 2).OR.(nsbc == 5) ) THEN
286
287         ztmpn(:,:,:)=0.0
288
289! x comp of wind stress (CI_1)
290! U point to F point
291         DO jj=1,jpjm1
292            DO ji=1,jpi
293               ztmp(ji,jj) = 0.5 * (  fr_iu(ji,jj) * utau(ji,jj)      &
294                                    + fr_iu(ji,jj+1) * utau(ji,jj+1) ) * fmask(ji,jj,1)
295            ENDDO
296         ENDDO
297         CALL nemo2cice(ztmp,strax,'F', -1. )
298
299! y comp of wind stress (CI_2)
300! V point to F point
301         DO jj=1,jpj
302            DO ji=1,jpim1
303               ztmp(ji,jj) = 0.5 * (  fr_iv(ji,jj) * vtau(ji,jj)      &
304                                    + fr_iv(ji+1,jj) * vtau(ji+1,jj) ) * fmask(ji,jj,1)
305            ENDDO
306         ENDDO
307         CALL nemo2cice(ztmp,stray,'F', -1. )
308
309! Surface downward latent heat flux (CI_5)
310         IF (nsbc == 2) THEN
311            DO jl=1,ncat
312               ztmpn(:,:,jl)=qla_ice(:,:,1)*a_i(:,:,jl)
313            ENDDO
314         ELSE
315! emp_ice is set in sbc_cpl_ice_flx as sublimation-snow
316            qla_ice(:,:,1)= - ( emp_ice(:,:)+sprecip(:,:) ) * Lsub
317! End of temporary code
318            DO jj=1,jpj
319               DO ji=1,jpi
320                  IF (fr_i(ji,jj).eq.0.0) THEN
321                     DO jl=1,ncat
322                        ztmpn(ji,jj,jl)=0.0
323                     ENDDO
324                     ! This will then be conserved in CICE
325                     ztmpn(ji,jj,1)=qla_ice(ji,jj,1)
326                  ELSE
327                     DO jl=1,ncat
328                        ztmpn(ji,jj,jl)=qla_ice(ji,jj,1)*a_i(ji,jj,jl)/fr_i(ji,jj)
329                     ENDDO
330                  ENDIF
331               ENDDO
332            ENDDO
333         ENDIF
334         DO jl=1,ncat
335            CALL nemo2cice(ztmpn(:,:,jl),flatn_f(:,:,jl,:),'T', 1. )
336
337! GBM conductive flux through ice (CI_6)
338!  Convert to GBM
339            IF (nsbc == 2) THEN
340               ztmp(:,:) = botmelt(:,:,jl)*a_i(:,:,jl)
341            ELSE
342               ztmp(:,:) = botmelt(:,:,jl)
343            ENDIF
344            CALL nemo2cice(ztmp,fcondtopn_f(:,:,jl,:),'T', 1. )
345
346! GBM surface heat flux (CI_7)
347!  Convert to GBM
348            IF (nsbc == 2) THEN
349               ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl))*a_i(:,:,jl) 
350            ELSE
351               ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl))
352            ENDIF
353            CALL nemo2cice(ztmp,fsurfn_f(:,:,jl,:),'T', 1. )
354         ENDDO
355
356      ELSE IF (nsbc == 4) THEN
357
358! Pass CORE forcing fields to CICE (which will calculate heat fluxes etc itself)
359! x comp and y comp of atmosphere surface wind (CICE expects on T points)
360         ztmp(:,:) = wndi_ice(:,:)
361         CALL nemo2cice(ztmp,uatm,'T', -1. )
362         ztmp(:,:) = wndj_ice(:,:)
363         CALL nemo2cice(ztmp,vatm,'T', -1. )
364         ztmp(:,:) = SQRT ( wndi_ice(:,:)**2 + wndj_ice(:,:)**2 )
365         CALL nemo2cice(ztmp,wind,'T', 1. )    ! Wind speed (m/s)
366         ztmp(:,:) = qsr_ice(:,:,1)
367         CALL nemo2cice(ztmp,fsw,'T', 1. )     ! Incoming short-wave (W/m^2)
368         ztmp(:,:) = qlw_ice(:,:,1)
369         CALL nemo2cice(ztmp,flw,'T', 1. )     ! Incoming long-wave (W/m^2)
370         ztmp(:,:) = tatm_ice(:,:)
371         CALL nemo2cice(ztmp,Tair,'T', 1. )    ! Air temperature (K)
372         CALL nemo2cice(ztmp,potT,'T', 1. )    ! Potential temp (K)
373! Following line uses MAX(....) to avoid problems if tatm_ice has unset halo rows 
374         ztmp(:,:) = 101000. / ( 287.04 * MAX(1.0,tatm_ice(:,:)) )   
375                                               ! Constant (101000.) atm pressure assumed
376         CALL nemo2cice(ztmp,rhoa,'T', 1. )    ! Air density (kg/m^3)
377         ztmp(:,:) = qatm_ice(:,:)
378         CALL nemo2cice(ztmp,Qa,'T', 1. )      ! Specific humidity (kg/kg)
379         ztmp(:,:)=10.0
380         CALL nemo2cice(ztmp,zlvl,'T', 1. )    ! Atmos level height (m)
381
382! May want to check all values are physically realistic (as in CICE routine
383! prepare_forcing)?
384
385! Divide shortwave into spectral bands (as in prepare_forcing)
386         ztmp(:,:)=qsr_ice(:,:,1)*frcvdr       ! visible direct
387         CALL nemo2cice(ztmp,swvdr,'T', 1. )             
388         ztmp(:,:)=qsr_ice(:,:,1)*frcvdf       ! visible diffuse
389         CALL nemo2cice(ztmp,swvdf,'T', 1. )             
390         ztmp(:,:)=qsr_ice(:,:,1)*frcidr       ! near IR direct
391         CALL nemo2cice(ztmp,swidr,'T', 1. )
392         ztmp(:,:)=qsr_ice(:,:,1)*frcidf       ! near IR diffuse
393         CALL nemo2cice(ztmp,swidf,'T', 1. )
394
395      ENDIF
396
397! Snowfall
398! Ensure fsnow is positive (as in CICE routine prepare_forcing) 
399      ztmp(:,:)=MAX(fr_i(:,:)*sprecip(:,:),0.0) 
400      CALL nemo2cice(ztmp,fsnow,'T', 1. ) 
401
402! Rainfall
403      ztmp(:,:)=fr_i(:,:)*(tprecip(:,:)-sprecip(:,:))
404      CALL nemo2cice(ztmp,frain,'T', 1. ) 
405
406! Freezing/melting potential
407! Calculated over NEMO leapfrog timestep (hence 2*dt)
408      nfrzmlt(:,:)=rau0*rcp*fse3t_m(:,:)*(Tocnfrz-sst_m(:,:))/(2.0*dt)
409
410      ztmp(:,:) = nfrzmlt(:,:)
411      CALL nemo2cice(ztmp,frzmlt,'T', 1. )
412
413! SST  and SSS
414
415      CALL nemo2cice(sst_m,sst,'T', 1. )
416      CALL nemo2cice(sss_m,sss,'T', 1. )
417
418! x comp and y comp of surface ocean current
419! U point to F point
420      DO jj=1,jpjm1
421         DO ji=1,jpi
422            ztmp(ji,jj)=0.5*(ssu_m(ji,jj)+ssu_m(ji,jj+1))*fmask(ji,jj,1)
423         ENDDO
424      ENDDO
425      CALL nemo2cice(ztmp,uocn,'F', -1. )
426
427! V point to F point
428      DO jj=1,jpj
429         DO ji=1,jpim1
430            ztmp(ji,jj)=0.5*(ssv_m(ji,jj)+ssv_m(ji+1,jj))*fmask(ji,jj,1)
431         ENDDO
432      ENDDO
433      CALL nemo2cice(ztmp,vocn,'F', -1. )
434
435      IF( nn_ice_embd == 2 ) THEN             !== embedded sea ice: compute representative ice top surface ==!
436          !
437          ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[n/nn_fsbc], n=0,nn_fsbc-1}
438          !                                               = (1/nn_fsbc)^2 * {SUM[n], n=0,nn_fsbc-1}
439         zintn = REAL( nn_fsbc - 1 ) / REAL( nn_fsbc ) * 0.5_wp
440          !
441          ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[1-n/nn_fsbc], n=0,nn_fsbc-1}
442          !                                               = (1/nn_fsbc)^2 * (nn_fsbc^2 - {SUM[n], n=0,nn_fsbc-1})
443         zintb = REAL( nn_fsbc + 1 ) / REAL( nn_fsbc ) * 0.5_wp
444          !
445         zpice(:,:) = ssh_m(:,:) + (  zintn * snwice_mass(:,:) +  zintb * snwice_mass_b(:,:)  ) * r1_rau0
446          !
447         !
448      ELSE                                    !== non-embedded sea ice: use ocean surface for slope calculation ==!
449         zpice(:,:) = ssh_m(:,:)
450      ENDIF
451
452! x comp and y comp of sea surface slope (on F points)
453! T point to F point
454      DO jj = 1, jpjm1
455         DO ji = 1, jpim1
456            ztmp(ji,jj)=0.5 * (  (zpice(ji+1,jj  )-zpice(ji,jj  )) * r1_e1u(ji,jj  )    &
457               &               + (zpice(ji+1,jj+1)-zpice(ji,jj+1)) * r1_e1u(ji,jj+1)  ) * fmask(ji,jj,1)
458         END DO
459      END DO
460      CALL nemo2cice( ztmp,ss_tltx,'F', -1. )
461
462! T point to F point
463      DO jj = 1, jpjm1
464         DO ji = 1, jpim1
465            ztmp(ji,jj)=0.5 * (  (zpice(ji  ,jj+1)-zpice(ji  ,jj)) * r1_e2v(ji  ,jj)    &
466               &               + (zpice(ji+1,jj+1)-zpice(ji+1,jj)) * r1_e2v(ji+1,jj)  ) *  fmask(ji,jj,1)
467         END DO
468      END DO
469      CALL nemo2cice(ztmp,ss_tlty,'F', -1. )
470
471      CALL wrk_dealloc( jpi,jpj, ztmp )
472      CALL wrk_dealloc( jpi,jpj,ncat, ztmpn )
473      !
474      IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_in')
475      !
476   END SUBROUTINE cice_sbc_in
477
478
479   SUBROUTINE cice_sbc_out (kt,nsbc)
480      !!---------------------------------------------------------------------
481      !!                    ***  ROUTINE cice_sbc_out  ***
482      !! ** Purpose: Get fields from CICE and set surface fields for NEMO
483      !!---------------------------------------------------------------------
484      INTEGER, INTENT( in  ) ::   kt   ! ocean time step
485      INTEGER, INTENT( in  ) ::   nsbc ! surface forcing type
486     
487      INTEGER  ::   ji, jj, jl                 ! dummy loop indices
488      REAL(wp), DIMENSION(:,:), POINTER :: ztmp1, ztmp2
489      !!---------------------------------------------------------------------
490
491      IF( nn_timing == 1 )  CALL timing_start('cice_sbc_out')
492      !
493      CALL wrk_alloc( jpi,jpj, ztmp1, ztmp2 )
494     
495      IF( kt == nit000 )  THEN
496         IF(lwp) WRITE(numout,*)'cice_sbc_out'
497      ENDIF
498     
499! x comp of ocean-ice stress
500      CALL cice2nemo(strocnx,ztmp1,'F', -1. )
501      ss_iou(:,:)=0.0
502! F point to U point
503      DO jj=2,jpjm1
504         DO ji=2,jpim1
505            ss_iou(ji,jj) = 0.5 * ( ztmp1(ji,jj-1) + ztmp1(ji,jj) ) * umask(ji,jj,1)
506         ENDDO
507      ENDDO
508      CALL lbc_lnk( ss_iou , 'U', -1. )
509
510! y comp of ocean-ice stress
511      CALL cice2nemo(strocny,ztmp1,'F', -1. )
512      ss_iov(:,:)=0.0
513! F point to V point
514
515      DO jj=1,jpjm1
516         DO ji=2,jpim1
517            ss_iov(ji,jj) = 0.5 * ( ztmp1(ji-1,jj) + ztmp1(ji,jj) ) * vmask(ji,jj,1)
518         ENDDO
519      ENDDO
520      CALL lbc_lnk( ss_iov , 'V', -1. )
521
522! x and y comps of surface stress
523! Combine wind stress and ocean-ice stress
524! [Note that fr_iu hasn't yet been updated, so still from start of CICE timestep]
525
526      utau(:,:) = ( 1.0 - fr_iu(:,:) ) * utau(:,:) - ss_iou(:,:)
527      vtau(:,:) = ( 1.0 - fr_iv(:,:) ) * vtau(:,:) - ss_iov(:,:)     
528
529! Freshwater fluxes
530
531      IF (nsbc == 2) THEN
532! Note that emp from the forcing files is evap*(1-aice)-(tprecip-aice*sprecip)
533! What we want here is evap*(1-aice)-tprecip*(1-aice) hence manipulation below
534! Not ideal since aice won't be the same as in the atmosphere. 
535! Better to use evap and tprecip? (but for now don't read in evap in this case)
536         emp(:,:)  = emp(:,:) + fr_i(:,:) * ( tprecip(:,:) - sprecip(:,:) )
537      ELSE IF (nsbc == 4) THEN
538         emp(:,:)  = ( 1.0 - fr_i(:,:) ) * emp(:,:)       
539      ELSE IF (nsbc ==5) THEN
540! emp_tot is set in sbc_cpl_ice_flx (called from cice_sbc_in above)
541! This is currently as required with the coupling fields from the UM atmosphere
542         emp(:,:) = emp_tot(:,:)+tprecip(:,:)*fr_i(:,:) 
543      ENDIF
544
545      CALL cice2nemo( fresh_gbm, ztmp1,'T', 1. )
546      CALL cice2nemo( fsalt_gbm, ztmp2,'T', 1. )
547
548! Check to avoid unphysical expression when ice is forming (ztmp1 negative)
549! Otherwise we are effectively allowing ice of higher salinity than the ocean to form
550! which has to be compensated for by the ocean salinity potentially going negative
551! This check breaks conservation but seems reasonable until we have prognostic ice salinity
552! Note the 1000.0 below is to convert from kg salt to g salt (needed for PSU)
553      WHERE (ztmp1(:,:) < 0._wp )   ztmp2(:,:) = MAX( ztmp2(:,:) , ztmp1(:,:)*sss_m(:,:)/1000._wp )
554      sfx(:,:) = ztmp2(:,:) * 1000.0
555      emp(:,:) = emp(:,:) - ztmp1(:,:)
556 
557      CALL lbc_lnk( emp , 'T', 1. )
558      CALL lbc_lnk( sfx , 'T', 1. )
559
560! Solar penetrative radiation and non solar surface heat flux
561
562! Scale qsr and qns according to ice fraction (bulk formulae only)
563
564      IF (nsbc == 4) THEN
565         qsr(:,:)=qsr(:,:)*(1.0-fr_i(:,:))
566         qns(:,:)=qns(:,:)*(1.0-fr_i(:,:))
567      ENDIF
568! Take into account snow melting except for fully coupled when already in qns_tot
569      IF (nsbc == 5) THEN
570         qsr(:,:)= qsr_tot(:,:)
571         qns(:,:)= qns_tot(:,:)
572      ELSE
573         qns(:,:)= qns(:,:)-sprecip(:,:)*Lfresh*(1.0-fr_i(:,:))
574      ENDIF
575
576! Now add in ice / snow related terms
577! [fswthru will be zero unless running with calc_Tsfc=T in CICE]
578      CALL cice2nemo( fswthru_gbm,ztmp1,'T', 1. )
579      qsr(:,:)=qsr(:,:)+ztmp1(:,:)
580      CALL lbc_lnk( qsr , 'T', 1. )
581
582      DO jj=1,jpj
583         DO ji=1,jpi
584            nfrzmlt(ji,jj) = MAX (nfrzmlt(ji,jj) , 0._wp )
585         ENDDO
586      ENDDO
587
588      CALL cice2nemo(fhocn_gbm,ztmp1,'T', 1. )
589      qns(:,:)=qns(:,:)+nfrzmlt(:,:)+ztmp1(:,:)
590
591      CALL lbc_lnk( qns , 'T', 1. )
592
593! Prepare for the following CICE time-step
594
595      CALL cice2nemo(aice,fr_i,'T', 1. )
596      IF ( (nsbc == 2).OR.(nsbc == 5) ) THEN
597         DO jl=1,ncat
598            CALL cice2nemo(aicen(:,:,jl,:),a_i(:,:,jl), 'T', 1. )
599         ENDDO
600      ENDIF
601
602! T point to U point
603! T point to V point
604      DO jj=1,jpjm1
605         DO ji=1,jpim1
606            fr_iu(ji,jj)=0.5*(fr_i(ji,jj)+fr_i(ji+1,jj))*umask(ji,jj,1)
607            fr_iv(ji,jj)=0.5*(fr_i(ji,jj)+fr_i(ji,jj+1))*vmask(ji,jj,1)
608         ENDDO
609      ENDDO
610
611      CALL lbc_lnk ( fr_iu , 'U', 1. )
612      CALL lbc_lnk ( fr_iv , 'V', 1. )
613
614      !                                      ! embedded sea ice
615      IF( nn_ice_embd /= 0 ) THEN            ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass
616         CALL cice2nemo(vsno(:,:,:),ztmp1,'T', 1. )
617         CALL cice2nemo(vice(:,:,:),ztmp2,'T', 1. )
618         snwice_mass  (:,:) = ( rhosn * ztmp1(:,:) + rhoic * ztmp2(:,:)  )
619         snwice_mass_b(:,:) = snwice_mass(:,:)
620         snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) / dt
621      ENDIF
622
623! Release work space
624
625      CALL wrk_dealloc( jpi,jpj, ztmp1, ztmp2 )
626      !
627      IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_out')
628      !
629   END SUBROUTINE cice_sbc_out
630
631
632#if defined key_oasis3 || defined key_oasis4
633   SUBROUTINE cice_sbc_hadgam( kt )
634      !!---------------------------------------------------------------------
635      !!                    ***  ROUTINE cice_sbc_hadgam  ***
636      !! ** Purpose: Prepare fields needed to pass to HadGAM3 atmosphere
637      !!
638      !!
639      INTEGER, INTENT( in  ) ::   kt   ! ocean time step
640      !!---------------------------------------------------------------------
641
642      INTEGER  ::   jl                        ! dummy loop index
643      INTEGER  ::   ierror
644
645      IF( nn_timing == 1 )  CALL timing_start('cice_sbc_hadgam')
646      !
647      IF( kt == nit000 )  THEN
648         IF(lwp) WRITE(numout,*)'cice_sbc_hadgam'
649         IF( sbc_cpl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' )
650      ENDIF
651
652      !                                         ! =========================== !
653      !                                         !   Prepare Coupling fields   !
654      !                                         ! =========================== !
655
656! x and y comp of ice velocity
657
658      CALL cice2nemo(uvel,u_ice,'F', -1. )
659      CALL cice2nemo(vvel,v_ice,'F', -1. )
660
661! Ice concentration (CO_1) = a_i calculated at end of cice_sbc_out 
662
663! Snow and ice thicknesses (CO_2 and CO_3)
664
665      DO jl = 1,ncat
666         CALL cice2nemo(vsnon(:,:,jl,:),ht_s(:,:,jl),'T', 1. )
667         CALL cice2nemo(vicen(:,:,jl,:),ht_i(:,:,jl),'T', 1. )
668      ENDDO
669      !
670      IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_hadgam')
671      !
672   END SUBROUTINE cice_sbc_hadgam
673
674#else
675   SUBROUTINE cice_sbc_hadgam( kt )    ! Dummy routine
676      INTEGER, INTENT( in  ) ::   kt   ! ocean time step
677      WRITE(*,*) 'cice_sbc_hadgam: You should not have seen this print! error?'
678   END SUBROUTINE cice_sbc_hadgam
679#endif
680
681   SUBROUTINE cice_sbc_final
682      !!---------------------------------------------------------------------
683      !!                    ***  ROUTINE cice_sbc_final  ***
684      !! ** Purpose: Finalize CICE
685      !!---------------------------------------------------------------------
686
687      IF(lwp) WRITE(numout,*)'cice_sbc_final'
688
689      CALL CICE_Finalize
690
691   END SUBROUTINE cice_sbc_final
692
693   SUBROUTINE cice_sbc_force (kt)
694      !!---------------------------------------------------------------------
695      !!                    ***  ROUTINE cice_sbc_force  ***
696      !! ** Purpose : Provide CICE forcing from files
697      !!
698      !!---------------------------------------------------------------------
699      !! ** Method  :   READ monthly flux file in NetCDF files
700      !!     
701      !!  snowfall   
702      !!  rainfall   
703      !!  sublimation rate   
704      !!  topmelt (category)
705      !!  botmelt (category)
706      !!
707      !! History :
708      !!----------------------------------------------------------------------
709      !! * Modules used
710      USE iom
711
712      !! * arguments
713      INTEGER, INTENT( in  ) ::   kt ! ocean time step
714
715      INTEGER  ::   ierror             ! return error code
716      INTEGER  ::   ifpr               ! dummy loop index
717      !!
718      CHARACTER(len=100) ::  cn_dir                            !   Root directory for location of CICE forcing files
719      TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i                 ! array of namelist informations on the fields to read
720      TYPE(FLD_N) ::   sn_snow, sn_rain, sn_sblm               ! informations about the fields to be read
721      TYPE(FLD_N) ::   sn_top1, sn_top2, sn_top3, sn_top4, sn_top5
722      TYPE(FLD_N) ::   sn_bot1, sn_bot2, sn_bot3, sn_bot4, sn_bot5 
723
724      !!
725      NAMELIST/namsbc_cice/ cn_dir, sn_snow, sn_rain, sn_sblm,   &
726         &                          sn_top1, sn_top2, sn_top3, sn_top4, sn_top5,  &
727         &                          sn_bot1, sn_bot2, sn_bot3, sn_bot4, sn_bot5
728      INTEGER :: ios
729      !!---------------------------------------------------------------------
730
731      !                                         ! ====================== !
732      IF( kt == nit000 ) THEN                   !  First call kt=nit000  !
733         !                                      ! ====================== !
734         REWIND( numnam_ref )              ! Namelist namsbc_cice in reference namelist :
735         READ  ( numnam_ref, namsbc_cice, IOSTAT = ios, ERR = 901)
736901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_cice in reference namelist', lwp )
737
738         REWIND( numnam_cfg )              ! Namelist namsbc_cice in configuration namelist : Parameters of the run
739         READ  ( numnam_cfg, namsbc_cice, IOSTAT = ios, ERR = 902 )
740902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_cice in configuration namelist', lwp )
741         WRITE ( numond, namsbc_cice )
742
743         ! store namelist information in an array
744         slf_i(jp_snow) = sn_snow   ;   slf_i(jp_rain) = sn_rain   ;   slf_i(jp_sblm) = sn_sblm
745         slf_i(jp_top1) = sn_top1   ;   slf_i(jp_top2) = sn_top2   ;   slf_i(jp_top3) = sn_top3
746         slf_i(jp_top4) = sn_top4   ;   slf_i(jp_top5) = sn_top5   ;   slf_i(jp_bot1) = sn_bot1
747         slf_i(jp_bot2) = sn_bot2   ;   slf_i(jp_bot3) = sn_bot3   ;   slf_i(jp_bot4) = sn_bot4
748         slf_i(jp_bot5) = sn_bot5
749         
750         ! set sf structure
751         ALLOCATE( sf(jpfld), STAT=ierror )
752         IF( ierror > 0 ) THEN
753            CALL ctl_stop( 'cice_sbc_force: unable to allocate sf structure' )   ;   RETURN
754         ENDIF
755
756         DO ifpr= 1, jpfld
757            ALLOCATE( sf(ifpr)%fnow(jpi,jpj,1) )
758            ALLOCATE( sf(ifpr)%fdta(jpi,jpj,1,2) )
759         END DO
760
761         ! fill sf with slf_i and control print
762         CALL fld_fill( sf, slf_i, cn_dir, 'cice_sbc_force', 'flux formulation for CICE', 'namsbc_cice' )
763         !
764      ENDIF
765
766      CALL fld_read( kt, nn_fsbc, sf )           ! Read input fields and provides the
767      !                                          ! input fields at the current time-step
768
769      ! set the fluxes from read fields
770      sprecip(:,:) = sf(jp_snow)%fnow(:,:,1)
771      tprecip(:,:) = sf(jp_snow)%fnow(:,:,1)+sf(jp_rain)%fnow(:,:,1)
772! May be better to do this conversion somewhere else
773      qla_ice(:,:,1) = -Lsub*sf(jp_sblm)%fnow(:,:,1)
774      topmelt(:,:,1) = sf(jp_top1)%fnow(:,:,1)
775      topmelt(:,:,2) = sf(jp_top2)%fnow(:,:,1)
776      topmelt(:,:,3) = sf(jp_top3)%fnow(:,:,1)
777      topmelt(:,:,4) = sf(jp_top4)%fnow(:,:,1)
778      topmelt(:,:,5) = sf(jp_top5)%fnow(:,:,1)
779      botmelt(:,:,1) = sf(jp_bot1)%fnow(:,:,1)
780      botmelt(:,:,2) = sf(jp_bot2)%fnow(:,:,1)
781      botmelt(:,:,3) = sf(jp_bot3)%fnow(:,:,1)
782      botmelt(:,:,4) = sf(jp_bot4)%fnow(:,:,1)
783      botmelt(:,:,5) = sf(jp_bot5)%fnow(:,:,1)
784
785      ! control print (if less than 100 time-step asked)
786      IF( nitend-nit000 <= 100 .AND. lwp ) THEN
787         WRITE(numout,*) 
788         WRITE(numout,*) '        read forcing fluxes for CICE OK'
789         CALL FLUSH(numout)
790      ENDIF
791
792   END SUBROUTINE cice_sbc_force
793
794   SUBROUTINE nemo2cice( pn, pc, cd_type, psgn)
795      !!---------------------------------------------------------------------
796      !!                    ***  ROUTINE nemo2cice  ***
797      !! ** Purpose :   Transfer field in NEMO array to field in CICE array. 
798#if defined key_nemocice_decomp
799      !!             
800      !!                NEMO and CICE PE sub domains are identical, hence
801      !!                there is no need to gather or scatter data from
802      !!                one PE configuration to another.
803#else
804      !!                Automatically gather/scatter between
805      !!                different processors and blocks
806      !! ** Method :    A. Ensure all haloes are filled in NEMO field (pn)
807      !!                B. Gather pn into global array (png)
808      !!                C. Map png into CICE global array (pcg)
809      !!                D. Scatter pcg to CICE blocks (pc) + update haloes 
810#endif
811      !!---------------------------------------------------------------------
812      CHARACTER(len=1), INTENT( in ) ::   cd_type   ! nature of pn grid-point (= T or F)
813      REAL(wp)        , INTENT( in ) ::   psgn      ! control of the sign change
814          !                                         !   =-1 , the sign is modified following the type of b.c. used
815          !                                         !   = 1 , no sign change
816      REAL(wp), DIMENSION(jpi,jpj) ::   pn          !!gm INTENT missing !!!!!
817#if !defined key_nemocice_decomp
818      REAL(wp), DIMENSION(jpiglo,jpjglo) :: png2
819      REAL (kind=dbl_kind), dimension(nx_global,ny_global) :: pcg
820#endif
821      REAL (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) ::   pc   !!gm INTENT missing !!!!
822      INTEGER (int_kind) :: &
823         field_type,        &! id for type of field (scalar, vector, angle)
824         grid_loc            ! id for location on horizontal grid
825                            !  (center, NEcorner, Nface, Eface)
826
827      INTEGER  ::   ji, jj, jn                      ! dummy loop indices
828      !!---------------------------------------------------------------------
829
830!     A. Ensure all haloes are filled in NEMO field (pn)
831
832      CALL lbc_lnk( pn , cd_type, psgn )
833
834#if defined key_nemocice_decomp
835
836      ! Copy local domain data from NEMO to CICE field
837      pc(:,:,1)=0.0
838      DO jj=2,ny_block-1
839         DO ji=2,nx_block-1
840            pc(ji,jj,1)=pn(ji-1+ji_off,jj-1+jj_off)
841         ENDDO
842      ENDDO
843
844#else
845
846!     B. Gather pn into global array (png)
847
848      IF ( jpnij > 1) THEN
849         CALL mppsync
850         CALL mppgather (pn,0,png) 
851         CALL mppsync
852      ELSE
853         png(:,:,1)=pn(:,:)
854      ENDIF
855
856!     C. Map png into CICE global array (pcg)
857
858! Need to make sure this is robust to changes in NEMO halo rows....
859! (may be OK but not 100% sure)
860
861      IF (nproc==0) THEN     
862!        pcg(:,:)=0.0
863         DO jn=1,jpnij
864            DO jj=nldjt(jn),nlejt(jn)
865               DO ji=nldit(jn),nleit(jn)
866                  png2(ji+nimppt(jn)-1,jj+njmppt(jn)-1)=png(ji,jj,jn)
867               ENDDO
868            ENDDO
869         ENDDO
870         DO jj=1,ny_global
871            DO ji=1,nx_global
872               pcg(ji,jj)=png2(ji+ji_off,jj+jj_off)
873            ENDDO
874         ENDDO
875      ENDIF
876
877#endif
878
879      SELECT CASE ( cd_type )
880         CASE ( 'T' )
881            grid_loc=field_loc_center
882         CASE ( 'F' )                             
883            grid_loc=field_loc_NEcorner
884      END SELECT
885
886      SELECT CASE ( NINT(psgn) )
887         CASE ( -1 )
888            field_type=field_type_vector
889         CASE ( 1 )                             
890            field_type=field_type_scalar
891      END SELECT
892
893#if defined key_nemocice_decomp
894      ! Ensure CICE halos are up to date
895      CALL ice_HaloUpdate (pc, halo_info, grid_loc, field_type)
896#else
897!     D. Scatter pcg to CICE blocks (pc) + update halos
898      CALL scatter_global(pc, pcg, 0, distrb_info, grid_loc, field_type)
899#endif
900
901   END SUBROUTINE nemo2cice
902
903   SUBROUTINE cice2nemo ( pc, pn, cd_type, psgn )
904      !!---------------------------------------------------------------------
905      !!                    ***  ROUTINE cice2nemo  ***
906      !! ** Purpose :   Transfer field in CICE array to field in NEMO array.
907#if defined key_nemocice_decomp
908      !!             
909      !!                NEMO and CICE PE sub domains are identical, hence
910      !!                there is no need to gather or scatter data from
911      !!                one PE configuration to another.
912#else 
913      !!                Automatically deal with scatter/gather between
914      !!                different processors and blocks
915      !! ** Method :    A. Gather CICE blocks (pc) into global array (pcg)
916      !!                B. Map pcg into NEMO global array (png)
917      !!                C. Scatter png into NEMO field (pn) for each processor
918      !!                D. Ensure all haloes are filled in pn
919#endif
920      !!---------------------------------------------------------------------
921
922      CHARACTER(len=1), INTENT( in ) ::   &
923          cd_type       ! nature of pn grid-point
924          !             !   = T or F gridpoints
925      REAL(wp), INTENT( in ) ::   &
926          psgn          ! control of the sign change
927          !             !   =-1 , the sign is modified following the type of b.c. used
928          !             !   = 1 , no sign change
929      REAL(wp), DIMENSION(jpi,jpj) :: pn
930
931#if defined key_nemocice_decomp
932      INTEGER (int_kind) :: &
933         field_type,        & ! id for type of field (scalar, vector, angle)
934         grid_loc             ! id for location on horizontal grid
935                              ! (center, NEcorner, Nface, Eface)
936#else
937      REAL (kind=dbl_kind), dimension(nx_global,ny_global) :: pcg
938#endif
939
940      REAL (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) :: pc
941
942      INTEGER  ::   ji, jj, jn                      ! dummy loop indices
943
944
945#if defined key_nemocice_decomp
946
947      SELECT CASE ( cd_type )
948         CASE ( 'T' )
949            grid_loc=field_loc_center
950         CASE ( 'F' )                             
951            grid_loc=field_loc_NEcorner
952      END SELECT
953
954      SELECT CASE ( NINT(psgn) )
955         CASE ( -1 )
956            field_type=field_type_vector
957         CASE ( 1 )                             
958            field_type=field_type_scalar
959      END SELECT
960
961      CALL ice_HaloUpdate (pc, halo_info, grid_loc, field_type)
962
963
964      pn(:,:)=0.0
965      DO jj=1,jpjm1
966         DO ji=1,jpim1
967            pn(ji,jj)=pc(ji+1-ji_off,jj+1-jj_off,1)
968         ENDDO
969      ENDDO
970
971#else
972
973!      A. Gather CICE blocks (pc) into global array (pcg)
974
975      CALL gather_global(pcg, pc, 0, distrb_info)
976
977!     B. Map pcg into NEMO global array (png)
978
979! Need to make sure this is robust to changes in NEMO halo rows....
980! (may be OK but not spent much time thinking about it)
981! Note that non-existent pcg elements may be used below, but
982! the lbclnk call on pn will replace these with sensible values
983
984      IF (nproc==0) THEN
985         png(:,:,:)=0.0
986         DO jn=1,jpnij
987            DO jj=nldjt(jn),nlejt(jn)
988               DO ji=nldit(jn),nleit(jn)
989                  png(ji,jj,jn)=pcg(ji+nimppt(jn)-1-ji_off,jj+njmppt(jn)-1-jj_off)
990               ENDDO
991            ENDDO
992         ENDDO
993      ENDIF
994
995!     C. Scatter png into NEMO field (pn) for each processor
996
997      IF ( jpnij > 1) THEN
998         CALL mppsync
999         CALL mppscatter (png,0,pn) 
1000         CALL mppsync
1001      ELSE
1002         pn(:,:)=png(:,:,1)
1003      ENDIF
1004
1005#endif
1006
1007!     D. Ensure all haloes are filled in pn
1008
1009      CALL lbc_lnk( pn , cd_type, psgn )
1010
1011   END SUBROUTINE cice2nemo
1012
1013#else
1014   !!----------------------------------------------------------------------
1015   !!   Default option           Dummy module         NO CICE sea-ice model
1016   !!----------------------------------------------------------------------
1017CONTAINS
1018
1019   SUBROUTINE sbc_ice_cice ( kt, nsbc )     ! Dummy routine
1020      WRITE(*,*) 'sbc_ice_cice: You should not have seen this print! error?', kt
1021   END SUBROUTINE sbc_ice_cice
1022
1023   SUBROUTINE cice_sbc_init (nsbc)    ! Dummy routine
1024      WRITE(*,*) 'cice_sbc_init: You should not have seen this print! error?'
1025   END SUBROUTINE cice_sbc_init
1026
1027   SUBROUTINE cice_sbc_final     ! Dummy routine
1028      WRITE(*,*) 'cice_sbc_final: You should not have seen this print! error?'
1029   END SUBROUTINE cice_sbc_final
1030
1031#endif
1032
1033   !!======================================================================
1034END MODULE sbcice_cice
Note: See TracBrowser for help on using the repository browser.