New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

sbcice_cice.F90 in NEMO/branches/2020/dev_r13923_Tiling_Cleanup_MPI3_LoopFusion/src/OCE/SBC – NEMO

Context Navigation

source: NEMO/branches/2020/dev_r13923_Tiling_Cleanup_MPI3_LoopFusion/src/OCE/SBC/sbcice_cice.F90 @ 13963

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

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

Download in other formats: