Resolução

Eu tenho matrizes e . é escasso e é com muito grande (pode ser da ordem de vários milhões.) é uma matriz de altura com bastante pequeno ( ) e cada coluna pode só tem uma única entrada com o resto sendo 's, de tal modo que . é enorme, por isso é realmente difícil de inverter, e eu posso resolver um sistema linear como iterativamente usando um método de subespaço de Krylov como , mas não tenho $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}$ explicitamente.

Eu quero resolver um sistema da forma: , onde e são vetores de comprimento . Uma maneira de fazer isso é usar um algoritmo iterativo dentro de um algoritmo iterativo para resolver para cada iteração do algoritmo iterativo externo. Isso seria extremamente caro computacionalmente, no entanto. Eu queria saber se existe uma maneira computacionalmente mais fácil de resolver esse problema. $(G^TA^{-1}G)x = b$ $x$ $b$ $m$ $A^{-1}$

linear-algebra linear-solver sparse-matrix krylov-method Costis
fonte

Acabei de adicionar à minha resposta uma observação sobre a exploração da estrutura 0-1.

Arnold Neumaier

Respostas:

Introduza o vetor $y:=-A^{-1}Gx$ e resolva o sistema acoplado grande $Ay+Gx=0$ , $G^Ty=-b$ para $(y,x)$ simultaneamente, usando um método iterativo. Se $A$ é simétrico (como parece provável, embora você não o indique explicitamente), o sistema é simétrico (mas indefinido, embora quase indefinido se $A$ é definitivo positivo), o que pode ajudá-lo a escolher um método apropriado. (palavras-chave relevantes: matriz KKT, matriz quaseidefinida).

Edit: Como é simétrico complexo, o mesmo ocorre com a matriz aumentada, mas não há quasefinibilidade. No entanto, você pode usar a rotina para calcular ; portanto, você pode adaptar um método como o QMR ftp://ftp.math.ucla.edu/pub/camreport/cam92-19.pdf (projetado para sistemas reais, mas você pode reescrevê-lo facilmente para sistemas complexos, usando o anexo em local da transposição) para resolver seu problema. $A$ $Ax$ $A^*x=\overline{A\overline x}$

Edit2: Na verdade, a estrutura (0,1) de significa que você pode eliminar e os componentes de simbolicamente, terminando assim com um sistema menor para resolver. Isso significa mexer com a estrutura de e paga apenas quando é fornecido explicitamente em formato esparso, e não como um operador linear. $G$ $x$ $G^Ty$ $A$ $A$

Arnold Neumaier
fonte

Obrigado! A é simétrico complexo. Existe razão para esperar que a condição da matriz aumentada seja pior do que a da matriz original

? Se m for pequeno, a matriz aumentada é apenas marginalmente maior em tamanho que A, então eu suspeitaria que resolver iterativamente esse sistema aumentado não deve ser muito mais difícil do que resolver um sistema com A?

A

$A$

Costis

O número de condição dos dois sistemas geralmente não é relacionado; depende muito do que

é. - Adicionei à minha resposta informações sobre como explorar simetria complexa.

G

$G$

Arnold Neumaier

Oi pessoal! Obrigado por todas as respostas; esse lugar é ótimo! Uma extensão da pergunta original: suponha agora que eu tenho

, onde G e A têm o mesmo significado que na pergunta original, mas B é uma matriz nxn com deficiência de classificação ( mesmo tamanho que A) e todo o

é a classificação completa. Como você resolveria o novo sistema, já que agora o inverso de B não existe, então você não pode ter

(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$ . Eu não acho que funcionaria simplesmente com o pseudoinverso de B também.

Costis

Introduza

e proceda de forma análoga ao caso resolvido. (Possivelmente você também precisa fator

em matrizes de posto completo e introduzir um vector intermediário adicional.)

y := A^{- 1} G x

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

z := A^{- H} B y

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

B

$B$

Arnold Neumaier

Oi Arnold. Obrigado, isso realmente funciona! Eu testei com alguns exemplos de teste muito pequenos e funciona muito bem. No entanto, meu solucionador iterativo está tendo problemas enormes ao inverter a matriz aumentada. Embora sejam necessárias apenas 80 iterações (alguns segundos) para resolver um sistema da forma

com a matriz A original, o sistema com a matriz aumentada (que é 2n + mx 2n + m ou 2n-mx 2n- m usando a abordagem de @ wolfgang-bangerth) leva dezenas de milhares de iterações (várias horas) para resolver um RHS. Existem estratégias para acelerar a convergência? talvez um pré-condicionador?

A x = b

$Ax=b$

Costis

Após a resposta de Arnold, há algo que você pode fazer para simplificar o problema. Reescreva especificamente o sistema como . Observe que, a partir da afirmação de que é alto e estreito e cada linha tem apenas 1 e zeros, a instrução significa que um subconjunto dos elementos de tem um valor fixo, a saber, os elementos de . $Ay+Gx=0, G^Ty=-b$ $G$ $G^Ty=-b$ $y$ $-b$

Digamos que, por simplicidade que tem colunas e linhas e que exatamente as primeiras linhas têm mais neles e que seja reordenar os elementos de que posso fazer com que tem o matriz identidade no topo e uma matriz zero na parte inferior. Então eu posso particionar em "restrito" e "livre" elementos para que $G$ $m$ $n$ $m$ $x$ $G$ $m \times m$ $n-m \times m$ $y=(y_c,y_f)$ $m$ $n-m$ . Também posso particionar para que . A partir da equação , obtive o seguinte: $y_c=-b$ $A$ $A=\begin{pmatrix} A_{cc} & A_{cf} \\ A_{fc} & A_{ff} \end{pmatrix}$ $Ay+Gx=0$ e usando o que sabemos sobre , temos a partir da segunda dessas equações e consequentemente Em outras palavras, a única matriz que você deve inverter é o subconjunto 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$ cujas linhas e colunas não são mencionadas em

(o espaço nulo de

). Isso você pode facilmente fazer: (i) calcular

; (ii) use qualquer solucionador que você precise resolver

; (iii) calcular

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$

Em outras palavras, dada a estrutura de , resolvendo o sistema linear você tem não é realmente mais difícil do que a solução de um único sistema linear com . $G$ $A$

Wolfgang Bangerth
fonte

$G$ $G^T$ $A$

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

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

$G^T G = I$ , then $G^T = G^{-1}$ , so $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$

Unless I've missed something, you don't need any iteration, or any solver to calculate x given $G$ , $A$ and $b$ .

Phil H
fonte

G^{T}

$G^T$ being a left inverse of

G

$G$ does not imply that it is also a right inverse. Consider

G = e_{1}

$G=e_1$ , where

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

$G^T=e_1^T$ is a left inverse, but

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

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

Jack Poulson

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

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

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

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

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

$GG^T\ne I_{n\times n}$ . Rather it's a projector on a subspace.

Deathbreath