[cig-commits] [commit] master: Create a postprocessor for complete material output. (05c01a5)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Fri Aug 1 19:23:13 PDT 2014
Repository : https://github.com/geodynamics/aspect
On branch : master
Link : https://github.com/geodynamics/aspect/compare/abae02cb1cb943692c91540d97c5066f3d74a54b...427ddaac5731732db91611fc4e41fc781a3f97ad
>---------------------------------------------------------------
commit 05c01a56e232f46bdac1e0a3c54f8fbfd16c05e3
Author: Rene Gassmoeller <R.Gassmoeller at mailbox.org>
Date: Mon Jul 28 19:55:29 2014 -0700
Create a postprocessor for complete material output.
>---------------------------------------------------------------
05c01a56e232f46bdac1e0a3c54f8fbfd16c05e3
...thermal_expansivity.h => material_properties.h} | 46 +++-
.../visualization/material_properties.cc | 239 +++++++++++++++++++++
2 files changed, 277 insertions(+), 8 deletions(-)
diff --git a/include/aspect/postprocess/visualization/thermal_expansivity.h b/include/aspect/postprocess/visualization/material_properties.h
similarity index 62%
copy from include/aspect/postprocess/visualization/thermal_expansivity.h
copy to include/aspect/postprocess/visualization/material_properties.h
index 50657f1..ed24bd7 100644
--- a/include/aspect/postprocess/visualization/thermal_expansivity.h
+++ b/include/aspect/postprocess/visualization/material_properties.h
@@ -1,5 +1,5 @@
/*
- Copyright (C) 2011, 2012 by the authors of the ASPECT code.
+ Copyright (C) 2014 by the authors of the ASPECT code.
This file is part of ASPECT.
@@ -19,8 +19,8 @@
*/
-#ifndef __aspect__postprocess_visualization_thermal_expansivity_h
-#define __aspect__postprocess_visualization_thermal_expansivity_h
+#ifndef __aspect__postprocess_visualization_material_properties_h
+#define __aspect__postprocess_visualization_material_properties_h
#include <aspect/postprocess/visualization.h>
#include <aspect/simulator_access.h>
@@ -36,20 +36,32 @@ namespace aspect
{
/**
* A class derived from DataPostprocessor that takes an output vector
- * and computes a variable that represents the thermal expansivity at
- * every point.
+ * and computes a variable that represents the reaction term for every
+ * compositional field at every point.
*
* The member functions are all implementations of those declared in the
* base class. See there for their meaning.
*/
template <int dim>
- class ThermalExpansivity
- : public DataPostprocessorScalar<dim>,
+ class MaterialProperties
+ : public DataPostprocessor<dim>,
public SimulatorAccess<dim>,
public Interface<dim>
{
public:
- ThermalExpansivity ();
+ MaterialProperties ();
+
+ virtual
+ std::vector<std::string>
+ get_names () const;
+
+ virtual
+ std::vector<DataComponentInterpretation::DataComponentInterpretation>
+ get_data_component_interpretation () const;
+
+ virtual
+ UpdateFlags
+ get_needed_update_flags () const;
virtual
void
@@ -59,6 +71,24 @@ namespace aspect
const std::vector<Point<dim> > &normals,
const std::vector<Point<dim> > &evaluation_points,
std::vector<Vector<double> > &computed_quantities) const;
+
+
+ /**
+ * Declare the parameters this class takes through input files.
+ */
+ static
+ void
+ declare_parameters (ParameterHandler &prm);
+
+ /**
+ * Read the parameters this class declares from the parameter file.
+ */
+ virtual
+ void
+ parse_parameters (ParameterHandler &prm);
+
+ private:
+ std::vector<std::string> property_names;
};
}
}
diff --git a/source/postprocess/visualization/material_properties.cc b/source/postprocess/visualization/material_properties.cc
new file mode 100644
index 0000000..b570393
--- /dev/null
+++ b/source/postprocess/visualization/material_properties.cc
@@ -0,0 +1,239 @@
+/*
+ Copyright (C) 2014 by the authors of the ASPECT code.
+
+ This file is part of ASPECT.
+
+ ASPECT is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2, or (at your option)
+ any later version.
+
+ ASPECT is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with ASPECT; see the file doc/COPYING. If not see
+ <http://www.gnu.org/licenses/>.
+*/
+
+
+#include <aspect/postprocess/visualization/material_properties.h>
+#include <aspect/simulator_access.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <algorithm>
+
+
+namespace aspect
+{
+ namespace Postprocess
+ {
+ namespace VisualizationPostprocessors
+ {
+ template <int dim>
+ MaterialProperties<dim>::
+ MaterialProperties ()
+ :
+ DataPostprocessor<dim> ()
+ {}
+
+ template <int dim>
+ std::vector<std::string>
+ MaterialProperties<dim>::
+ get_names () const
+ {
+ std::vector<std::string> solution_names;
+
+ for (unsigned int i=0; i<property_names.size(); ++i)
+ if (property_names[i] == "reaction terms")
+ {
+ for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+ solution_names.push_back (this->introspection().name_for_compositional_index(c) + "_change");
+ }
+ else
+ {
+ solution_names.push_back(property_names[i]);
+ std::replace(solution_names.back().begin(),solution_names.back().end(),' ', '_');
+ }
+
+ return solution_names;
+ }
+
+ template <int dim>
+ std::vector<DataComponentInterpretation::DataComponentInterpretation>
+ MaterialProperties<dim>::
+ get_data_component_interpretation () const
+ {
+ std::vector<DataComponentInterpretation::DataComponentInterpretation> interpretation;
+ for (unsigned int i=0; i<property_names.size(); ++i)
+ {
+ if (property_names[i] == "reaction terms")
+ {
+ for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+ interpretation.push_back (DataComponentInterpretation::component_is_scalar);
+ }
+ else
+ interpretation.push_back (DataComponentInterpretation::component_is_scalar);
+ }
+
+ return interpretation;
+ }
+
+ template <int dim>
+ UpdateFlags
+ MaterialProperties<dim>::
+ get_needed_update_flags () const
+ {
+ return update_gradients | update_values | update_q_points;
+ }
+
+ template <int dim>
+ void
+ MaterialProperties<dim>::
+ compute_derived_quantities_vector (const std::vector<Vector<double> > &uh,
+ const std::vector<std::vector<Tensor<1,dim> > > &duh,
+ const std::vector<std::vector<Tensor<2,dim> > > &dduh,
+ const std::vector<Point<dim> > &normals,
+ const std::vector<Point<dim> > &evaluation_points,
+ std::vector<Vector<double> > &computed_quantities) const
+ {
+ const unsigned int n_quadrature_points = uh.size();
+ Assert (computed_quantities.size() == n_quadrature_points, ExcInternalError());
+ Assert (uh[0].size() == this->introspection().n_components, ExcInternalError());
+
+ typename MaterialModel::Interface<dim>::MaterialModelInputs in(n_quadrature_points,
+ this->n_compositional_fields());
+ typename MaterialModel::Interface<dim>::MaterialModelOutputs out(n_quadrature_points,
+ this->n_compositional_fields());
+
+ in.position = evaluation_points;
+ for (unsigned int q=0; q<n_quadrature_points; ++q)
+ {
+ Tensor<2,dim> grad_u;
+ for (unsigned int d=0; d<dim; ++d)
+ grad_u[d] = duh[q][d];
+ in.strain_rate[q] = symmetrize (grad_u);
+
+ in.pressure[q]=uh[q][this->introspection().component_indices.pressure];
+ in.temperature[q]=uh[q][this->introspection().component_indices.temperature];
+
+ for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+ in.composition[q][c] = uh[q][this->introspection().component_indices.compositional_fields[c]];
+ }
+
+ this->get_material_model().evaluate(in, out);
+
+ for (unsigned int q=0; q<n_quadrature_points; ++q)
+ {
+ unsigned k = 0;
+ for (unsigned int i=0; i<property_names.size(); ++i)
+ {
+ if (property_names[i] == "viscosity")
+ computed_quantities[q][i+k] = out.viscosities[q];
+
+ else if (property_names[i] == "density")
+ computed_quantities[q][i+k] = out.densities[q];
+
+ else if (property_names[i] == "thermal expansivity")
+ computed_quantities[q][i+k] = out.thermal_expansion_coefficients[q];
+
+ else if (property_names[i] == "specific heat")
+ computed_quantities[q][i+k] = out.specific_heat[q];
+
+ else if (property_names[i] == "thermal conductivity")
+ computed_quantities[q][i+k] = out.thermal_conductivities[q];
+
+ else if (property_names[i] == "compressibility")
+ computed_quantities[q][i+k] = out.compressibilities[q];
+
+ else if (property_names[i] == "entropy derivative pressure")
+ computed_quantities[q][i+k] = out.entropy_derivative_pressure[q];
+
+ else if (property_names[i] == "entropy derivative temperature")
+ computed_quantities[q][i+k] = out.entropy_derivative_temperature[q];
+
+ else if (property_names[i] == "reaction terms")
+ {
+ for (k=0; k<this->n_compositional_fields(); ++k)
+ {
+ computed_quantities[q][i+k] = out.reaction_terms[q][k];
+ }
+ --k;
+ }
+ }
+ }
+ }
+
+ template <int dim>
+ void
+ MaterialProperties<dim>::declare_parameters (ParameterHandler &prm)
+ {
+ prm.enter_subsection("Postprocess");
+ {
+ prm.enter_subsection("Visualization");
+ {
+ prm.enter_subsection("Material properties");
+ {
+ const std::string pattern_of_names
+ = "viscosity|density|thermal expansivity|specific heat|"
+ "thermal conductivity|compressibility|entropy derivative temperature|"
+ "entropy derivative pressure|reaction terms";
+
+ prm.declare_entry("List of material properties",
+ "density,thermal expansivity,specific heat,viscosity",
+ Patterns::MultipleSelection(pattern_of_names),
+ "A comma separated list of material properties that should be "
+ "written whenever writing graphical output. By default, the "
+ "material properties will always contain the density, thermal "
+ "expansivity, specific heat and viscosity. "
+ "The following material properties are available:\n\n"
+ +
+ pattern_of_names);
+ }
+ prm.leave_subsection();
+ }
+ prm.leave_subsection();
+ }
+ prm.leave_subsection();
+ }
+
+
+ template <int dim>
+ void
+ MaterialProperties<dim>::parse_parameters (ParameterHandler &prm)
+ {
+ prm.enter_subsection("Postprocess");
+ {
+ prm.enter_subsection("Visualization");
+ {
+ prm.enter_subsection("Material properties");
+ {
+ property_names = Utilities::split_string_list(prm.get ("List of material properties"));
+ }
+ prm.leave_subsection();
+ }
+ prm.leave_subsection();
+ }
+ prm.leave_subsection();
+ }
+ }
+ }
+}
+
+
+// explicit instantiations
+namespace aspect
+{
+ namespace Postprocess
+ {
+ namespace VisualizationPostprocessors
+ {
+ ASPECT_REGISTER_VISUALIZATION_POSTPROCESSOR(MaterialProperties,
+ "material properties",
+ "A visualization output object that generates output "
+ "for the material properties given by the material model.")
+ }
+ }
+}
More information about the CIG-COMMITS
mailing list