Se voglio risolvere un sistema triangolare superiore completo, posso chiamare linsolve(A,b,'UT')
. Tuttavia questo al momento non è supportato per matrici sparse. Come posso superare questo?Risolvi * sparse * sistema triangolare superiore
risposta
Modifica Dal momento che quello che vi serve è una procedura di risolvere triangolare, chiamato anche indietro/sostituzione in avanti, è possibile utilizzare ordinaria MATLAB backslash \
operatore che:
x = U\b
Come accennato nella risposta originale, MATLAB riconoscerà il fatto che la matrice è triangolare. Per essere sicuri di ciò, è possibile confrontare la prestazione con la procedura cs_usolve
che si trova in SuiteSparse. È una funzione di messaggistica implementata in C che risolve la soluzione triangolare sparsa per la matrice sparsa triangolare superiore (vi sono anche funzioni simili: cs_lsolve
, cs_utsolve
e cs_ltsolve
).
È possibile dare un'occhiata a performance comparison di MATLAB nativo e cs_l(t)solve
nel contesto di fattorizzazione sparsa di Cholesky. In sostanza, le prestazioni di MATLAB sono buone. L'unica insidia è che se si vuole risolvere un sistema di trasposto
x = U'\b
MATLAB non non riconoscere che e crea in modo esplicito una trasposizione di U
. In tal caso, chiamare esplicitamente cs_utsolve
.
risposta originale Se il sistema è simmetrica e si memorizza solo la parte matrice triangolare superiore (è così che ho capito piena in questione), e se la decomposizione di Cholesky è adatto per voi, chol gestisce matrici simmetriche, se la tua matrice è definita positiva. Per le matrici indefinite è possibile utilizzare ldl. Entrambi gestiscono la memorizzazione sparsa e lavorano sulle parti della matrice simmetrica.
Le nuove versioni di MATLAB utilizzano cholmod and suitesparse per quello. Questa è la fattorizzazione Cholesky di gran lunga migliore che io conosca. In matlab è anche parallelizzato in BALS paralleli.
Il fattore si ottiene dalle funzioni sopra è superiore matrice triangolare L tale che
A=LL'
Basta fare ora è eseguire l'avanzamento e sostituzione all'indietro, che sia semplice ed economico. In MATLAB questo viene fatto automaticamente operatore tha backslash
x=L'\(L\b)
la matrice può essere sparsa, e MATLAB riconoscerà che è superiore/triangolare inferiore. Utilizzerai anche questa chiamata insieme alla sostituzione diretta per i fattori ottenuti utilizzando la fattorizzazione di Cholesky.
Penso che significhi 'A = triu (...)' (completo) vs 'A = sparse (triu (...))' (sparse) –
@RodyOldenhuis Oh, ora l'ho letto di nuovo Penso che tu abbia ragione . Ma la mia risposta include comunque informazioni sui triangolari risolti (sostituzione avanti/indietro) - alla fine questo è quello che fai dopo aver fattorizzato la tua matrice comunque :) – angainor
è possibile utilizzare MLDIVIDE (\) o MRDIVIDE (/) gli operatori sulle matrici sparse ...
I sistemi UT e LT sono tra i sistemi più facili da risolvere. Avere una lettura on the wiki su di esso. Sapendo questo, è facile scrivere il proprio UT o LT risolutore:
%# some example data
A = sparse(triu(rand(100)));
b = rand(100,1);
%# solve UT system by back substitution
x = zeros(size(b));
for n = size(A,1):-1:1
x(n) = (b(n) - A(n,n+1:end)*x(n+1:end))/A(n,n);
end
La procedura è molto simile per i sistemi LT.
Detto questo, è generalmente molto più facile e veloce da usare di Matlab operatore backslash:
x = A\b
che funziona anche per pezzi di ricambio matrici, come nate già indicato.
Si noti che questo operatore risolve anche i sistemi UT che hanno A
non quadrato o se A
ha alcuni elementi uguali a zero (o < eps
) sulla diagonale principale. Risolve questi casi in un senso dei minimi quadrati, che potrebbe essere o non essere desiderabile per te. Si potrebbe verificare la presenza di questi casi, prima di effettuare la soluzione:
if size(A,1)==size(A,2) && all(abs(diag(A)) > eps)
x = A\b;
else
%# error, warning, whatever you want
end
Per saperne di più sul (posteriore) taglio su operatore digitando
>> help \
o
>> help slash
sul prompt dei comandi di Matlab .
Naturalmente posso implementare quella sostituzione posteriore da solo, ho pensato che fosse ovvio :) il problema è che i loop for sono generalmente evitati in MATLAB poiché sono molto lenti. L'operazione di barra è garantita per utilizzare la sostituzione posteriore per le matrici triangolari? – olamundo
@noam: dai un'occhiata [qui] (http://scicomp.stackexchange.com/questions/1001/how-does-the-matlab-backslash-operator-solve-ax-b-for-square-matrices). –
\ è una barra rovesciata, o una divisione sinistra ('mldivide'), mentre'/'è una barra o divisione destra (' mrdivide'). – angainor
Utilizzare ['completo'] (http://www.mathworks.com/help/matlab/ref/full.html)? – chaohuang
@chaohuang Una pessima idea. Usa 'sparse' per un motivo. – angainor
Dai un'occhiata alla mia risposta aggiornata. – angainor