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

Context Navigation

source: branches/UKMO/dev_r6501_GO6_package_trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 6507

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

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

Download in other formats: