La solution d'un système d'équations linéaire peut-elle être approchée uniquement pour les premières variables?

J'ai un système linéaire d'équations de taille mxm, où m est grand. Cependant, les variables qui m'intéressent ne sont que les n premières variables (n est petit par rapport à m). Existe-t-il un moyen d'approximer la solution pour les premières valeurs m sans avoir à résoudre l'ensemble du système? Dans l'affirmative, cette approximation serait-elle plus rapide que la résolution du système linéaire complet?

Comme d'autres l'ont souligné, cela est difficile à faire avec un solveur direct. Cela dit, ce n'est pas si difficile à faire avec des solveurs itératifs. À cette fin, notez que la plupart des solveurs itératifs minimisent d'une manière ou d'une autre l'erreur par rapport à une norme. Souvent, cette norme est soit induite par la matrice elle-même, mais parfois c'est aussi juste la norme du vecteur l2. Mais cela ne doit pas être le cas: vous pouvez choisir la norme dans laquelle vous souhaitez minimiser l'erreur (ou le résiduel), et vous pouvez, par exemple, choisir une norme dans laquelle vous pesez les composants qui vous intéressent avec 1 et tous les autres avec 1e-12, par exemple quelque chose comme (1e-24) ∑ N i = 6 x 2 i et produit scalaire correspondant. Ensuite, écrivez toutes les étapes du solveur itératif par rapport à cette norme et produit scalaire, et vous obtenez un solveur itératif qui accorde beaucoup plus d'attention aux éléments vectoriels qui vous intéressent qu'aux autres.$|| x ||^2 = \sum_{i=1}^5 x_i^2 +$$\sum_{i=6}^N x_i^2$

La question est bien sûr de savoir si vous avez besoin de moins d'itérations qu'avec le produit standard / scalaire qui pèse tous les composants de manière égale. Mais cela devrait en effet être le cas: disons que vous ne vous souciez que des cinq premiers éléments vectoriels. Ensuite, vous devez avoir au plus cinq itérations pour réduire l'erreur d'un facteur 1e12, car cinq itérations sont nécessaires pour le système 5x5 qui les décrit. Ce n'est pas une preuve, mais je suis assez certain que vous devriez en effet vous en tirer avec un nombre d'itérations bien plus petit si le poids dans la norme (1e-12 ci-dessus) est plus petit que la tolérance avec laquelle vous voulez résoudre le système linéaire de manière itérative .

Former le complément Schur

Supposons que vous ayez permuté et partitionné votre matrice sous la forme

de telle sorte que contient vos degrés de liberté d'intérêt et est beaucoup plus petit que A 11 , alors on peut former le complément Schur$A_{22}$$A_{11}$

soit par une factorisation LU partielle droite, soit par la formule explicite, puis peut être compris dans le sens suivant:$S_{22}$

où représente la partie «sans intérêt» de la solution. Ainsi, à condition qu'un côté droit ne soit pas nul dans les degrés de liberté du complément Schur S 22 , il suffit de résoudre contre S 22 pour obtenir la portion de la solution correspondant à ces degrés de liberté.$\star$$S_{22}$$S_{22}$

Complexité informatique dans un cas dense non structuré

Réglage de à la hauteur de A et n à la hauteur de A 22 , la méthode de référence pour le calcul de S 22 est d'abord facteur L 11 U 11 : = A 11 (Ignorons pivotant pour le moment) à peu près deux / trois ( N - n ) 3 travaux, puis former$N$$A$$n$$A_{22}$$S_{22}$$L_{11} U_{11} := A_{11}$$2/3 (N-n)^3$

en utilisant deux résolutions de triangle nécessitant travaux chacune, puis effectuer la mise à jour vers A 22 en 2 n 2 ( N - n ) travaux.$n(N-n)^2$$A_{22}$$2n^2 (N-n)$

Ainsi, le travail total est d' environ . Quand n est très faible, N - n ≈ N , de sorte que le coût peut être considérée comme à peu près 2 / 3 N 3 , qui est le coût d'une factorisation complète.$2/3 (N-n)^3 + 2n(N-n)^2 + 2n^2 (N-n)$$n$$N-n \approx N$$2/3 N^3$

L'avantage est que, s'il y a un très grand nombre de côtés droits à résoudre avec le même système d'équations, alors pourrait potentiellement être réutilisé un grand nombre de fois, où chaque résolution ne nécessiterait que 2 n 2 de travail (plutôt que 2 N 2 ) si S 22 est pris en compte.$S_{22}$$2n^2$$2N^2$$S_{22}$

Complexité informatique dans le cas rare (typique)

Si votre système clairsemé est né d'un certain type d'approximation par différences finies ou éléments finis, les solveurs directs clairsemés seront presque certainement en mesure d'exploiter une partie de la structure; 2d systèmes peuvent être résolus avec travail et O ( N log N ) de stockage, tandis que les systèmes 3D peuvent être résolus avec O ( N 2 ) de travail et O ( N 4 / 3 ) de stockage. Les systèmes factorisés peuvent alors être résolus avec la même quantité de travail que les besoins de stockage.$O(N^{3/2})$$O(N \log N)$$O(N^2)$$O(N^{4/3})$

Le point de faire apparaître les complexités de calcul est que, si et vous avez un système 2d, alors comme le complément de Schur sera probablement dense, la complexité de résolution étant donné le complément de Schur factorisé seraO(n2)=O(N), qui ne manque qu'un facteur logarithmique par rapport à la résolution complète système! En 3D, il fautO(N)travaillieu deO(N 4 / 3 ).$n \approx \sqrt{N}$$O(n^2) = O(N)$$O(N)$$O(N^{4/3})$

Il est donc important de garder à l’esprit que, dans votre cas où , il n'y aura des économies importantes que si vous travaillez dans plusieurs dimensions et que vous avez de nombreux problèmes à résoudre.$n=\sqrt{N}$

L'approche de réduction du modèle

Depuis que Paul a demandé, je vais parler de ce qui se passe si vous utilisez des méthodes de réduction de modèle basées sur la projection sur ce problème. Supposons que vous puissiez trouver un projecteur tel que la plage de P , notée R ( P ) , contienne la solution de votre système linéaire A x = b , et ait la dimension k , où k est le nombre d'inconnues pour lesquelles vous souhaitent résoudre dans un système linéaire.$\mathbf{P}$$\mathbf{P}$$\mathcal{R}(\mathbf{P})$$\mathbf{Ax} = \mathbf{b}$$k$$k$

Une décomposition en valeurs singulières de donnera la matrice partitionnée suivante:$\mathbf{P}$

Les matrices obscurcies par les étoiles comptent pour d'autres choses (comme l'erreur d'estimation, etc.), mais pour l'instant, nous éviterons de traiter des détails étrangers. Il s'ensuit que

est une décomposition de plein rang de .$\mathbf{P}$

Essentiellement, vous allez résoudre le système

d'une manière intelligente, parce que et W ont aussi la propriété que W T V = I . En multipliant les deux côtés de P A x = P b par W T et en laissant y = V x est une approximation de x rendements$\mathbf{V}$$\mathbf{W}$$\mathbf{W}^{T}\mathbf{V} = \mathbf{I}$$\mathbf{PAx} = \mathbf{Pb}$$\mathbf{W}^{T}$$\mathbf{y} = \mathbf{V}\widehat{\mathbf{x}}$$\mathbf{x}$

Résoudre pour x , Prémultiplier par V , et vous avez y , votre approximation pour x .$\widehat{\mathbf{x}}$$\mathbf{V}$$\mathbf{y}$$\mathbf{x}$

Pourquoi l'approche du complément Schur est probablement meilleure

Pour commencer, vous devez choisir manière ou d'une autre. Si la solution de A x = b est dans R ( P ) , alors y = x , et y n'est pas une approximation. Sinon, y ≠ x , et vous introduisez une erreur d'approximation. Cette approche ne tire pas vraiment parti de la structure que vous avez mentionnée vouloir exploiter. Si nous choisissons P tel que sa plage soit la base d'unité standard dans les coordonnées de x que vous souhaitez calculer, les coordonnées correspondantes de y auront des erreurs. Vous ne savez pas comment choisir$\mathbf{P}$$\mathbf{Ax} = \mathbf{b}$$\mathcal{R}(\mathbf{P})$$\mathbf{y} = \mathbf{x}$$\mathbf{y}$$\mathbf{y} \neq \mathbf{x}$$\mathbf{P}$$\mathbf{x}$$\mathbf{y}$ . Vous pouvez utiliser un SVD de A , par exemple, et sélectionner P pour être le produit des k premiersvecteurs singuliers gauches de A et l'adjoint des k premiersvecteurs singuliers droits de A , en supposant que les vecteurs singuliers sont disposés par ordre décroissant de valeur singulière. Ce choix de projecteur équivaudrait à effectuer une décomposition orthogonale appropriée sur A , et il minimiserait l'erreurL 2 dans la solution approximative.$\mathbf{P}$$\mathbf{A}$$\mathbf{P}$$k$$\mathbf{A}$$k$$\mathbf{A}$$\mathbf{A}$$_{2}$

$\mathbf{V}$$\mathbf{W}$$\mathbf{P}$

Les inconvénients ressemblent beaucoup à l'approche de JackPoulson, sauf que vous ne tirez pas tout à fait parti de la structure que vous avez mentionnée.

La réponse longue est ... en quelque sorte.

$k$

$k$

$n-k$$n$

De plus, gardez à l'esprit que restreindre l'ordre dans lequel vous allez effectuer une rétro-substitution peut restreindre la forme de la matrice (cela enlève la possibilité d'échanger des colonnes), ce qui pourrait éventuellement conduire à un système mal conditionné, mais je ne suis pas sûr de cela - juste quelque chose à garder à l'esprit.