[cig-commits] (no subject)
dealii.demon at gmail.com
dealii.demon at gmail.com
Wed Jul 24 13:15:00 PDT 2013
-e To: cig-commits at geodynamics.org
Subject: commit by buerg to /var/svn/dealii/aspect
Revision 1819
Split DoFHandler.
U trunk/aspire/cantera.prm
U trunk/aspire/include/aspect/introspection.h
U trunk/aspire/include/aspect/simulator.h
U trunk/aspire/include/aspect/simulator_access.h
U trunk/aspire/source/mesh_refinement/composition.cc
U trunk/aspire/source/mesh_refinement/density.cc
U trunk/aspire/source/mesh_refinement/temperature.cc
U trunk/aspire/source/mesh_refinement/thermal_energy_density.cc
U trunk/aspire/source/mesh_refinement/velocity.cc
U trunk/aspire/source/postprocess/composition_statistics.cc
U trunk/aspire/source/postprocess/heat_flux_statistics.cc
U trunk/aspire/source/postprocess/temperature_statistics.cc
U trunk/aspire/source/postprocess/velocity_statistics.cc
U trunk/aspire/source/postprocess/visualization.cc
U trunk/aspire/source/simulator/assembly.cc
U trunk/aspire/source/simulator/checkpoint_restart.cc
U trunk/aspire/source/simulator/core.cc
U trunk/aspire/source/simulator/helper_functions.cc
U trunk/aspire/source/simulator/initial_conditions.cc
U trunk/aspire/source/simulator/introspection.cc
U trunk/aspire/source/simulator/parameters.cc
U trunk/aspire/source/simulator/simulator_access.cc
U trunk/aspire/source/simulator/solver.cc
U trunk/aspire/source/termination_criteria/steady_rms_velocity.cc
http://www.dealii.org/websvn/revision.php?repname=Aspect+Repository&path=%2F&rev=1819&peg=1819
Diff:
Modified: trunk/aspire/cantera.prm
===================================================================
--- trunk/aspire/cantera.prm 2013-07-16 23:37:10 UTC (rev 1818)
+++ trunk/aspire/cantera.prm 2013-07-24 20:14:34 UTC (rev 1819)
@@ -1,9 +1,6 @@
set End time = 100
-set Time step = 1e-4
-set Iterate viscosity = false
-set Nonlinear solver scheme = fixed point iteration
-set Nonlinear solver tolerance = 1e-3
-set Max nonlinear iterations = 20
+set Time step = 1e-12
+set Nonlinear solver scheme = IMPES
set Number of cheap Stokes solver steps = 0
subsection Discretization
@@ -45,12 +42,12 @@
end
subsection Geometry model
-# set Model name = small axisymmetric jetflow apparatus
- set Model name = file
- subsection File
- set format = msh
- set name = jet.msh
- end
+ set Model name = axisymmetric jetflow apparatus
+#set Model name = file
+# subsection File
+# set format = msh
+# set name = jet.msh
+# end
end
subsection Gravity model
@@ -66,12 +63,12 @@
end
subsection Mesh refinement
- set Initial global refinement = 1
- set Initial adaptive refinement = 3
- set Time steps between mesh refinement = 4
+ set Initial global refinement = 0
+ set Initial adaptive refinement = 0
+ set Time steps between mesh refinement = 1000000000
set Strategy = density
- set Coarsening fraction = 0.01
- set Refinement fraction = 0.3
+ set Coarsening fraction = 0.0
+ set Refinement fraction = 0.0
end
@@ -85,11 +82,6 @@
end
end
-subsection Mesh refinement
- set Initial global refinement = 1
- set Initial adaptive refinement = 0
-end
-
subsection Model settings
set Prescribed velocity boundary indicators = 1: axisymmetric jetflow apparatus
set Tangential velocity boundary indicators = 0,2
Modified: trunk/aspire/include/aspect/introspection.h
===================================================================
--- trunk/aspire/include/aspect/introspection.h 2013-07-16 23:37:10 UTC (rev 1818)
+++ trunk/aspire/include/aspect/introspection.h 2013-07-24 20:14:34 UTC (rev 1819)
@@ -68,34 +68,15 @@
* @{
*/
/**
- * The number of vector components used by the finite element
- * description of this problem. It equals $d+2+n_c$ where
- * $d$ is the dimension and equals the number of velocity components,
- * and $n_c$ is the number of advected (compositional) fields. The
- * remaining components are the scalar pressure and temperature
- * fields.
- */
- const unsigned int n_components;
-
- /**
- * The number of vector blocks. This equals $3+n_c$ where,
- * in comparison to the n_components field, the velocity components
- * form a single block.
- */
- const unsigned int n_blocks;
-
- /**
* A structure that contains FEValues extractors for every block
* of the finite element used in the overall description.
*/
struct Extractors
{
- Extractors (const unsigned int n_compositional_fields);
+ Extractors ();
- const FEValuesExtractors::Vector velocities;
- const FEValuesExtractors::Scalar pressure;
- const FEValuesExtractors::Scalar temperature;
- const std::vector<FEValuesExtractors::Scalar> compositional_fields;
+ const FEValuesExtractors::Vector velocities;
+ const FEValuesExtractors::Scalar pressure;
};
/**
* A variable that contains extractors for every block
@@ -110,12 +91,8 @@
*/
struct ComponentIndices
{
- ComponentIndices (const unsigned int n_compositional_fields);
-
- static const unsigned int velocities[dim];
- static const unsigned int pressure = dim;
- static const unsigned int temperature = dim+1;
- const std::vector<unsigned int> compositional_fields;
+ static const unsigned int velocities[dim];
+ static const unsigned int pressure = dim;
};
/**
* A variable that enumerates the vector components of the finite
@@ -129,12 +106,8 @@
*/
struct BlockIndices
{
- BlockIndices (const unsigned int n_compositional_fields);
-
- static const unsigned int velocities = 0;
- static const unsigned int pressure = 1;
- static const unsigned int temperature = 2;
- const std::vector<unsigned int> compositional_fields;
+ static const unsigned int velocities = 0;
+ static const unsigned int pressure = 1;
};
/**
* A variable that enumerates the vector blocks of the finite
@@ -152,9 +125,6 @@
{
ComponentMask velocities;
ComponentMask pressure;
- ComponentMask temperature;
- std::vector<ComponentMask> compositional_fields;
- ComponentMask all_compositional_fields;
};
/**
* A variable that contains component masks for each of the variables
@@ -167,7 +137,7 @@
* A variable that describes for each vector component which vector
* block it corresponds to.
*/
- const std::vector<unsigned int> components_to_blocks;
+ const std::vector<unsigned int> velocity_components_to_blocks;
/**
* @}
@@ -178,11 +148,11 @@
* @{
*/
/**
- * A variable that describes how many of the degrees of freedom
- * on the current mesh belong to each of the n_blocks blocks
- * of the finite element.
+ * A variable that describes how many of the velocity degrees of
+ * freedom on the current mesh belong to each of the n_blocks
+ * blocks of the finite element.
*/
- std::vector<types::global_dof_index> system_dofs_per_block;
+ std::vector<types::global_dof_index> velocity_dofs_per_block;
/**
* A structure that contains index sets describing which of the
@@ -192,30 +162,61 @@
struct IndexSets
{
/**
- * An index set that indicates which (among all) degrees of freedom
+ * An index set that indicates which (among all) velocity degrees of freedom
* are relevant to the current processor. See the deal.II documentation
* for the definition of the term "locally relevant degrees of freedom".
*/
- IndexSet system_relevant_set;
+ IndexSet velocity_relevant_set;
/**
+ * An index set that indicates which (among all) temperature degrees of freedom
+ * are relevant to the current processor. See the deal.II documentation
+ * for the definition of the term "locally relevant degrees of freedom".
+ */
+ IndexSet temperature_relevant_set;
+
+ /**
+ * An index set that indicates which (among all) composition degrees of freedom
+ * are relevant to the current processor. See the deal.II documentation
+ * for the definition of the term "locally relevant degrees of freedom".
+ */
+ IndexSet composition_relevant_set;
+
+ /**
* A collection of index sets that for each of the vector blocks
* of this finite element represents the global indices of the
- * degrees of freedom owned by this processor. The n_blocks
+ * velocity degrees of freedom owned by this processor. The n_blocks
* elements of this array form a mutually exclusive decomposition
- * of the index set containing all locally owned degrees of freedom.
+ * of the index set containing all locally owned velocity degrees
+ * of freedom.
*/
- std::vector<IndexSet> system_partitioning;
+ std::vector<IndexSet> velocity_partitioning;
/**
+ * An index set that represents the global indices of the temperature
+ * degrees of freedom owned by this processor. The n_blocks elements
+ * of this array form a mutually exclusive decomposition of the index
+ * set containing all locally owned temperature degrees of freedom.
+ */
+ IndexSet temperature_partitioning;
+
+ /**
+ * An index set that represents the global indices of the composition
+ * degrees of freedom owned by this processor. The n_blocks elements
+ * of this array form a mutually exclusive decomposition of the index
+ * set containing all locally owned composition degrees of freedom.
+ */
+ IndexSet composition_partitioning;
+
+ /**
* A collection of index sets that for each of the vector blocks
* of this finite element represents the global indices of the
- * degrees of freedom are relevant to this processor. The n_blocks
+ * velocity degrees of freedom are relevant to this processor. The n_blocks
* elements of this array form a mutually exclusive decomposition
* of the index set containing all locally relevant degrees of freedom,
- * i.e., of the system_relevant_set index set.
+ * i.e., of the velocity_relevant_set index set.
*/
- std::vector<IndexSet> system_relevant_partitioning;
+ std::vector<IndexSet> velocity_relevant_partitioning;
};
/**
* A variable that contains index sets describing which of the
Modified: trunk/aspire/include/aspect/simulator.h
===================================================================
--- trunk/aspire/include/aspect/simulator.h 2013-07-16 23:37:10 UTC (rev 1818)
+++ trunk/aspire/include/aspect/simulator.h 2013-07-24 20:14:34 UTC (rev 1819)
@@ -27,7 +27,10 @@
#include <deal.II/base/parameter_handler.h>
#include <deal.II/base/conditional_ostream.h>
#include <deal.II/base/symmetric_tensor.h>
+#include <deal.II/base/synchronous_iterator.h>
+#include <deal.II/grid/filtered_iterator.h>
+
#include <deal.II/lac/trilinos_block_vector.h>
#include <deal.II/lac/trilinos_block_sparse_matrix.h>
#include <deal.II/lac/trilinos_precondition.h>
@@ -67,14 +70,16 @@
{
template <int dim> struct StokesPreconditioner;
template <int dim> struct StokesSystem;
- template <int dim> struct AdvectionSystem;
+ template <int dim> struct TemperatureSystem;
+ template <int dim> struct CompositionSystem;
}
namespace CopyData
{
template <int dim> struct StokesPreconditioner;
template <int dim> struct StokesSystem;
- template <int dim> struct AdvectionSystem;
+ template <int dim> struct TemperatureSystem;
+ template <int dim> struct CompositionSystem;
}
}
@@ -98,10 +103,7 @@
{
enum Kind
{
- fixed_point_iteration,
- IMPES,
- iterated_IMPES,
- iterated_Stokes
+ IMPES
};
};
@@ -169,11 +171,8 @@
unsigned int timing_output_frequency;
double linear_solver_tolerance;
unsigned int n_cheap_stokes_solver_steps;
- double nonlinear_solver_tolerance;
- unsigned int max_nonlinear_iterations;
double temperature_solver_tolerance;
double composition_solver_tolerance;
- bool iterate_viscosity;
/**
* @}
*/
@@ -307,73 +306,6 @@
private:
/**
- * A structure that is used as an argument to functions that can
- * work on both the temperature and the compositional variables
- * and that need to be told which one of the two, as well as on
- * which of the compositional variables.
- */
- struct TemperatureOrComposition
- {
- /**
- * An enum indicating whether the identified variable is the
- * temperature or one of the compositional fields.
- */
- enum FieldType { temperature_field, compositional_field };
-
- /**
- * A variable indicating whether the identified variable is the
- * temperature or one of the compositional fields.
- */
- const FieldType field_type;
-
- /**
- * A variable identifying which of the compositional fields is
- * selected. This variable is meaningless if the temperature
- * is selected.
- */
- const unsigned int compositional_variable;
-
- /**
- * Constructor.
- * @param field_type Determines whether this variable should select
- * the temperature field or a compositional field.
- * @param compositional_variable The number of the compositional field
- * if the first argument in fact chooses a compositional variable.
- * Meaningless if the first argument equals temperature.
- *
- * This function is implemented in
- * <code>source/simulator/helper_functions.cc</code>.
- */
- TemperatureOrComposition (const FieldType field_type,
- const unsigned int compositional_variable = numbers::invalid_unsigned_int);
-
- /**
- * A static function that creates an object identifying the temperature.
- *
- * This function is implemented in
- * <code>source/simulator/helper_functions.cc</code>.
- */
- static
- TemperatureOrComposition temperature ();
-
- /**
- * A static function that creates an object identifying given
- * compositional field.
- *
- * This function is implemented in
- * <code>source/simulator/helper_functions.cc</code>.
- */
- static
- TemperatureOrComposition composition (const unsigned int compositional_variable);
-
- /**
- * Return whether this object refers to the temperature field.
- */
- bool
- is_temperature () const;
- };
-
- /**
* @name Top-level functions in the overall flow of the numerical algorithm
* @{
*/
@@ -399,8 +331,8 @@
void setup_introspection ();
/**
- * A function that is responsible for initializing the temperature/compositional
- * field before the first time step. The temperature field then serves as the
+ * A function that is responsible for initializing the temperature field
+ * before the first time step. The temperature field then serves as the
* temperature from which the velocity is computed during the first time
* step, and is subsequently overwritten by the temperature field one gets
* by advancing by one time step.
@@ -408,9 +340,18 @@
* This function is implemented in
* <code>source/simulator/initial_conditions.cc</code>.
*/
- void set_initial_temperature_and_compositional_fields ();
+ void set_initial_temperature ();
/**
+ * A function that is responsible for initializing the compositional field
+ * before the first time step.
+ *
+ * This function is implemented in
+ * <code>source/simulator/initial_conditions.cc</code>.
+ */
+ void set_initial_compositional_fields ();
+
+ /**
* Do some housekeeping at the beginning of each time step. This includes
* generating some screen output, adding some information to the statistics
* file, and interpolating time-dependent boundary conditions specific to
@@ -443,14 +384,22 @@
void build_stokes_preconditioner ();
/**
- * Initialize the preconditioner for the advection equation of
+ * Initialize the preconditioner for the temperature equation.
+ *
+ * This function is implemented in
+ * <code>source/simulator/assembly.cc</code>.
+ */
+ void build_temperature_preconditioner (std_cxx1x::shared_ptr<aspect::LinearAlgebra::PreconditionBase> &preconditioner);
+
+ /**
+ * Initialize the preconditioner for the composition equation of
* field index.
*
* This function is implemented in
* <code>source/simulator/assembly.cc</code>.
*/
- void build_advection_preconditioner (const TemperatureOrComposition &temperature_or_composition,
- std_cxx1x::shared_ptr<aspect::LinearAlgebra::PreconditionBase> &preconditioner);
+ void build_composition_preconditioner (const unsigned int c,
+ std_cxx1x::shared_ptr<aspect::LinearAlgebra::PreconditionBase> &preconditioner);
/**
* Initiate the assembly of the Stokes matrix and right hand side.
@@ -461,43 +410,55 @@
void assemble_stokes_system ();
/**
- * Initiate the assembly of one advection matrix and right hand side
+ * Initiate the assembly of the temperature matrix and right hand side
+ * and build a preconditioner for this matrix.
+ *
+ * This function is implemented in
+ * <code>source/simulator/assembly.cc</code>.
+ */
+ void assemble_temperature_system ();
+
+ /**
+ * Initiate the assembly of one composition matrix and right hand side
* and build a preconditioner for the matrix.
*
* This function is implemented in
* <code>source/simulator/assembly.cc</code>.
*/
- void assemble_advection_system (const TemperatureOrComposition &temperature_or_composition);
+ void assemble_composition_system (const unsigned int c);
/**
- * Solve one block of the the temperature/composition linear system. Return the initial
- * nonlinear residual, i.e., if the linear system to be solved is $Ax=b$, then
- * return $\|Ax_0-b\|$ where $x_0$ is the initial guess for the solution variable
- * and is taken from the current_linearization_point member variable.
+ * Solve the temperature linear system. Return the initial nonlinear residual, i.e.,
+ * if the linear system to be solved is $Ax=b$, then return $\|Ax_0-b\|$ where $x_0$
+ * is the initial guess for the solution variable and is taken from the
+ * current_linearization_point_temperature member variable.
*
* This function is implemented in
* <code>source/simulator/solver.cc</code>.
*/
- double solve_advection (const TemperatureOrComposition &temperature_or_composition);
+ double solve_temperature ();
/**
- * Solve the Stokes linear system. Return the initial nonlinear residual,
- * i.e., if the linear system to be solved is $Ax=b$, then return $\|Ax_0-b\|$
- * where $x_0$ is the initial guess for the solution variable and is taken from
- * the current_linearization_point member variable.
+ * Solve one block of the the composition linear system. Return the initial nonlinear
+ * residual, i.e., if the linear system to be solved is $Ax=b$, then return
+ * $\|Ax_0-b\|$ where $x_0$ is the initial guess for the solution variable and is
+ * taken from the current_linearization_point member variable.
*
* This function is implemented in
* <code>source/simulator/solver.cc</code>.
*/
- double solve_stokes ();
+ double solve_composition (const unsigned int c);
/**
- * Evolve the compositional variables by one time step.
+ * Solve the Stokes linear system. Return the initial nonlinear residual,
+ * i.e., if the linear system to be solved is $Ax=b$, then return $\|Ax_0-b\|$
+ * where $x_0$ is the initial guess for the solution variable and is taken from
+ * the current_linearization_point member variable.
*
* This function is implemented in
* <code>source/simulator/solver.cc</code>.
*/
- void evolve_composition ();
+ double solve_stokes ();
/**
* This function is called at the end of every time step. It
@@ -578,15 +539,33 @@
*/
/**
* Set up the size and structure of the matrix used to store the
- * elements of the linear system.
+ * elements of the linear system for the velocity.
*
* This function is implemented in
* <code>source/simulator/core.cc</code>.
*/
- void setup_system_matrix (const std::vector<IndexSet> &system_partitioning);
+ void setup_system_matrix_velocity (const std::vector<IndexSet> &system_partitioning);
/**
* Set up the size and structure of the matrix used to store the
+ * elements of the linear system for the temperature.
+ *
+ * This function is implemented in
+ * <code>source/simulator/core.cc</code>.
+ */
+ void setup_system_matrix_temperature (const IndexSet &system_partitioning);
+
+ /**
+ * Set up the size and structure of the matrix used to store the
+ * elements of the linear system for the composition.
+ *
+ * This function is implemented in
+ * <code>source/simulator/core.cc</code>.
+ */
+ void setup_system_matrix_composition (const IndexSet &system_partitioning);
+
+ /**
+ * Set up the size and structure of the matrix used to store the
* elements of the matrix that is used to build the preconditioner for
* the system.
*
@@ -619,7 +598,9 @@
* <code>source/simulator/assembly.cc</code>.
*/
void
- local_assemble_stokes_preconditioner (const typename DoFHandler<dim>::active_cell_iterator &cell,
+ local_assemble_stokes_preconditioner (const SynchronousIterators<std_cxx1x::tuple<FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator> > > &cell,
internal::Assembly::Scratch::StokesPreconditioner<dim> &scratch,
internal::Assembly::CopyData::StokesPreconditioner<dim> &data);
@@ -641,7 +622,9 @@
* <code>source/simulator/assembly.cc</code>.
*/
void
- local_assemble_stokes_system (const typename DoFHandler<dim>::active_cell_iterator &cell,
+ local_assemble_stokes_system (const SynchronousIterators<std_cxx1x::tuple<FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator> > > &cell,
internal::Assembly::Scratch::StokesSystem<dim> &scratch,
internal::Assembly::CopyData::StokesSystem<dim> &data);
@@ -656,46 +639,62 @@
copy_local_to_global_stokes_system (const internal::Assembly::CopyData::StokesSystem<dim> &data);
/**
- * Compute the integrals for one advection matrix and right hand side
+ * Compute the integrals for the temperature matrix and right hand side
* on a single cell.
*
* This function is implemented in
* <code>source/simulator/assembly.cc</code>.
*/
void
- local_assemble_advection_system (const TemperatureOrComposition &temperature_or_composition,
- const std::pair<double,double> global_field_range,
- const double global_max_velocity,
- const double global_entropy_variation,
- const typename DoFHandler<dim>::active_cell_iterator &cell,
- internal::Assembly::Scratch::AdvectionSystem<dim> &scratch,
- internal::Assembly::CopyData::AdvectionSystem<dim> &data);
+ local_assemble_temperature_system (const std::pair<double,double> global_field_range,
+ const double global_max_velocity,
+ const double global_entropy_variation,
+ const SynchronousIterators<std_cxx1x::tuple<FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator> > > &cell,
+ internal::Assembly::Scratch::TemperatureSystem<dim> &scratch,
+ internal::Assembly::CopyData::TemperatureSystem<dim> &data);
/**
- * Compute the heating term for the advection system index.
- * Currently the heating term is 0
- * for compositional fields, but this can be changed in the
- * future to allow for interactions between compositional fields.
+ * Compute the integrals for one composition matrix and right hand side
+ * on a single cell.
+ *
+ * This function is implemented in
+ * <code>source/simulator/assembly.cc</code>.
+ */
+ void
+ local_assemble_composition_system (const unsigned int c,
+ const std::pair<double,double> global_field_range,
+ const double global_max_velocity,
+ const double global_entropy_variation,
+ const SynchronousIterators<std_cxx1x::tuple<FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator>,
+ FilteredIterator<typename DoFHandler<dim>::active_cell_iterator> > > &cell,
+ internal::Assembly::Scratch::CompositionSystem<dim> &scratch,
+ internal::Assembly::CopyData::CompositionSystem<dim> &data);
+
+
+ /**
+ * Copy the contribution to the temperature system
+ * from a single cell into the global matrix that stores these elements.
*
* This function is implemented in
* <code>source/simulator/assembly.cc</code>.
*/
- double compute_heating_term(const internal::Assembly::Scratch::AdvectionSystem<dim> &scratch,
- typename MaterialModel::Interface<dim>::MaterialModelInputs &material_model_inputs,
- typename MaterialModel::Interface<dim>::MaterialModelOutputs &material_model_outputs,
- const TemperatureOrComposition &temperature_or_composition,
- const unsigned int q) const;
+ void
+ copy_local_to_global_temperature_system (const internal::Assembly::CopyData::TemperatureSystem<dim> &data);
/**
- * Copy the contribution to the advection system
+ * Copy the contribution to the composition system
* from a single cell into the global matrix that stores these elements.
*
* This function is implemented in
* <code>source/simulator/assembly.cc</code>.
*/
void
- copy_local_to_global_advection_system (const internal::Assembly::CopyData::AdvectionSystem<dim> &data);
+ copy_local_to_global_composition_system (const unsigned int c,
+ const internal::Assembly::CopyData::CompositionSystem<dim> &data);
/**
* @}
@@ -768,70 +767,44 @@
* This function is implemented in
* <code>source/simulator/assembly.cc</code>.
*/
- double get_entropy_variation (const double average_value,
- const TemperatureOrComposition &temperature_or_composition) const;
+ double get_temperature_variation (const double average_value) const;
/**
- * Compute the minimal and maximal temperature througout the domain from a
- * solution vector extrapolated from the previous time steps. This is needed
- * to compute the artificial diffusion stabilization terms.
+ * Compute the variation (i.e., the difference between maximal and
+ * minimal value) of the entropy $(T-ar T)^2$ where $ar T$ is the
+ * average temperature throughout the domain given as argument to
+ * this function.
*
+ * This function is used in computing the artificial diffusion
+ * stabilization term.
+ *
* This function is implemented in
- * <code>source/simulator/helper_functions.cc</code>.
+ * <code>source/simulator/assembly.cc</code>.
*/
- std::pair<double,double>
- get_extrapolated_temperature_or_composition_range (const TemperatureOrComposition &temperature_or_composition) const;
+ double get_composition_variation (const double average_value,
+ const unsigned int c) const;
/**
- * Compute the size of the next time step from the mesh size and
- * the velocity on each cell. The computed time step has to satisfy
- * the CFL number chosen in the input parameter file on each cell
- * of the mesh. If specified in the parameter file, the time step
- * will be the minimum of the convection *and* conduction time
- * steps. Also returns whether the timestep is dominated by
- * convection (true) or conduction (false).
+ * Compute the minimal and maximal composition througout the domain from a
+ * solution vector extrapolated from the previous time steps. This is needed
+ * to compute the artificial diffusion stabilization terms.
*
* This function is implemented in
* <code>source/simulator/helper_functions.cc</code>.
*/
- std::pair<double,bool> compute_time_step () const;
+ std::pair<double,double> get_extrapolated_composition_range (const unsigned int n) const;
/**
- * Compute the artificial diffusion coefficient value on a cell
- * given the values and gradients of the solution passed as
- * arguments.
+ * Compute the minimal and maximal temperature througout the domain from a
+ * solution vector extrapolated from the previous time steps. This is needed
+ * to compute the artificial diffusion stabilization terms.
*
* This function is implemented in
- * <code>source/simulator/assembly.cc</code>.
+ * <code>source/simulator/helper_functions.cc</code>.
*/
- double
- compute_viscosity(internal::Assembly::Scratch::AdvectionSystem<dim> &scratch,
- const double global_u_infty,
- const double global_T_variation,
- const double average_temperature,
- const double global_entropy_variation,
- const double cell_diameter,
- const TemperatureOrComposition &temperature_or_composition) const;
+ std::pair<double,double> get_extrapolated_temperature_range () const;
/**
- * Compute the residual of one advection equation to be used
- * for the artificial diffusion coefficient value on a cell
- * given the values and gradients of the solution
- * passed as arguments.
- *
- * This function is implemented in
- * <code>source/simulator/assembly.cc</code>.
- */
- void
- compute_advection_system_residual(internal::Assembly::Scratch::AdvectionSystem<dim> &scratch,
- const double average_field,
- const TemperatureOrComposition &temperature_or_composition,
- double &max_residual,
- double &max_velocity,
- double &max_density,
- double &max_specific_heat) const;
-
- /**
* Extract the values of temperature, pressure, composition and
* optional strain rate for the current linearization point.
* These values are stored as input arguments for the material
@@ -848,29 +821,15 @@
* <code>source/simulator/assembly.cc</code>.
*/
void
- compute_material_model_input_values (const TrilinosWrappers::MPI::BlockVector &input_solution,
- const FEValues<dim,dim> &input_finite_element_values,
+ compute_material_model_input_values (const TrilinosWrappers::MPI::BlockVector &input_solution_velocity,
+ const TrilinosWrappers::MPI::Vector &input_solution_temperature,
+ const TrilinosWrappers::MPI::BlockVector &input_solution_composition,
+ const FEValues<dim,dim> &input_fe_values_velocity,
+ const FEValues<dim,dim> &input_fe_values_temperature,
+ const FEValues<dim,dim> &input_fe_values_composition,
typename MaterialModel::Interface<dim>::MaterialModelInputs &material_model_inputs) const;
/**
- * Return whether the Stokes matrix depends on the values of the
- * solution at the previous time step. This is the case is
- * the coefficients that appear in the matrix (i.e., the
- * viscosity and, in the case of a compressible model, the
- * density) depend on the solution.
- *
- * This function exists to ensure that the Stokes matrix is
- * rebuilt in time steps where it may have changed, while we
- * want to save the effort of rebuilding it whenever we don't
- * need to.
- *
- * This function is implemented in
- * <code>source/simulator/helper_functions.cc</code>.
- */
- bool
- stokes_matrix_depends_on_solution () const;
-
- /**
* Generate and output some statistics like timing information and memory consumption.
* Whether this function does anything or not is controlled through the
* variable aspect::output_parallel_statistics.
@@ -981,17 +940,22 @@
const MappingQ<dim> mapping;
- const FESystem<dim> finite_element;
+ const FESystem<dim> fe_velocity;
+ const FE_Q<dim> fe_temperature;
+ const FE_Q<dim> fe_composition;
- DoFHandler<dim> dof_handler;
+ DoFHandler<dim> dof_handler_velocity;
+ DoFHandler<dim> dof_handler_temperature;
+ DoFHandler<dim> dof_handler_composition;
- ConstraintMatrix constraints;
- ConstraintMatrix current_constraints;
+ ConstraintMatrix constraints_velocity;
+ ConstraintMatrix constraints_temperature;
+ ConstraintMatrix current_constraints_velocity;
+ std::vector<ConstraintMatrix> current_constraints_composition;
double pressure_adjustment;
double pressure_scaling;
double pressure_scaling_old;
- double viscosity_scaling;
/**
* @}
@@ -1002,15 +966,27 @@
* @name Variables that describe the linear systems and solution vectors
* @{
*/
- LinearAlgebra::BlockSparseMatrix system_matrix;
+ LinearAlgebra::BlockSparseMatrix system_matrix_velocity;
+ LinearAlgebra::SparseMatrix system_matrix_temperature;
+ LinearAlgebra::BlockSparseMatrix system_matrix_composition;
LinearAlgebra::BlockSparseMatrix system_preconditioner_matrix;
- LinearAlgebra::BlockVector solution;
- LinearAlgebra::BlockVector old_solution;
- LinearAlgebra::BlockVector old_old_solution;
- LinearAlgebra::BlockVector system_rhs;
+ LinearAlgebra::BlockVector solution_velocity;
+ LinearAlgebra::Vector solution_temperature;
+ LinearAlgebra::BlockVector solution_composition;
+ LinearAlgebra::BlockVector old_solution_velocity;
+ LinearAlgebra::Vector old_solution_temperature;
+ LinearAlgebra::BlockVector old_solution_composition;
+ LinearAlgebra::BlockVector old_old_solution_velocity;
+ LinearAlgebra::Vector old_old_solution_temperature;
+ LinearAlgebra::BlockVector old_old_solution_composition;
+ LinearAlgebra::BlockVector system_rhs_velocity;
+ LinearAlgebra::Vector system_rhs_temperature;
+ LinearAlgebra::BlockVector system_rhs_composition;
- LinearAlgebra::BlockVector current_linearization_point;
+ LinearAlgebra::BlockVector current_linearization_point_velocity;
+ LinearAlgebra::Vector current_linearization_point_temperature;
+ LinearAlgebra::BlockVector current_linearization_point_composition;
// only used if is_compressible()
LinearAlgebra::BlockVector pressure_shape_function_integrals;
Modified: trunk/aspire/include/aspect/simulator_access.h
===================================================================
--- trunk/aspire/include/aspect/simulator_access.h 2013-07-16 23:37:10 UTC (rev 1818)
+++ trunk/aspire/include/aspect/simulator_access.h 2013-07-24 20:14:34 UTC (rev 1819)
@@ -187,49 +187,87 @@
/** @name Accessing variables that identify the solution of the problem */
/** @{ */
-
/**
* Return a reference to the vector that has the current
- * solution of the entire system, i.e. the velocity and
- * pressure variables as well as the temperature. This vector
- * is associated with the DoFHandler object returned by
- * get_dof_handler().
+ * solution of the velocity system, i.e. the velocity and
+ * pressure variables. This vector is associated with the
+ * DoFHandler object returned by get_dof_handler_velocity ().
*
* @note In general the vector is a distributed vector; however, it
* contains ghost elements for all locally relevant degrees of freedom.
*/
const LinearAlgebra::BlockVector &
- get_solution () const;
+ get_solution_velocity () const;
/**
- * Return a reference to the vector that has the solution
- * of the entire system at the previous time step.
- * This vector is associated with the DoFHandler object returned by
- * get_stokes_dof_handler().
+ * Return a reference to the vector that has the current
+ * solution of the temperature system. This vector is
+ * associated with the DoFHandler object returned by
+ * get_dof_handler_temperature ().
*
* @note In general the vector is a distributed vector; however, it
* contains ghost elements for all locally relevant degrees of freedom.
*/
+ const LinearAlgebra::Vector &
+ get_solution_temperature () const;
+
+ /**
+ * Return a reference to the vector that has the current
+ * solution of the composition system. This vector is
+ * associated with the DoFHandler object returned by
+ * get_dof_handler_composition ().
+ *
+ * @note In general the vector is a distributed vector; however, it
+ * contains ghost elements for all locally relevant degrees of freedom.
+ */
const LinearAlgebra::BlockVector &
- get_old_solution () const;
+ get_solution_composition () const;
/**
* Return a reference to the DoFHandler that is used to discretize
- * the variables at the current time step.
+ * the velocity variables at the current time step.
*/
const DoFHandler<dim> &
- get_dof_handler () const;
+ get_dof_handler_velocity () const;
/**
+ * Return a reference to the DoFHandler that is used to discretize
+ * the temperature at the current time step.
+ */
+ const DoFHandler<dim> &
+ get_dof_handler_temperature () const;
+
+ /**
+ * Return a reference to the DoFHandler that is used to discretize
+ * the composition variables at the current time step.
+ */
+ const DoFHandler<dim> &
+ get_dof_handler_composition () const;
+
+ /**
* Return a reference to the finite element that the DoFHandler that
- * is used to discretize the variables at the current time step is built
- * on. This is the finite element for the entire, couple problem, i.e.,
- * it contains sub-elements for velocity, pressure, temperature
- * and all other variables in this problem (e.g., compositional variables,
- * if used in this simulation).
+ * is used to discretize the velocity variables at the current time step
+ * is built on. This is the finite element for the velocity problem, i.e.,
+ * it contains sub-elements for velocity and pressure.
*/
const FiniteElement<dim> &
- get_fe () const;
+ get_fe_velocity () const;
+
+ /**
+ * Return a reference to the finite element that the DoFHandler that
+ * is used to discretize the temperature at the current time step is built
+ * on.
+ */
+ const FiniteElement<dim> &
+ get_fe_temperature () const;
+
+ /**
+ * Return a reference to the finite element that the DoFHandler that
+ * is used to discretize the composition at the current time step is built
+ * on.
+ */
+ const FiniteElement<dim> &
+ get_fe_composition () const;
/** @} */
Modified: trunk/aspire/source/mesh_refinement/composition.cc
===================================================================
--- trunk/aspire/source/mesh_refinement/composition.cc 2013-07-16 23:37:10 UTC (rev 1818)
+++ trunk/aspire/source/mesh_refinement/composition.cc 2013-07-24 20:14:34 UTC (rev 1819)
@@ -42,16 +42,13 @@
{
Vector<float> this_indicator (indicators.size());
- KellyErrorEstimator<dim>::estimate (this->get_dof_handler(),
-//TODO: Replace the 2 by something reasonable, adjusted to the polynomial degree
- QGauss<dim-1>(2),
+ KellyErrorEstimator<dim>::estimate (this->get_mapping (),
+ this->get_dof_handler_composition (),
+ QGauss<dim - 1> (this->get_fe_composition ().degree + 1),
typename FunctionMap<dim>::type(),
- this->get_solution(),
- this_indicator,
- this->introspection().component_masks.compositional_fields[c],
- 0,
- 0,
- this->get_triangulation().locally_owned_subdomain());
+ this->get_solution_composition ().block (c),
+ this_indicator, ComponentMask (), 0, 0,
+ this->get_triangulation ().locally_owned_subdomain ());
indicators += this_indicator;
}
}
Modified: trunk/aspire/source/mesh_refinement/density.cc
===================================================================
--- trunk/aspire/source/mesh_refinement/density.cc 2013-07-16 23:37:10 UTC (rev 1818)
+++ trunk/aspire/source/mesh_refinement/density.cc 2013-07-24 20:14:34 UTC (rev 1819)
@@ -41,7 +41,7 @@
// then we can get away with simply interpolating it spatially
- // create a vector in which we set the temperature block to
+ // create a vector in which we set the temperature vector to
// be a finite element interpolation of the density.
// we do so by setting up a quadrature formula with the
// temperature unit support points, then looping over these
@@ -50,16 +50,23 @@
// (because quadrature points and temperature dofs are,
// by design of the quadrature formula, numbered in the
// same way)
- LinearAlgebra::BlockVector vec_distributed (this->introspection().index_sets.system_partitioning,
- this->get_mpi_communicator());
+ LinearAlgebra::Vector vec_distributed (this->introspection ().index_sets.temperature_partitioning,
+ this->get_mpi_communicator ());
- const Quadrature<dim> quadrature(this->get_fe().base_element(2).get_unit_support_points());
- std::vector<types::global_dof_index> local_dof_indices (this->get_fe().dofs_per_cell);
- FEValues<dim> fe_values (this->get_mapping(),
- this->get_fe(),
- quadrature,
More information about the CIG-COMMITS
mailing list