Résolution

J'ai des matrices et . est clairsemé et est avec très grand (peut être de l'ordre de plusieurs millions.) est une matrice haut avec plutôt petit ( ) et chaque colonne peut seulement un seul entrée avec le reste étant « s, de telle sorte que . est énorme, il est donc très difficile à inverser, et je peux résoudre un système linéaire tel que manière itérative en utilisant une méthode de sous-espace de Krylov telle que , mais je n'ai pas $A$ $G$ $A$ $n\times n$ $n$ $G$ $n\times m$ $m$ $1 \lt m \lt 1000$ $1$ $0$ $G^TG = I$ $A$ $Ax = b$ $\mathrm{BiCGStab}(l)$ $A^{-1}$ explicitement.

Je veux résoudre un système de la forme: , où et sont des vecteurs de longueur . Une façon de le faire est d'utiliser un algorithme itératif dans un algorithme itératif pour résoudre pour chaque itération de l'algorithme itératif externe. Ce serait cependant extrêmement coûteux en calcul. Je me demandais s'il y avait un moyen plus simple sur le plan informatique de résoudre ce problème. $(G^TA^{-1}G)x = b$ $x$ $b$ $m$ $A^{-1}$

linear-algebra linear-solver sparse-matrix krylov-method Costis
la source

Je viens d'ajouter à ma réponse une remarque sur l'exploitation de la structure 0-1.

Arnold Neumaier

Réponses:

Introduisez le vecteur et résolvez le grand système couplé , pour simultanément, en utilisant une méthode itérative. Si est symétrique (comme il semble probable que vous ne le déclariez pas explicitement), alors le système est symétrique (mais indéfini, bien que quasi-défini si $y:=-A^{-1}Gx$ $Ay+Gx=0$ $G^Ty=-b$ $(y,x)$ $A$ $A$ est définie positive), ce qui pourrait vous aider à choisir une méthode appropriée. (mots clés pertinents: matrice KKT, matrice quasi-définie).

Edit: Comme est symétrique complexe, la matrice augmentée l'est aussi, mais il n'y a pas de quasi-infinité. Vous pouvez cependant utiliser la routine pour calculer ; vous pouvez donc adapter une méthode telle que QMR ftp://ftp.math.ucla.edu/pub/camreport/cam92-19.pdf (conçue pour les systèmes réels, mais vous pouvez facilement la réécrire pour les systèmes complexes, en utilisant l'adjoint dans lieu de la transposition) pour résoudre votre problème. $A$ $Ax$ $A^*x=\overline{A\overline x}$

Edit2: En fait, la structure (0,1) de signifie que vous pouvez éliminer et les composants de symboliquement, se retrouvant ainsi avec un système plus petit à résoudre. Cela signifie jouer avec la structure de , et ne paie que lorsque est donné explicitement dans un format clairsemé plutôt que comme un opérateur linéaire. $G$ $x$ $G^Ty$ $A$ $A$

Arnold Neumaier
la source

Merci! A est symétrique complexe. Y a-t-il des raisons de s'attendre à ce que l'état de la matrice augmentée soit pire que celui de la matrice

? Si m est petit, la matrice augmentée n'est que légèrement plus grande que A, donc je soupçonne que la résolution itérative de ce système augmenté ne devrait pas être beaucoup plus difficile que la résolution d'un système avec A?

A

$A$

Costis

Le numéro de condition des deux systèmes est généralement assez indépendant; cela dépend beaucoup de ce que

est. - J'ai ajouté à ma réponse des informations sur la façon d'exploiter la symétrie complexe.

G

$G$

Arnold Neumaier

Salut les gars! Merci pour toutes vos réponses; cet endroit est super! Une extension de la question d'origine: supposons maintenant que j'ai

, où G et A ont la même signification que dans la question d'origine mais B est une matrice nxn déficiente en rang ( même taille que A) et l'ensemble

est de rang complet. Comment allez-vous résoudre le nouveau système, puisque maintenant l'inverse de B n'existe pas, vous ne pouvez donc pas avoir

(G^{T} A^{- H} B A^{- 1} G) x = b

$(G^TA^{-H}BA^{-1}G)x=b$

G^{T} A^{- H} B A^{- 1} G

$G^TA^{-H}BA^{-1}G$

A B^{- 1} A^{H}

$AB^{-1}A^H$ . Je ne pense pas non plus que cela fonctionnerait simplement avec le pseudoinverse de B.

Costis

Introduisez

, et procédez par analogie avec le cas élaboré. (Probablement, vous devez également factoriser

dans les matrices de rang complet et introduire un vecteur intermédiaire supplémentaire.)

y := A^{- 1} G x

$y:=A^{-1}Gx$

z := A^{- H} B y

$z:=A^{-H}By$

B

$B$

Arnold Neumaier

Salut Arnold. Merci, cela fonctionne vraiment! Je l'ai testé avec de très petits exemples de test, et cela fonctionne très bien. Cependant, mon solveur itératif a d'énormes problèmes pour inverser la matrice augmentée. Alors qu'il ne faut que 80 itérations environ (quelques secondes) pour résoudre un système de la forme

avec la matrice A d'origine, le système avec la matrice augmentée (qui est 2n + mx 2n + m ou 2n-mx 2n- m en utilisant l'approche de @ wolfgang-bangerth) prend des dizaines de milliers d'itérations (plusieurs heures) à résoudre pour un RHS. Existe-t-il des stratégies pour accélérer la convergence? peut-être un préconditionneur?

A x = b

$Ax=b$

Costis

Après la réponse d'Arnold, vous pouvez faire quelque chose pour simplifier le problème. Plus précisément, réécrivez le système sous la forme . Notez ensuite que d'après l'énoncé que est grand et étroit et que chaque ligne n'a qu'un 1 et des zéros sinon, l'énoncé signifie qu'un sous-ensemble des éléments de a une valeur fixe, à savoir les éléments de . $Ay+Gx=0, G^Ty=-b$ $G$ $G^Ty=-b$ $y$ $-b$

Disons que pour simplifier que a colonnes et lignes et qu'exactement les premières lignes en contiennent et que réorganiser les éléments de je puisse faire en sorte que ait la matrice d'identité en haut et une matrice zéro en bas. Ensuite, je peux partitionner en "contraints" et "libres" afin que $G$ $m$ $n$ $m$ $x$ $G$ $m \times m$ $n-m \times m$ $y=(y_c,y_f)$ $m$ $n-m$ . Je peux également partitionner pour que . De l'équation j'obtiens alors: $y_c=-b$ $A$ $A=\begin{pmatrix} A_{cc} & A_{cf} \\ A_{fc} & A_{ff} \end{pmatrix}$ $Ay+Gx=0$ et en utilisant ce que nous savons sur nous avons de la seconde de ces équations et par conséquent En d'autres termes, la seule matrice que vous devez inverser est le sous-ensemble de

A_{c c} y_{c} + A_{c f} y_{f} + x = 0, A_{f c} y_{c} + A_{f f} y_{f} = 0

$A_{cc} y_c + A_{cf} y_f + x = 0, \\ A_{fc} y_c + A_{ff} y_f = 0$

y_{c}

$y_c$

A_{f f} y_{f} = A_{f c} b

$A_{ff} y_f = A_{fc} b$

x = A_{c c} b - A_{c f} A_{f f}^{- 1} A_{f c} b .

$x = A_{cc} b - A_{cf} A_{ff}^{-1} A_{fc} b.$

A

$A$ dont les lignes et les colonnes ne sont pas mentionnées dans

(l'espace nul de

). Vous pouvez facilement le faire: (i) calculer

; (ii) utilisez le solveur que vous avez pour résoudre

; (iii) calculer

G

$G$

G

$G$

z = A_{f c} b

$z=A_{fc} b$

A_{f f} h = z

$A_{ff} h = z$

x = A_{c c} b - A_{c f} h

$x = A_{cc} b - A_{cf} h$

En d' autres termes, étant donné la structure de , la résolution du système linéaire que vous avez est vraiment pas plus difficile que la résolution d' un seul système linéaire avec . $G$ $A$

Wolfgang Bangerth
la source

Mais nous connaissons , et , donc $G$ $G^T$ $A$

$G^T A^{-1} G x = b$

$G G^T A^{-1} G x = G b$

Puisque , alors , donc : $G^T G = I$ $G^T = G^{-1}$ $G G^T=I$

$A^{-1} G x = G b$

$A A^{-1} G x = A G b$

$G x = A G b$

$G^T G x = G^T A G b$

$x = G^T A G b$

Sauf si j'ai raté quelque chose, vous n'avez besoin d'aucune itération, ni d'aucun solveur pour calculer x étant donné , et . $G$ $A$ $b$

Phil H
la source

étant un inverse gauche de

n'implique pas qu'il soit également un inverse droit. Considérons

, où

est un inversegauche, mais

G^{T}

$G^T$

G

$G$

G = e_{1}

$G=e_1$

G^{T} = e_{1}^{T}

$G^T=e_1^T$

G G^{T} = e_{1} e_{1}^{T} \neq I

$GG^T=e_1e_1^T\neq I$

Jack Poulson

, donc

, mais

. C'est plutôt un projecteur sur un sous-espace.

G \in C^{n} \otimes C^{m}

$G\in \mathbb{C}^n\otimes \mathbb{C}^m$

G^{T} G = I_{m \times m}

$G^TG=I_{m\times m}$

G G^{T} \neq I_{n \times n}

$GG^T\ne I_{n\times n}$

Deathbreath