I did some testing to reproduce Miro's results. Here are times for
a symmetric Poisson-like problem on an N x N x N mesh of the unit cube
with CG1 elements. If I use the MUMPS solver, then setting
solver.parameters['symmetric'] to True, and so using Cholesky,
results in a modest memory saving, and about 25%-30% faster solve.
If I use PETSc instead of MUMPS, then setting symmetric equal to True
results in a modest *increase* in memory and a *huge increase* in
solver time. For example, on a 40 x 40 x 40 mesh the PETSc LU solver
uses 37 seconds while the PETSc Cholesky solver uses 198 seconds.
PETSc takes 5 times longer than MUMPS for the LU solve and a whopping
39 times longer for the Cholesky solve.
As Garth indicated, if you want to solve large systems by direct solve
don't use PETSc, (and, if you must, don't tell PETSc that your
matrix is symmetric).
Here are the actual timings:
Solution of Poisson problem, CG1 on N x N x N mesh of cube
symmetric=False symmetric=True
N time memory time memory
MUMPS
20 mumps 0.3 942976 0.2 937212
22 mumps 0.4 967744 0.4 947736
24 mumps 0.6 992196 0.5 964740
26 mumps 1.0 1032872 0.8 992284
28 mumps 1.3 1078404 1.1 1025716
30 mumps 1.9 1141312 1.4 1065532
32 mumps 2.4 1186708 1.8 1095176
34 mumps 3.1 1262252 2.4 1141552
36 mumps 4.8 1374788 3.4 1221876
38 mumps 5.4 1456412 3.9 1269188
40 mumps 7.3 1582028 5.1 1351760
PETSC
20 petsc 0.6 950972 2.1 950656
22 petsc 1.1 976664 3.8 979100
24 petsc 1.8 1007896 6.4 1015736
26 petsc 2.9 1098964 11.0 1067084
28 petsc 4.5 1149148 17.8 1136260
30 petsc 6.7 1242468 28.0 1222572
32 petsc 9.8 1322416 43.0 1322780
34 petsc 13.9 1333408 64.8 1457924
36 petsc 19.7 1440184 95.9 1621312
38 petsc 27.3 1669508 139.4 1826692
40 petsc 37.3 1723864 198.2 2076156
and here is the program that generated them:
# sample call: python ch.py 50 'mumps' True
from dolfin import *
import sys
N = int(sys.argv[1])
method = sys.argv[2]
symmetric = (sys.argv[3] == 'True')
mesh = UnitCubeMesh(N, N, N)
V = FunctionSpace(mesh, 'CG', 1)
u = TrialFunction(V)
v = TestFunction(V)
a = ( dot(grad(u), grad(v)) + u*v ) *dx
L = v*dx
u = Function(V)
problem = LinearVariationalProblem(a, L, u)
solver = LinearVariationalSolver(problem)
solver.parameters['linear_solver'] = method
solver.parameters['symmetric'] = symmetric
timer = Timer('solver')
timer.start()
solver.solve()
tim = timer.stop()
mem = memory_usage(as_string=False)
print "{:6} {:8} {:1} {:7.1f} {:9}".\
format(N, method, symmetric, tim, mem[1])
On 03/01/2015 08:58 AM, Garth N. Wells wrote:
The PETSc built-in direct solver is slow for large systems. It really
there for cases where LU is needed as part of another algorithm, e.g.
the coarse level is multigrid.
If you want to solve large systems, use one of the specialised direct
solvers.
Garth
On Sunday, 1 March 2015, Miro Kuchta <[email protected]
<mailto:[email protected]>> wrote:
Hi,
please consider the attached script. Following this
<https://bitbucket.org/fenics-project/dolfin/pull-request/2/use-cholesky-rather-than-lu-decomposition/diff#chg-dolfin/fem/LinearVariationalSolver.cpp>
discussion, if method is mumps, petsc or pastix and
we have symmetric=True the linear system is solved with Cholesky
factorization (it this so?). While testing different method/symmetry
combinations I noticed that PETSc's own symmetric solver is easily
10 times slower then mumps (I don't have pastix to compare against).
Can anyone else reproduce
this? Thanks.
Regards, Miro
--
Garth N. Wells
Department of Engineering, University of Cambridge
http://www.eng.cam.ac.uk/~gnw20
_______________________________________________
fenics mailing list
[email protected]
http://fenicsproject.org/mailman/listinfo/fenics
_______________________________________________
fenics mailing list
[email protected]
http://fenicsproject.org/mailman/listinfo/fenics
