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

Context Navigation

source: branches/UKMO/dev_r5107_hadgem3_cplseq/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 5482

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

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

Download in other formats: