Context Navigation

sbcice_cice.F90 @ 9828

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

Download in other formats: