Context Navigation

sbcice_cice.F90 @ 8731

Last change on this file since 8731 was 8731, checked in by dancopsey, 6 years ago
Merged in NEMO3.6 version of this branch (dev_r5518_GC3_couple_pkg)
File size: 44.9 KB

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

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

Context Navigation

source: branches/UKMO/dev_r8183_GC_couple_pkg/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 8731

Download in other formats: