PETScLUSolver constructor, cpp

Asked by Paul Constantine on 2013-02-01

The following code snippet -- adapted from the Poisson demo won't compile.

====
  // Create mesh and function space
  UnitSquareMesh mesh(1024, 1024);
  Poisson::FunctionSpace V(mesh);

  // Define boundary condition
  Constant u0(0.0);
  DirichletBoundary boundary;
  DirichletBC bc(V, u0, boundary);

  // Define variational forms
  Poisson::BilinearForm a(V, V);
  Poisson::LinearForm L(V);
  Source f;
  dUdN g;
  L.f = f;
  L.g = g;

  // Compute solution
  Function u(V);

  PETScMatrix A;
  PETScVector b;
  assemble_system(A,b,a,L,bc);

  PETScLUSolver lu(A,"mumps");
  lu.solve(*u.vector(),b);
====

Here's the error:
====
$ make
[100%] Building CXX object CMakeFiles/demo_poisson.dir/main.cpp.o
/home/paulcon/FEniCS/demo/pde/poisson/cpp/main.cpp(96): error: no instance of constructor "dolfin::PETScLUSolver::PETScLUSolver" matches the argument list
            argument types are: (dolfin::PETScMatrix, const char [6])
    PETScLUSolver lu(A,"mumps");
                     ^

compilation aborted for /home/paulcon/FEniCS/demo/pde/poisson/cpp/main.cpp (code 2)
make[2]: *** [CMakeFiles/demo_poisson.dir/main.cpp.o] Error 2
make[1]: *** [CMakeFiles/demo_poisson.dir/all] Error 2
make: *** [all] Error 2
====

What am I missing? The documentation online is unclear, and there aren't any examples of PETScLUSolver in the demos.

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
Paul Constantine
Solved:
2013-02-01
Last query:
2013-02-01
Last reply:
2013-02-01

This question was reopened

Okay, it works with

PETScLUSolver lu(A.copy(),"mumps");

But that's strange, no? Can't I pass a pointer?

Sorry. Actually, the above did not work. I had to do:

PETScLUSolver lu("mumps");
lu.set_operator(A.copy());

Still odd, right?

Garth Wells (garth-wells) said : #3

On 1 February 2013 03:50, Paul Constantine
<email address hidden> wrote:
> Question #220709 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/220709
>
> Status: Solved => Open
>
> Paul Constantine is still having a problem:
> Sorry. Actually, the above did not work. I had to do:
>
> PETScLUSolver lu("mumps");
> lu.set_operator(A.copy());
>
> Still odd, right?
>

I wrote the sample code form memory rather than testing, . . .

It should be

  boost::shared_ptr<PETScMatrix> A(new PETScMatrix);
  PETScVector b;
  assemble_system(*A,b,a,L,bc);

  PETScLUSolver lu(A,"mumps");
  lu.solve(*u.vector(),b);

(note the shared pointer for A). I tested the above.

We use a shared pointer so the solver objects can be passed around and
used even if A goes out of scope.

Garth

> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.

Thanks!