$BC_2$ type multivariable matrix functions and matrix spherical functions
Matrix spherical functions associated to the compact symmetric pair $(\mathrm{SU}(m+2), \mathrm{S}(\mathrm{U}(2)\times \mathrm{U}(m))$, having reduced root system of type $\mathrm{BC}_2$, are studied. We consider an irreducible $K$-representation $(π,V)$ arising from the $\mathrm{U}(2)$-part of $K$, and the induced representation $\mathrm{Ind}_K^G π$ splits multiplicity free. The corresponding spherical functions, i.e. $Φ\colon G \to \mathrm{End}(V)$ satisfying $Φ(k_1gk_2)=π(k_1)Φ(g)π(k_2)$ for all $g\in G$, $k_1,k_2\in K$, are studied by studying certain leading coefficients which involve hypergeometric functions. This is done explicitly using the action of the radial part of the Casimir operator on these functions and their leading coefficients. To suitably grouped matrix spherical functions we associate two-variable matrix orthogonal polynomials giving a matrix analogue of Koornwinder's 1970s two-variable orthogonal polynomials, which are Heckman-Opdam polynomials for $\mathrm{BC}_2$. In particular, we find explicit orthogonality relations and the matrix polynomials being eigenfunctions to an explicit second order matrix partial differential operator. The scalar part of the matrix weight is less general than Koornwinder's weight.