Risolvi sistema triangolare superiore spoglio

Question

Sto cercando di capire come risolvere un sistema triangolare sparse in modo efficiente, Au * x = b in scarse sparse.

Ad esempio, si può costruire una rada matrice triangolare superiore, Au, e un lato b destro con:

import scipy.sparse as sp 
import scipy.sparse.linalg as sla 
import numpy as np 

n = 2000 
A = sp.rand(n, n, density=0.4) + sp.eye(n) 
Au = sp.triu(A).tocsr() 
b = np.random.normal(size=(n)) 

Possiamo ottenere una soluzione al problema utilizzando spsolve, tuttavia è chiaro che il la struttura triangolare non viene sfruttata. Questo può essere dimostrato cronometrando la soluzione e confrontandola con il metodo di risoluzione in splu. (Qui usando la magia tempo% del ipython)

%time x1 = sla.spsolve(Au,b) 
CPU times: user 3.63 s, sys: 79.1 ms, total: 3.71 s 
Wall time: 1.1 s 

%time Au_lu = sla.splu(Au) 
CPU times: user 3.61 s, sys: 62.2 ms, total: 3.67 s 
Wall time: 1.08 s 

%time x2 = Au_lu.solve(b) 
CPU times: user 25 ms, sys: 332 µs, total: 25.4 ms 
Wall time: 7.01 ms 

Come Au è già superiore triangolare la chiamata a splu in realtà non dovrebbe fare molto di qualche cosa, tuttavia, come n diventa grande questa chiamata diventa proibitivo (come fa l'uso di spsolve), mentre il tempo di risoluzione rimane piccolo.

C'è un modo per utilizzare il risolutore triangolare di superLU senza prima chiamare splu? O c'è un modo migliore per farlo del tutto?

Answer 1

Ho paura che questo non sia terribilmente istruttivo, ma hai provato a cambiare la permutazione della colonna? Quando uso "NATURAL", ottengo enormi velocità.

%time x1 = sla.spsolve(Au, b, permc_spec="NATURAL") 
CPU time: user 46.7 ms, sys: 0 ns, total: 46.7 ms 
Wall time: 49 ms 

Per me, non è abbastanza veloce come con i metodi della funzione di uscita splu risolvere, ma sembra per ottenere ragionevolmente vicino (ed evita di chiamare splu). Forse sarà sufficiente? Scipy Docs