Hello,

I want to use Libmesh in a multiphysics context where
different solvers, including solvers not written in Libmesh, share information.

The interface between solvers is very simple: they send
and receive constant or writeable references to vectors of doubles. I
want to send both actual solutions and gradients.

For example, in the case of a first order electrostatic
(Poisson) solver , I want to generate :

a vector <double> V
( containing the containing the nodal values of the potential and numbered the same
indices as the nodes of the mesh )

three vector<double>
Ei ( containing the containing the values of the electric field components in
the centroids and numbered the same indices as the elements of the mesh )

Now here are my questions:

1) I
can get the values ordered per DOF via:

std::vector<double> v;

equation_systems.build_discontinuous_solution_vector
(v);

The indexing of the
values in the vector v is not exactly the way I want it to be , but I figured
out that by looping through the elements I can recover the nodal values for
each element. This works, but it is not very elegant ( or efficient) as each
node is being treated several times.

Question : Do you have
any suggestions how to do this better ? Ideally, I want to loop through the nodes
and find the corresponding value in v.

2) How
do I recover the gradients ? I can go through the whole exercise of building
again an FE object, and recuperate the derivatives of the shape functions in the
qpoints, which should allow to calculate the gradients ( along the lines of the
construction of grad_u_old in ex.9). That seems a bit overkill to me, but I
cannot think of a simpler way.

Question: What is the simplest
way to do this? If going through building FE is the only way, can you please write
a simple example ? There is plenty of opportunity to make errors while
converting between the different indices, and I ‘d feel more comfortable
if I get an example.

All the best,

Jacques Kools