Context Navigation

sbcice_cice.F90 @ 8879

Last change on this file since 8879 was 8879, checked in by frrh, 6 years ago
Merge in http://fcm3/projects/NEMO.xm/log/branches/UKMO/dev_r8183_ICEMODEL_svn_removed revisions 8738:8847 inclusive.
File size: 44.4 KB

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

Download in other formats: