The Relationship Between Brightest Cluster Galaxy Star Formation and the Intracluster Medium in CLASH

TL;DR

This study reveals a strong link between the thermodynamic state of the intracluster medium and star formation activity in brightest cluster galaxies, highlighting the role of AGN feedback in regulating galaxy growth.

Contribution

It establishes a tight correlation between BCG star formation rates and the ICM cooling time to free-fall time ratio, providing observational evidence for AGN-induced condensation processes.

Findings

BCG SFR correlates with t_cool/t_ff ratio.

Starburst durations relate to ICM cooling times.

BCGs with high SFRs have dust masses similar to extreme starbursts.

Abstract

We study the nature of feedback mechanisms in the 11 CLASH brightest cluster galaxies (BCGs) that exhibit extended ultraviolet and nebular line emission features. We estimate star formation rates (SFRs), dust masses, and starburst durations using a Bayesian photometry-fitting technique that accounts for both stellar and dust emission from the UV through far IR. By comparing these quantities to intracluster medium (ICM) cooling times and free-fall times derived from X-ray observations and lensing estimates of the cluster mass distribution, we discover a tight relationship between the BCG SFR and the ICM cooling time to free-fall time ratio, $t_{cool}/t_{ff}$, with an upper limit on the intrinsic scatter of 0.15 dex. Furthermore, starburst durations may correlate with ICM cooling times at a radius of $0.025R_{500}$, and the two quantities converge upon reaching the Gyr regime. Our results provide a direct observational link between the thermodynamical state of the ICM and the intensity and duration of BCG star formation activity, and appear consistent with a scenario where active galactic nuclei (AGN) induce condensation of thermally unstable ICM overdensities that fuel long-duration ($>$ 1 Gyr) BCG starbursts. This scenario can explain (a) how gas with a low cooling time is depleted without causing a cooling flow and (b) the scaling relationship between SFR and $t_{cool}/t_{ff}$. We also find that the scaling relation between SFR and dust mass in BCGs with SFRs $<100$ M$_{\odot}$ yr$^{-1}$ is similar to star-forming field galaxies; BCGs with large ($>100$ M$_{\odot}$ yr$^{-1}$) SFRs have dust masses comparable to extreme starbursts.

Tables

Table 1. Table 1 Herschel Observations

	Instrument	Observation ID	Exposure Time
BCG			(seconds)
Abell 383	PACS	1342189151	7704
		1342189152	7704
		1342189153	7185
		1342189154	7185
	SPIRE	1342189503	5803
		1342201147	3172
MACS0329.7$-$0211	PACS	1342249280	7704
		1342249281	7704
	SPIRE	1342214564	169
		1342239844	1411
MACS0429.6$-$0253	PACS	1342250641	7704
		1342250836	7704
	SPIRE	1342239932	169
		1342241124	1411
MACS1115.9+0219	PACS	1342247672	7704
		1342247691	7704
	SPIRE	1342223226	169
		1342256866	1411
MACS1423.8+2404	PACS	1342188215	9850
		1342188216	9850
	SPIRE	1342188159	6636
MACS1720.3+3536	PACS	1342243800	7704
		1342243801	7704
	SPIRE	1342229601	169
		1342239976	1411
MACS1931.8$-$2653	PACS	1342241619	7704
		1342241681	7704
	SPIRE	1342215993	169
		1342254639	1411
MS2137$-$2353	PACS	1342187803	9850
		1342187804	9850
	SPIRE	1342195938	5786
RXJ1347.5$-$1145	PACS	1342213836	9420
		1342213837	9420
	SPIRE	1342201256	859
		1342201257	859
		1342201258	859
		1342201259	859
		1342201260	859
		1342201261	859
		1342201262	859
		1342201263	859
		1342247859	1952
		1342247860	1952
		1342247861	1952
RXJ1532.9+3021	PACS	1342258435	7704
		1342258435	7704
	SPIRE	1342234776	169
		1342261681	1411
RXJ2129.7+0005	PACS	1342187256	9569
		1342187257	9569
		1342196791	7185
		1342196792	7185
	SPIRE	1342188167	6636
		1342209312	3172

Table 2. Table 2 SED Fitting Parameters

Stellar Population Model
Synthetic Stellar Population	Bruzual & Charlot (2003)
Initial Mass Function	Salpeter (1955)
Attenuation Law	Calzetti et al. (2000)
Dust Emission	Draine & Li (2007)
	Minimum	Maximum	Sampling
Model Parameter Space Constraints	Value	Value	Interval
Old Stellar Population
Age, $t$	6 Gyr	$z_{Age}^a$	Linear
Decay Rate, $\tau$	$0.05t$	$0.2t$	Linear
Metallicity	$3\times {10}^{-2}{Z}_{\odot }$	$1.5Z_{\odot }$	Linear
Burst Population
Burst Duration, $\Delta t_b$	${10}^{-2}$ Gyr	5.0 Gyr	Logarithmicb
Burst Decay Percentage	0.01	0.99	Linear
Relative Burst Mass, fburstc	0.0016	6.4	Logarithmic
Dust Parameters
Attenuation AV	0	2	Linear
PAH Abundance Index $q_{PAH}$	0.10	4.58	Lineard
${\gamma }^e$	0.0	1.0	Linear
$U_{min}^e$	0.10	25.0	Logarithmic
$U_{max}^e$	${10}^3$	${10}^7$	Logarithmic

Table 3. Table 3 BCG Stellar and Dust Parameters

	log10 SFR	log10 $\Delta t_b$	log10 $\Delta t_b$	log10 $M_d$	A${}_{\text{V}}$	$R_S$a	$R_S$	$R_S$
			restr. fburstb			SFR$\times \Delta t_b$	SFR$\times$ $M_d$	M${}_{\text{d}}{}^{}\times \Delta t_b$
BCG	(M⊙ yr-1)	(Gyr)	(Gyr)	(M⊙)	(mag)
Abell 383	${0.18}_{-0.25}^{+0.24}$	$>-{0.49}^d$	${0.19}_{-0.34}^{+0.35}$	${7.47}_{-0.87}^{+0.84}$	${0.48}_{-0.24}^{+0.22}$	-0.16	0.36	-0.03
MACS0329.7-0211	${1.6}_{-0.2}^{+0.18}$	$>-0.73$	${0.01}_{-0.40}^{+0.37}$	${8.4}_{-0.73}^{+0.71}$	${0.56}_{-0.19}^{+0.18}$	-0.06	0.07	0.07
MACS0429.6-0253	${1.53}_{-0.24}^{+0.23}$	$>-0.36$	${0.09}_{-0.37}^{+0.36}$	${8.49}_{-0.68}^{+0.7}$	${0.75}_{-0.24}^{+0.23}$	0.03	0.20	0.07
MACS1115.9+0219	${0.85}_{-0.28}^{+0.29}$	$>-0.40$	${0.19}_{-0.36}^{+0.35}$	${7.6}_{-0.81}^{+0.85}$	${0.44}_{-0.3}^{+0.3}$	-0.03	0.44	-0.01
MACS1423.8+2404	${1.41}_{-0.18}^{+0.18}$	$>-0.65$	$>-0.03$	${8.47}_{-0.72}^{+0.69}$	${0.42}_{-0.18}^{+0.18}$	-0.11	0.25	0.02
MACS1720.3+3536	${0.19}_{-0.26}^{+0.27}$	$>-0.34$	$>-0.01$	${7.6}_{-0.71}^{+0.67}$	${0.6}_{-0.26}^{+0.26}$	-0.10	0.10	-0.0
MACS1931.8-2653	${2.42}_{-0.18}^{+0.22}$	$-{0.97}_{-0.56}^{+0.53}$	$-{1.01}_{-0.35}^{+0.34}$	${8.88}_{-0.4}^{+0.39}$	${0.87}_{-0.21}^{+0.21}$	-0.52	-0.23	0.19
MS2137-2353	${0.25}_{-0.29}^{+0.29}$	$>-0.43$	$>0.08$	${7.41}_{-0.87}^{+0.86}$	${0.43}_{-0.31}^{+0.31}$	-0.23	0.47	-0.21
RXJ1347.5-1145	${1.07}_{-0.23}^{+0.23}$	$>-0.39$	${0.13}_{-0.35}^{+0.36}$	${7.76}_{-0.9}^{+0.96}$	${0.31}_{-0.23}^{+0.22}$	0.03	0.47	-0.10
RXJ1532.9+3021	${1.99}_{-0.19}^{+0.2}$	$-{0.29}_{-0.59}^{+0.54}$	$-{0.39}_{-0.41}^{+0.38}$	${8.77}_{-0.51}^{+0.47}$	${0.9}_{-0.19}^{+0.19}$	-0.40	-0.16	0.11
RXJ2129.7+0005	$-{0.5}_{-0.23}^{+0.22}$	${0.29}_{-0.31}^{+0.28}$	$>0.25$	${6.81}_{-0.71}^{+0.87}$	${0.28}_{-0.17}^{+0.14}$	-0.19	0.12	0.01

Table 4. Table 4 Best Fit χ 2 superscript 𝜒 2 \chi^{2} Values

	Calzetti	Modified Calzetti	Witt Clumpy SMC	Chabrier IMF
BCG	${\chi }_{\nu }^2$	${\chi }_{\nu }^2$	${\chi }_{\nu }^2$	${\chi }_{\nu }^2$
Abell 383	0.99	1.33	1.54	1.41
MACS0329.7$-$0211	0.57	0.76	0.60	0.64
MACS0429.6$-$0253	1.53	1.53	2.25	1.39
MACS1115.9+0219	0.59	0.53	0.46	0.52
MACS1423.8+2404	0.80	1.44	1.24	0.80
MACS1720.3+3536	1.41	1.48	1.02	1.31
MACS1931.8$-$2653	0.66	1.39	1.22	1.17
MS2137$-$2353	0.64	0.70	0.68	0.67
RXJ1347.5$-$1145	0.72	0.72	0.66	0.78
RXJ1532.9+3021	1.63	1.51	1.28	1.40
RXJ2129.7+0005	0.70	0.73	0.79	0.73

Table 5. Table 5 BCG Stellar and Dust Parameters assuming a Chabrier IMF and Calzetti attenuation law

	log10 SFR	log10 $\Delta t_b$	log10 $M_d$	A${}_{\text{V}}$	$R_S$a	$R_S$	$R_S$
					SFR$\times \Delta t_b$	SFR$\times$ $M_d$	M${}_{\text{d}}{}^{}\times \Delta t_b$
BCG	(M⊙ yr-1)	(Gyr)	(M⊙)	(mag)
Abell 383	$-{0.02}_{-0.24}^{+0.23}$	$>-0.47$	${7.39}_{-0.82}^{+0.78}$	${0.48}_{-0.24}^{+0.22}$	-0.11	0.30	0.04
MACS0329.7-0211	${1.39}_{-0.22}^{+0.19}$	$>-0.68$	${8.39}_{-0.64}^{+0.66}$	${0.53}_{-0.19}^{+0.18}$	0.07	0.01	-0.07
MACS0429.6-0253	${1.32}_{-0.22}^{+0.23}$	$>-0.73$	${8.49}_{-0.65}^{+0.67}$	${0.74}_{-0.22}^{+0.24}$	0.10	0.17	0.09
MACS1115.9+0219	${0.65}_{-0.26}^{+0.27}$	$>-0.75$	${7.54}_{-0.88}^{+0.96}$	${0.43}_{-0.29}^{+0.29}$	-0.08	0.45	-0.04
MACS1423.8+2404	${1.21}_{-0.2}^{+0.19}$	$>-0.71$	${8.42}_{-0.75}^{+0.74}$	${0.41}_{-0.18}^{+0.17}$	-0.10	0.18	0.04
MACS1720.3+3536	${0.01}_{-0.23}^{+0.26}$	${0.13}_{-0.41}^{+0.4}$	${7.61}_{-0.68}^{+0.74}$	${0.62}_{-0.26}^{+0.26}$	-0.19	0.17	-0.0
MACS1931.8-2653	${2.21}_{-0.19}^{+0.19}$	$-{0.93}_{-0.51}^{+0.49}$	${8.86}_{-0.39}^{+0.38}$	${0.88}_{-0.21}^{+0.19}$	-0.43	-0.27	0.11
MS2137-2353	${0.05}_{-0.29}^{+0.27}$	${0.15}_{-0.44}^{+0.4}$	${7.43}_{-0.86}^{+0.8}$	${0.44}_{-0.3}^{+0.32}$	-0.27	0.36	-0.15
RXJ1347.5-1145	${0.87}_{-0.22}^{+0.22}$	$>-0.66$	${7.74}_{-0.81}^{+0.91}$	${0.31}_{-0.22}^{+0.23}$	-0.11	0.45	0.04
RXJ1532.9+3021	${1.8}_{-0.18}^{+0.2}$	$-{0.2}_{-0.55}^{+0.51}$	${8.79}_{-0.48}^{+0.47}$	${0.91}_{-0.21}^{+0.21}$	-0.34	-0.10	0.07
RXJ2129.7+0005	$-{0.7}_{-0.21}^{+0.24}$	$>-0.01$	${6.71}_{-0.73}^{+0.81}$	${0.28}_{-0.14}^{+0.14}$	-0.30	0.27	-0.19

Table 6. Table 6 BCG Stellar and Dust Parameters assuming a modified Calzetti attenuation law

	log10 SFR	log10 $\Delta t_b$	log10 $M_d$	A${}_{\text{V}}$	$R_S$a	$R_S$	$R_S$
					SFR$\times \Delta t_b$	SFR$\times$ $M_d$	M${}_{\text{d}}{}^{}\times \Delta t_b$
BCG	(M⊙ yr-1)	(Gyr)	(M⊙)	(mag)
Abell 383	${0.18}_{-0.26}^{+0.24}$	$>-{0.54}^b$	${7.41}_{-0.78}^{+0.83}$	${0.47}_{-0.26}^{+0.25}$	-0.13	0.35	0.05
MACS0329.7-0211	${1.6}_{-0.21}^{+0.19}$	$>-0.31$	${8.39}_{-0.69}^{+0.72}$	${0.53}_{-0.18}^{+0.19}$	-0.11	0.08	0.04
MACS0429.6-0253	${1.54}_{-0.25}^{+0.24}$	$>-0.80$	${8.5}_{-0.75}^{+0.75}$	${0.76}_{-0.24}^{+0.24}$	0.02	0.26	-0.01
MACS1115.9+0219	${0.85}_{-0.28}^{+0.28}$	$>-0.38$	${7.61}_{-0.82}^{+0.87}$	${0.43}_{-0.28}^{+0.28}$	-0.0	0.42	0.03
MACS1423.8+2404	${1.4}_{-0.23}^{+0.21}$	$>-0.76$	${8.37}_{-0.78}^{+0.77}$	${0.42}_{-0.2}^{+0.19}$	0.02	0.21	-0.02
MACS1720.3+3536	${0.19}_{-0.24}^{+0.24}$	$>-0.37$	${7.64}_{-0.74}^{+0.68}$	${0.59}_{-0.24}^{+0.23}$	-0.17	0.27	-0.04
MACS1931.8-2653	${2.41}_{-0.19}^{+0.21}$	$-{0.95}_{-0.58}^{+0.52}$	${8.89}_{-0.37}^{+0.37}$	${0.85}_{-0.2}^{+0.2}$	-0.53	-0.25	0.15
MS2137-2353	${0.27}_{-0.27}^{+0.3}$	${0.14}_{-0.52}^{+0.44}$	${7.43}_{-0.78}^{+0.75}$	${0.44}_{-0.34}^{+0.33}$	-0.21	0.36	0.0
RXJ1347.5-1145	${1.07}_{-0.21}^{+0.23}$	$-{0.11}_{-0.61}^{+0.59}$	${7.71}_{-0.87}^{+0.92}$	${0.31}_{-0.22}^{+0.23}$	0.05	0.37	0.09
RXJ1532.9+3021	${1.99}_{-0.19}^{+0.2}$	$-{0.26}_{-0.55}^{+0.57}$	${8.77}_{-0.49}^{+0.51}$	${0.88}_{-0.19}^{+0.2}$	-0.43	-0.12	0.03
RXJ2129.7+0005	$-{0.53}_{-0.24}^{+0.24}$	$>0.06$	${6.65}_{-0.73}^{+0.83}$	${0.26}_{-0.17}^{+0.18}$	-0.14	0.04	0.18

Table 7. Table 7 BCG Stellar and Dust Parameters assuming Witt clumpy SMC-like dust in a shell geometry

	log10 SFR	log10 $\Delta t_b$	log10 $M_d$	A${}_{\text{V}}$	$R_S$a	$R_S$	$R_S$
					SFR$\times \Delta t_b$	SFR$\times$ $M_d$	M${}_{\text{d}}{}^{}\times \Delta t_b$
BCG	(M⊙ yr-1)	(Gyr)	(M⊙)	(mag)
Abell 383	${0.36}_{-0.24}^{+0.25}$	$>-0.39$	${7.6}_{-0.78}^{+0.75}$	${0.9}_{-0.56}^{+0.61}$	-0.23	0.27	-0.06
MACS0329.7-0211	${1.61}_{-0.19}^{+0.19}$	$-{0.31}_{-0.83}^{+0.74}$	${8.41}_{-0.69}^{+0.7}$	${0.58}_{-0.42}^{+0.39}$	-0.19	-0.08	-0.01
MACS0429.6-0253	${1.57}_{-0.2}^{+0.21}$	$>-{1.24}^b$	${8.54}_{-0.64}^{+0.65}$	${1.05}_{-0.57}^{+0.56}$	-0.10	0.17	0.01
MACS1115.9+0219	${0.98}_{-0.26}^{+0.26}$	$>-1.06$	${7.84}_{-0.81}^{+0.82}$	${0.67}_{-0.52}^{+0.57}$	-0.15	0.27	-0.04
MACS1423.8+2404	${1.51}_{-0.19}^{+0.19}$	$-{0.27}_{-0.77}^{+0.72}$	${8.49}_{-0.75}^{+0.73}$	${0.56}_{-0.42}^{+0.4}$	-0.15	0.05	-0.01
MACS1720.3+3536	${0.43}_{-0.26}^{+0.26}$	$-{0.1}_{-0.62}^{+0.57}$	${7.66}_{-0.71}^{+0.7}$	${1.12}_{-0.62}^{+0.6}$	-0.19	0.11	-0.01
MACS1931.8-2653	${2.36}_{-0.18}^{+0.18}$	$-{1.12}_{-0.56}^{+0.54}$	${8.88}_{-0.37}^{+0.37}$	${1.09}_{-0.52}^{+0.51}$	-0.33	-0.23	0.03
MS2137-2353	${0.55}_{-0.31}^{+0.31}$	$-{0.07}_{-0.61}^{+0.58}$	${7.66}_{-0.82}^{+0.79}$	${1.01}_{-0.68}^{+0.68}$	-0.21	0.31	-0.08
RXJ1347.5-1145	${1.21}_{-0.21}^{+0.19}$	$-{0.26}_{-0.76}^{+0.68}$	${8.08}_{-0.78}^{+0.78}$	${0.43}_{-0.33}^{+0.32}$	-0.09	0.23	-0.04
RXJ1532.9+3021	${1.98}_{-0.18}^{+0.18}$	$-{0.49}_{-0.81}^{+0.67}$	${8.79}_{-0.5}^{+0.5}$	${1.18}_{-0.5}^{+0.46}$	-0.36	-0.13	0.04
RXJ2129.7+0005	$-{0.24}_{-0.24}^{+0.27}$	$>-0.32$	${6.82}_{-0.71}^{+0.75}$	${0.83}_{-0.62}^{+0.67}$	-0.31	0.07	0.04

Equations

$${\chi }^2=\sum _{i=1}^N\frac{{\left(F_i-\mathcal{A}{M}_i\right)}^2}{{\sigma }_i^2},$$

$$\begin{array}{rl}{t}_{cool}\left(r\right)& =\frac{3}{2}\frac{n\left(r\right)k_{B}T\left(r\right)}{n_e\left(r\right)n_H\left(r\right)\Lambda \left[T\left(r\right),Z\left(r\right)\right]}\\ & =\frac{6.9}{2}\frac{k_{B}T\left(r\right)}{n_e\left(r\right)\Lambda \left[T\left(r\right),Z\left(r\right)\right]},\end{array}$$

$$M_{enc}\left(r\right)=4\pi {\rho }_0{r}_s^3\left[\ln \left(\frac{r_s+r}{r_s}\right)-\frac{r}{r_s+r}\right],$$

$$t_{ff}\left(r\right)=\sqrt{\frac{2r^3}{G{M}_{enc}\left(r\right)}}.$$

$$\frac{M_{\text{burst}}}{M_{\text{aperture}}}=\left(\frac{M_{\text{burst}}}{M_{\text{BCG}}}\right)\times \left(\frac{L_{F160,\text{BCG}}}{L_{F160,\text{aperture}}}\right),$$

$${\log }_{10}\frac{t_{cool}}{t_{ff}}=\left(1.6\pm 0.4\right)-\left(0.15\pm 0.03\right){\log }_{10}\frac{\text{SFR}}{{\text{M}}_{\odot }{\text{ yr}}^{-1}},$$

$${\log }_{10}\frac{M_d}{{\text{M}}_{\odot }}=\left(7.1_{-0.3}^{+0.4}\right)+\left(0.97_{0.24}^{+0.34}\right){\log }_{10}\frac{\text{SFR}}{{\text{M}}_{\odot }{\text{ yr}}^{-1}}$$

$$\begin{array}{rl}{ϵ}_{MD}& ={\int }_{{\delta }_c\left(\frac{t_{cool}}{t_{ff}}\right)}^{\infty }\frac{1}{\sqrt{2\pi }\sigma \rho }{e}^{-\frac{{\left(\ln \rho \right)}^2}{2{\sigma }^2}}d\rho \\ & =\frac{1}{2}Erfc\left(\frac{\ln {\delta }_c\left(\frac{t_{cool}}{t_{ff}}\right)}{\sqrt{2}\sigma }\right),\end{array}$$

$${\log }_{10}{ϵ}_{MD}\sim {\log }_{10}\text{SFR}-{\log }_{10}{\dot{M}}_g=\frac{t_{cool}/t_{ff}-10}{-17.1_{-2.3}^{+1.9}},$$

Full text

††slugcomment: Accepted for publication in the Astrophysical Journal

1. Introduction

The brightest cluster galaxy (BCG) in the center of the Perseus cluster, NGC1275, was long known to be an emission-line system with ionized hydrogen gas spanning $\sim 100$ kpc (Lynds, 1970). While BCGs are typically quiescent systems, the correlation of emission-line BCGs inside so-called “cooling flow” clusters was noted by the early studies of X-ray clusters and the optical spectra of their BCGs (e.g. Heckman, 1981; Cowie et al., 1983; Hu et al., 1985; Johnstone et al., 1987). Such clusters, more recently termed “cool core” clusters since their central gas X-ray temperatures are somewhat cooler than their outskirts, have highly peaked central X-ray surface brightness profiles arising from relatively high central gas densities. Initially, it was thought that since the cooling time for this gas (the ratio of the thermal energy content to the radiative loss rate) was short compared to the Hubble time, that such gas would cool, lose pressure, and gradually allow more gas to settle into the center of such a cluster: a “cooling flow”. The rates inferred were 100-1000 M*⊙* yr*-1* in some cases (see, e.g. Fabian (1994) for a review). Subsequent high-resolution X-ray spectroscopy from XMM-Newton convincingly showed that this simple model was incorrect (e.g. Peterson et al., 2003; Peterson & Fabian, 2006). Nevertheless, observations of the BCGs in samples of X-ray luminous galaxy clusters have revealed that in up to 70$\%$ of cool-core clusters, the otherwise quiescent elliptical BCG shows signs of ultraviolet and nebular line emission. This emission, in some cases, is consistent with star formation rates (SFRs) on the order of 100 M*⊙* yr*-1* (e.g. McNamara & O’Connell, 1989; Crawford et al., 1999; Edwards et al., 2007; Donahue et al., 2010). This activity is observed to be related to the presence of low entropy intracluster gas (the intracluster medium, or ICM) in the cluster core. Since low-entropy ICM gas has a short cooling time, this correlation suggests that the observed star formation is fueled by cold gas that has condensed from a hot gas reservoir, however at a far gentler rate than the simple cooling flow model predicted (Rafferty et al., 2006, 2008; Cavagnolo et al., 2008, 2009; McDonald et al., 2010; Hoffer et al., 2012). One viable source of heat replacing energy lost via radiative cooling is energy injected by active galactic nuclei (AGN) into the ICM. The mechanical work indicated by the size and ICM gas pressure surrounding X-ray cavities filled by radio-emitting plasma demonstrates that there is sufficient energy to counter radiative cooling (Rafferty et al., 2006; McNamara & Nulsen, 2007; Fabian, 2012; McNamara & Nulsen, 2012).

There is a growing body of evidence supporting AGN feedback-driven precipitation and condensation as being responsible for balancing heating and cooling in cool-core galaxy clusters. AGN jet feedback is thought to inject energy into the ICM, thus offsetting cooling while also triggering condensation of thermodynamically unstable volumes of ICM plasma, which then precipitate onto the BCG and fuel star formation and further AGN activity (Voit et al., 2015; Voit & Donahue, 2015; Li & Bryan, 2014; Li et al., 2015; Gaspari et al., 2012, 2013, 2015; Fogarty et al., 2015; Gaspari et al., 2016; Prasad et al., 2015; Yang & Reynolds, 2016; Meece et al., 2016). Recently, the precipitation aspect of this model was shown to be consistent with observations of molecular gas accretion onto the AGN in the Abell 2597 BCG (Tremblay et al., 2016). The morphologies of UV and nebular emission structures in observed BCGs in cool core clusters (e.g. Donahue et al., 2015; Tremblay et al., 2015) are reproduced in simulations of AGN feedback-regulated condensation and precipitation in low-entropy ICM gas (Li et al., 2015; Gaspari et al., 2017). Furthermore, submillimeter observations of cool core clusters reveal reservoirs of as much as $\sim 10^{11}$ M*⊙* of molecular gas (Edge et al., 2002; O’Dea et al., 2008; Russell et al., 2016). The picture that is emerging is one of a complex feedback-driven interaction between the BCG and ICM in the cluster core which produces substantial gas condensation but suppresses runaway ICM cooling.

The condensation model follows thermal instabilities that are both triggered and regulated by AGN feedback. In this model, the hot ICM gas is near hydrostatic equilibrium with the gravitational potential, and the gas entropy at large cluster-centric distances is governed by cosmological processes. At very small radii close to the central AGN, the gas can develop a nearly flat entropy profile that allows thermal instabilities to rapidly grow. However, Meece et al. (2016) showed that even a modest fraction of AGN feedback in the form of mechanical jets can transport energy beyond the region local to the AGN. Rather than resulting in unphysical catastrophic cooling, such systems can self-regulate. The resulting entropy profile within 5-10 kpc of the simulated AGN with jet feedback is somewhat shallower than in the outer parts of the cluster, but condensation of low-entropy inhomogeneities which have been uplifted to greater altitudes can produce multi-phase structure at radii larger than 10 kpc (Meece et al., 2016; Voit et al., 2017). Key to this model is the prediction that ICM plasma with a cooling-to-freefall time ratio ($t_{cool}/t_{ff}$) $\lesssim 10$ is sufficiently unstable that jets trigger condensation. Supporting this prediction are observations indicating that $t_{cool}/t_{ff}\leq 20$ in the center of a galaxy cluster is a good predictor of BCG activity (Voit & Donahue, 2015; Voit et al., 2017). In simulations, material will condense out of the ambient ICM when $t_{cool}/t_{ff}\lesssim 10$ locally (Li & Bryan, 2014; Li et al., 2015).

Alternatively, McNamara et al. (2016) posit that most cold molecular gas condenses in the vicinity of BCGs when AGN jets uplift low entropy plasma from within the BCG to high enough radii for the plasma to condense. In this scenario, jet uplifting can cause core plasma to condense if drag prevents it from sinking back to its original altitude in a time shorter than the cooling time of the gas. Observations of massive (up to $10^{10}$ M*⊙) flows of molecular gas with velocities of several hundred km s-1* in cool core clusters such as Abell 1835 and Abell 1664 suggest that this mode of condensation and precipitation plays an important role in AGN-regulated feedback (McNamara et al., 2014; Russell et al., 2014). Both condensation due to uplifting and condensation in the ambient ICM in a cluster core are studied in models of cluster-scale feedback (Li et al., 2015; Gaspari et al., 2016; Voit et al., 2017).

Molendi et al. (2016) propose a possible scenario where star formation in BCGs lags behind cooling by up to a Gyr, owing to molecular gas in the environs of the BCG forming stars at only a few percent efficiency, as one possible explanation for the presence of star forming BCGs in the absence of cooling flows. This scenario can explain relationships between BCG SFRs and the thermodynamical state of the ICM, but would imply that the density and temperature of the hot ICM in the vicinity of the BCG is related to the amount of gas that had cooled in a previous cooling phase.

In this paper, we analyze the spectral energy distributions of the UV-luminous BCGs in the full X-ray selected cluster sample from the CLASH111Cluster Lensing And Supernova survey with Hubble program. Our aim is to explore physical connections between BCG star formation activity and the properties of the surrounding ICM and constrain models for possible mechanisms of condensation in the ICM. We take full advantage deep, multi-band imaging of the relatively uniform and massive sample of galaxy clusters observed by CLASH and we are able to reveal relationships that may have previously been confounded by analyses of less deeply observed and/or less uniform samples. Specifically, by studying the BCGs of cool-core clusters in CLASH, we explore how cooling and thermal instability timescales ($t_{cool}$ and $t_{cool}/t_{ff}$) in the ICM may relate to BCG star formation.

Our paper is structured as follows: in Section 2 we describe the observational data set used to construct the spectral energy distributions (SEDs) for the 11 active CLASH BCGs studied in this sample. In Section 3, we describe constructing and fitting UV-through-FIR SEDs, as well as the X-ray derived parameters we use to study the ICM. In Section 4 we present our results, which we discuss in Section 5. We summarize our conclusions in Section 6. Throughout our analysis we adopt a $\Lambda$CDM cosmology with $\Omega_{m}=0.3$, $\Omega_{\Lambda}=0.7$, $H_{0}$ = 70 km/s/Mpc, and $h=0.7$.

2. Observations

The near UV through far IR SEDs in this paper are based on photometry from the CLASH HST data set in combination with mid- and far-IR data from Spitzer and Herschel. X-ray data used for measuring the temperature, density and metallicity profiles of the ICM were taken from the Chandra archive.

2.1. HST: UV through near IR photometry

A detailed summary of the science-level data products for CLASH may be found in Postman et al. (2012). For the 11 X-ray selected clusters with evidence of BCG star formation activity, we used the CLASH photometric data set covering 16 filters spanning from $\sim$2000-17000 Å in the observer frame. We used drizzled mosaics with a 65 milliarcsecond pixel scale, the same image data used in Fogarty et al. (2015). These data are publicly available via MAST HLSP 222https://archive.stsci.edu/prepds/clash/. Drizzled mosaics were constructed using the MosaicDrizzle pipeline (Koekemoer et al., 2011). In keeping with our previous work, we calculated a single Milky Way foreground reddening correction for each BCG in each filter using the Schlegel et al. (1998) dust maps. Our drizzled images are background corrected using an iterative 3-sigma clipping technique. Owing to additional uncertainty in the flat-fielding of WFC3/UVIS photometry over large spatial scales, discussed in Fogarty et al. (2015), the median flux measured in an annulus around each BCG was subtracted from each of the UV filters as well.

2.2. Spitzer: mid-IR photometry

Spitzer/IRAC 3.6 $\mu$m and 4.5 $\mu$m mosaicked observations and catalogs are available for all of the BCGs studied in this paper, and 5.7 $\mu$m and 7.9 $\mu$m observations are available for Abell 383, MACS1423.8+2404, and RXJ1347.5-1145. Fluxes for Spitzer IR sources in the CLASH fields were taken from the publicly available CLASH/Spitzer catalog 333http://irsa.ipac.caltech.edu/data/SPITZER/CLASH/. We describe aperture selection and our method for correcting for crowded fields in the Spitzer photometry in Section 3.1.

Photometry for each channel was measured on mosaic images generated using the MOPEX software package. The default MOPEX settings for the catalog use the 3.6 $\mu$m channel, and use the fiducial image frame for this channel to generate the mosaics for longer wavelength IRAC channels for each CLASH cluster. The Spitzer mosaic images have a pixel scale of 0.6 arcseconds. Flux values were obtained with Source Extractor in double-image mode using MOPEX weights (Bertin & Arnouts, 1996). The Source Extractor-generated catalogs consist of aperture photometry for sources in fixed-diameter apertures ranging from 2 to 40 pixels. Full details about the parameters used with MOPEX and Source Extractor to generate the CLASH catalog are given in the online documentation444http://irsa.ipac.caltech.edu/data/SPITZER/CLASH/docs/

README.CLASHSpitzer.

2.3. Herschel: far-IR photometry

We used archival data from Herschel/PACS (100 $\mu$m, 160 $\mu$m) and Herschel/SPIRE (250 $\mu$m, 350 $\mu$m, and 500 $\mu$m) to extend BCG photometry into the far-infrared. Table 1 details the observations used and their exposure times. We obtained level 2 archival data using the Herschel Science Archive reduction pipeline for SPIRE observations and the HSA MADMAP reduction pipeline for PACS observations. Observations were co-added and photometric parameters were measured using the HIPE software package. We describe Herschel aperture selection and background subtraction in Section 3.1.

2.4. Chandra: X-ray Observations

Chandra data exists for all CLASH clusters, with exposure times for individual observations ranging from 19.3 ks (for Abell 383) to 115.1 ks (for MACS1423.8+2404). These data were assembled in Donahue et al. (2014) to construct the parameter profiles that we use to analyze the thermodynamics of the ICM in the environs of each BCG. We also examine the archival data for evidence of X-ray loud AGN, in order to determine whether it is important to consider AGN emission models in our SED fits. We found evidence for only one X-ray loud AGN, in the cluster MACS1931.8-2635. The AGN classification is based on the presence of a Chandra point source, with a spectrum showing evidence for hard X-ray ($>$5 keV) emission over that expected from hot gas with $kT\sim 5-7$ keV.

3. Methods

In Fogarty et al. (2015), we identified 11 BCGs in the CLASH X-ray selected sample with evidence of ongoing star formation. In each of these BCGs, we detected extended UV emission in WFC3/UVIS photometry and extended H$\alpha$ + [N II] emission by differencing ACS images. The extended UV and line emission features were shown in Donahue et al. (2015) and Fogarty et al. (2015) to be the site of nebular emission and star formation. These 11 BCGs form the star-forming BCG sample we study in this paper.

3.1. BCG Photometry

We constructed multi-instrument UV-through-IR SEDs for the UV-bright filamentary features in each UV-luminous CLASH BCG. Since we are interested in measuring the star formation properties of the BCGs, for the Hubble UV through near-IR photometry we extracted flux from within the apertures described in Fogarty et al. (2015) (apertures for each BCG are shown in Figure 3 in that paper). These apertures were selected to encompass the region in each BCG exhibiting a UV luminosity of $7.14\times 10^{24}$ erg s*-1* Hz*-1* pix*-2*, corresponding to a star formation surface density of $\Sigma_{\rm SFR}$ $\geq 0.001$ M*⊙* yr*-1* pix*-2*, after accounting for dust reddening. Measuring fluxes in these apertures maximizes the contribution made by the recently formed stellar population to the SED, while it minimizes the contribution made by the passive stellar population in the bulk of the BCGs. Our procedure for measuring fluxes minimizes the risk of underestimating dust attenuation in the star forming regions of the BCGs since we are not averaging dust attenuation over the dusty star-forming and relatively dust-free quiescent parts of the galaxy. We do not match the Hubble point spread functions (PSF) since the photometric aperture sizes we use are much larger than the sizes of the PSFs for the HST passbands.

Estimating mid-IR fluxes using Spitzer was more complicated, since star-forming features are not well resolved spatially in Spitzer photometry. We measured IR fluxes from Spitzer in apertures in the CLASH/Spitzer catalogs that were selected to encompass the BCG while minimizing inclusion of satellite galaxies. Both the old stellar population and dust emission (primarily in the form of PAH features) contribute significantly to the flux in these filters. Since the angular resolution of Spitzer ($>$ 1.45 arcsec) is insufficient to resolve the star forming features we wish to study, we needed to subtract the contribution made by old stellar light from outside our Hubble apertures to the Spitzer IR fluxes. We accomplished this by first measuring the HST/WFC3 F160W flux in the region covered by the Spitzer apertures but not the Hubble apertures. We assumed this flux is due entirely to old stellar light, and estimated (IRAC-F160W) color by modeling old stellar populations experiencing dust attenuation $\leq 0.5$ A${}_{\textrm{V}}$ using the Bayesian SED fitting algorithm iSEDfit (see Section 3.2 for details) for each BCG. We used these colors to scale the F160W fluxes corresponding to the old stellar population outside the Hubble apertures to Spitzer IRAC fluxes and subtracted these scaled fluxes from our Spitzer photometry. The resulting Spitzer IRAC photometry corresponds to the fluxes in the Hubble apertures. An example pair of Spitzer and Hubble apertures is shown in Figure 1.

We note that H2 vibrational modes and other near-IR emission lines have been detected between 5-25$\mu$m in star-forming BCGs (Donahue et al., 2011). However, contamination from emission lines similar to those observed in Donahue et al. (2011) would only contribute a few percent to the mid-IR fluxes derived for the stellar and dust emission in star forming regions, while the uncertainties on these fluxes after subtracting the excess old stellar light component is $\gtrsim 20\%$. Therefore, we do not attempt to estimate the contribution of these lines to the Spitzer mid-IR photometry.

We extracted Herschel PACS and SPIRE photometry from archival data using the HIPE software package. Since Herschel bands are dominated by dust re-emission of UV flux, we assumed that all of the Herschel flux in each BCG comes from the star forming features. We measured photometry in circular apertures large enough to cover the PSF of each Herschel filter. Aperture radii were chosen to be 8$"$, 12$"$ , 18$"$ , 24$"$ , and 36$"$ for the 100 $\mu$m, 160 $\mu$m, 250 $\mu$m, 350 $\mu$m, and 500 $\mu$m filters respectively. To obtain fluxes we took the mean of pairs of cross scans, while we estimated uncertainties using the difference between scans. For each Herschel band, we found typical uncertainties on the order of $\lesssim 2\%$ for the PACS filters and $\sim 10\%$ for the SPIRE filters. Since the BCGs occupy crowded fields, we measured crowded source backgrounds in annuli centered on the BCG with inner radii of 16$"$, 24$"$, 36$"$, 48$"$, 72$"$ for the respective PACS and SPIRE filters and outer radii of 180$"$.

3.1.1 Photometric Errors

Error budgets for the fluxes in our SEDs include uncertainty from counting statistics, systematic uncertainties from weight maps (produced as output from Drizzling in the case of the Hubble data and provided by the level 2 data pipeline in the case of the Herschel PACS and SPIRE data), and the absolute calibration uncertainty for each instrument. In the case of the Spitzer IRAC data, the total counting and formal systematic uncertainties are available in the CLASH/Spitzer catalogs from the fluxes and MOPEX weight maps (see Section 2.2). Calculating uncertainties for the mid-IR fluxes obtained using Spitzer IRAC photometry required propagating the additional uncertainty involved in estimating and subtracting the excess old stellar light component from the mid-IR photometry.

Source confusion noise is an important component of the error budget for SPIRE detections (Nguyen et al., 2010). We incorporate confusion noise terms of 5.8 $\mu$Jy, 6.3 $\mu$Jy, and 6.8 $\mu$Jy from Nguyen et al. (2010) into our estimates of the uncertainty on SPIRE 250 $\mu$m, 350 $\mu$m, and 500 $\mu$m fluxes respectively. In bands where we do not detect significant flux from the BCG, we use the confusion and instrument noise estimates computed in Nguyen et al. (2010) to estimate 3$\sigma$ upper limits.

For HST ACS and WFC3 data, we included a 5$\%$ total absolute and relative calibration uncertainty in our error budget. For the ACS and WFC3/IR filters, the absolute calibration is the dominant term in the overall error budget. Since several of the early CLASH WFC3/UVIS observations were affected by non-uniform flat-fielding at large (hundreds to thousands of pixels) scales, we calculated the scatter in identical apertures placed in empty patches of sky for each UV observation. This extra uncertainty component in the WFC3/UVIS error budget is in most cases similar to the absolute calibration uncertainty, although for several filters with only faintly detected BCG UV emission it is the main source of uncertainty.

We adopted an additional 3$\%$ uncertainty for the Spitzer IRAC fluxes, and 10$\%$ uncertainty for Herschel PACS and SPIRE fluxes. The 3$\%$ figure accounts for the absolute flux calibration uncertainty in Spitzer (Reach et al., 2005). For Herschel, we incorporate the absolute calibration uncertainty ($\sim 5\%$), the relative calibration uncertainty ($\sim 2\%$), and allow for an additional factor of $\sim 2$ applied to the systematic uncertainty to account for the fact that we are measuring the fluxes of marginally extended sources and that our method of estimating the uncertainty in the Herschel images may underestimate the uncertainty (Griffin et al., 2013; Balog et al., 2014).

For most of the BCGs, the error budget for the Hubble photometry is dominated by the absolute calibration uncertainty, the Spitzer error budget is dominated by the uncertainty from modeling and subtracting the old stellar excess, the Herschel PACS uncertainties are dominated either by counting statistics or the calibration, and the Herschel SPIRE uncertainties are dominated by the calibration uncertainty or the confusion noise.

3.2. SED Fitting

We fit SEDs for the CLASH BCGs using the Bayesian fitting code iSEDfit (Moustakas et al., 2013). In order to take full advantage of the CLASH UV-FIR observations, we incorporated modifications to iSEDfit to account for emission from dust. The iSEDfit package fits SEDs by first generating a grid of parameters that describe models of dust–attenuated stellar emission by a composite stellar population. Grids are produced by randomly sampling a bounded section of parameter space. The sampling is weighted with either a uniform prior or a logarithmic scale prior. With these parameters, iSEDfit constructs a grid of model spectra given a choice of synthetic stellar population (SSP), initial mass function (IMF), and dust attenuation model. The package then computes synthetic photometry by convolving these model spectra with the filter responses of the CLASH SEDs, and uses this grid of synthetic SEDs to sample the posterior probability distribution that the model stellar population and dust parameters produce the observed SED.

The iSEDfit package uses a Bayesian Monte Carlo approach to estimate probability distributions for model parameters. For each BCG, we constructed a model grid consisting of 20000 models. The likelihood and mass scaling $\mathcal{A}$ of each model given the observed SED was calculated by minimizing $\chi^{2}$ for each model,

[TABLE]

where $F_{i}$ and $\sigma_{i}$ are the SED fluxes and uncertainties, and $M_{i}$ are the model fluxes. The minimum $\chi^{2}$ and $\mathcal{A}$ are found by solving for $\partial\chi^{2}/\partial\mathcal{A}=0$. We obtained a posterior probability distribution for the models using a weighted random sampling of the model grid, where model weights were determined by the likelihoods. Posterior probability distributions for individual physical parameters were obtained by taking the distribution of parameters of the sampled models. For the SFR, the posterior consists of the distribution of the instantaneous normalized SFR determined by the parameterized star formation history, multiplied by $\mathcal{A}$. For a detailed discussion of iSEDfit, see Moustakas et al. (2013).

We model the star formation history of each BCG with an exponentially decaying curve (the early-type population) and a super-imposed exponentially decaying burst at recent times. The initial exponentially decaying curve is parameterized by the age of the BCG $t$ and the decay rate of the curve $\tau$, while the exponentially decaying burst is parameterized by the duration of the starburst $\Delta t_{b}$, the decay rate of the burst, and the mass of the burst relative to the old stellar population (see Figure 2). The bounds on parameter space for the entire model and assumptions we made for the BCG stellar populations are given in Table 2.

We modified the SED–fitting routine to incorporate dust emission models from Draine & Li (2007) into the parameter grid and synthetic spectra used in iSEDfit. The dust parameter space was sampled using bounds and priors given in Table 2. For each synthetic spectrum, the total luminosity of the dust spectrum was normalized by the difference between the un-attenuated and attenuated stellar spectrum. Figure 3 shows an example of the synthetic attenuated stellar plus dust spectrum, and the Appendix shows best-fit synthetic spectrum for each BCG.

With our modified version of iSEDfit, we were able to test multiple stellar population and dust attenuation models in order to examine the dependence of our results on the choice of SSP, IMF, and attenuation law. We adopt a Bruzual & Charlot (2003) SSP, the Salpeter (1955) IMF defined over the interval 0.1-100 M*⊙, and the Calzetti et al. (2000) attenuation law. We chose to adopt a Salpeter (1955) IMF in order to produce SFRs consistent with the SFRs estimated in our previous paper using the Kennicutt (1998) relationship. In order to examine the effect of an $A_{\textrm{V}}$-dependent attenuation curve, we also performed fits using a modified Calzetti attenuation law, where we multiplied the Calzetti curve with an attenuation-dependent slope that matches the attenuation dependence of the curves published in Witt & Gordon (2000). To test if our results have any dependence on our choice to adopt a Calzetti attenuation law, we also ran iSEDfit assuming clumpy SMC-like dust in a shell geometry. We found our results to be largely consistent with results using the regular Calzetti curve, as well as the attenuation-dependent version of the Calzetti curve. This is because, in the typical CLASH BCG, the attenuation is $A_{\textrm{V}}\lesssim 1.0$, where overall attenuation dependence on the shape of the curve has a only a modest effect. As a further test of the model dependence of our results, we performed SED fits assuming a Chabrier (2003) IMF. Altering the IMF shifts the SFRs downwards by 0.25 dex, but otherwise does not significantly alter the stellar and dust parameters we seek to measure ($\Delta t_{b}$, Md*, A${}_{\textrm{V}}$).

3.3. Cooling and Freefall Time Profiles

We calculated radial profiles of the cooling time, defined to be the ratio of the thermal energy density to the rate of radiative energy density loss for an optically thin plasma undergoing Bremsstrahlung emission, for each cluster using

[TABLE]

where $k_{B}$ is the Boltzmann constant, $n\left(r\right)$ is the total number density profile of the plasma, $n_{H}\left(r\right)$ is the H number density, $n_{e}\left(r\right)$ is the electron number density, $T\left(r\right)$ is the temperature profile, $Z\left(r\right)$ the metallicity profile, and $\Lambda\left(T,Z\right)$ is the cooling function. We obtained values of the cooling function for specific temperatures and metallicities by interpolating the cooling function values in Sutherland & Dopita (1993), and assumed $n\approx 2.3n_{H}$ (Cavagnolo et al., 2009). The ICM density, temperature, and metallicity profiles used in this study are available in Donahue et al. (2014). We used the non-parametric Joint Analysis of Cluster Observations (JACO) profiles (see Mahdavi et al. (2007) and Mahdavi et al. (2013) for a description of the JACO algorithm), which are reported in concentric shells spaced so that each annulus contains at least 1500 X-ray counts.

We measured $t_{cool}$ at specific radii by interpolating on $n_{e}\left(r\right)$, $T\left(r\right)$, and $Z\left(r\right)$ at the desired radius and solving Equation 2. In order to determine the uncertainty on $t_{cool}$ we produced an ensemble of 1000 $n_{e}$, $T$, and $Z$ profiles, where the values in each profile were obtained by drawing from normal distributions defined by the observed values and uncertainties in each profile. By calculating the distribution of $t_{cool}$ given the interpolated values of $n_{e}$, $T$, and $Z$ for the ensemble of profiles, we sampled the probability distribution of $t_{cool}$ and thus estimated uncertainties. We constrained the ensemble of metallicity profiles to only allow values of $Z$ in the range $0.15-1.5Z_{\odot}$, since ICM metallicities are typically $0.3Z_{\odot}$ and tend to be $\gtrsim 0.6-0.8Z_{\odot}$ in the centers of cool core clusters (De Grandi et al., 2004).

We calculated freefall times by assuming cluster masses obey an NFW profile (Navarro et al., 1997). NFW mass concentration parameters and values of M200 were obtained from Merten et al. (2015). We computed the enclosed mass for each cluster as a function of radius,

[TABLE]

where $\rho_{0}$ and $r_{s}$ are NFW scale factors determined by the mass and concentration parameter of the galaxy cluster, as well as the critical density at the cluster redshift. We then calculated freefall time profiles,

[TABLE]

We estimated the contribution of BCG stellar mass to the free-fall time using the Cooke et al. (2016) stellar mass estimates of CLASH BCGs. We did not estimate BCG stellar masses with our SED fits since our photometry does not cover the entire extent of the BCG. Stellar mass density profiles were estimated from the stellar mass of each BCG assuming a Hernquist (1990) profile and the Shen et al. (2003) mass-size relationship (Ruszkowski & Springel, 2009; Laporte et al., 2013). The inclusion of stellar mass in our estimates of $M_{enc}\left(r\right)$ does not significantly affect the freefall time profiles at radii $\gtrsim 25$ kpc (or similarly $\gtrsim 0.025R_{500}$), which is not surprising given that BCG stellar mass dominates the cluster density profile only within the central $\sim 10$ kpc of a cluster (Newman et al., 2013; Monna et al., 2017; Caminha et al., 2017). Still, we include the BCG stellar mass in our estimation of the free-fall time profiles and find that the BCG stellar masses alter the free fall times by $\lesssim 4\%$ at $\sim 10$ kpc and $\lesssim 1\%$ at $\sim 20$ kpc.

Comparing the characteristics of BCG starbursts to ICM thermodynamic parameters requires choosing a radius in the ICM profile at which to measure these parameters. We measure ICM parameters at $0.025R_{500}$, which is the smallest fraction of $R_{500}$ that does not require extrapolation of the X-ray profiles for the 11 CLASH clusters studied. We chose this radius since we expect feedback effects to be strongest near the centers of cool cores, and negligible outside the core. At radii outside the cool core the thermodynamical state of the ICM is only weakly tied to feedback and cooling, and therefore should be weakly related to the BCG (Cowie et al., 1983; Hu et al., 1985).

We also base our expectation of the radial dependence between the properties of the ICM and feedback on BCG simulations in Li et al. (2015). At large radii ($\gtrsim 100$ kpc), the variation in $t_{cool}/t_{ff}$ with time in the simulation (and with the level of star formation and AGN-driven feedback) is weak, so a relationship between SFR and cluster dynamical state would be difficult to detect (see their Figure 9). On the other hand, $t_{cool}/t_{ff}$ varies wildly at small radii ($\lesssim$ 10 kpc), so in this case too, relationships involving $t_{cool}/t_{ff}$ and other parameters related to cooling and feedback in the cluster would be difficult to interpret (although the variation in simulations may be due to the idealized AGN entering periods of complete quiescence as the fuel supply reaches zero). However, since $0.025R_{500}$ is greater than 10 kpc for all the CLASH clusters, even if this effect is physically realistic we would not expect to observe it.

4. Results

Table 3 lists our estimates for the best-fit SFR, the duration of the exponentially decaying starburst $\Delta t_{b}$, and dust mass $M_{d}$ for each active BCG, expressed as the mean of the marginalized posterior probability distribution with uncertainties given as the 68.3% confidence interval. We also list the best-fit intrinsic optical attenuation, $A_{\textrm{V}}$. The SFRs in the CLASH BCGs span the range from $\sim 0.5$ M*⊙* yr*-1* to $\sim 250$ M*⊙* yr*-1*. The BCGs have dust masses ranging from $10^{7}$ M*⊙* in the case of RXJ2129.7+0005 to $10^{9}$ $M_{\odot}$ in the case of MACS1931.8-2653. Our multiband photometry and the corresponding best fit SEDs for each star forming CLASH BCG are shown in the appendix in Figure 18.

Burst durations ($\Delta t_{b}$) range from $\lesssim 100$ Myr to several Gyr with very large uncertainties for any given galaxy. Estimating the posterior probability distribution of $\Delta t_{b}$ is more complicated than estimating the distribution of the SFR or $M_{d}$, owing to the strong dependence of $\Delta t_{b}$ on fburst, the ratio of the stellar mass of the burst population to the stellar mass of the old population in the apertures we use. We find the mean Spearman correlation coefficient between $\Delta t_{b}$ and fburst is $R_{S}=0.72$ (see Figure 15 in the Appendix).

Since the uncertainties on $\Delta t_{b}$ are large and the probability distribution of $\Delta t_{b}$ is highly correlated with fburst, we can increase the precision of our $\Delta t_{b}$ estimates by imposing physically motivated constraints on fburst. We do this by using literature estimates of the fraction of stellar mass contributed to BCGs by star formation at a redshift $<1$ as the range of plausible burst mass percentages (referring back to Table 2, the percentage of stellar mass due to the starburst in our models is fburst/(1 + fburst)). This range is on the order of $1\%-10\%$ (e.g. McIntosh et al., 2008; Inagaki et al., 2015; Cooke et al., in prep.). By restricting the burst mass percentage (and therefore fburst) to a single order of magnitude, we reduce the uncertainty in our estimates of $\Delta t_{b}$ by slightly more than a factor of 1.5, and can provide estimates of the best-fit $\Delta t_{b}$ in some cases where we otherwise report lower limits.

Restricting fburst for each SED is slightly complicated by the fact that we measure our SEDs in apertures that do not cover the entire surface area of the BCG. While the range of burst mass percentages we cite above correspond to $M_{\textrm{burst}}/M_{\textrm{BCG}}$, the stellar mass contributed by the starburst over the stellar mass of the BCG, the burst mass percentage in our SED models is $M_{\textrm{burst}}/M_{\textrm{aperture}}$, the stellar mass contributed by the starburst over the stellar mass in the aperture. However,

[TABLE]

where we assume that the F160W near-IR luminosity, $L_{F160}$, is a reasonable proxy for stellar mass. We measured $L_{F160,\textrm{BCG}}/L_{F160,\textrm{aperture}}$ using the 2D image model of the BCG from the isophotal fitting procedure described in $\S$ 2.3.2 in Lauer et al. (2014). The 2D image models are free of contamination by light from surrounding galaxies. By multiplying the $1\%-10\%$ burst mass percentage range by $L_{F160,\textrm{BCG}}/L_{F160,\textrm{aperture}}$ for each BCG, we were able to find appropriate ranges of $M_{\textrm{burst}}/M_{\textrm{aperture}}$, and therefore fburst, for each SED. For example, $L_{F160,\textrm{BCG}}/L_{F160,\textrm{aperture}}$ is $\sim 2.5-2.7$ in the case of RXJ1532.9+3021, implying that $2.5\%\lesssim M_{\textrm{burst}}/M_{\textrm{aperture}}\lesssim 25\%$. This implies a range of fburst of 0.025-0.4. Since MACS1931.8-2635 does not have a 2D image model available, we used the same range for fburst, but otherwise calculated the restricted range of fburst for each BCG individually.

The best fit values for $\Delta t_{b}$ assuming restricted ranges of fburst are reported in Table 3. By restricting the range of fburst to values that correspond to modest contributions to the BCG stellar mass from recent star formation, we produce a more precise set of $\Delta t_{b}$ results compared to fits with no restrictions on fburst. We adopt these values in the rest of our analysis.

The trends we report below between quantities measured with SED fits are not caused by correlations in the posterior probability distributions of individual fits. The posterior probability distributions of SFR, $\Delta t_{b}$, and $M_{d}$ for individual BCGs are either uncorrelated or only weakly correlated. We report Spearman correlation coefficients in Table 3 and two-dimensional posterior probability distributions in Figure 12 in the Appendix. Across the sample of CLASH BCGs, the posterior probability distributions of SFR and $M_{d}$ are the most correlated, but only with a mean Spearman correlation coefficient of 0.19, while $M_{d}$ and $\Delta t_{b}$ are the least correlated, with a mean Spearman correlation coefficient of 0.01.

Our modified version of iSEDfit produces fits with a mean $\chi_{\nu}^{2}$ of 0.93, and a range between $\chi_{\nu}^{2}=0.57$ and $\chi_{\nu}^{2}=1.63$. The $\chi_{\nu}^{2}$ statistic for the best fitting model in the grid for each BCG is listed in Table 4. Table 4 lists $\chi_{\nu}^{2}$ values we obtain for each dust attenuation law we try using, coupled with a Salpeter (1955) IMF, and the values we obtain using a Calzetti et al. (2000) attenuation law coupled with a Chabrier (2003) IMF. The choice of the adopted attenuation law and IMF can affect the quality of the fits to each individual SED but does not affect the average $\chi_{\nu}^{2}$ for the full sample nor does it affect our qualitative conclusions. When the Calzetti et al. (2000) attenuation law is replaced by Witt & Gordon (2000) clumpy SMC attenuation, the mean $\chi_{\nu}^{2}$ increases by 16% to 1.08; when we use a modified Calzetti et al. (2000) law instead the mean is 1.10.

The best-fit stellar and dust parameters are similarly minimally affected by choice of attenuation law. Adopting a Chabrier (2003) IMF systematically shifts the SFRs $\sim$ 0.25 dex downwards but does not otherwise significantly alter the probability distributions or correlations between parameters, and does not strongly affect values of $\chi_{\nu}^{2}$. The fits using the Chabrier IMF have a mean $\chi_{\nu}^{2}$ of 0.98, and any offsets in SED fit parameters are comparable to, or smaller than, our reported uncertainties. We report the parameters obtained for each model used in the appendix, but in our analyses of the results we only discuss the SED fits obtained using the Calzetti et al. (2000) law and Salpeter (1955) IMF, unless otherwise noted.

4.1. The Starburst - ICM Connection

For BCGs with detectable star-formation activity, we find a tight correlation between the BCG SFR and the ratio $t_{cool}/t_{ff}$ measured at $0.025R_{500}$. The trend and the data are shown in Figure 4. The best fit line between these two parameters is

[TABLE]

and is shown as the black curve in Figure 4. The data are fit by this trend with no detectable intrinsic scatter ($\sigma_{i}<0.15$ dex at 3$\sigma$). We estimated the best fit lines and intrinsic scatters using the least squares method in Hogg et al. (2010) for fitting data with uncertainties in two dimensions and intrinsic scatter. The data have a Spearman correlation coefficient $R_{S}=-0.98$ and a Pearson correlation coefficient $R=-0.95$.

We estimated both the strength of the relationship we observe, and the strength of the claim that the relationship has little intrinsic scatter. We first investigated the possibility that there is no relationship between SFR and $t_{cool}/t_{ff}$, by calculating the marginal probability distribution of the slope using Equation 35 of Hogg et al. (2010), and find that the slope of the relationship is less than 0 at 99.993% confidence.

We then examined the possibility that a relationship obeying the trend in Equation 6 but with a large intrinsic scatter produced the results we observe. We created an ensemble with $10^{4}$ synthetic datasets, consisting of 11 SFR values randomly sampled in logarithmic units from 0.1 M*⊙* yr*-1* to 1000 M*⊙* yr*-1*, and calculated the corresponding predicted values of log10 $t_{cool}/t_{ff}$. Assuming a given intrinsic scatter, $\sigma_{i}$, we generated an offset from the trend line for each of the 11 points by sampling a normal distribution of standard deviation $\sigma_{i}$ to obtain the amplitude of the offset and determined a direction of the offset by sampling a uniform distribution between 0 and $2\pi$. We added an additional offset for each point in the x-direction and in the y-direction by sampling normal distributions with standard deviations equal to the mean uncertainty of the observed SFRs and the observed $t_{cool}/t_{ff}$ values, respectively. The probability that data more tightly correlated than our dataset ($\|R\|=0.95$) would be produced given $\sigma_{i}$ was calculated for the ensemble of datasets. Since we used log-normally distributed variables to generate our synthetic data, we used the Pearson coefficient $R$ to measure the strength of the correlation. We repeated this procedure for $\sigma_{i}$ ranging from 0.0 to 0.6, and found that for $\sigma_{i}\gtrsim 0.22$ dex, the probability of drawing a dataset with a trend tighter than the one we observe is $<0.3\%$. The results of the above correlation test are presented in Figure 5.

The ratio $t_{cool}/t_{ff}$ is thought to be a proxy for the thermal instability of ICM gas – when the gas can cool quickly relative to the time it takes for it to infall, it can more readily become thermally unstable and collapse (Gaspari et al., 2012; Li & Bryan, 2014). Assuming ICM instability increases with decreasing $t_{cool}/t_{ff}$, our results suggest that BCG SFRs will increase with the thermal instability of the surrounding ICM.

The relationship between the star formation duration $\Delta t_{b}$ and cooling time at $0.025R_{500}$, $t_{\textrm{cool}}$, is shown in Figure 6. The relationship between these two parameters is difficult to quantify owing to the large uncertainties on $\Delta t_{b}$. In order to estimate the strength of this correlation, we calculated the Spearman coefficients both including and excluding MACS1931.8-2653. In performing these computations, we also need to address the fact that some of the $\Delta t_{b}$ values are only lower limits. If we use the lower limit values to compute the $\Delta t_{b}$ ranking, then $\Delta t_{b}$ and $t_{cool}$ are strongly correlated, with a Spearman coefficient of 0.86 including MACS1931.8-2653, and a coefficient of 0.81 excluding it. If we assume that the lower limits in $\Delta t_{b}$ are tied for the highest rank, then the Spearman coefficient including MACS1931.8-2653 is 0.57, implying a weak correlation ($p<0.065$), while if MACS1931.8-2653 is excluded, the data are not significantly correlated ($R_{S}=0.427$, $p<0.22$). In the limit of the worst possible $\Delta t_{b}$ ranking for the lower limit data (e.g., the data point with the highest ranked $t_{cool}$ amongst the ‘lower limit’ points has the lowest relative $\Delta t_{b}$ rank), the two datasets are not correlated ($R_{S}=0.433$, $p<0.18$ with MACS1931.8-2653, $R_{S}=0.24$, $p<0.50$ without). Therefore, cooling times and burst durations are likely weakly positively correlated (although this observation may be driven by MACS1931.8-2653), and $\Delta t_{b}$ and $t_{\textrm{cool}}$ become comparable in amplitude when $\Delta t_{b}$ reaches Gyr timescales.

All the BCGs either lie to the left of the line where $\Delta t_{b}=t_{\textrm{cool}}$ or are consistent with lying to the left of this line. This is the region where starbursts occur on timescales shorter than the timescale for the ICM to cool. In clusters where $\Delta t_{b}$ and $t_{cool}$ fall to the left of this line, star formation has been ongoing for a shorter duration than the time it would take for the ICM at $0.025R_{500}$ to radiatively cool. Therefore, assuming $t_{cool}$ profiles typically increase monotonically with radius, if the cold gas reservoir fueling star formation is the result of the ICM radiatively cooling, then the body of low-cooling time gas inside $0.025R_{500}$ that was present at the onset of star formation will be depleted. Alternatively, in clusters where $\Delta t_{b}$ and $t_{cool}$ fall to the right of this line, either higher cooling time gas will have radiatively cooled and replenished the gas inside $0.025R_{500}$ in order to continue forming the cold gas fueling star formation, or radiative cooling has been arrested and $t_{cool}$ is locally static.

4.2. Star Formation and Dust Parameters

We find that both the dust mass and burst duration of CLASH BCGs are correlated with their SFRs. The relationship between $\Delta t_{b}$ and SFR is shown in Figure 7. Large SFRs are consistent with star formation episodes that have recently begun, and as the bursts persist to $\sim$ Gyr timescales, the SFRs diminish by several orders of magnitude.

Figure 8 shows the relationship between $M_{d}$ and SFR. In order to be consistent with studies of SFR and Md conducted by da Cunha et al. (2010) and Hjorth et al. (2014), when comparing these two quantities we analyze SFRs and dust masses obtained with a Chabrier (2003) IMF. We continue to use Salpeter (1955) when discussing results in the rest of our paper.

Dust masses and SFRs are correlated, with a best-fit trendline of

[TABLE]

with an intrinsic scatter of $<0.83$ dex ($3\sigma$ limit). The analysis of da Cunha et al. (2010) derived dust masses and SFRs for field galaxies using a UV-IR SED fitting technique that was similar to ours. We overlay their best-fit trendline, which has a best-fit slope of $1.11\pm 0.01$ and intercept of $7.1\pm 0.01$ in Figure 8. In order to constrain the behavior of the trend at the large SFR end, we include the dust mass and SFR of the Phoenix Cluster BCG (e.g. McDonald et al., 2013). We adopt the Mittal et al. (2017) estimate of 454-494 M*⊙* yr*-1* for the Phoenix SFR and fit the far-IR SED of Phoenix to estimate $M_{d}$. This extended dataset shows a flattening slope at the high-SFR end of the relationship. When we overlay a trend drawn from Hjorth et al. (2014) that takes into account the evolution of star formation and dust for starbursting galaxies with large ($\gtrsim 1000$ M*⊙* yr*-1*) SFRs, we find excellent correspondence to our data across the range of SFRs studied. While our confidence that we observe a change in the SFR-Md trend at the high SFR-end is limited by the size of our sample, our results are fully consistent with starforming BCGs producing dust reservoirs like starbursts in the field.

4.3. Testing Impact of AGN Emission on SED Results on MACS1931.8-2653

We wish to understand the possible impact of AGN emission on our SED fitting results. Therefore, we analyzed the impact of AGN emission on the SED fit to MACS1931.8-2653, since this BCG shows the strongest evidence for AGN emission. If the AGN component has little impact on the SED fit to MACS1931.8-2653, we do not believe there will be substantial AGN contamination of our fits to the other starbursting BCGs in CLASH, which host weaker or quiescent AGN. The X-ray point source in this BCG is well-fit by an AGN power spectrum with a 2-10 keV luminosity of $5.32^{+0.40}_{-0.37}\times 10^{43}$ ergs s*-1* cm*-2* (Hlavacek-Larrondo et al., 2013). Moreover, the IR template-fitting results of Santos et al. (2016) imply a substantial AGN contribution may exist in the case of MACS1931.8-2653, so it is important to investigate the effect of adding an AGN component to the model used to fit the UV-IR SED of this X-ray loud BCG.

As Santos et al. (2016) postulate a potentially large AGN contribution to the UV-IR SED of MACS1931.8-2653, we ran a separate SED fitting analysis for this BCG wherein we include the IR AGN model of Siebenmorgen et al. (2015). We allowed iSEDfit to sample the full range of parameters that describe the Siebenmorgen et al. (2015) AGN emission model library, which for our purposes are nuisance parameters. The contribution of the AGN to the total UV-IR luminosity of the galaxy was allowed to range between $0.001\times$ and $10\times$ the stellar contribution.

We find that the effect of AGN emission in MACS1931.8-2653 is marginal. The best fit $\chi^{2}$ degrades slightly, to 1.09, most likely owing to the addition of AGN model nuisance parameters. After incorporating AGN emission, we find the log10 SFR = $2.36\pm 0.19$ M*⊙* yr*-1*, log10 $\Delta t_{b}$ = $-0.69^{+0.46}_{-0.80}$ Gyr, and that log10 Md = $8.86^{+0.34}_{-0.38}$ M*⊙. Compared to our AGN-free model, we find the SFR changes by $\sim 0.3\sigma$, $\Delta t_{b}$ by $\sim 0.5\sigma$, and Md* by $\sim 0.05\sigma$. This is not surprising, since the AGN contribution to the UV-FIR luminosity relative to the stellar contribution is log10 fAGN = $-1.65\pm 0.87$.

One possible explanation for the discrepancy between these findings and the results published in Santos et al. (2016) for MACS1931.8-2653 is our decision not to incorporate WISE photometry in our SED fitting. WISE W3 and W4 filters cover the region of the MIR spectrum ($\sim 10$ $\mu$m) sensitive to the contribution of AGN flux. Since for our purposes, the AGN flux is a contaminant, it does not make sense to include these filters in our SED. The difference may also be explained, in part, by the fitting technique– the theoretical modelling we employ to fit the data allows us greater flexibility to fit the data than the empirical templates used in Santos et al. (2016).

5. Discussion

CLASH galaxy clusters provide a data-rich sample for the study of feedback in the environments of BCGs. Just over 50$\%$ of the CLASH X-ray selected clusters host BCGs that exhibit signs of vigorous feedback and their relatively low redshifts allow us to study their properties in detail. We find a strong observational relationship between BCG star formation and the ratio $t_{cool}/t_{ff}$, which is a proxy for thermal instability in the ICM. As we will discuss below, this relationship appears to strongly support the AGN-regulated cooling mechanism in the condensation and precipitation model advocated in Voit et al. (2015) and Voit et al. (2017). However, several of the implications of our findings have not been anticipated by models of cooling and feedback in the ICM, and raise interesting questions about the dynamics of condensation and heating in the BCG and its environs.

5.1. $t_{cool}/t_{ff}$ As A Proxy For Thermal Instability and ICM Condensation

Models of AGN-driven condensation and precipitation that involve $t_{cool}/t_{ff}$ as a proxy for thermal instability in the ICM provide a natural foundation for interpreting our observed correlation between between $t_{cool}/t_{ff}$ and the SFR. If $t_{cool}/t_{ff}$ is related to the rate of molecular gas production around an active BCG, it should be correlated with the SFR. Specifically, recent simulations show that $t_{cool}/t_{ff}$ determines the critical overdensity at which ICM density perturbations can condense (Singh & Sharma, 2015). We hypothesize this critical overdensity, in turn, determines the mass deposition efficiency $\epsilon_{MD}$, defined to be the mass fraction in a region of the ICM that cools into molecular gas. Assuming that density perturbations of the ICM follow a log-normal distribution in the core of the cluster, $\epsilon_{MD}$ at a radius $r$ ought to be

[TABLE]

where $\delta_{c}\left(\frac{t_{cool}}{t_{ff}}\right)$ is the critical overdensity for ICM condensation as a function of $t_{cool}/t_{ff}$ and $\sigma$ is the width of the ICM density perturbation distribution. The observed relation between SFR and $t_{cool}/t_{ff}$ would therefore imply that SFR scales with $\epsilon_{MD}$.

We can use the $t_{cool}/t_{ff}$-SFR relationship to infer properties of the relationship between $t_{cool}/t_{ff}$ and $\epsilon_{MD}$, and therefore constrain models of feedback-regulated cooling. As a first-order approximation, we assume that a condition close to equilibrium exists between cooling and star-formation for most of the duration of an episode of feedback-regulated cooling, so that SFR $\sim$ $\dot{M}_{g,real}$, where $\dot{M}_{g,real}$ is the actual ICM cooling rate. By measuring $\dot{M}_{g}\equiv M_{g}/t_{cool}$ within $0.025R_{500}$ with the X-ray parameter profiles used in this paper, we find that

[TABLE]

where for the purposes of this simple model we have set $\log_{10}\epsilon_{MD}=0.0$ when $t_{cool}/t_{ff}=10.0$, which is a ‘critical’ value for mechanical feedback-triggered condensation (Voit & Donahue, 2015). In a relatively uniform sample of clusters like CLASH, where core gas masses occupy a narrow range ($\sim$0.4-1.4$\times 10^{13}$ M*⊙*), $\dot{M}_{g}$ does plays a relatively minor role in Equation 9, resulting in a tight relationship between SFR and $t_{cool}/t_{ff}$. While it is important to bear in mind that we have made a simple estimate of $\epsilon_{MD}$ which may not be the actual mass deposition efficiency, our analysis demonstrates how the relationship between $t_{cool}/t_{ff}$ and star formation may be used to constrain the processes governing the cooling and condensation of gas in the ICM.

In Fogarty et al. (2015), we attempted to estimate the cooling rate of the ICM by measuring $\dot{M}_{g}$ for gas that was at a radius $<35$ kpc or for gas that had an average $t_{cool}/t_{ff}$ ratio below 70 555The cooling times quoted in Fogarty et al. (2015) were taken from the ACCEPT website (http://www.pa.msu.edu/astro/MC2/accept/). ACCEPT website cooling times are incorrect and need to be multiplied by a factor of 6.9/2. Values of $t_{cool}$ were calculated using ACCEPT profiles and assuming fixed metallicity.. We found that while reddening-corrected UV photometric SFRs scale with $\dot{M}_{g}$ in both cases, star formation appears to be increasingly ‘inefficient’ as the SFR decreases. The lowest SFRs we observed were $\sim 0.1-1\%$ of $\dot{M}_{g}$, while at the other extreme the two quantities were comparable. Such a trend would be expected if, as Equation  9 indicates, the $\epsilon_{MD}$ scales with the SFR.

5.1.1 The Role of $t_{ff}$ in the SFR-$t_{cool}/t_{ff}$ relationship

Recent work examining the critical condition for the onset of BCG activity in a cool core cluster suggests that this activity is driven by $t_{cool}$, not $t_{cool}/t_{ff}$ (Hogan et al., 2017). These results may imply that we ought to observe a relationship between SFR-$t_{cool}$, without any significant contribution from $t_{ff}$. An SFR-$t_{cool}$ relationship may be observed if the dominant driver of condensation is the AGN jet propelling low cooling time gas from lower to higher altitudes, where it can condense. Since the cooling time of the uplifted gas determines the radius where the gas becomes thermally unstable, it also determines how much work must be done by the jet to lift the gas to a radius where condensation can occur. Therefore, if the uplift of plasma is driving condensation, we expect to see a SFR-mass condensation efficiency relationship in the form of a SFR-$t_{cool}$ relationship. Such a relationship would appear similar to an SFR-$t_{cool}/t_{ff}$ relationship in a sample such as CLASH with a narrow range of cluster masses, since all of the clusters in CLASH have similar free-fall times in their cores.

We thus investigated the possibility that the underlying relationship we observe in Figure 5 is primarily between SFR and $t_{cool}$. The range of $t_{ff}$ at $0.025R_{500}$ in CLASH clusters with star-forming BCGs is smaller than the range of $t_{cool}$ – while $t_{ff}$ varies by a factor of $\sim$ 1.5 across the sample, $t_{cool}$ varies by a factor of $\sim$ 3. Hence, the effect of $t_{ff}$ on the relationship with SFR is modest, and our ability to distinguish between an underlying SFR-$t_{cool}/t_{ff}$ relationship vs. a SFR-$t_{cool}$ relationship is limited. Several lines of reasoning provide evidence for a non-negligible contribution from $t_{ff}$, although we cannot rule out the interpretation that the underlying relationship we observe is solely between SFR and $t_{cool}$.

SFR and $t_{cool}$ are correlated as shown in Figure 9, with a Spearman correlation coefficient of -0.90 (and a Pearson coefficient of -0.90). We calculate an intrinsic scatter between the two quantities of $0.05^{+0.03}_{-0.01}$ dex ($<0.18$ dex at $3\sigma$). The Spearman correlation coefficient in the SFR-$t_{cool}/t_{ff}$ and SFR-$t_{cool}$ relationships as a function of sampling radius are presented in Figure 10. We do not find a substantial difference between $t_{cool}$ and $t_{cool}/t_{ff}$ in terms of how tightly these quantities relate to the SFR. However, the distinction between $t_{cool}$ in the star-forming and non-starforming CLASH clusters is less clear than the distinction between $t_{cool}/t_{ff}$ in these two populations, and $t_{cool}$ for two clusters with non-starforming BCGs is comparable to $t_{cool}$ in the clusters with BCGs exhibiting $0.1-1$ M*⊙* yr*-1* of star formation.

The situation becomes clearer when we examine the relationship between $t_{cool}$, $t_{ff}$, and the residuals in the fits to the CLASH SFR-$t_{cool}/t_{ff}$ and SFR-$t_{cool}$ datasets. In Figure 11, we show $t_{ff}$ vs. the residuals when we fit a log-log relationship to the SFR-$t_{cool}$ dataset at a radius of $0.025R_{500}$ and at $0.075R_{500}$. We also show $t_{ff}$ vs. the residuals of SFR-$t_{cool}/t_{ff}$ at $0.025R_{500}$ and $t_{cool}$ vs. the residuals of SFR-$t_{cool}/t_{ff}$ at $0.025R_{500}$. We find that at $0.025R_{500}$, $t_{ff}/\langle t_{ff}\rangle$, where $\langle t_{ff}\rangle$ is the mean $t_{ff}$ for the sample, is correlated with the residuals in the fit to the SFR-$t_{cool}$ dataset, $\frac{t_{cool}}{Predicted\textrm{ }t_{cool}}$, where $Predicted\textrm{ }t_{cool}$ is $t_{cool}$ predicted by the SFR-$t_{cool}$ relationship for a given SFR. These two quantities are consistent with $\frac{t_{cool}}{Predicted\textrm{ }t_{cool}}=t_{ff}/\langle t_{ff}\rangle$, implying that dividing $t_{cool}$ by $t_{ff}$ will offset the residuals in the SFR-$t_{cool}$ relationship. The relationship between $t_{ff}$ and $\frac{t_{cool}}{Predicted\textrm{ }t_{cool}}$ at $0.025R_{500}$ has a positive slope with $\sim 98\%$ confidence. The plot in Figure 11 does not have obvious outliers, so the scatter reduction seen at $0.025R_{500}$ by dividing $t_{cool}$ by $t_{ff}$ is not attributable to reducing the residuals in an extreme outlier. There is also no evidence that $t_{ff}$ or $t_{cool}$ are correlated with the other residuals we examined.

5.2. BCG Activity and Low Cooling Time Gas

We find a relationship between $\Delta t_{b}$ and $t_{cool}$ that shows a possible positive correlation between these two quantities with $\Delta t_{b}$ approaching $t_{cool}$ at Gyr timescales. Since $\Delta t_{b}$ is potentially correlated with the SFR, and individual measurement uncertainties are relatively large ($\sim$0.3 – 0.4 dex) compared to the range of $\Delta t_{b}$ values measured in the sample (1.5 dex), the observed relationship has limited power to constrain models of cluster-scale feedback evolution. However, our measurements of the starburst durations have implications for understanding how AGN-regulated feedback progresses over time.

Four of the 11 BCG starburst durations are reported as lower limits. These BCGs may be undergoing continuous star formation. However, given the association between longer $\Delta t_{b}$ and lower SFR, and shorter $\Delta t_{b}$ and higher SFR, we suspect that star formation in these BCGs is slowly decaying. Alternatively, BCGs undergoing feedback may initially exhibit large starbursts before settling down to a relatively steady state with anywhere between $\sim$ 1-10 M*⊙* yr*-1* of star formation. Star formation may also ‘flicker’ on timescales that are short relative to the values of $\Delta t_{b}$ we measure.

Along with recent spectroscopic observations in the far-UV, the long burst durations we observe suggest that star formation and the thermodynamical state of the ICM in the cool-core is temporally decoupled from the gas condensation rate fueling AGN feedback, which may feature irregular spikes and dips over time. Recent observations of RXJ1532.9+3021 with Cosmic Origins Spectrograph (COS) (Green et al., 2012) reveal that star formation in this BCG exceeds the upper limit on gas cooling measured by N V and O VI by a factor of $\sim$10, suggesting the initial build up of molecular gas in this system has run its course (Donahue et al., 2016). Meanwhile similar spectroscopy in Abell 1795 and the Phoenix cluster reveal that gas cooling outstrips the SFR (McDonald et al., 2014a, 2015). Delayed consumption of molecular gas may allow the star formation history of BCGs to ‘smooth over’ these intermittent spikes. Slowly decaying BCG starbursts may also trace a transition from an initial mode of rapid condensation (with $\sim$100-1000 M*⊙* yr*-1* of gas condensing) to a long-duration mode of more modest condensation (with $\sim$1-10 M*⊙* yr*-1* of condensation). If star formation lags behind cooling in making this transition, the resulting observables would be consistent both with our results and the offset between SFRs and cooling rates seen in Molendi et al. (2016) and Donahue et al. (2016).

A positive correlation between $\Delta t_{b}$ and $t_{cool}$ and the existence of long-lived ($>1$ Gyr) starbursts are consistent with the mass condensation efficiency interpretation of the $t_{cool}/t_{ff}$-SFR scaling relationship discussed in Section 5.1. Assuming that the duration of the BCG starburst is a proxy for the duration of AGN feedback, then longer periods of star formation in a BCG imply more energy injection into the surrounding ICM, which in turn raises the cooling time of the surrounding ICM. This in turn diminishes the rate at which molecular gas condenses, and gradually shuts off both the starburst and feedback. The duration $\Delta t_{b}$ and $t_{cool}$ ought to converge in the limit of a long starburst, although feedback mechanisms could run out of fuel and shut down before the two quantities converge depending on how inefficiently mass condenses out of the ICM. Alternatively, a quasi-steady state may be reached with effectively continuous star formation in the BCG coupled with ICM plasma with $t_{cool}$ of a few Gyr being replenished by hotter gas at about the same rate it condenses. In this case, star formation and feedback may decay very slowly, or not at all. Our results do not rule this latter scenario out, leaving the following possible options: (1) BCG activity is cyclical, (2) BCG activity approaches a slowly decaying quasi-steady state, or (3) BCG activity is cyclical but the duty cycle is comparable to the age of the cool core of the cluster.

Finally, it has been noted that there are few examples of BCGs with post-starburst spectra in most BCG samples (Liu et al., 2012; Loubser et al., 2016). Indeed, none of the SOAR or SDSS spectra of CLASH BCGs in Fogarty et al. (2015) show post–starburst features. The lack of post-starburst BCGs is consistent with our estimates of $\sim$ Gyr duration episodes of star formation. Wild et al. (2009) notes that for a galaxy with an exponentially decaying starburst to feature a post-starburst spectrum, the burst’s decay timescale would have to be $\lesssim 0.1$ Gyr and the burst would have to account for at least 5-10$\%$ of the galaxy’s stellar mass. Given the nature of long-duration star formation in our observations and models of feedback-regulated cooling, we would not expect to see a substantial population of post-starburst BCGs.

5.3. Characteristics of Large Starbursts in BCGs

Large BCG starbursts (with SFRs $\gtrsim 100$ M*⊙* yr*-1*) may differ from their more modest counterparts in several ways. Firstly, the SFRs and dust masses in the CLASH BCGs are consistent with using Hjorth et al. (2014) to describe the relationship between BCG star formation and dust mass. Specifically, MACS1931.8-2653 and the Phoenix cluster match the flattening slope and eventual turnover in the Hjorth et al. (2014) relation closely, raising the possibility that these BCGs harbor starbursts similar to the starbursts in the submillimeter-detected population of galaxies studied by Hjorth et al. (2014). The Hjorth et al. sample extends the study of dust and star formation conducted in da Cunha et al. (2010) to cover starbursts to galaxies forming stars at $\gtrsim 1000$ M*⊙* yr*-1*. Since the largest BCG starbursts in the CLASH sample also began forming stars more recently we hypothesize that, like these massive field galaxy starbursts, BCGs with large star formation rates are either forming or building up their dust reservoirs.

Secondly, both MACS1931.8-2653 and RXJ1532.9 +3021 boast prodigious SFRs and distinctive X-ray cavities. Furthermore, MACS1931.8-2653 exhibits an X-ray loud AGN, a feature which is not obvious in the other CLASH clusters. These features are noteworthy because the starbursts in RXJ1532.9+3021 and MACS1931.8-2653 are the youngest in our sample at log10 $\Delta t_{b}$ = $-0.39^{+0.38}_{-0.41}$ and $-1.01^{+0.35}_{-0.34}$ Gyr, respectively. Our findings suggest that stronger X-ray features may be associated with younger, larger BCG starbursts, such as the extreme example of BCG star formation present in the Phoenix Cluster (McDonald et al., 2012, 2013, 2014b).

6. Conclusions

We use multi-wavelength observations from HST, Spitzer and Herschel to derive SFRs, starburst durations and dust masses from fitting the SEDs of the 11 CLASH BCGs with extended UV and nebular line emission features. The SFRs and dust masses span nearly three orders of magnitude, with SFRs ranging from $\sim$0.3 to $\sim$ 250 M*⊙* yr*-1*, and dust masses ranging from $\sim 10^{6}$ to $\sim 10^{9}$ M*⊙*. BCG starbursts are $\lesssim$ 100 Myr to several Gyr old.

We find compelling evidence for a direct link between the thermodynamic state of the ICM and BCG star formation. Specifically, we observe a tight SFR-$t_{cool}/t_{ff}$ relationship and a relationship with an intrinsic scatter $\sigma_{i}<0.15$ dex and a slope of $-0.15\pm 0.03$. These results strongly suggest that thermally unstable ICM plasma with a low cooling time is the source of material that forms the reservoir of cool gas fueling star formation in the CLASH BCGs and that BCG star formation and feedback either exhausts the supply of this material on Gyr timescales or settles into a state with relatively modest ($\sim 1-10$ M*⊙* yr*-1*) continuous star formation.

Even if the trend we observe between SFR and $t_{cool}/t_{ff}$ is due to an underlying trend between SFR and $t_{cool}$ we would still find the AGN-driven condensation and precipitation model to provide a compelling explanation. However, since we find a $\sim 2\sigma$ detection of a correlation between $t_{ff}$ and the residuals from the best-fit line SFR-$t_{cool}$ relationship, we suspect that $t_{ff}$ plays a role in the physics governing the relationship between SFR and the ICM. With larger datasets spanning a larger range of $t_{ff}$ that include both measurements of $t_{cool}$ and $t_{ff}$, it will be possible to constrain the role of $t_{ff}$ in cluster core dynamics with greater confidence.

While our results are not a direct observation of feedback-induced condensation, the condensation model provides predictions that are consistent with our observations and provides a framework for understanding why the star formation rate in BCGs would scale with $t_{cool}/t_{ff}$. If $t_{cool}/t_{ff}$ scales with the critical density for ICM perturbation collapse, $t_{cool}/t_{ff}$ measures the efficiency of condensation, and our findings can be interpreted as a relationship between SFR and the condensation efficiency of the ICM.

Our study also raises several questions about the life-cycle of cooling and feedback. We present evidence that star formation episodes in BCGs with larger SFRs are younger relative to BCGs with more modest SFRs, suggesting that BCG starbursts decay over time. We also present evidence that BCG star formation can persist over $\gtrsim$ Gyr timescales. However, it is not clear whether the starbursts in BCGs are cyclical, or slowly decaying single events that may eventually settle into a low level of persistent star formation. It is also possible we are observing the superposition of several shorter-lived events.

The dust and star formation in BCGs is consistent with the SFR-Md relationship described in Hjorth et al. (2014). This consistency holds even when we include the SFR and dust mass measured for the Phoenix cluster BCG, which forms stars at a rate of $\sim 500$ M*⊙* yr*-1*. Our results lead us to hypothesize that, while uncommon in BCGs, large starbursts like those in MACS1931.8-2635 and the Phoenix cluster may have properties in common with young, violent starbursts in field galaxies.

Our work shows a direct link between the amount of star formation occurring in a BCG and the thermodynamical state of the surrounding ICM. The quality of the data we used in combination with the sample selected allowed us to study specifically the interaction between BCG and the ICM in a uniformly-selected sample of cool-core clusters. These clusters and their BCGs have properties that are consistent with a process of cooling and feedback and the results presented herein bring us closer to a complete understanding of feedback in galaxy clusters. A clear extension of this study would be to examine how the BCG-ICM relationships evolve with redshift and cluster mass by analyzing additional deep, multi-wavelength galaxy cluster surveys.

Acknowledgments

We acknowledge insightful discussions with Mark Voit and Anton Koekemoer that helped enhance our interpretation of the observations. We thank Michael McDonald for providing the infrared SED of the Phoenix Cluster BCG. We also thank the anonymous referee for insightful comments that significantly improved the quality of this paper. This research is supported, in part, by NASA grants HSTGO-12065.01-A and HSTGO-13367. The CLASH Multi-Cycle Treasury Program is based on observations made with the NASA/ESA Hubble Space Telescope and which is operated by the Space Telescope Science Institute. This work is also based, in part, on observations made with the Spitzer Space Telescope and the Herschel Space Observatory. The Spitzer Space Telescope is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. Herschel is an ESA space observatory with science instruments provided by European-led Principal Investigator consortia and with important participation from NASA.