Topography on lower surface (e.g. CMB) #1519

jaustermann · 2017-05-08T18:46:34Z

I've worked with @ian-r-rose on this extension to the existing dynamic topography processor. It now also calculates topography at the lower surface of the model domain. I've only added this to the visualization postprocessor for now. The equivalent extension has to be done in the dynamic topography postprocessor that outputs the values in a file. I will do that once I have your comments on this file. This pull request will also break all existing dynamic topography tests and I can update those once the tester has run.

gassmoeller

Looks good except for the one variable name and a few small typos. Feel free to modify the other postprocessor as well, either in this PR or in a separate one. Also this change definitely deserves a file in doc/modules/changes.

gassmoeller · 2017-05-08T20:24:40Z

source/postprocess/visualization/dynamic_topography.cc

+                  {
+                    const double depth_face_center = this->get_geometry_model().depth (cell->face(f)->center());
+                    const double upper_depth_cutoff = cell->face(f)->minimum_vertex_distance()/3;
+                    const double lower_depth_cutoff = this->get_geometry_model().maximal_depth() -  cell->face(f)->minimum_vertex_distance()/3;


there is an additional space in the formula after the minus

gassmoeller · 2017-05-08T20:26:50Z

source/postprocess/visualization/dynamic_topography.cc

-                    if (cell->at_boundary(f) && this->get_geometry_model().depth (cell->face(f)->center()) < cell->face(f)->minimum_vertex_distance()/3)
+                // see if the cell is at the *top* or *bottom* boundary, not just any boundary
+                unsigned int face_idx = numbers::invalid_unsigned_int;
+                bool at_upper_surface;


Please add a comment here that the variable is set to false if the face is at the lower boundary. Maybe a name like top_or_bottom_boundary_face or lower_or_upper_boundary_face would make that clearer?

Actually I think at_upper_surface might make more sense because then true has a clear meaning, which isn't the case if the name is top_or_bottom_boundary_face. But I agree that I should add a comment that false corresponds to the lower boundary. On it!

You are right, just add the comment and it should be fine.

gassmoeller · 2017-05-08T20:27:47Z

source/postprocess/visualization/dynamic_topography.cc

+                {
+                  const double depth_face_center = this->get_geometry_model().depth (cell->face(f)->center());
+                  const double upper_depth_cutoff = cell->face(f)->minimum_vertex_distance()/3;
+                  const double lower_depth_cutoff = this->get_geometry_model().maximal_depth() -  cell->face(f)->minimum_vertex_distance()/3;


same double space as above

gassmoeller · 2017-05-08T20:28:14Z

source/postprocess/visualization/dynamic_topography.cc

@@ -215,7 +253,8 @@ namespace aspect
              prm.declare_entry ("Subtract mean of dynamic topography", "false",
                                 Patterns::Bool (),
                                 "Option to remove the mean dynamic topography "
-                                 "in the outputted data file (not visualization). "
+                                 "in the visualization. The mean dynamic topography is calculated "
+                                 "and subtracted separately for the upper and lower boundary. "


Nice catch!

gassmoeller · 2017-05-08T20:28:37Z

source/postprocess/visualization/dynamic_topography.cc

+              prm.declare_entry ("Density below","9900",
+                                 Patterns::Double (0),
+                                 "Dynamic topography is calculated as the excess or lack of mass that is supported by mantle flow. "
+                                 "This value depends on the density of material that is moved up or down, i.e. mantle above CMB,  and the "


double space

tjhei · 2017-05-09T00:30:30Z

source/postprocess/visualization/dynamic_topography.cc

-                    if (cell->at_boundary(f) && this->get_geometry_model().depth (cell->face(f)->center()) < cell->face(f)->minimum_vertex_distance()/3)
+                // see if the cell is at the *top* or *bottom* boundary, not just any boundary
+                unsigned int face_idx = numbers::invalid_unsigned_int;
+                bool at_upper_surface;


can it happen that this variable stays undefined? (neither one of the if's below is true)

Additionally, can a cell be on the "upper" and the "lower" surface at the same time? If not, can this be checked with an Assert()?

@tjhei If at_upper_surface is undefined face_idx is also not assigned and we check for that in the next if loop
if (face_idx == numbers::invalid_unsigned_int)
continue
If it's not assigned we continue to the next cell so I don't think there's an issue.
I will add an Assert statement that addresses your second comment.

tjhei · 2017-05-09T00:30:46Z

source/postprocess/visualization/dynamic_topography.cc

+                for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+                  {
+                    const double depth_face_center = this->get_geometry_model().depth (cell->face(f)->center());
+                    const double upper_depth_cutoff = cell->face(f)->minimum_vertex_distance()/3;


can you write 3.0 instead of 3 here and at the other places?

tjhei · 2017-05-09T00:32:17Z

source/postprocess/visualization/dynamic_topography.cc

@@ -215,7 +253,8 @@ namespace aspect
              prm.declare_entry ("Subtract mean of dynamic topography", "false",
                                 Patterns::Bool (),
                                 "Option to remove the mean dynamic topography "
-                                 "in the outputted data file (not visualization). "
+                                 "in the visualization. The mean dynamic topography is calculated "
+                                 "and subtracted separately for the upper and lower boundary. "
                                 "'true' subtracts the mean, 'false' leaves "
                                 "the calculated dynamic topography as is. ");


remove space at the end

tjhei · 2017-05-09T00:32:26Z

source/postprocess/visualization/dynamic_topography.cc

+                                 "density of the material that is displaced (generally outer core material). While the density of mantle rock "
+                                 "is part of the material model, this parameter `Density below' allows the user to specify the density "
+                                 "value of material that is displaced below the solid surface. By default this material is assumed to "
+                                 "be outer core material with a density of 9900. "


remove space at the end

jaustermann · 2017-05-09T02:47:10Z

@gassmoeller @tjhei I took another stab at this and also included it in the postprocessor/dynamic_topography.cc. In the latter the CMB topography is outputted as a separate file and only if the parameter 'Output lower surface topography' is set to true. If it's not set to true it's not even calculated.

gassmoeller · 2017-05-09T03:09:03Z

/run-tests

jaustermann · 2017-05-09T14:45:20Z

@gassmoeller: I rebased to master and updated existing dynamic topography tests to also output / test the lower boundary topography

jaustermann · 2017-05-09T15:54:14Z

@gassmoeller There seem to be three tests related to the stokes solver that don't fail because they compile. I'm not sure what the reason is for that. I checked and fixed all the topography related warnings and tests.

ian-r-rose · 2017-05-09T21:36:37Z

source/postprocess/dynamic_topography.cc

+                  if (depth_face_center < upper_depth_cutoff && depth_face_center > lower_depth_cutoff)
+                    AssertThrow(false, ExcMessage("Your geometry model so small that the upper and lower boundary of "
+                                                  "the domain are bordered by the same cell. "
+                                                  "Consider using a higher mesh resolution.") );


I like this error message.

bangerth · 2017-05-09T21:45:41Z

source/postprocess/visualization/dynamic_topography.cc

+
+                    // Check if cell is at upper and lower surface at the same time
+                    if (depth_face_center < upper_depth_cutoff && depth_face_center > lower_depth_cutoff)
+                      AssertThrow(false, ExcMessage("Your geometry model so small that the upper and lower boundary of "


this sentence no verb :-)

gassmoeller · 2017-05-09T21:54:33Z

@jaustermann These tests crash with your PR, and they use the dynamic topography postprocessor. Please run:

ctest -R postprocess_iterations_Stokes_only
cd tests
../aspect postprocess_iterations_Stokes_only.x.prm

take a look at the error message and think about why that parameter file causes a floating point exception. Let me know if you want to discuss how to debug those errors.

jaustermann · 2017-05-10T15:35:46Z

@gassmoeller This might be ready to be looked at again whenever you have a chance.

gassmoeller

Ok, 4 small documentation issues remaining. Ready to merge once you address them.

gassmoeller · 2017-05-10T18:55:18Z

source/postprocess/dynamic_topography.cc

-              {
-                for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
-                  if (cell->at_boundary(f) && this->get_geometry_model().depth (cell->face(f)->center()) < cell->face(f)->minimum_vertex_distance()/3)
+// see if the cell is at the *top* or *bottom* boundary, not just any boundary


gassmoeller · 2017-05-10T18:55:39Z

include/aspect/postprocess/dynamic_topography.h

@@ -70,9 +70,19 @@ namespace aspect
        bool subtract_mean_dyn_topography;

        /**
+         * A parameter that specifies whether the lower topography is outputted.
+               */


gassmoeller · 2017-05-10T18:56:42Z

source/postprocess/dynamic_topography.cc

          prm.declare_entry ("Density above","0",
                             Patterns::Double (0),
                             "Dynamic topography is calculated as the excess or lack of mass that is supported by mantle flow. "
                             "This value depends on the density of material that is moved up or down, i.e. crustal rock, and the "
                             "density of the material that is displaced (generally water or air). While the density of crustal rock "
                             "is part of the material model, this parameter `Density above' allows the user to specify the density "
                             "value of material that is displaced above the solid surface. By default this material is assumed to "
-                             "be air, with a density of 0. "
+                             "be air, with a density of 0."


Keep the space at the end of the line, there is another line below.

gassmoeller · 2017-05-10T18:56:56Z

source/postprocess/dynamic_topography.cc

+                             "density of the material that is displaced (generally outer core material). While the density of mantle rock "
+                             "is part of the material model, this parameter `Density below' allows the user to specify the density "
+                             "value of material that is displaced below the solid surface. By default this material is assumed to "
+                             "be outer core material with a density of 9900."


also here, missing space at end of line

…opography at the lower surface; Adjust tests

gassmoeller · 2017-05-10T20:37:33Z

All tests pass, and all comments have been addressed.

gassmoeller reviewed May 8, 2017

View reviewed changes

tjhei requested changes May 9, 2017

View reviewed changes

jaustermann force-pushed the CMB_topography branch 2 times, most recently from 57043f3 to eb7b98b Compare May 9, 2017 02:45

jaustermann force-pushed the CMB_topography branch from eb7b98b to 14fb7a3 Compare May 9, 2017 14:44

jaustermann force-pushed the CMB_topography branch from 14fb7a3 to 6c7c731 Compare May 9, 2017 15:52

jaustermann force-pushed the CMB_topography branch 2 times, most recently from 8a88c41 to dafbb28 Compare May 9, 2017 19:24

ian-r-rose reviewed May 9, 2017

View reviewed changes

bangerth requested changes May 9, 2017

View reviewed changes

jaustermann force-pushed the CMB_topography branch 2 times, most recently from 354dcfe to e5706e9 Compare May 10, 2017 01:49

gassmoeller approved these changes May 10, 2017

View reviewed changes

jaustermann force-pushed the CMB_topography branch from e5706e9 to 5e21482 Compare May 10, 2017 19:02

Extension to dynamic topography postprocessor to also calculate the t…

33c6c09

…opography at the lower surface; Adjust tests

jaustermann force-pushed the CMB_topography branch from 5e21482 to 33c6c09 Compare May 10, 2017 19:05

gassmoeller merged commit 1e20e7b into geodynamics:master May 10, 2017

jaustermann deleted the CMB_topography branch June 22, 2018 21:09

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Topography on lower surface (e.g. CMB) #1519

Topography on lower surface (e.g. CMB) #1519

jaustermann commented May 8, 2017

gassmoeller left a comment

gassmoeller May 8, 2017

gassmoeller May 8, 2017

jaustermann May 8, 2017

gassmoeller May 9, 2017

gassmoeller May 8, 2017

gassmoeller May 8, 2017

gassmoeller May 8, 2017

tjhei May 9, 2017

tjhei May 9, 2017

jaustermann May 9, 2017

tjhei May 9, 2017

tjhei May 9, 2017

tjhei May 9, 2017

jaustermann commented May 9, 2017

gassmoeller commented May 9, 2017

jaustermann commented May 9, 2017 •

edited

jaustermann commented May 9, 2017

ian-r-rose May 9, 2017

bangerth May 9, 2017

gassmoeller commented May 9, 2017

jaustermann commented May 10, 2017

gassmoeller left a comment

gassmoeller May 10, 2017

gassmoeller May 10, 2017

gassmoeller May 10, 2017

gassmoeller May 10, 2017

gassmoeller commented May 10, 2017

Topography on lower surface (e.g. CMB) #1519

Topography on lower surface (e.g. CMB) #1519

Conversation

jaustermann commented May 8, 2017

gassmoeller left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jaustermann commented May 9, 2017

gassmoeller commented May 9, 2017

jaustermann commented May 9, 2017 • edited

jaustermann commented May 9, 2017

Choose a reason for hiding this comment

Choose a reason for hiding this comment

gassmoeller commented May 9, 2017

jaustermann commented May 10, 2017

gassmoeller left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

gassmoeller commented May 10, 2017

jaustermann commented May 9, 2017 •

edited