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 NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/SBC – NEMO

source: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/SBC/sbcice_cice.F90 @ 12656

Last change on this file since 12656 was 12656, checked in by techene, 4 years ago

change key_QCO into key_qco, stepLF & traatf: add substitute, dynatfQCO: remove pe3 arguments

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