A Lipschitz continuous function $f : [-1,1] \to \mathbb{C}$ has a unique representation as a series in terms of the Chebyshev polynomials $T_k$, $$ f(x) = \sum_{k = 0}^\infty a_k \, T_k(x) \qquad \forall x \in [-1,1] . $$ See e.g. Theorem 3.1 in Trefethen's "Approximation Theory and Approximation Practice".

Is there an analogue of this result for arbitrary compact sets $S \subset \mathbb{C}$? That is, given such an $S$, is there a sequence of polynomials $T_k^{(S)}$ such that every Lipschitz-continuous function $f : S \to \mathbb{C}$ can be written as $$ f(x) = \sum_{k = 0}^\infty a_k \, T_k^{(S)}(x) \qquad \forall x \in S ? $$ I'm actually only interested in the case of $f$ being analytic in a neighbourhood of $S$, so if that makes a difference then never mind the Lipschitz-continuous case.

What I have so far:

Given $f$ analytic in a neighbourhood of $S$, potential theory tells us that polynomial approximation converges exponentialy, i.e. $$ \min_{p_k \in \mathcal{P}_k} \|f - p_k\|_{\infty,S} \leq C \, \exp(-\gamma \, k) $$ for some $C, \gamma > 0$. Given this, it is tempting to take any polynomial basis $p_k$ and look at the limit $$ \DeclareMathOperator*{\argmin}{arg\,min} a = \lim_{n \to \infty} \argmin_{a^{(n)} \in \mathbb{C}^n} \left\|f - \sum_{k = 0}^n a_k^{(n)} p_k \right\|_{\infty,S}, $$ but this doesn't even exist in general, as outlined in my previous question here.

This example shows that we must be careful about the choice of the basis $p_k$. The Fekete polynomials as defined in Saff's "Logarithmic Potential Theory with Applications to Approximation Theory" look like a candidate because they generalise Chebyshev polynomials in the sense that they are the monic polynomials with asymptotically the smallest extremum. I don't see how this relates to the above problem though.

Why I am interested in this:

Given a matrix $A$ and a function $f$, I would like to estimate the size of the elements of $f(A)$. If I had the desired result, I could write $$ f(A)_{ij} = \sum_{k = 0}^\infty a_k \, T_k^{(S)}(A_{ij}) $$ and would be done if I can get both a bound on $a_k$ (probably from a polynomial approximation estimate) and a bound on $p_k(A_{ij})$ (by assuming certain decay properties of $A$). If $f$ is defined on a single interval, I can do exactly that using the theory on Chebyshev functions, but I am now interested in the case when the domain of $f$ consists of two or more intervals.