Yes. We have a necessary and sufficient characterization for when the Cauchy-Schwarz inequality becomes equality.
For this post, $U,V,W$ shall denote finite dimensional complex Hilbert spaces. Suppose that $A_1,\dots,A_n:U\rightarrow U,B_1,\dots,B_n:V\rightarrow V$ are linear operators. Then define an operator $\Gamma(A_1,\dots,A_n;B_1,\dots,B_n):L(V,U)\rightarrow L(V,U)$ by letting $$\Gamma(A_1,\dots,A_n;B_1,\dots,B_n)(X)=A_1XB_1^*+\dots+A_nXB_n^*.$$
Observe that $\Gamma(A_1,\dots,A_n;B_1,\dots,B_n)$ is similar to $A_1\otimes \overline{B_1}+\dots+A_n\otimes\overline{B_n}$ (Here, $\overline{C}=(C^{T})^{*}=(C^*)^T$, so $\overline{C}$ is the matrix obtained by replacing every entry in $C$ with its complex conjugate), so the spectral radius Cauchy-Schwarz inequality
$$\rho(\Gamma(A_1,\dots,A_n;B_1,\dots,B_n))
\leq\rho(\Phi(A_1,\dots,A_n))^{1/2}\cdot\rho(\Phi(B_1,\dots,B_n))^{1/2}$$
hold in all cases. It is not too hard to prove the spectral radius Cauchy-Schwarz inequality by using the conventional Cauchy-Schwarz inequality and the characterization of $\rho(\Phi(A_1,\dots,A_n))^{1/2}$ given in the 1998 paper The $p$-norm joint spectral radius for even integers by Ding-Xuan Zhou.
Observe that if there is a $\lambda$ and an invertible $B$ with $B_j=\lambda BA_jB^{-1}$ for $1\leq j\leq n$, then
$$\rho(\Gamma(A_1,\dots,A_n;B_1,\dots,B_n))
=\rho(\Phi(A_1,\dots,A_n))^{1/2}\cdot\rho(\Phi(B_1,\dots,B_n))^{1/2}.$$
I claim that when $(A_1,\dots,A_n)$ and $(B_1,\dots,B_n)$ have no invariant subspaces, this is the only way in which the spectral radius Cauchy-Schwarz inequality becomes equality. By decomposing $(A_1,\dots,A_n)$ and $(B_1,\dots,B_n)$ according to their invariant subspaces, we obtain necessary and sufficient conditions for when the Cauchy-Schwarz inequality is actually an equality.
Recall that a linear operator $\mathcal{E}:L(U)\rightarrow L(V)$ is called positive if $\mathcal{E}(P)$ is positive semidefinite whenever $P$ is positive semidefinite. A linear operator $\mathcal{E}:L(U)\rightarrow L(V)$ is called completely positive if
$\mathcal{E}\otimes 1_W:L(U\otimes W)\rightarrow L(V\otimes W)$ is positive whenever $W$ is a finite dimensional Hilbert space. We say that an operator $\mathcal{E}:L(U)\rightarrow L(V)$ is trace preserving if $\text{Tr}(\mathcal{E}(X))=\text{Tr}(X)$ for all $X\in L(U)$. An operator $\mathcal{E}:L(U)\rightarrow L(V)$ is said to be a quantum channel if $\mathcal{E}$ is both trace preserving and completely positive.
The completely positive mappings from $L(V)$ to $L(V)$ are precisely the mappings of the form $\Phi(A_1,\dots,A_n)$. Observe that any linear mapping from $L(V,U)$ to $L(V,U)$ is of the form $\Gamma(A_1,\dots,A_n;B_1,\dots,B_n)$. It is not too hard to show using the Frobenius inner product that the mapping $\Gamma(A_1,\dots,A_n;B_1,\dots,B_n)$ is trace preserving if and only if $A_1^*B_1+\dots+A_n^*B_n=1_V$. In particular, the quantum channels are precisely the mappings of the form $\Phi(A_1,\dots,A_n)$ where $A_1^*A_1+\dots+A_n^*A_n=1_V$.
Observe that if $\mathcal{E}$ is a quantum channel, then $\rho(\mathcal{E})=1$.
Let $E_{U,n}$ be the collection of tuples $(A_1,\dots,A_n)\in L(U)^n$ such that there is some complex number $\lambda$ and invertible $B$ such that if $B_j=\lambda BA_jB^{-1}$ for $1\leq j\leq n$, then $\Phi(B_1,\dots,B_n)$ is a quantum channel. By this answer, $L(U)^n\setminus E_{U,n}$ is a quite small set whenever $n>1$. In particular, if $(A_1,\dots,A_n)$ has no-invariant subspace, then $(A_1,\dots,A_n)\in E_{U,n}$.
Theorem: Suppose that $(A_1,\dots,A_n)\in L(U)^n,(B_1,\dots,B_n)\in L(V)^n$ have no invariant subspace. Then $$\rho(\Gamma(A_1,\dots,A_n;B_1,\dots,B_n))=\rho(\Phi(A_1,\dots,A_n))^{1/2}\cdot\rho(\Phi(B_1,\dots,B_n))^{1/2}$$ if and only if there is some invertible matrix $C$ along with some complex number $\eta\neq 0$ where
$A_j=\eta CB_jC^{-1}$ for $1\leq j\leq n$.
Proof: Since $(A_1,\dots,A_n),(B_1,\dots,B_n)$ have no invariant subspace, there are non-zero complex numbers $\mu,\nu$ along with invertible matrices $A,B$ where if we set $R_j=\mu AA_jA^{-1},S_j=\nu BB_jB^{-1}$ for $1\leq j\leq n$, then $\mathcal{E}=\Phi(R_1,\dots,R_n),\mathcal{F}=\Phi(S_1,\dots,S_n)$ are quantum channels. Let $W$ be a complex inner product space with orthonormal basis $(e_1,\dots,e_n).$ Now, let $R,S\in L(V,V\otimes W)$ be the operators defined by letting
$R=\sum_{j=1}^{n}R_j\otimes e_j,S=\sum_{j=1}^{n}S_j\otimes e_j$. Then
$$\mathcal{E}(X)=\text{Tr}_{W}(RXR^*),\mathcal{F}(X)=\text{Tr}_{W}(SXS^*).$$
Define a mapping $\mathcal{G}$ by setting $$\mathcal{G}(X)=\text{Tr}_{W}(RXS^*)=\sum_{k=1}^nR_kXS_k^*.$$
Since $R_1^*R_1+\dots+R_n^*R_n=S_1^*S_1+\dots+S_n^*S_n=1_V$, the mappings $R,S$ are isometries.
Now, let $\lambda$ be an eigenvalue of $\mathcal{G}$. Then there is some
eigenvector $X$ with $\mathcal{G}(X)=\lambda X$. Now perform a polar decomposition of $X$ to write $X=HP$ where $H$ is an isometry and $P$ is a positive semidefinite matrix. Therefore, we have
$$\lambda HP=\lambda X=\mathcal{G}(X)=\text{Tr}_{W}(RXS^*)=\text{Tr}_W(RHPS^*).$$ Now, set $T=(H^*\otimes 1_W)RH$. Then
$$\lambda P=\lambda H^*HP=H^*\text{Tr}_W(RHPS^*)=\text{Tr}_W((H^*\otimes 1_W)RHPS^*)=\text{Tr}_W(TPS^*).$$ Now let $P=\sum_{k=1}^{n}\sigma_ke_ke_k^*$. Then $\text{Tr}(\lambda P)=\lambda\cdot\sum_{k=1}^n\sigma_k$. Therefore, we have
$$\text{Tr}(\lambda P)=\text{Tr}(\text{Tr}_W(TPS^*))=\text{Tr}(TPS^*)=\sum_{k=1}^{n}\sigma_k\text{Tr}(Te_ke_k^*S^*)=\sum_{k=1}^n\sigma_k\text{Tr}(Te_k(Se_k)^*)=\sum_{k=1}^n\sigma_k\langle Te_k,Se_k\rangle.$$
Therefore, since $$\lambda\cdot\sum_{k=1}^n\sigma_k=\sum_{k=1}^n\sigma_k\langle Te_k,Se_k\rangle,$$ we know that $|\lambda|\leq 1$. Furthermore, if $|\lambda|=1$, then we know that $Te_k=\lambda Se_k$ whenever $\sigma_k>0$.
Therefore, we have $T|_{\text{Im}(P)}=\lambda\cdot S|_{\text{Im}(P)}$. In this case, we have $\lambda P=\text{Tr}_W(TPS^*)=\text{Tr}_W(\lambda SPS^*)$, so $P=\text{Tr}_W(SPS^*)$. Since
$P=\text{Tr}_W(SPS^*)$ and since $(S_1,\dots,S_n)$ has no invariant subspace, we know that $P$ is positive semidefinite. Therefore, $(H^*\otimes 1_W)RH=T=\lambda S$.
Now, if $v\in V$, then $\|RH v|=\|v\|=\|\lambda Sv\|=\|(H^*\otimes 1_W)RHv\|$. Therefore, if $v\in V$, then $RHv\in\text{Im}(H)\otimes 1_W$. Therefore, since $\text{Im}(H)$ is a non-trivial invariance subspace of $U$, we know that $\text{Im}(H)=U$. Therefore, $H$ is a unitary operator.
Thus, $H^*R_jH=\lambda S_j$ for $1\leq j\leq n$. Thus, $H^*\mu AA_jA^{-1}H=\lambda \nu BB_jB^{-1}$. We conclude that
$$A_j=\mu^{-1}\lambda\nu A^{-1}HBB_jB^{-1}H^{-1}A=\mu^{-1}\lambda\nu A^{-1}HBB_j(A^{-1}HB)^{-1}.$$
Q.E.D.
Theorem: Suppose that $A_1,\dots,A_n:U\rightarrow U,B_1,\dots,B_n:V\rightarrow V$ be linear operators. Then assign $U,V$ bases so that
$$A_j=\begin{bmatrix}
A_{j,1,1}&\dots&A_{j,1,u}\\
\vdots&\ddots&\vdots\\
A_{j,u,1}&\dots&A_{j,u,u}
\end{bmatrix}$$
and
$$B_j=\begin{bmatrix}
B_{j,1,1}&\dots&B_{j,1,v}\\
\vdots&\ddots&\vdots\\
B_{j,v,1}&\dots&B_{j,v,v}
\end{bmatrix}$$
and where for each $\alpha,\beta$ each of the $n$ matrices $A_{1,\alpha,\beta},\dots,A_{n,\alpha,\beta}$ have the same dimensions,
for each $\alpha,\beta$, each of the $n$ matrices $B_{1,\alpha,\beta},\dots,B_{n,\alpha,\beta}$ have the same dimensions, and where $A_{j,\alpha,\beta}=0$ whenever $\alpha>\beta$, and where $B_{j,\alpha,\beta}=0$ whenever $\alpha>\beta$, and where for $1\leq\alpha\leq u$, the matrices
$(A_{1,\alpha,\alpha},\dots,A_{n,\alpha,\alpha})$ have no non-trivial invariant subspace, and where if $1\leq\beta\leq v$, the matrices
$(B_{1,\beta,\beta},\dots,B_{n,\beta,\beta})$ have non-trivial no invariant subspace either.
Then $\rho(\Gamma(A_1,\dots,A_n;B_1,\dots,B_n))=\rho(\Phi(A_1,\dots,A_n))^{1/2}\cdot\rho(\Phi(B_1,\dots,B_n))^{1/2}$ if and only if there are $\alpha,\beta$ with $1\leq\alpha\leq u,1\leq\beta\leq v$ and where
$\rho(\Phi(A_1,\dots,A_n))=\rho(\Phi(A_{1,\alpha,\alpha},\dots,A_{n,\alpha,\alpha}))$
$\rho(\Phi(B_1,\dots,B_n))=\rho(\Phi(A_{1,\beta,\beta},\dots,A_{n,\beta,\beta}))$, and
There is an invertible matrix $C$ and some complex number $\lambda\neq 0$ such that $A_{j,\alpha,\beta}=\lambda CB_{j,\alpha,\beta}C^{-1}$ for $1\leq j\leq n$.
A few observations:
In this answer, we actually have a couple of different proofs of the spectral radius Cauchy-Schwarz inequality. To prove the Cauchy-Schwarz inequality, we need to prove that $\rho(\Gamma(R_1,\dots,R_n;S_1,\dots,S_n))\leq 1$ whenever $\Phi(R_1,\dots,R_n),\Phi(S_1,\dots,S_n)$ are quantum channels, and the general case will follow from the fact that $E_{V,n}$ is dense in $L(V)^n$ whenever $n>1$. We have already shown that $\rho(\Gamma(R_1,\dots,R_n;S_1,\dots,S_n))\leq 1$, but there is another way of showing this using the induced trace norm.
If $\mathcal{H}:L(V)\rightarrow L(V)$, then define the induced trace norm of $\mathcal{H}$ to be $\|\mathcal{H}\|_1=\max\{\|\mathcal{H}(X)\|_1:\|X\|_1\leq 1\}.$ Recall that if $R,S$ are isometries and $A$ is a complex matrix, then $\|RAS^*\|_1=\|A\|_1$ whenever the matrix multiplication exists. Also, recall that $\|\text{Tr}_{V}(X)\|_1\leq\|X\|_1$ whenever this inequality makes sense.
We have
$$\|\Gamma(R_1,\dots,R_n;S_1,\dots,S_n)(X)\|_1=\|\text{Tr}_{W}(RXS^*)\|_1\leq\|RXS^*\|_1=\|X\|_1.$$
Therefore, $\|\Gamma(R_1,\dots,R_n;S_1,\dots,S_n)\|_1\leq 1$, so
$\rho(\Gamma(R_1,\dots,R_n;S_1,\dots,S_n))\leq 1$ since the induced trace norm is submultiplicative.
Empirical verification
I have empirically verified using computer calculations that the conclusions that we have made are reasonable. Define a fitness function $F:M_d(\mathbb{C})^n\times M_d(\mathbb{C})^n\rightarrow\mathbb{R}$ by letting
$$F(A_1,\dots,A_n;B_1,\dots,B_n)=\frac{\rho(\Gamma(A_1,\dots,A_n;B_1,\dots,B_n))}{\rho(\Phi(A_1,\dots,A_n))^{1/2}\rho(\Phi(B_1,\dots,B_n))^{1/2}}.$$
One can maximize the value of $F(A_1,\dots,A_n;B_1,\dots,B_n)$ using gradient ascent to obtain examples of tuples $(A_1,\dots,A_n;B_1,\dots,B_n)$ with
$F(A_1,\dots,A_n;B_1,\dots,B_n)\approx 1$, but in each of these examples, we always have a complex number $\lambda$ along with some invertible $B$ where
$B_j\approx \lambda BA_jB^{-1}$ for $1\leq j\leq n$.
Yes. The equality can in fact be reached. Our strategy will be to produce a compact set $K_{d,r}\subseteq M_{d}(\mathbb{C})^{r}$ such that if $(X_1,\dots,X_r)\in K_{d,r}$, then $\rho(\Phi(X_1,\dots,X_r))=1$, and where
$$\rho_{2,d}(A_1,\dots,A_r)=\max\{\rho(A_1\otimes X_1+\dots+A_r\otimes X_r)\mid(X_1,\dots,X_r)\in K_{d,r}\}.$$
Let us go over a few definitions and facts to give some context to our construction. These facts can easily be found in John Watrous's 2018 book called The Theory of Quantum Information.
A mapping $\mathcal{E}:L(V)\rightarrow L(W)$ is said to be positive if whenever $X$ is positive semidefinite, $\mathcal{E}(X)$ is also positive semidefinite. A mapping $\mathcal{E}:L(V)\rightarrow L(W)$ is said to be completely positive if the mapping $\mathcal{E}\otimes 1_{U}:L(V\otimes U)\rightarrow L(W\otimes U)$ is positive for each finite dimensional complex Hilbert space $U$.
A mapping $\mathcal{E}:L(V)\rightarrow L(W)$ is said to be trace preserving if $\text{Tr}(\mathcal{E}(X))=\text{Tr}(X)$ for each $X\in L(V)$.
Proposition: Let $\mathcal{E}:L(V)\rightarrow L(W)$ be a linear operator.
$\mathcal{E}$ is completely positive if and only if there are $A_1,\dots,A_r$ where $\mathcal{E}=\Phi(X_1,\dots,X_r)$.
If $\mathcal{E}$ is defined by letting $\mathcal{E}(X)=A_1XB_1^*+\dots+A_rXB_r^*$, then $\mathcal{E}$ is trace preserving if and only if $A_1^*B_1+\dots+A_r^*B_r=1_V$.
A quantum channel is a completely positive trace preserving operator $\mathcal{E}:L(V)\rightarrow L(W)$. In quantum information theory, the quantum channels are the main morphisms between quantum states.
Proposition: If $\mathcal{E}$ is a quantum channel, then $\rho(\mathcal{E})=1$.
Let $K_{d,r}$ be the collection of all tuples $(X_1,\dots,X_r)\in M_{d}(\mathbb{C})^r$ such that $X_1^*X_1+\dots+X_r^*X_r=1_d$. Said differently, $K_{d,r}$ is the collection of all tuples $(X_1,\dots,X_r)\in M_{d}(\mathbb{C})^r$ such that
$\Phi(X_1,\dots,X_r)$ is a quantum channel. The set $K_{d,r}$ is clearly a closed set, and $K_{d,r}$ is bounded since $d=\text{Tr}(X_1^*X_1+\dots+X_r^*X_r)=\|X_1\|_2^2+\dots+\|X_r\|_2^2$ where $\|\cdot\|_2$ denotes the Frobenius norm, so $K_{d,r}$ is compact.
Lemma: Let $A_1,\dots,A_r\in L(V)$. Suppose that there is no $x\in V\setminus\{0\}$ with $A_1x=\dots=A_rx=0$. Furthermore, suppose that there is no subspace $W\subseteq V$ with $W\neq\{0\},W\neq V$, and $W=A_1[W]+\dots+A_r[W]$. Then there is a $\lambda>0$ along with a positive definite $P$ with $\Phi(A_1,\dots,A_r)(P)=\lambda P$.
Proof: Now, let $\mathcal{Q}$ be the collection of all positive semidefinite matrices in $L(V)$ with trace $1$, and let $F:\mathcal{Q}\rightarrow\mathcal{Q}$ be the mapping defined by letting $$F(P)=\frac{\Phi(A_1,\dots,A_r)(P)}{\text{Tr}\big(\Phi(A_1,\dots,A_r)(P)\big)}.$$ Then $\mathcal{Q}$ is convex, and $F$ is a continuous bijection, so by the Brouwer fixed point theorem, there is some $P\in\mathcal{Q}$ with
$F(P)=P$. Therefore, we have $\Phi(A_1,\dots,A_r)(P)=\lambda P$ for some positive $\lambda$.
Now,
$$\text{Im}(P)=\text{Im}(\lambda P)=\text{Im}(A_1PA_1^*+\dots+A_rPA_r^*)$$
$$=\text{Im}(A_1PA_1^*)+\dots+\text{Im}(A_rPA_r^*)=\text{Im}(A_1P)+\dots+\text{Im}(A_rP)$$
$$=A_1[\text{Im}(P)]+\dots+A_r[\text{Im}(P)].$$
This is only possible if $\text{Im}(P)=V$. $\square$
Let $O_{d,r}$ be the collection of all $(X_1,\dots,X_r)\in M_{d}(\mathbb{C})^{r}$ where $\Phi(X_1,\dots,X_r)$ is not nilpotent. Let $E_{d,r}$ be the collection of all $(X_1,\dots,X_r)\in M_d(\mathbb{C})^r$ where there is no subspace $W\subseteq\mathbb{C}^d$ with
$W=X_1^*[W]+\dots+X_r^*[W]$ and $W\neq\{0\},W\neq\mathbb{C}^d$ and where there is no $x$ with $X_j^*x=0$ for $1\leq j\leq r$.
The sets $O_{d,r},E_{d,r}$ are dense subsets of $M_{d}(\mathbb{C})^{r}$ with $E_{d,r}\subseteq O_{d,r}$.
Proposition: Suppose that $(X_1,\dots,X_r)\in E_{d,r}$. Then there is an invertible matrix $B$ and some positive number $\lambda$ with $(\lambda BX_1B^{-1},\dots,\lambda BX_rB^{-1})\in K_{d,r}$.
Proof: Suppose that $B$ is invertible and $\lambda$ is positive. Then $(\lambda BX_1B^{-1},\dots,\lambda BX_rB^{-1})\in K_{d,r}$ if and only if
$$\sum_{k=1}^{r}\lambda^2 (B^{-1})^*X_k^*B^*BX_k^*B^{-1}=I$$ if and only if
$$\sum_{k=1}^{r}\lambda^2 X_k^*B^*BX_k=B^*B$$.
Therefore, there are $\lambda,B$ with $(\lambda BX_1B^{-1},\dots,\lambda BX_rB^{-1})\in K_{d,r}$ if and only if there is a positive $\mu$ and a positive definite $P$ with $$\mu\Phi(X_1^*,\dots,X_k^*)P=P.$$
On the other hand, the existence of such a positive definite $P$ and positive $\mu$ is guaranteed by the above lemma. $\square$
Now, define a mapping $G_{A_1,\dots,A_r}:O_{d,r}\rightarrow\mathbb{R}$ by letting
$$G_{A_1,\dots,A_r}(X_1,\dots,X_r)=\frac{\rho(A_1\otimes X_1+\dots+A_r\otimes X_r)}{\rho(\Phi(X_1,\dots,X_r))^{1/2}}.$$ Since the mapping $G$ is continuous, we have
$$\rho_{2,d}(A_1,\dots,A_r)=\sup\{G_{A_1,\dots,A_r}(X_1,\dots,X_r)\mid (X_1,\dots,X_r)\in E_{d,r}\}.$$
Since $$G_{A_1,\dots,A_r}(X_1,\dots,X_r)=G_{A_1,\dots,A_r}(\lambda BX_1B^{-1},\dots,\lambda BX_rB^{-1})$$ whenever $B$ is invertible and $\lambda$ is a non-zero complex number, by the above proposition, we know that
$$\{G_{A_1,\dots,A_r}(X_1,\dots,X_r)\mid (X_1,\dots,X_r)\in E_{d,r}\}$$
$$=\{G_{A_1,\dots,A_r}(X_1,\dots,X_r)\mid (X_1,\dots,X_r)\in K_{d,r}\}.$$
Therefore, since $K_{d,r}$ is compact, there is some $(Z_1,\dots,Z_r)\in K_{d,r}$ with
$$G_{A_1,\dots,A_r}(Z_1,\dots,Z_r)=\max\{G_{A_1,\dots,A_r}(X_1,\dots,X_r)\mid (X_1,\dots,X_r)\in K_{d,r}\}$$
$$=\rho_{2,d}(A_1,\dots,A_r).$$
I have ran computations that maximize $G_{A_1,\dots,A_r}$, and in these computations the maximum seems to actually be reached.
Best Answer
I claim that there is a somewhat abstract notion of a multi-spectral radius and that there is probably an abstract theory behind this abstract notion. I will try to justify this abstract multi-spectral radius by showing that it captures the specific examples of multi-spectral radii that I have mentioned in the question and that the simplest examples of these multi-spectral radii are reasonable mathematical objects. With that being said, there are notions of a multi-spectral radius that I have not shown to fit within this framework, so more research on this topic is needed.
Suppose that $A$ is a complex Banach algebra. We say that a function $\rho:A^r\rightarrow[0,\infty)$ is a multi-spectral radius if there is an isometric embedding $\iota:A\rightarrow B$ along with a bounded subset $\mathcal{C}\subseteq B^r$ where
$x_j\iota(a)=\iota(a)x_j$ whenever $a\in A,(x_1,\dots,x_r)\in\mathcal{C}$,
if $(x_1,\dots,x_r)\in\mathcal{C}$, then $(\lambda_1x_1,\dots,\lambda_rx_r)\in\mathcal{C}$ whenever $|\lambda_j|=1$ for $1\leq j\leq r$, and
$$\rho(a_1,\dots,a_r)=\rho_{\iota,\mathcal{C}}(a_1,\dots,a_r)=\sup_{(x_1,\dots,x_r)\in\mathcal{C}}\rho(x_1\iota(a_1)+\dots+x_r\iota(a_r))$$ whenever $a_1,\dots,a_r\in A$.
One may also want to require that if $(x_1,\dots,x_r)\in\mathcal{C}$, then $(\lambda_1x_1,\dots,\lambda_rx_r)\in\mathcal{C}$ whenever $|\lambda_j|=1$ for $1\leq j\leq r$, but this condition is not necessary. One can show that if $a,b\in B$, then the mapping from $\mathbb{C}$ to $[-\infty,\infty)$ defined by $\lambda\mapsto \ln(\rho(a+\lambda b))$ is subharmonic, so by the maximum principle,
$\max\{\rho(\lambda_1x_1\iota(a_1)+\dots+\lambda_rx_r\iota(a_r)):|\lambda_1|=\dots=|\lambda_r|=1\}=\max\{\rho(\lambda_1x_1\iota(a_1)+\dots+\lambda_rx_r\iota(a_r)):|\lambda_1|\leq 1,\dots,|\lambda_r|\leq 1\}.$
The $L_1$-spectral radius can be characterized in terms of our framework.
Theorem: $\rho_1(a_1,\dots,a_r)$ is the maximum value of $\rho(x_1\iota(a_1)+\dots+x_r\iota(a_r))$ where $\iota:A\rightarrow B$ is an isometric embedding of Banach algebras, and $\|x_j\|\leq 1$ for $1\leq j\leq r$.
The proof of the above result is not too hard, and I have given a proof of the above result in this answer.
We say that a multi-spectral radius $\rho$ is unitary invariant if $\rho(a_1,\dots,a_r)=\rho(b_1,\dots,b_r)$ whenever there is an $n\times n$-unitary matrix $(u_{i,j})_{i,j}$ where $b_j=\sum_{i=1}^ru_{i,j}a_i$ for $1\leq j\leq r$. The following lemma is a standard result from quantum information theory.
Lemma: Suppose that $A_1,\dots,A_r,B_1,\dots,B_r\in M_n(\mathbb{C})$. Then $\Phi(A_1,\dots,A_r)=\Phi(B_1,\dots,B_r)$ if and only if there is an $r\times r$-unitary matrix $(u_{i,j})_{i,j}$ where $B_j=\sum_{i=1}^ru_{i,j}A_i$ for $1\leq j\leq r$.
Therefore, a multi-spectral radius $\rho:M_n(\mathbb{C})^r\rightarrow[0,\infty)$ is unitary invariant if and only if $\rho(A_1,\dots,A_r)=\rho(B_1,\dots,B_r)$ whenever $\Phi(A_1,\dots,A_r)=\Phi(B_1,\dots,B_r)$. The continuous unitary invariant multi-spectral radii are completely determined by the mapping $\Phi(A_1,\dots,A_r)\mapsto\rho(A_1,\dots,A_r)$ where $\Phi(A_1,\dots,A_r)$ is completely positive and trace preserving (a completely positive trace preserving map is known as a quantum channel).
The following easy lemmas show that how we can always upgrade a multi-spectral radius to a unitary invariant multi-spectral radius.
Lemma: Let $A$ be an algebra over a field $K$. Suppose that $(a_1,\dots,a_r),(b_1,\dots,b_r),(x_1,\dots,x_r),(y_1,\dots,y_r)\in A^r$. Let $(u_{i,j})_{i,j},(v_{i,j})_{i,j}\in M_r(K)$ be inverse matrices. Suppose that $a_k=\sum_{i=1}^ru_{i,k}b_i$ and $x_k=\sum_{j=1}^rv_{k,j}y_j$ for $1\leq k\leq r$. Then $$\sum_{k=1}^ra_kx_k=\sum_{i=1}^rb_iy_j.$$
Lemma: Suppose that $K$ is a field and $A$ is an algebra over $K$. Let $(u_{i,k})_{i,k}\in M_r(K)$. Suppose furthermore that $a_1,\dots,a_r,b_1,\dots,b_r,x_1,\dots,x_r,y_1,\dots,y_r\in A$ and $a_k=\sum_{i=1}^ru_{i,k}b_i$ for $1\leq k\leq r$ and $x_i=\sum_{k=1}^ru_{i,k}y_k$ for $1\leq k\leq r$. Then $\sum_{k=1}^ra_ky_k=\sum_{i=1}^rb_ix_i$.
Proposition: Let $A,B$ be Banach algebras. Let $\iota:A\rightarrow B$ be an isometric embedding. Suppose that $\mathcal{C}\subseteq B^r$ is a bounded subset with $x_j\iota(a)=\iota(a)x_j$ for $1\leq j\leq r$. Let $\mathcal{D}$ be the collection of all tuples $(y_1,\dots,y_r)$ where there is some $r\times r$-unitary matrix $(u_{i,j})_{i,j}$ and $(x_1,\dots,x_r)\in\mathcal{C}$ where $y_j=\sum_{i=1}^ru_{i,j}x_i$ for $1\leq i\leq r$. Then $$\rho_{\iota,\mathcal{D}}(x_1,\dots,x_r)=\sup\{\rho_{\iota,\mathcal{C}}(\sum_ju_{1,j}x_j,\dots,\sum_ju_{r,j}x_j)\mid (u_{i,j})_{i,j}\in U(r)\}.$$
By using the following version of Holder's inequality that can be proven using the classical Holder's inequality, we can show that the $L_2$-spectral radius is a multi-spectral radius.
Theorem: $\rho(A_1\otimes B_1+\dots+A_r\otimes B_r)\leq \rho_p(A_1,\dots,A_r)\cdot\rho_q(B_1,\dots,B_r)$ whenever $p,q\in(1,\infty)$ and $\frac{1}{p}+\frac{1}{q}=1$.
As a consequence, if $d\geq n$ and $A_1,\dots,A_r\in M_n(\mathbb{C})$, then $$\rho_2(A_1,\dots,A_r)=\max_{(X_1,\dots,X_r)\in M_d(\mathbb{C})}\frac{\rho(A_1\otimes X_1+\dots+A_r\otimes X_r)}{\rho_2(X_1,\dots,X_r)}.$$
We have another construction that allows us to show that the $L_p$-spectral radius is a multi-spectral radius for $1\leq p<\infty$. Suppose now that $1\leq p<\infty$. Now, let $A$ be a Banach algebra. Let $x_1,\dots,x_r$ be non-commutating variables. Let $B$ be the collection of all sums of the form $\sum_{k=0}^n\sum_{i_1,\dots,i_k\in\{1,\dots,r\}}a_{i_1,\dots,i_k}x_{i_1}\dots x_{i_k}$. We observe that for $p>1$ the Banach space $\ell^p$ indexed with the natural numbers cannot be endowed with a convolution operation since $(1/n)_{n=1}^{\infty}*(1/n)_{n=1}^{\infty}=(+\infty)_{n=1}^\infty$. We can give $B$ a norm that combines the $\ell^p$ and the $\ell^1$ norms that makes the completion of $B$ into a Banach algebra.
Then give $B$ the norm $$\|\sum_{k=0}^n\sum_{i_1,\dots,i_k\in\{1,\dots,r\}}a_{i_1,\dots,i_k}x_{i_1}\dots x_{i_k}\|=\sum_{k=0}^n\|(a_{i_1,\dots,i_k})_{i_1,\dots,i_k}\|_p.$$ Give $B$ the multiplication defined by bilinearity along with the condition that $$(a\cdot x_{i_1}\dots x_{i_m})\cdot (b\cdot x_{j_1}\dots x_{j_n})= ab\cdot x_{i_1}\dots x_{i_m}x_{j_1}\dots x_{j_n}.$$ In other words, each element in $A$ commutes with each variable $x_j$, but we do not impose any other version of commutativity.
$B$ is submultiplicative: Let $$u=\sum_{j=0}^\infty\sum_{i_1,\dots,i_j\in\{1,\dots,r\}}a_{i_1,\dots,i_j}x_{i_1}\dots x_{i_j}$$ and let $$v=\sum_{j=0}^\infty\sum_{i_1,\dots,i_j\in\{1,\dots,r\}}a_{i_1,\dots,i_j}x_{i_1}\dots x_{i_j}$$ where only finitely many terms of these 'non-commutative polynomials' are non-zero.
Then
$$\|u\cdot v\|$$ $$=\|(\sum_{k=0}^\infty\sum_{i_1,\dots,i_k\in\{1,\dots,r\}}a_{i_1,\dots,i_k}x_{i_1}\dots x_{i_k})\cdot (\sum_{k=0}^\infty\sum_{i_1,\dots,i_k\in\{1,\dots,r\}}b_{i_1,\dots,i_k}x_{i_1}\dots x_{i_k})\|$$
$$=\|\sum_{k=0}^{\infty}\sum_{j=0}^k\sum_{i_1,\dots,i_j\in\{1,\dots,r\}}\sum_{i_{j+1},\dots,i_k}a_{i_1,\dots,i_j}b_{i_{j+1},\dots,i_k}x_{i_1}\dots x_{i_k}\|$$
$$\leq\sum_{k=0}^{\infty}\sum_{j=0}^k\|\sum_{i_1,\dots,i_k\{1,\dots,r\}}a_{i_1,\dots,i_j}b_{i_{j+1},\dots,i_k}x_{i_1}\dots x_{i_k}\|$$
$$=\sum_{k=0}^{\infty}\sum_{j=0}^k\|(a_{i_1,\dots,i_j}\cdot b_{i_{j+1},\dots, i_k})_{i_1,\dots,i_k\in\{1,\dots,r\}}\|_p$$ $$\leq\sum_{k=0}^{\infty}\sum_{j=0}^k\|(a_{i_1,\dots,i_j})_{i_1,\dots,i_j\in\{1,\dots,r\}}\|_p\cdot \|(b_{i_{j+1},\dots,i_k})_{i_{j+1},\dots,i_k\in\{1,\dots,r\}}\|_p$$ $$=\sum_{j=0}^\infty\|(a_{i_1,\dots,i_j})_{i_1,\dots,i_j\in\{1,\dots,r\}}\|_p\cdot\sum_{k=0}^\infty\|(b_{i_1,\dots,i_k})_{i_1,\dots,i_k\in\{1,\dots,r\}}\|_p=\|u\|\cdot\|v\|.$$
Therefore, the completion $\overline{B}$ of $B$ is a Banach algebra, and the original Banach algebra $A$ embeds into $\overline{B}$. In this case, we simply have $\rho_p(a_1,\dots,a_r)=\rho(a_1x_1+\dots+a_rx_r)$.
One should be able to generalize the above construction to most sensible notions of a multi-spectral radius.
Other examples:
In order for our notion of a multi-spectral radius to be sensible, one would expect that the functions $\rho_{\iota,\mathcal{C}}$ would be coherent and interesting for the simplest possible cases of $\iota,\mathcal{C}$. For example, if $\iota:A\rightarrow A$ is the identity function and $1\leq p\leq\infty$, and $\mathcal{C}$ is the unit ball in $\mathbb{C}^r$ with respect to the $p$-norm, then one should expect for $\rho_{\iota,\mathcal{C}}$ to be about as reasonable of a function as the $L_p$-spectral radii, and experimental computations indicate that this is indeed the case.
Define a mapping $F_{\iota,\mathcal{C},a_1,\dots,a_r}:\mathcal{C}\rightarrow[0,\infty)$ by $F_{\iota,\mathcal{C},a_1,\dots,a_r}(x_1,\dots,x_r)= \rho(x_1\iota(a_1)+\dots+x_r\iota(a_r))$. Experimental computations suggest that the local maxima $(x_1,\dots,x_r)$ for the function $F_{\iota,\mathcal{C},a_1,\dots,a_r}$ tend to resemble a sort of conjugate of $a_1,\dots,a_r$.
Let $\iota_n:M_n(\mathbb{C})\rightarrow M_n(\mathbb{C})$ be the identity mapping, and let $\mathcal{L}_{r;p}=\{(\lambda_1,\dots,\lambda_r)\in \mathbb{C}^r:\|(\lambda_1,\dots,\lambda_r)\|_p=1\}$.
In some of my experiments with $A_1,\dots,A_r\in M_n(\mathbb{R})$ and in all of my experiments with $A_1,\dots,A_r\in M_n(\mathbb{C})$ that are Hermitian or real symmetric, when I computed $\lambda_1,\dots,\lambda_r$ locally maximizes $F_{\iota_n,S_1^r,A_1,\dots,A_r}$, then one can find a $\lambda\in S_1$ and $e_1,\dots,e_r\in\{-1,1\}$ where $\lambda_j=\lambda\cdot e_j$ for $1\leq j\leq r$. A similar phenomenon holds when I locally maximized $F_{\iota_n,\mathcal{L}_{r;p},A_1,\dots,A_r}$ for $1\leq p\leq\infty$ even though this phenomenon seems to break down as $p$ gets close to $1$ and it holds better for Hermitian matrices than it does for non-symmetric real matrices.
My computer experiments indicate that if we locally maximize $F_{\iota,\mathcal{C},a_1,\dots,a_r}$, then as $\mathcal{C}$ better approximates $A$, the local maxima $(x_1,\dots,x_r)$ will become more and more similar to a conjugate version of $(a_1,\dots,a_r)$. On the other hand, if $\mathcal{C}$ is too complicated and has too much room to work with, then the local maxima $(x_1,\dots,x_r)$ will again poorly represent the elements in $A$. Therefore, in order to best represent the conjugates of the elements in $A$, it is best if $\mathcal{C}$ is a little bit simpler than $A$.
Let $\iota_{r;n,d}:M_n(\mathbb{C})\rightarrow M_{n\times d}(\mathbb{C})$ be the algebra homomorphism defined by $\iota_{r;n,d}(A)=A\otimes I_d$. Let $\mathcal{C}_{r;n,d}$ be the collection of all tuples $(I_n\otimes X_1,\dots,I_n\otimes X_r)$ where $\rho_2(X_1,\dots,X_r)=1$.
Theorem: Suppose that $A_1,\dots,A_r,B_1,\dots,B_r$ are $n\times n$-complex matrices where $A_1,\dots,A_r$ do not have a common invariant subspace. Suppose furthermore that $\rho_2(A_1,\dots,A_r)>0,\rho_2(B_1,\dots,B_r)>0$. Then $\rho(A_1\otimes B_1+\dots+A_r\otimes B_r)=\rho_2(A_1,\dots,A_r)\rho_2(B_1,\dots,B_r)$ if and only if there is some $\lambda$ and invertible $C$ where $B_j=\overline{\lambda\cdot C\cdot A_j\cdot C^{-1}}$ for $1\leq j\leq r$.
See this answer or this link for proofs that I gave of the above result.
From the above result, we see that if $I_n\otimes\overline{X_1},\dots,I_n\otimes\overline{X_r}\in M_{n\times n}(\mathbb{C})$ globally maximizes $F_{\iota_{r;n,n},\mathcal{C}_{r;n,n},A_1,\dots,A_r}$ and $A_1,\dots,A_r$ have no common invariant subspace, then there are $C,\lambda$ where $X_j=\lambda CA_jC^{-1}$ for $1\leq j\leq r.$
If $A_1,\dots,A_r\in M_n(\mathbb{C})$ does not have a common invariant subspace and $I_n\otimes\overline{X_1},\dots,I_n\otimes\overline{X_r}\in M_{n\times d}(\mathbb{C})$ locally maximizes $F_{\iota_{r;n,d},\mathcal{C}_{r;n,d},A_1,\dots,A_r}$, then the matrices $X_1,\dots,X_r$ will (up-to-similarity and a constant factor) resemble $A_1,\dots,A_r$. For example, if $A_1,\dots,A_r$ are all real, complex symmetric, real symmetric, Hermitian, real positive semidefinite, complex positive semidefinite, quaternionic, rank $\leq k$, etc, and $I_n\otimes\overline{X_1},\dots,I_n\otimes\overline{X_r}$ locally maximizes $F_{\iota_{r;n,d},\mathcal{C}_{r;n,d},A_1,\dots,A_r}$, then one will often be able to find a constant $\lambda$ and invertible matrix $C$ where $Y_j=\lambda CX_rC^{-1}$ satisfy those properties respectively. Furthermore, one will often be able to find matrices $R,S$ where $Y_j=RA_jS$ for $1\leq j\leq r$. In this case, $RS=I_d$ and $P=SR$ will be a (non-orthogonal) projection matrix. Define linear operators $F,G:M_n(\mathbb{C})\rightarrow M_n(\mathbb{C})$ by setting $F(X)=\sum_{k=1}^rA_kX(PA_kP)^*$ and $G(X)=\sum_{k=1}^rA_k^*XPA_kP$ (here $F=G^*$). Define $U_0=I_n,V_0=I_n$ and set $U_{n+1}=F(U_n)/\|F(U_n)\|,V_{n+1}=G(V_n)/\|G(V_n)\|$ for $n\geq 0$. Then $U_n,V_n$ experimentally converge to positive semidefinite matrices $U,V$, the dominant eigenvectors of $F$ and $G$. It seems like the strategy that the optimization algorithm chose for locally maximizing the spectral radius was to make the dominant eigenvalues of $F,G$ positive semidefinite matrices of rank $d$, but the best way to retain the positive semidefiniteness of the dominant eigenvectors of $F,G$ is to make the operators $PA_kP$ closely related to the operators $A_k$. It seems like the reason this strategy works is that in order for a spectral radius of a matrix $A$ to be large, the matrix $A$ should be designed to maximize a particular eigenvalue, and by making the operators $PA_kP$ related to $A_k$, we can maximize the spectral radius of $F,G$. Since the local maximum values of $F_{\iota_{r;n,d},\mathcal{C}_{r;n,d},A_1,\dots,A_r}$ are closely related to the tuples $(A_1,\dots,A_r)$ themselves, I would regard the multi-spectral radius $\rho_{\iota_{r;n,d},\mathcal{C}_{r;n,d}}$ as a legitimate generalization of the notion of the spectral radius to multiple operators which I call the $L_{2,d}$-spectral radius $\rho_{2,d}$.
Other multi-spectral radii $\rho_{\iota,\mathcal{C}}$ are probably reasonably well-behaved, but more computer experiments are needed to verify whether other multi-spectral radii $\rho_{\iota,\mathcal{C}}$ behave nearly as well as $\rho_{\iota_{r;n,d},\mathcal{C}_{r;n,d}}$.
One can find more details on $\rho_{\iota_{r;n,d},\mathcal{C}_{r;n,d}}$ at my site here, and here is another page where I apply $\rho_{\iota_{r;n,d},\mathcal{C}_{r;n,d}}$ to evaluate cryptographic algorithms. I also gave some experimental observations of $\rho_{\iota_{r;n,d},\mathcal{C}_{r;n,d}}$ right here.
Multi-spectrum:
There seems to be a somewhat reasonable definition of a multi-spectrum of a collection of operators.
Suppose that $\rho_{\iota,\mathcal{C}}$ is a multi-spectral radius. If $(x_1,\dots,x_r)\in\mathcal{C}$ and $\rho(x_1\iota(a_1)+\dots+x_r\iota(a_r))=\rho_{\iota,\mathcal{C}}(a_1,\dots,a_r)$, then we say that the spectrum of $x_1\iota(a_1)+\dots+x_r\iota(a_r)$ is a multi-spectrum of $a_1,\dots,a_r$ with respect to the embedding $\iota$ and set $\mathcal{C}$.
Suppose that $(x_1,\dots,x_r)\in\mathcal{C}$ and for every neighborhood $U$ of $(x_1,\dots,x_r)$ with respect to the topology induced by the norm on $A$, then whenever $(y_1,\dots,y_r)\in U\cap\mathcal{C}$, we have $\rho(x_1\iota(a_1)+\dots+x_r\iota(a_r))\geq\rho(y_1\iota(a_1)+\dots+y_r\iota(a_r))$; then we say that the spectrum of $x_1\iota(a_1)+\dots+x_r\iota(a_r)$ is a local multi-spectrum of $(a_1,\dots,a_r)$ with respect to the embedding $\iota$ and the set $\mathcal{C}$.
This notion of a multi-spectrum depends on the choice of $\iota,\mathcal{C}$ and not only on the multi-spectral radius $\rho_{\iota,\mathcal{C}}$. For example, if $A_1,\dots,A_r$ are complex matrices with no common invariant subspace, then the multi-spectrum of $A_1,\dots,A_r$ with respect to $\iota_{r;n,n}$ and $\mathcal{C}_{r;n,n}$ is simply the spectrum of $A_1\otimes\overline{A_1}+\dots+A_r\otimes\overline{A_r}$. On the other hand, suppose that $x_1,\dots,x_r$ are the non-commuting variables in $\overline{B}$ where $\rho_p(A_1,\dots,A_r)=\rho(x_1\iota(A_1)+\dots+x_r\iota(A_r))$ for all matrices $A_1,\dots,A_r$ and suppose that $\mathcal{C}=\{(\lambda_1 x_1,\dots,\lambda_r x_r)\mid \lambda_1,\dots,\lambda_r\in S_1\}$. Then a multi-spectrum of $(A_1,\dots,A_r)$ with respect to $\iota$ and $\mathcal{C}$ is the spectrum of $x_1A_1+\dots+x_rA_r$. If $\lambda$ is complex number with $|\lambda|=1$, then there is an automorphism $\phi$ of the Banach algebra $B$ with $\phi(x_1A_1+\dots+x_rA_r)=\lambda(x_1A_1+\dots+x_rA_r)$. Therefore, since $x_1A_1+\dots+x_rA_r$ has the same spectrum as $\lambda(x_1A_1+\dots+x_rA_r)$, there is some compact set $C\subseteq[0,\infty)$ with $\sigma(x_1A_1+\dots+x_rA_r)=\{\lambda t:|\lambda|=1,t\in C\}$.
Unresolved properties
If the set $\mathcal{C}$ is compact, then the function $\rho_{\iota,\mathcal{C}}$ is automatically upper-semicontinuous. If $\rho_{\iota,\mathcal{C}}$ is not upper-semicontinuous, then we can just take the upper-semicontinuous regularization of $\rho_{\iota,\mathcal{C}}$, but the possible lack of upper-semicontinuity is a potential problem with the abstract theory that I am proposing. I do not know if we should require $\mathcal{C}$ to always be compact in order to make $\rho_{\iota,\mathcal{C}}$ always upper-semicontinuous.
So far, I have mainly experimental results about multi-spectral radii, but I would like for there to be more theorems about multi-spectral radii.