Functor material properties

[lindsayad]

@GiudGiud @rwcarlsen @joe61vette

As per https://github.com/idaholab/moose/commit/bec76f305fcc6d6f12f8bb34bcf5b6ab24276633 I think I should make a FunctorMaterial or something that doesn't get added in triplicate like our current Materials do

To-do:

• [x] Make all functor material property lambdas block restricted
• [x] Demonstrate feasibility of caching
• [x] Allow functions to be called like our other Functors. This would allow a FunctorMaterial version of, for instance, GenericFunctionMaterial

Closes #16809

[moosebuild]

Job Documentation on 0653f11 wanted to post the following:

View the site here

This comment will be updated on new commits.

[rwcarlsen]

This has some really good ideas. And a few things that I'm uncomfortable with. I love that we get some of that on-demand property calcs stuff I've been wanting. But I am concerned at the limitations of your Elem * arg to the lambda. What about qps, on nodes, etc.? If we are going to make the effort to start supporting a new type of material property, we should try to make it work for more use cases than just the particular (FV) itch being scratched here. Along similar lines, the lack of support for old/older/stateful props is an additional inconsistency. I think the API approach we take for a new type of material property should be carefully thought out and discussed so we can be sure we know what we want before we commit to bifurcating the APIs.

[lindsayad]

But I am concerned at the limitations of your Elem * arg to the lambda

You can add new geometric overloads at any time. For instance today I will add a FaceInfo version because I'm tired of the material property ghosting. Node can be similarly added. qp is a little less obvious for an on-demand calculation. I think a better argument for an on-demand calculation would be something like a std::pair<const Elem *, Point>

[GiudGiud]

looking at the commit rather than the message, wrt to INSFV you have not changed the behavior?

not sure how we would evaluate rho * v directly (eg not rho_f v_f) on the face unless we make it the variable instead of v

[lindsayad]

wrt to INSFV you have not changed the behavior?

The behavior should be changed. In the interior we would interpolate rhou together before, e.g. rhou_avg = (rho_elem u_elem + rho_neighbor u_neighbor) / 2. Now rho and u are interpolated separately and multiplied, e.g. rhou_avg = (rho_elem + rho_neighbor) / 2. * (u_elem + u_neighbor) / 2.

not sure how we would evaluate rho * v directly (eg not rho_f v_f) on the face unless we make it the variable instead of v

Correct. Right now with the functor implementation this would be the only way to "interpolate them together"

[GiudGiud]

Oh rho_u used to be the advected quantity, true true

[GiudGiud]

btw I was thinking you were going to evaluate rho(T_f, P_f) at the face, not really (rho_elem + rho_neighbor) / 2. The former requires a full reimplementation of fluid_properties right?

[lindsayad]

btw I was thinking you were going to evaluate rho(T_f, P_f) at the face

Yes any call to evaluate rho at a face will result in a query of fluid properties for rho_from_T_p(T_f, p_f). However a user could still ask for rho_elem and then they would get rho_from_T_p(T_elem, p_elem)

[lindsayad]

@joe61vette thinks are looking promising. Here are for a couple different Rayleigh numbers

[joe61vette]

I would say it looks more than a little promising. I also prefer keeping to the FVMatAdvection scheme whenever possible. It is so easy (using the AD system) to build up compound material properties (eg, alfarhou) that can be combinations of nonlinear variables and fluid properties to be advected.

[lindsayad]

It is still easy to build compound properties with the functor property scheme and it halves the necessary data members in residual objects. But if being able to do interpolations together, in the sense that (rhou)_avg != rho_avg * u_avg, is important, then that is a technical issue that would require more thought.

[joe61vette]

I didn't make my comment very clear (sorry). For the interpolations, I would want to interpolate between adv_quantity_elem and adv_quantity_neighbor. Where adv_quantity_elem would be a compound material property (eg, alfa rho u / eps). Where alfa & u would be non-linear variables (volume fraction and superficial velocity), rho would be a fluid property (fn(P,T)), and eps would be an AuxVar for the porosity. I believe this is what you are proposing.

[lindsayad]

So, functionally speaking, you can do this very easily with the current functor implementation, but your result would be this:

interp(alfa rho u / eps) = interp(alfa) interp(rho) interp(u) / interp(eps)

If this is not desirable to you, then I can think harder about how to do aggregated interpolations with the functor system

[joe61vette]

It has been a long time since I have looked into the proper way to go from the PDEs to the finite volume formulation. But I believe that interpolating each quantity separately could lead to problems. For example, let's consider a compound material property say "momentum_x", for our conditions momentum_x_elem < momentum_x_neighbor. If we interpolate each quantity separately, I don't believe there is a guarantee that:

momentum_x_elem < momentum_x_face < momentum_x_neighbor

I admit that I am not sure about this. Perhaps @GiudGiud can weigh in as well. At any rate, I thought this was the direction you were headed based on the above comment:

"The behavior should be changed. In the interior we would interpolate rhou together before, e.g. rhou_avg = (rho_elem u_elem + rho_neighbor u_neighbor) / 2. Now rho and u are interpolated separately and multiplied, e.g. rhou_avg = (rho_elem + rho_neighbor) / 2. * (u_elem + u_neighbor) / 2."

[GiudGiud]

I would rather we do not lose momentum and energy conservation over this. There could be an option to interpolate the combined property over evaluating each of its components. It's probably much cheaper to do too (but as often, might not show up in profiling)

[lindsayad]

It's probably much cheaper to do to

How is that? In my head there's the same number of operations

[lindsayad]

My hatred of performing material property evaluation "at neighbor" stems from the fact that it requires material property ghosting which has raised a whole slate of issues for us both in terms of code maintenance and in terms of user input file construction. If we can train developers to think about what they do at domain boundaries when using FVFluxKernels then some-day we could remove the need to do material property ghosting

[lindsayad]

I just don't like ghosting period. If we are at a domain boundary, we should evaluate a quantity based on Dirichlet information or with getExtrapolatedBoundaryFaceValue

[GiudGiud]

I guess I'm assuming the neighbor property has already been evaluated. But we could get rid of that entirely. An alternative to ghosting is using auxvariables for these quantities.

[lindsayad]

Ok this restores the aggregate interpolation in FVMatAdvection

[lindsayad]

With the most recent commit, we have the ability to either do a combination of atomic interpolations or a single interpolation. So we now have to decide what we want. @GiudGiud and I have been having some discussion about this on Slack. Moukalled's text does combo-atomic in some places and single interpolation in others. One thorn is that if you were to adopt the stance that all things should be interpolated together, then you run into an immediate road-block with Rhie-Chow interpolation of velocity. That Rhie-Chow interpolated velocity is then used as a piece of data which helps to determine other interpolations, such as upwind. So it's impossible to move ahead with a stance that all things should be interpolated together.

But one can reasonably ask whether the advected quantity should always be interpolated as a unit even though we are using this unique interpolation scheme for the advector.

I don't really have the answer. Intuitively, I feel like it should be interpolated as a unit, but as @GiudGiud will point out equation 16.6 in Moukalled does a multiplication of atomic interpolations. And then to further muddy the waters, a lot of quantities in section 15.4 are interpolated as a unit. The D coefficient, which we compute in our code, is equal to V / a, and both currently in our code (e.g. MOOSE master/devel) and in Moukalled's text we perform interpolation of the V / a quotient as opposed to performing division of atomic V and a interpolations.

So in summary I don't know what to do. We have the code capability with this PR to do either. I would like something that is as close to uniformly consistent as possible. I'll reiterate that my intuitive feeling, and seemingly @joe61vette 's as well, is to interpolate things as a unit (except of course for the Rhie-Chow velocity) but a literature reference would be preferable to intuition.

[lindsayad]

Maybe @tanoret has a thought on whether for a weakly compressible implementation, we should interpolate say rhou for momentum advection as rho_f u_f or as (rhou)_f ? (With full momentum flux given by rhou v_RC.) For a pure upwind scheme the two are equivalent, but they differ the moment we introduce any degree of central differencing

[lindsayad]

Mmm I don't think so? If a Material is not block restricted, then I query for the mesh subdomain IDs from the _mesh and then add key-value pairs for each block ID and the provided lambda.

Something that is true right now is that I do not do error checking for whether a functor material property is added multiply by different materials on the same block. I started out with error checking for that but then was reminded that we add materials in triplicate (volume, face, and neighbor). So I think at some point in this PR I will create a FunctorMaterial that is only singleton added. (the whole idea of these functor materials is that we don't need these degenerate materials)

[lindsayad]

In OpenFOAM's SIMPLE application, the mass flow rate at the face is the product of the linear interpolation of (density * velocity) dotted with the surface vector, but that is only in the constructor

[rwcarlsen]

Joe is right - part of the reason I didn't create an interpolation system that does more for you magically is because the different problems/methods in FV need to interpolate different aggregate quantities - so I opted to just leave it all up to the kernel devs to implement as necessary. A system that only allows us to do it one way (interpolate things all together or interpolate them individually and then combine) will always fail to work properly for certain use cases.

[dschwen]

Well, a customizable interpolation system would be very helpful to support stateful material properties with adaptivity. Currently we do neatest neighbor, which gives pronounced artifacts when running inelastic models.

On Mon, Aug 2, 2021, 10:56 AM Robert Carlsen @.***> wrote:

Joe is right - part of the reason I didn't create an interpolation system that doesn't do more for you magically is because the different problems/methods in FV need to interpolate different aggregate quantities - so I opted to just leave it all up to the kernel devs to implement as necessary. A system that only allows us to do it one way (interpolate things all together or interpolate them individually and then combine) will always fail to work properly for certain use cases.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/idaholab/moose/pull/18395#issuecomment-891179477, or unsubscribe https://github.com/notifications/unsubscribe-auth/AABRMPSXR5FKF5JO6HJQJYDT23E3NANCNFSM5AWO3ZAQ .

[lindsayad]

Joe is right - part of the reason I didn't create an interpolation system that does more for you magically is because the different problems/methods in FV need to interpolate different aggregate quantities - so I opted to just leave it all up to the kernel devs to implement as necessary. A system that only allows us to do it one way (interpolate things all together or interpolate them individually and then combine) will always fail to work properly for certain use cases.

We have not removed any of the global Moose::FV functions that allow you to do aggregate interpolations.

There are now overloads on FunctorInterface that allow you to ask for a quantity at a

• face. If you are calling this overload with a material property that is a function of multiple variables, then the result will be a function of the variables (individually) evaluated at the face
• element, but with an attached face information object as well. If the element has no dofs associated with it, then we assume that this is a ghost element, in which case we can determine what the ghost value should be

With the addition of the latter overload, the functor property is a complete substitute for regular properties in FVFluxKernels because you can (safely) ask for the functor property at both the element and neighbor and then do whatever interpolation you like to the face using our global Moose::FV functionality that we've always had

[lindsayad]

Ok it is my understanding that after Monday's call with the Griffin team (@YaqiWang) there is interest in using this new system in order to do material property caching, raised as an issue in #18202. If that is the case, then that is good motivation to split out the navier_stokes module modifications from this branch (which probably should be done anyway).

[lindsayad]

It would be good if someone pointed me to a problem whose cost is known to be heavy in material property evaluation and which from a Newton solve standpoint can have some or all of the property evaluations moved off of linear. I believe that @vincentlaboure had concocted such an example IIRC

[vincentlaboure]

@lindsayad This Griffin branch has such examples.

[lindsayad]

@vincentlaboure can you grant me access to your Griffin fork please?

[vincentlaboure]

@lindsayad Sorry for the delay, you should be good to go now

[lindsayad]

This is ready for some review

[rwcarlsen]

I'll give it some love on monday.

[GiudGiud]

quick look, I ll come back to it. Powerful stuff for sure

[lindsayad]

Ok if I run my caching test with n=20, I get the following timings:

• always clear cache (always evaluate): 172.3 seconds
• clear cache on linear iteration begin: 26.78 seconds
• clear cache on nonlinear iteration begin: 7.806 seconds
• clear cache on timestep begin: 6.829 seconds

[moosebuild]

Job Create Mac Stack on 624b162 : invalidated by @lindsayad

[lindsayad]

Want to add WCNSFVTimeDerivative kernels since that's basically the only thing missing for using WCNSFV?

I would like to wait on more modules-specific stuff for a follow-on PR if possible. This obviously already has some but its more just updating the code to use the new system, and then as a happy coincidence I can setup a WCNSFV steady test without adding any module code, just changing it.

[lindsayad]

I think this is ready for further review

[GiudGiud]

once this passes Precheck I m good to approve. @rwcarlsen want to review more? @friedmud I know you have opinions on the material system

[moosebuild]

Job Code coverage on c50f0ce wanted to post the following:

Generated code coverage:

moose
chemical_reactions
combined
contact
external_petsc_solver
fluid_properties
fsi
functional_expansion_tools
geochemistry
heat_conduction
level_set
misc
navier_stokes
peridynamics
phase_field
porous_flow
ray_tracing
rdg
richards
stochastic_tools
tensor_mechanics
xfem

This comment will be updated on new commits.

[lindsayad]

Ok I'm hoping I've addressed all conceptual concerns with user-facing design. I'm going to be looking at code coverage and making sure we're hitting everything

[lindsayad]

Ok I've completed coverage of FunctorInterface and demonstrated creation of FunctorMaterialProperty sub-classes. Have to start writing reports and applying this to get good results for NRC by December so if there's any major design overhauls requested this will probably sit in a branch for quite a while 😄

[GiudGiud]

where can I see the coverage?

[lindsayad]

@loganharbour and I are still getting everything sorted out. We will need https://github.com/idaholab/moose/pull/18778 to go in before we can see accurate coverage reports from CIVET

[loganharbour]

Yeah. Once #18778 gets into next and then devel, we should get accurate reporting.

[moosebuild]

Job Python 3.6 Tests on f8ee33a : invalidated by @lindsayad

[GiudGiud]

so for the Rhie Chow user object I m thinking we'll want to have neighbor overloads for the element, so we can compute the friction term on the neighbor completely (using the friction coef neighbor value, and the variable value there too). This system can do that right? (just need to add a new overload right). I m not asking for this to be added now btw

Then if we have that it's trivial to get some finite difference gradient going or some Green Gauss one

[GiudGiud]

oh yeah bad question, just need to fill the elem argument one with the neighbor and we re good. My bad!

[GiudGiud]

Thoughts about adding a ParsedFunctorMaterial, which creates the same functor for all overloads based on the parsed function we provide. I d be happy to do it