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 NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/SBC – NEMO

Context Navigation

source: NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/SBC/sbcice_cice.F90 @ 12624

Last change on this file since 12624 was 12624, checked in by techene, 4 years ago
all: add e3 substitute and limit precompiled files lines to about 130 character, change key_LF into key_QCO, change module name (dynatfQCO, traatfQCO, stepLF)
Property svn:keywords set to `Id`
File size: 42.9 KB

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

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

Download in other formats: