Existe-t-il un algorithme efficace pour les fractions continues à valeurs matricielles?

Supposons que j'ai une équation matricielle définie récursivement comme

A[n] = inverse([1 - b[n]A[n+1]]) * a[n]

L'équation pour A [1] ressemble alors à une fraction continue, pour laquelle il existe des méthodes très efficaces qui évitent un recalcul fastidieux (voir «Recettes numériques» pour quelques exemples).

Cependant, je me demande s'il existe des méthodes analogues qui permettent aux coefficients b [n] et a [n] d'être des matrices, avec la seule contrainte que b [n] A [n + 1] soit une matrice carrée de sorte que la matrice

1 - b[n]A[n+1]

est en fait inversible.

algorithms Lagerbaer
la source

Telle est la question que vous avez posée en mathématiques.SE quelques mois auparavant, non? Est carré ou rectangulaire?

A

$A$

Je me souviens que quelqu'un dans les commentaires de math.SE a suggéré que je pose cette question ici une fois la version bêta en ligne :) Dans mon cas spécial, A est rectangulaire. Les équations récursives correspondent à un ensemble hiérarchique d'équations, et le nombre de quantités croît avec . Dans mon cas, la dimension de A [n] est nx (n-1)

n

$n$

Lagerbaer

Juste curieux, pour quelle application voulez-vous l'utiliser?

Hjulle

Très brièvement, en utilisant l'identité de Dyson pour un hamiltonien particulier génère des fonctions de Green que je peux marquer avec un certain indice

N

$N$ . La collecte de toutes les fonctions avec le même indice dans un vecteur

V_{N}

$V_N$ me permet d'écrire

V_{N} = α_{N} V_{N - 1} + β_{N} V_{N + 1}

$V_N = \alpha_N V_{N-1} + \beta_N V_{N+1}$ utilisant l'identité de Dyson et une approximation appropriée. Utiliser un seuil de coupure pour que

V_{N} = 0

$V_N = 0$ pour tout

n \geq N

$n \ge N$ me permette de trouver des matrices

A_{n}

$A_n$ pour que

et ces matrices sont données par mon équation de style à fraction continue. Cette technique peut, par exemple, calculer les fonctions du réseau de Green pour les modèles à liaison étroite.

V_{n} = A_{n} V_{n - 1}

$V_n = A_n V_{n-1}$

Lagerbaer

Ce n'est pas mon domaine, mais j'étais il y a quelque temps à un séminaire où quelque chose concernant ce problème a été présenté. [Ici] [1] est la seule trace que j'ai pu en trouver en ligne. Je ne sais vraiment pas si ça aide. [1]: mh2009.ulb.ac.be/ResActiv.pdf

user189035

Réponses:

Les deux méthodes suivantes sont données dans Fonctions des matrices: théorie et calcul par Nicholas Higham, à la page 81. Ces formules évaluent

oùest une matrice carrée.

r (X) = b_{0} + \frac{{une}_{1} X}{b_{1} + \frac{{une}_{2} X}{b_{2} + \dots + \frac{{une}_{2 m - 1} X}{b_{2 m - 1} + \frac{{une}_{2 m} X}{b_{2 m}}}}}

$r(X) = b_0 + \frac{a_1 X}{b_1+\frac{a_2 X}{b_2+\cdots + \frac{a_{2m-1} X}{b_{2m-1} + \frac{a_{2m}X}{b_{2m}}}}}$

X

$X$

Méthode descendante:

$P_{−1}=I, Q_{−1}=0,P_0 =b_0 I,Q_0 =I$

pour j = 1: 2m

$P_j =b_j P_{j−1}+a_jXP_{j−2}$

$Q_j =b_j Q_{j−1}+a_j X Q_{j−2}$

fin

$r_m =P_{2m} Q_{2m}^{-1}$

Méthode ascendante:

$Y_{2m} = (a_{2m}/b_{2m})X$

pour j = 2m − 1: −1: 1

Résoudre pour . $(b_jI+Y_{j+1})Y_j =a_jX$ $Y_j$

fin

$r_m=b_0I+Y_1$

La question demande une évaluation de la forme plus générale

b_{0} + \frac{a_{1} X_{1}}{b_{1} + \frac{a_{2} X_{2}}{b_{2} + \dots + \frac{a_{2 m - 1} X_{2 m - 1}}{b_{2 m - 1} + \frac{a_{2 m} X_{2 m}}{b_{2 m}}}}}

$b_0 + \frac{a_1 X_1}{b_1+\frac{a_2 X_2}{b_2+\cdots + \frac{a_{2m-1} X_{2m-1}}{b_{2m-1} + \frac{a_{2m}X_{2m}}{b_{2m}}}}}$

Cela peut être évalué par une simple généralisation des formules ci-dessus; par exemple, la méthode ascendante devient

$Y_{2m} = (a_{2m}/b_{2m})X_{2m}$

pour j = 2m − 1: −1: 1

Résoudre pour . $(b_jI+Y_{j+1})Y_j =a_j X_j$ $Y_j$

fin

. $r_m=b_0I+Y_1$

David Ketcheson
la source

Cela semble très intéressant. Je vais voir si je peux l'appliquer à mon problème spécifique mais cela répond à la question puisque mon b [n] * A [n + 1] est une matrice carrée

Lagerbaer

Ah, mais je viens de remarquer que la matrice

est la même partout dans votre solution, mais pas la mienne.

X

$X$

Lagerbaer

D'accord, je l'ai généralisé.

David Ketcheson

Je sais que cette réponse fait beaucoup d'hypothèses, mais elle généralise au moins votre algorithme:

Supposons que , et la matrice d'ensemencement, , forment tous une famille de navettage de matrices normales, où les décompositions de valeurs propres de et sont connues a priori, disons , et $\{A_n\}$ $\{B_n\}$ $V_N$ $\{A_n\}$ $\{B_n\}$ $U' V_N U = \Lambda_N$ $U' A_n U = \Omega_n$ $U' B_n U = \Delta_n$ , où est unitaire et , et sont des matrices diagonales à valeurs complexes. $U$ $\Lambda_N$ $\{\Omega_n\}$ $\{\Delta_n\}$

Une fois que nous avons dit la décomposition, par induction,

V_{n} = (I - B_{n} V_{n + 1})^{- 1} A_{n} = (I - U Δ_{n} U^{'} U Λ_{n + 1} U^{'})^{- 1} U Ω_{n} U^{'},

$V_n = (I - B_n V_{n+1})^{-1} A_n = (I - U \Delta_n U' U\Lambda_{n+1} U')^{-1} U \Omega_n U',$

qui peut être réorganisé sous la forme

V_{n} = U (I - Δ_{n} Λ_{n + 1})^{- 1} Ω_{n} U^{'} \equiv U Λ_{n} U^{'},

$V_n = U (I - \Delta_n \Lambda_{n+1})^{-1} \Omega_n U' \equiv U \Lambda_n U',$

où est bien sûr toujours diagonal, donc toute la famille va nécessairement commuter avec les autres opérateurs, et nous avons montré que les valeurs diagonales de chaque sont découplées, de sorte que votre formule de récursion scalaire rapide peut être appliquée indépendamment sur les valeurs propres de et les matrices de coefficients. $\Lambda_n$ $\{V_n\}$ $\Lambda_n$ $V_N$

$A_n \equiv \alpha_n I$ $B_n \equiv \beta_n I$ $V_N$

Jack Poulson
la source