Solve for two RHSs

Question

Dear all,

I have a simple problem that needs some optimization. My problem is that I have two linear systems sharing the same form: A u = p and A v = q, with the same BCs (neumann and dirichlet).

Right now what I do is simple, I assemble the system separately, and solve them:

PETScMatrix A;
PETScVector b;

assemble_system(A, b, a, L, bcs);

PETScLUSolver solver;

solver.solve(A, *(du1->vector()), b);

PETScVector b_t;

assemble_system(A, b_t, a, L_t, bcs);

solver.solve(A, *(du2->vector()), b_t);

Since in essence I have just one matrix, and one dirichlet BCs vector, and that won't change in computing the two solutions, what can I do to optimize my code?

I am not bound to LU solvers, just in case another one would be better.

Thanks!

Comments

So... is it possible to do something?

Answer 1

I believe this kind of functionality is available in Trilinos, with their MultiVector.

Comments

Yes, it is, in the past I've used Epetra_MultiVector.

I have two problems with Epetra or Trilinos. First, as far as I remember, they're not part of Fenics, and second, I'd very much like to avoid having two assemble() calls for the RHSs.

How do you suggest I could use Epetra?

PS. The people I'm working with are not really computer-savvy, so the Fenics application is the perfect solution (we're on Mac computers).