diff --git a/doc/modules/changes/20170609_gassmoeller b/doc/modules/changes/20170609_gassmoeller new file mode 100644 index 00000000000..44e74f158cb --- /dev/null +++ b/doc/modules/changes/20170609_gassmoeller @@ -0,0 +1,5 @@ +New: A new initial temperature condition was added that adds random +noise of user prescribed magnitude to the initial temperature field. +The same type of plugin was added for initial compositions. +
+(Rene Gassmoeller, 2017/06/19) diff --git a/include/aspect/initial_composition/random_perturbation.h b/include/aspect/initial_composition/random_perturbation.h new file mode 100644 index 00000000000..7a0c113e8aa --- /dev/null +++ b/include/aspect/initial_composition/random_perturbation.h @@ -0,0 +1,96 @@ +/* + Copyright (C) 2017 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 + . +*/ + + +#ifndef _aspect_initial_composition_random_perturbation_h +#define _aspect_initial_composition_random_perturbation_h + +#include +#include + +DEAL_II_DISABLE_EXTRA_DIAGNOSTICS +#include +DEAL_II_ENABLE_EXTRA_DIAGNOSTICS + +namespace aspect +{ + namespace InitialComposition + { + using namespace dealii; + + /** + * A class that describes a perturbation for the initial composition fields + * for any geometry model or dimension. The perturbation follows + * a random noise of a prescribed magnitude. + * + * @ingroup InitialCompositions + */ + template + class RandomPerturbation : public Interface, public aspect::SimulatorAccess + + { + public: + void initialize () override; + /** + * Return the initial composition as a function of position. + */ + + double initial_composition (const Point &position, const unsigned int composition_index) const override; + + /** + * 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. + */ + void + parse_parameters (ParameterHandler &prm) override; + + + private: + /** + * The maximal magnitude of the random noise. + */ + double magnitude; + + /** + * Whether to use a random seed for the random + * number generator. This parameter controls whether + * this plugin generates different or identical + * perturbations for subsequent model runs of + * the same setup. + */ + bool use_random_seed; + + /** + * A random number generator used by this class to get + * random composition perturbations. Should return the same + * random numbers every time since it is seeded with one. + */ + mutable std::mt19937 random_number_generator; + }; + } +} + +#endif diff --git a/include/aspect/initial_temperature/random_perturbation.h b/include/aspect/initial_temperature/random_perturbation.h new file mode 100644 index 00000000000..d1fdf0b1cdd --- /dev/null +++ b/include/aspect/initial_temperature/random_perturbation.h @@ -0,0 +1,95 @@ +/* + Copyright (C) 2017 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 + . +*/ + + +#ifndef _aspect_initial_temperature_random_perturbation_h +#define _aspect_initial_temperature_random_perturbation_h + +#include +#include + +DEAL_II_DISABLE_EXTRA_DIAGNOSTICS +#include +DEAL_II_ENABLE_EXTRA_DIAGNOSTICS + +namespace aspect +{ + namespace InitialTemperature + { + using namespace dealii; + + /** + * A class that describes a perturbation for the initial temperature field + * for any geometry model or dimension. The perturbation follows + * a random noise of a prescribed magnitude. + * + * @ingroup InitialTemperatures + */ + template + class RandomPerturbation : public Interface, public aspect::SimulatorAccess + { + public: + void initialize () override; + /** + * Return the initial temperature as a function of position. + */ + + double initial_temperature (const Point &position) const override; + + /** + * 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. + */ + void + parse_parameters (ParameterHandler &prm) override; + + private: + /** + * The maximal magnitude of the random noise. + */ + double magnitude; + + /** + * Whether to use a random seed for the random + * number generator. This parameter controls whether + * this plugin generates different or identical + * perturbations for subsequent model runs of + * the same setup. + */ + bool use_random_seed; + + /** + * A random number generator used by this class to get + * random temperature perturbations. Should return the same + * random numbers every time since it is seeded with one. + */ + mutable std::mt19937 random_number_generator; + + }; + } +} + +#endif diff --git a/source/initial_composition/random_perturbation.cc b/source/initial_composition/random_perturbation.cc new file mode 100644 index 00000000000..425dff106df --- /dev/null +++ b/source/initial_composition/random_perturbation.cc @@ -0,0 +1,109 @@ +/* + Copyright (C) 2017 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 + . +*/ + + +#include + +// #include +#include +#include + +namespace aspect +{ + namespace InitialComposition + { + template + void + RandomPerturbation:: + initialize() + { + const unsigned int my_rank = Utilities::MPI::this_mpi_process(this->get_mpi_communicator()); + if (use_random_seed) + random_number_generator.seed(time(NULL)+my_rank); + else + random_number_generator.seed(9088+my_rank); + } + + template + double + RandomPerturbation:: + initial_composition (const Point &position, const unsigned int composition_index) const + { + // Uniform distribution on the interval [-magnitude,magnitude). This + // will be used to generate the random temperature perturbation. + std::uniform_real_distribution uniform_distribution(-magnitude,magnitude); + return uniform_distribution(random_number_generator); + } + + template + void + RandomPerturbation::declare_parameters (ParameterHandler &prm) + { + prm.enter_subsection ("Initial composition model"); + { + prm.enter_subsection("Random perturbation"); + { + prm.declare_entry ("Magnitude", "1.0", + Patterns::Double (0), + "The magnitude of the random perturbation."); + prm.declare_entry ("Use random seed", "false", + Patterns::Bool(), + "Whether to use a random seed for the random " + "number generator. This parameter controls whether " + "this plugin generates different or identical " + "perturbations for subsequent model runs of" + "the same setup."); + } + prm.leave_subsection (); + } + prm.leave_subsection (); + } + + + template + void + RandomPerturbation::parse_parameters (ParameterHandler &prm) + { + prm.enter_subsection ("Initial composition model"); + { + prm.enter_subsection("Random perturbation"); + { + magnitude = prm.get_double ("Magnitude"); + use_random_seed = prm.get_bool("Use random seed"); + } + prm.leave_subsection (); + } + prm.leave_subsection (); + } + } +} + +// explicit instantiations +namespace aspect +{ + namespace InitialComposition + { + ASPECT_REGISTER_INITIAL_COMPOSITION_MODEL(RandomPerturbation, + "random perturbation", + "An initial composition anomaly that perturbs the composition " + "following a random noise with uniform distribution and user " + "specified magnitude.") + } +} diff --git a/source/initial_temperature/random_perturbation.cc b/source/initial_temperature/random_perturbation.cc new file mode 100644 index 00000000000..24f186af38f --- /dev/null +++ b/source/initial_temperature/random_perturbation.cc @@ -0,0 +1,108 @@ +/* + Copyright (C) 2017 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 + . +*/ + + +#include + +// #include +#include +#include + +namespace aspect +{ + namespace InitialTemperature + { + template + void + RandomPerturbation:: + initialize() + { + const unsigned int my_rank = Utilities::MPI::this_mpi_process(this->get_mpi_communicator()); + if (use_random_seed) + random_number_generator.seed(time(NULL)+my_rank); + else + random_number_generator.seed(9088+my_rank); + } + + template + double + RandomPerturbation:: + initial_temperature (const Point &position) const + { + // Uniform distribution on the interval [-magnitude,magnitude). This + // will be used to generate the random temperature perturbation. + std::uniform_real_distribution uniform_distribution(-magnitude,magnitude); + return uniform_distribution(random_number_generator); + } + + template + void + RandomPerturbation::declare_parameters (ParameterHandler &prm) + { + prm.enter_subsection ("Initial temperature model"); + { + prm.enter_subsection("Random perturbation"); + { + prm.declare_entry ("Magnitude", "1.0", + Patterns::Double (0), + "The magnitude of the random perturbation."); + prm.declare_entry ("Use random seed", "false", + Patterns::Bool(), + "Whether to use a random seed for the random " + "number generator. This parameter controls whether " + "this plugin generates different or identical " + "perturbations for different model runs."); + } + prm.leave_subsection (); + } + prm.leave_subsection (); + } + + + template + void + RandomPerturbation::parse_parameters (ParameterHandler &prm) + { + prm.enter_subsection ("Initial temperature model"); + { + prm.enter_subsection("Random perturbation"); + { + magnitude = prm.get_double ("Magnitude"); + use_random_seed = prm.get_bool("Use random seed"); + } + prm.leave_subsection (); + } + prm.leave_subsection (); + } + } +} + +// explicit instantiations +namespace aspect +{ + namespace InitialTemperature + { + ASPECT_REGISTER_INITIAL_TEMPERATURE_MODEL(RandomPerturbation, + "random perturbation", + "An initial temperature anomaly that perturbs the temperature " + "following a random noise with uniform distribution and user " + "specified magnitude.") + } +} diff --git a/tests/random_composition_perturbation.prm b/tests/random_composition_perturbation.prm new file mode 100644 index 00000000000..bd82910303d --- /dev/null +++ b/tests/random_composition_perturbation.prm @@ -0,0 +1,55 @@ +# Test to check the status of the random initial composition perturbation +# plugin. + +subsection Boundary velocity model + set Tangential velocity boundary indicators = left, right, bottom, top +end + +subsection Mesh refinement + set Initial global refinement = 3 +end + +subsection Postprocess + set List of postprocessors = composition statistics +end + +subsection Compositional fields + set Number of fields = 2 +end + +subsection Material model + set Model name = simpler +end + +set End time = 0 +subsection Geometry model + set Model name = box +end + +subsection Gravity model + set Model name = vertical +end + +subsection Initial temperature model + set List of model names = function + subsection Function + set Function expression = 0.5 + end + +end + +subsection Initial composition model + set List of model names = function, random perturbation + subsection Random perturbation + set Magnitude = 0.2 + end + + subsection Function + set Function expression = 0.5;0.0 + end + +end + +subsection Boundary temperature model + set List of model names = box +end diff --git a/tests/random_composition_perturbation/screen-output b/tests/random_composition_perturbation/screen-output new file mode 100644 index 00000000000..1a746fdca3b --- /dev/null +++ b/tests/random_composition_perturbation/screen-output @@ -0,0 +1,21 @@ + +Number of active cells: 64 (on 4 levels) +Number of degrees of freedom: 1,526 (578+81+289+289+289) + +*** Timestep 0: t=0 years + Solving temperature system... 0 iterations. + Solving C_1 system ... 0 iterations. + Solving C_2 system ... 0 iterations. + Rebuilding Stokes preconditioner... + Solving Stokes system... 2+0 iterations. + + Postprocessing: + Compositions min/max/mass: 0.3002/0.6997/0.5069 // -0.199/0.1999/-0.007404 + +Termination requested by criterion: end time + + ++---------------------------------------------+------------+------------+ ++---------------------------------+-----------+------------+------------+ ++---------------------------------+-----------+------------+------------+ + diff --git a/tests/random_composition_perturbation/statistics b/tests/random_composition_perturbation/statistics new file mode 100644 index 00000000000..9743ff8bf2b --- /dev/null +++ b/tests/random_composition_perturbation/statistics @@ -0,0 +1,20 @@ +# 1: Time step number +# 2: Time (years) +# 3: Time step size (years) +# 4: Number of mesh cells +# 5: Number of Stokes degrees of freedom +# 6: Number of temperature degrees of freedom +# 7: Number of degrees of freedom for all compositions +# 8: Iterations for temperature solver +# 9: Iterations for composition solver 1 +# 10: Iterations for composition solver 2 +# 11: Iterations for Stokes solver +# 12: Velocity iterations in Stokes preconditioner +# 13: Schur complement iterations in Stokes preconditioner +# 14: Minimal value for composition C_1 +# 15: Maximal value for composition C_1 +# 16: Global mass for composition C_1 +# 17: Minimal value for composition C_2 +# 18: Maximal value for composition C_2 +# 19: Global mass for composition C_2 +0 0.000000000000e+00 0.000000000000e+00 64 659 289 578 0 0 0 1 3 2 3.00222459e-01 6.99744700e-01 5.06886136e-01 -1.99012748e-01 1.99867409e-01 -7.40357014e-03 diff --git a/tests/random_temperature_perturbation.prm b/tests/random_temperature_perturbation.prm new file mode 100644 index 00000000000..9d8eeded384 --- /dev/null +++ b/tests/random_temperature_perturbation.prm @@ -0,0 +1,43 @@ +# Test to check the status of the random initial temperature perturbation +# plugin. + +subsection Boundary velocity model + set Tangential velocity boundary indicators = left, right, bottom, top +end + +subsection Mesh refinement + set Initial global refinement = 3 +end + +subsection Postprocess + set List of postprocessors = temperature statistics +end + +subsection Material model + set Model name = simpler +end + +set End time = 0 +subsection Geometry model + set Model name = box +end + +subsection Gravity model + set Model name = vertical +end + +subsection Initial temperature model + set List of model names = function, random perturbation + subsection Random perturbation + set Magnitude = 0.2 + end + + subsection Function + set Function expression = 0.5 + end + +end + +subsection Boundary temperature model + set List of model names = box +end diff --git a/tests/random_temperature_perturbation/screen-output b/tests/random_temperature_perturbation/screen-output new file mode 100644 index 00000000000..944cd96a0d8 --- /dev/null +++ b/tests/random_temperature_perturbation/screen-output @@ -0,0 +1,19 @@ + +Number of active cells: 64 (on 4 levels) +Number of degrees of freedom: 948 (578+81+289) + +*** Timestep 0: t=0 years + Solving temperature system... 0 iterations. + Rebuilding Stokes preconditioner... + Solving Stokes system... 14+0 iterations. + + Postprocessing: + Temperature min/avg/max: 0.302 K, 0.4957 K, 0.6998 K + +Termination requested by criterion: end time + + ++---------------------------------------------+------------+------------+ ++---------------------------------+-----------+------------+------------+ ++---------------------------------+-----------+------------+------------+ + diff --git a/tests/random_temperature_perturbation/statistics b/tests/random_temperature_perturbation/statistics new file mode 100644 index 00000000000..b5bfb68e845 --- /dev/null +++ b/tests/random_temperature_perturbation/statistics @@ -0,0 +1,14 @@ +# 1: Time step number +# 2: Time (years) +# 3: Time step size (years) +# 4: Number of mesh cells +# 5: Number of Stokes degrees of freedom +# 6: Number of temperature degrees of freedom +# 7: Iterations for temperature solver +# 8: Iterations for Stokes solver +# 9: Velocity iterations in Stokes preconditioner +# 10: Schur complement iterations in Stokes preconditioner +# 11: Minimal temperature (K) +# 12: Average temperature (K) +# 13: Maximal temperature (K) +0 0.000000000000e+00 0.000000000000e+00 64 659 289 0 13 15 14 3.02047032e-01 4.95676273e-01 6.99758663e-01