Context Navigation

sbcice_cice.F90 @ 8882

Last change on this file since 8882 was 8882, checked in by flavoni, 6 years ago
dev_CNRS_2017 branch: merged dev_r7881_ENHANCE09_RK3 with trunk r8864
Property svn:keywords set to `Id`
File size: 44.5 KB

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

Note: See TracBrowser for help on using the repository browser.

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

source: branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 8882

Download in other formats: