On Sun, 1 Mar 2015 17:48:42 +0000
Jan Blechta <[email protected]> wrote:

> For the sake of objectivity, it should be noted that PETSc LU and
> Cholesky perform so bad in DOLFIN because reordering is not employed,
> see
> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCFactorSetMatOrderingType.html
> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatOrderingType.html
> 
> Also note that in DOLFIN we somehow manipulate MUMPS and SuperLU_dist
> ordering, see
> https://bitbucket.org/fenics-project/dolfin/src/21ba7f92a2acb2f45351d47eca6658a24e7affc2/dolfin/common/SubSystemsManager.cpp?at=master#cl-191
> to workaround some problems with SCOTCH and/or ParMETIS and/or 64-bit
> PetscInt. MUMPS automatic choice ICNTL(7) == 7 may perform better.
> This hack is little bit unexpectable because the code is only

I mean                    ^^^^^^^^^^^^
                          unpredictable

Jan

> executed if DOLFIN takes care of PETSc initialization.
> 
> Jan
> 
> 
> On Sun, 01 Mar 2015 11:15:19 -0600
> Douglas N Arnold <[email protected]> wrote:
> 
> > 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
> 
> _______________________________________________
> fenics mailing list
> [email protected]
> http://fenicsproject.org/mailman/listinfo/fenics

_______________________________________________
fenics mailing list
[email protected]
http://fenicsproject.org/mailman/listinfo/fenics

Reply via email to