[WIP] force vector

[tjhei]

This is one out of several PRs to get Stokes solver improvements (with and without melt) with @rrgrove6 into mainline.

This adds an optional material model output force vector to be applied to the RHS of the Stokes system. We need this to be able to construct manufactured solutions.

Questions:

• Should I also add this to the non-melt assembly? It is not strictly needed, because one can abuse the "rho*g" in the RHS (like the burstedde test does). I feel like this approach would be much cleaner though.
• Is the creation of the instance inside evaluate() (see make_shared) the best way to do this?

[bangerth]

then throw an exception

[bangerth]

no need for the empty line

[bangerth]

It's not quite clear to me why this should be part of a material model. What were the design choices here, and why did you go this way?

[tjhei]

It's not quite clear to me why this should be part of a material model. What were the design choices here, and why did you go this way?

All other variables/terms in the equation come from the material model, so that was the natural place to put it. Also, AdditionalOutputs provided an easy way to put it in. Is your suggestion to create something like a GravityModel or GeometryModel? That would make this feature leak into many different parts of the code even though that is not desirable (I assume this feature is used very rarely).

[bangerth]

Well, but that's not true that all other terms come from the material model. Dirichlet as well as traction boundary conditions come from their own set of plugins. So do prescribed velocities in the interior, and as you point out the gravity vector also comes from somewhere else.

I think the most natural way would really be to have it somewhere separate, e.g., in a system of plugins. That also allows to neatly separate concerns -- for example, one could describe the forcing function as a function expression in the .prm file, and adjust it by hand to whatever material model is chosen. The way you implement it, one is essentially forced to re-implement the forcing term for (i) each material model you want to test, (ii) for each testcase you want to test.

[tjhei]

I think the most natural way would really be to have it somewhere separate

Okay. I am fine with making that change. My goal when implementing it was to make it as invisible as possible because it is not a commonly used feature. But that is why we do WIP PRs. ;-)

[tjhei]

Options:

1. Write an assembler within plugin for each problem that needs a force vector (see "alternative" above)
2. AdditionalOutputs (see this PR)
3. Create a base class for 1. so the implementation is much smaller.
4. Implement a new plugin type for ForceVector and implement Function and/or new plugins.

[bangerth]

I'm for 1, 3, 4, 2, in this order.

[jdannberg]

I'm for 2, 3, 4, 1 (in that order).
If we could also somehow fit elasticity into this structure (it also needs a new RHS, but only for the velocity) that would be a bonus.

[bangerth]

I guess we will need to come up with a rank-based ordering system if we get multiple ballots with different orders :-)

I'll just restate that I think that 2) is a bad idea because the body force is not related to the material model. I know that it would be convenient to put it there, but I think we will come to regret doing that in the long run.

[tjhei]

@gassmoeller ? ;-)

[jdannberg]

I don't like 1) because we will end up copying a lot of code each time we make a new model that uses the force vector, and I think it is a useful feature and could be used for almost every benchmark with a manufactured solution.

[bangerth]

@jdannberg 's argument is good. I change my mind to 3, 4, 1, 2.

[bangerth]

Superseded by #1547.