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

Context Navigation

source: branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 10394

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

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

Download in other formats: