Inverse of symmetric tridiagonal block Toeplitz matrix There is a triagonal block matrix $M$ of form:
$$
M = \begin{bmatrix}
A & B^T & 0 & 0 & \cdots & 0 & 0 \\
B & A & B^T & 0 & \cdots & 0 & 0 \\
0 & B & A & B^T & \cdots & 0 & 0 \\
\vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\
0 & 0 & 0 & 0 & \cdots & B & A 
\end{bmatrix}
$$
where $A, B$ are real-valued square matrices of the same size. 
Also, $A$ is positive definite and symmetric. Later makes $M$ symmetric as well. 
My interest is in closed-form solution for elements of $M^{-1}$. 
From "Explicit inverses of some tridiagonal matrices" C.M. da Fonseca, J. Petronilho, I am aware of the closed-form solution for tridiagonal toeplitz matrix of form:
$$
T = \begin{bmatrix}
a & b & 0 & 0 & \cdots & 0 & 0 \\
c & a & b & 0 & \cdots & 0 & 0 \\
0 & c & a & b & \cdots & 0 & 0 \\
\vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\
0 & 0 & 0 & 0 & \cdots & c & a 
\end{bmatrix}
$$
where $a, b, c$ - scalars.
$$
(T^{-1})_{ij} = \begin{cases}
(-1)^{i+j}\frac{b^{j-i}}{\left(\sqrt{bc}\right)^{j-i+1}}\frac{U_{i-1}(d)U_{n-j}(d)}{U_{n}(d)} \quad \text{if} \quad i \le j \\
(-1)^{i+j}\frac{c^{i-j}}{\left(\sqrt{bc}\right)^{i-j+1}}\frac{U_{j-1}(d)U_{n-i}(d)}{U_{n}(d)} \quad \text{if} \quad i \gt j
\end{cases}
$$
where $d = \frac{a}{2\sqrt{bc}}$, $U_{k}(x)$ - Chebyshev polynomials of the second kind. 
Don't see a way to extend it to block matrices though. Does anybody know how it can be done? or any alternative way? 
 A: After some extended search, found an article "Numerically Stable Algorithms for Inversion of Block Tridiagonal and Banded Matrices" that gives direct algorithm for even more general problem. 
$$
M = \begin{bmatrix}
A_1 & -B_1 & 0 & 0 & \cdots & 0\\
-B_1^T & A_2 & -B_2 & 0 & \cdots & 0\\
0 & -B_2^T & A_3 & -B_3 & \cdots & 0\\
\vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\
0 & 0 & 0 & 0 & \cdots & A_N 
\end{bmatrix}
$$
where all $A_i$ are symmetric, all $B_i$ are non-singular. 
Then, symmetry of $M$, as a whole, is exploited to get decomposition $(M^{-1})_{ij}=U_iV_j^T,\quad j \ge i$. 
Note, when inverse of a symmetric matrix exists (and we assume it does) it should be symmetric as well.      
$U$ and $V$ are given recursively:
$$
U_1 = I, \quad U_2=B_1^{-1}A_1\\
U_{i+1}=B_i^{-1}(A_iU_i-B^T_{i-1}U_{i-1}),\quad i=2,...,N-1\\
V_N^T=(A_NU_N-B_{N-1}^TU_{N-1}^T)^{-1}, \quad V_{N-1}^T=V_N^TA_NB_{N-1}^{-1}\\
V_i^T=(V_{i+1}^TA_{i+1}-V_{i+2}^TB_{i+1}^T)B_i^{-1}, \quad i=N-2,...,1
$$
Under constraints, that all $A_i$ are equal, all $B_i$ are equal, and $B=B^T$ (which I missed when posting original question!) equations become even simpler:
$$
U_1 = I, \quad U_2 = B^{-1}A \\
U_{i+1} = B^{-1}AU_i - U_{i-1},\quad i=2,...,N-1\\
V_N^T=(AU_N-BU_{N-1}^T)^{-1}, \quad V_{N-1}^T=V_N^TAB^{-1}\\
V_i^T=V_{i+1}^TAB^{-1}-V_{i+2}^T, \quad i=N-2,...,1
$$ 
