queso-0.57.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
LibMeshFunction.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // QUESO - a library to support the Quantification of Uncertainty
5 // for Estimation, Simulation and Optimization
6 //
7 // Copyright (C) 2008-2017 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 #include <queso/Defines.h>
26 
27 #ifdef QUESO_HAVE_LIBMESH_SLEPC
28 
29 #include <iostream>
30 #include <queso/LibMeshFunction.h>
31 #include <queso/FunctionOperatorBuilder.h>
32 #include <libmesh/libmesh.h>
33 #include <libmesh/mesh.h>
34 #include <libmesh/mesh_generation.h>
35 #include <libmesh/equation_systems.h>
36 #include <libmesh/explicit_system.h>
37 #include <libmesh/system_norm.h>
38 #include <libmesh/exodusII_io.h>
39 
40 // Define the Finite Element object.
41 #include <libmesh/fe.h>
42 
43 // Define Gauss quadrature rules.
44 #include <libmesh/quadrature_gauss.h>
45 
46 // Define useful datatypes for finite element
47 // matrix and vector components.
48 #include <libmesh/sparse_matrix.h>
49 #include <libmesh/numeric_vector.h>
50 #include <libmesh/dense_matrix.h>
51 #include <libmesh/dense_vector.h>
52 #include <libmesh/elem.h>
53 
54 // Define the DofMap, which handles degree of freedom
55 // indexing.
56 #include <libmesh/dof_map.h>
57 #include <libmesh/utility.h>
58 #include <libmesh/string_to_enum.h>
59 #include <libmesh/enum_order.h>
60 #include <libmesh/enum_fe_family.h>
61 
62 namespace QUESO {
63 
65  const FunctionOperatorBuilder & builder, libMesh::MeshBase & m)
66  : FunctionBase(),
67  builder(builder),
68  equation_systems(new libMesh::EquationSystems(m))
69 {
70  this->equation_systems->add_system<libMesh::ExplicitSystem>("Function");
71  this->equation_systems->get_system("Function").add_variable("u",
72  libMesh::Utility::string_to_enum<libMeshEnums::Order>(this->builder.order),
73  libMesh::Utility::string_to_enum<libMeshEnums::FEFamily>(this->builder.family));
74  this->equation_systems->init();
75 }
76 
78 {
79 }
80 
82 {
83  // Print information about the mesh to the screen.
84  this->equation_systems->get_mesh().print_info(std::cerr);
85 }
86 
87 void LibMeshFunction::save_function(const std::string & filename, double time) const
88 {
89  // The "1" is hardcoded for the number of time steps because the ExodusII
90  // manual states that it should be the number of timesteps within the file.
91  // Here, we are explicitly only doing one timestep per file.
92  libMesh::ExodusII_IO(this->equation_systems->get_mesh()).write_timestep(
93  filename, *this->equation_systems, 1, time);
94 }
95 
96 void LibMeshFunction::add(double scale, const FunctionBase & rhs) {
97  // We know we're dealing with a derived class type, so cast
98  const LibMeshFunction & rhs_derived = dynamic_cast<
99  const LibMeshFunction &>(rhs);
100 
101  this->equation_systems->get_system<libMesh::ExplicitSystem>("Function").solution->add(
102  scale, *(rhs_derived.equation_systems->get_system<libMesh::ExplicitSystem>(
103  "Function").solution));
104 }
105 
107  const FunctionBase & f2)
108 {
109  const LibMeshFunction & f1_derived = static_cast<
110  const LibMeshFunction &>(f1);
111  const LibMeshFunction & f2_derived = static_cast<
112  const LibMeshFunction &>(f2);
113 
114  this->equation_systems->get_system<libMesh::ExplicitSystem>("Function").solution->pointwise_mult(
115  *(f1_derived.equation_systems->get_system<libMesh::ExplicitSystem>("Function").solution),
116  *(f2_derived.equation_systems->get_system<libMesh::ExplicitSystem>("Function").solution));
117 }
118 
119 void LibMeshFunction::scale(double scale) {
120  this->equation_systems->get_system<libMesh::ExplicitSystem>(
121  "Function").solution->scale(scale);
122 }
123 
125  this->equation_systems->get_system<libMesh::ExplicitSystem>(
126  "Function").solution->zero();
127 }
128 
129 double LibMeshFunction::L2_norm() const {
130  libMesh::ExplicitSystem & system =
131  this->equation_systems->get_system<libMesh::ExplicitSystem>("Function");
132 
133  double norm = system.calculate_norm(*system.solution,
134  libMesh::SystemNorm(libMeshEnums::L2));
135  return norm;
136 }
137 
139 {
140  LibMeshFunction * clone = new LibMeshFunction(this->builder,
141  this->equation_systems->get_mesh());
142  clone->equation_systems->get_system<libMesh::ExplicitSystem>(
143  "Function").solution->zero();
144 
146  return ptr;
147 }
148 
151  return this->equation_systems;
152 }
153 
154 } // End namespace QUESO
155 
156 #endif // QUESO_HAVE_LIBMESH_SLEPC
std::shared_ptr< T > Type
Definition: SharedPtr.h:69
std::string order
String to store the polynomial order to use. Default is &quot;FIRST&quot;.
virtual double L2_norm() const
Return the L2-norm of this.
SharedPtr< libMesh::EquationSystems >::Type equation_systems
const FunctionOperatorBuilder & builder
virtual SharedPtr< FunctionBase >::Type zero_clone() const
Create a zero function copy of this and return pointer to it.
virtual SharedPtr< libMesh::EquationSystems >::Type get_equation_systems() const
Return the internal libmesh equation systems object.
void print_info() const
Will print mesh-related libMesh foo to std::cerr.
~LibMeshFunction()
Destructor.
virtual void zero()
Set this to zero everywhere.
virtual void pointwise_mult(const FunctionBase &f1, const FunctionBase &f2)
Pointwise multiply f1 and f2 and store the result in *this.
virtual void save_function(const std::string &filename, double time) const
Save the current function to an Exodus file called filename. time is the time to attach to the functi...
Function objects using libMesh for the backend.
virtual void scale(double scale)
Execute this *= scale.
virtual void add(double scale, const FunctionBase &rhs)
Execute this += scale * rhs.
LibMeshFunction(const FunctionOperatorBuilder &builder, libMesh::MeshBase &m)
Construct a function with a user-provided mesh m and builder.
Abstract base class for function objects.
Definition: FunctionBase.h:47

Generated on Tue Jun 5 2018 19:48:55 for queso-0.57.1 by  doxygen 1.8.5