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

Context Navigation

source: branches/2015/dev_r5021_UKMO1_CICE_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 5445

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

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

Download in other formats: