Prescribed Solution Post-processor by lhy11009 · Pull Request #6882 · geodynamics/aspect

lhy11009 · 2026-03-05T15:15:35Z

Pull Request Checklist. Please read and check each box with an X. Delete any part not applicable. Ask on the forum if you need help with any step.

Add Post-processor to the Prescribed Solution Plugins

Before your first pull request:

I have read the guidelines in our CONTRIBUTING.md document.

For all pull requests:

I have followed the instructions for indenting my code.

For new features/models or changes of existing features:

I have tested my new feature locally to ensure it is correct.
I have created a testcase for the new feature/benchmark in the tests/ directory.
I have added a changelog entry in the doc/modules/changes directory that will inform other users of my change.

lhy11009 · 2026-03-05T15:25:51Z

Hi @mibillen,

I figured it would be better for me to work on this myself for now since the structure involved is probably a bit too complicated for you and Sophia at this point. I can also reuse it in my current project.

In this unfinished PR, I have already addressed the following:

Included the prescribed solution in the simulator access interface.

Added a postprocess class in the same interface.cc file of the prescribed solution. I plan to move it later to a separate file under the postprocessors folder. This class is already registered and currently outputs a data array with trivial values for the prescribed conditions and the indicators with the test prm I include.

At this point, I have a specific question that I think @gassmoeller might be able to answer.

If you take a look at the evaluate_vector_field function in this new class: how do you think we should propagate the results from computing constraints on cells (which I already copied from the Manager class) to assigning them to the global vector? My understanding is that this is different from the current_constraints object used in the Manager class.

Also, please feel free to comment on any other parts of the code.

lhy11009 · 2026-03-05T15:27:21Z

+        };
+
+        void
+        evaluate_vector_field(const DataPostprocessorInputs::Vector<dim> &input_data,


@gassmoeller Here

mibillen · 2026-03-05T16:41:45Z

Hi Haoyuan, Thank you for taking the lead on this. I really appreciate it. Please let Sophia and I know if we can help with testing or anything else. Magali From: lhy11009 ***@***.***> Reply-To: geodynamics/aspect ***@***.***> Date: Thursday, March 5, 2026 at 7:26 AM To: geodynamics/aspect ***@***.***> Cc: Magali Isabelle Billen ***@***.***>, Mention ***@***.***> Subject: Re: [geodynamics/aspect] Prescribed Solution Post-processor (PR #6882) [https://avatars.githubusercontent.com/u/28665365?s=20&v=4]lhy11009 left a comment (geodynamics/aspect#6882)<#6882 (comment)> Hi @mibillen<https://github.com/mibillen>, I figured it would be better for me to work on this myself for now since the structure involved is probably a bit too complicated for you and Sophia at this point. I can also reuse it in my current project. In this unfinished PR, I have already addressed the following: Included the prescribed solution in the simulator access interface. Added a postprocess class in the same interface.cc file of the prescribed solution. I plan to move it later to a separate file under the postprocessors folder. This class is already registered and currently outputs a data array with trivial values for the prescribed conditions and the indicators with the test prm I include. At this point, I have a specific question that I think @gassmoeller<https://github.com/gassmoeller> might be able to answer. If you take a look at the evaluate_vector_field function in this new class: how do you think we should propagate the results from computing constraints on cells (which I already copied from the Manager class) to assigning them to the global vector? My understanding is that this is different from the current_constraints object used in the Manager class. Also, please feel free to comment on any other parts of the code. — Reply to this email directly, view it on GitHub<#6882 (comment)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/ACS64UUFVEHCG7ZB2XLAZJT4PGMBNAVCNFSM6AAAAACWIMKN2KVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHM2DAMBVHA3TAMZSGU>. You are receiving this because you were mentioned.Message ID: ***@***.***>

gassmoeller

I will wait with more comments until the code is in the right places, but generally speaking this is close to what I had imagined.

Nice work figuring out the details!

gassmoeller · 2026-03-06T15:04:17Z

+          std::vector<DataComponentInterpretation::DataComponentInterpretation> interpretation;
+          for (auto &model_name : prescribed_solution_manager.get_active_plugin_names())
+            {
+              // for values
+              if (model_name == "velocity function")
+                interpretation.push_back (DataComponentInterpretation::component_is_part_of_vector);
+              else
+                interpretation.push_back (DataComponentInterpretation::component_is_scalar);
+
+              // for indicators
+              interpretation.push_back (DataComponentInterpretation::component_is_scalar);
+            }


I see you wanted to make this generally useful. I would recommend to choose one of two simplifications:

either make this a scalar postprocessor that only outputs for a single (possibly user chosen) component where it is prescribed

or let it always create as many output fields as there are components (introspection.n_components), and then output for every field where it is prescribed (just outputting 0 for all fields that are not prescribed). This saves you the problems of figuring out which prescribed solution plugins are active.

I took the first approach. In the visualization section, I let the user select the prescribed solution model they want to visualize. But I didn't restrict it to scalar fields, because I still include the velocity function as an option, and that is a vector field.

gassmoeller · 2026-03-06T15:09:13Z

+          for (unsigned int q=0; q<n_quadrature_points; ++q)
+            {
+              for (unsigned i = 0; i < plugin_names.size(); ++i)
+                {
+                  computed_quantities[q][i] = 1.0; // for values
+
+                  ++i;
+                  computed_quantities[q][i] = 1.0; // for indicators
+                }
+            }


Do you intend to have different output fields for indicators and values? I would suggest to have just set of fields, and set it to 0 where should_be_constrained is false, and set it to the value where should_be_constrained is true.

I still have both indicators and values as output, and they are laid out as "temperature function" and "temperature function indicator". My reasoning is that the prescribed conditions are not something that users could directly inspect from other outputs, so outputting both the indicator and the values give user much better control.

For velocity, although there are 3 different directions, I only write one indicator, and that is the maximum of the three. This is similar to what we did when the spherical unit vectors were used to prescribe the velocities.

lhy11009 · 2026-03-07T18:14:01Z

@gassmoeller. This is ready to be fully reviewed. And I'll then mark some places where my code needs some suggestions to be improved. Otherwise, this is serving our purposes. Now, it can separately output indicators and values. I could also handle fields related to both the "temperature function" and "velocity function" in parallel.
For example, these are the visualization results of the temperature at timestep 0.

gassmoeller

There are still some smaller comments about minor issues I need to make, but it is late on a Friday, so I will get back to you when I can.

But I had a larger comment that I wanted to leave now already, to check what you think about it.

gassmoeller · 2026-03-13T14:59:12Z

+              const std::string pattern_of_names
+                = "temperature function|velocity function";
+
+              prm.declare_entry("List of model names",
+                                "temperature function",
+                                Patterns::MultipleSelection(pattern_of_names),
+                                "A comma-separated list of prescribed solution models that "
+                                "will be written as both indicators and values"
+                                "The following prescribed solution models are available for output:\n\n"
+                                +
+                                pattern_of_names);


I dont think the approach of selecting the output by the prescribed solution plugin is useful. It has at least two problems:

selecting a plugin by name is hard (see above), and your pattern_of_names is hard-coded and can get out of sync with the actual plugin list

What do you do if one plugin is written to prescribe the velocity and compositions 2,4 and 5, but not the other ones? Then you dont know how many components this postprocessor should have.

Like I mentioned in my earlier review, I think you really only have the two options of:

Only output one component, temperature, one compositional field, one velocity component and let the user specify which component (this is similar to what we do in the artificial_viscosity_composition postprocessor)

or: output all components and let the postprocessor essentially output two fields for every solution component

Both is possible, and both is simpler to implement than what you are trying at the moment. The only point I could see in selecting prescribed solution plugins would be to check them individually, but you could achieve the same result by only activating a single plugin at a time.

In the end I think of this postprocessor as a debugging tool, something that will mostly be used to set up a model, and disabled for production runs. It would be best if kept as simple as possible (i.e. without complicated logic about which plugin to query).

@gassmoeller I might have misunderstood this point from the start. For the first option, it seems that this is not referring to the type of plugin (e.g., temperature function vs. initial function), but rather to applying it to all prescribed temperature plugins, with the user specifying which component to constrain. By "let the user specify which component", do you mean the user should explicitly indicate components such as "temperature" or "velocity"?

For the second point, do you mean that we should loop over all solution components (e.g., temperature, velocity, etc.) rather than iterating over each specific plugin?

I'll make sure I get this right before proceeding.

I changed the behavior to output all components and let the postprocessor essentially output two fields for every solution component. This indeed simplifies the code structure a lot.

lhy11009 · 2026-03-25T06:27:39Z

I have changed the implementation and adjusted the tests. But let's wait for the other PR (#6889) to be merged so that I can use it to test the composition output.

lhy11009 · 2026-03-25T06:28:52Z

In addition, below are the velocity indicator and velocity magnitude at timestep 1, along with the temperature outputs shown before.

lhy11009 commented Mar 5, 2026

View reviewed changes

lhy11009 mentioned this pull request Mar 5, 2026

Work plan to merge my 3-D free subduction model as cookbook #6838

Open

10 tasks

lhy11009 changed the title ~~Prescribed Solution Post-processor~~ Prescribed Solution Post-processor (WIP) Mar 6, 2026

gassmoeller reviewed Mar 6, 2026

View reviewed changes

lhy11009 force-pushed the prescribed_solution_postp branch 2 times, most recently from 1ec7c01 to cf62e6d Compare March 7, 2026 04:11

lhy11009 force-pushed the prescribed_solution_postp branch from cf62e6d to 978940b Compare March 7, 2026 18:16