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 @ 5033

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

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

Download in other formats: