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/UKMO/nemo_v3_6_STABLE_pkg/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

Context Navigation

source: branches/UKMO/nemo_v3_6_STABLE_pkg/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 6249

Last change on this file since 6249 was 6249, checked in by frrh, 8 years ago

Merge in

rev 5662:5668 from branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8
since we HAVE to do this before we can merge in
branches/UKMO/dev_r5518_coupling_GSI7_GSI8_landice
because that is a branch of a branch and directly tryingto merge
that will make a complete mess of things.

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

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

Download in other formats: