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

Context Navigation

source: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90 @ 4616

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

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

Download in other formats: