Discontinuous finite element for advection equations #3766
Unanswered
nanguaxiaofendui
asked this question in
Q&A
Replies: 1 comment
-
I've seen my mistake, thank you! |
Beta Was this translation helpful? Give feedback.
0 replies
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Uh oh!
There was an error while loading. Please reload this page.
Uh oh!
There was an error while loading. Please reload this page.
-
Discontinuous finite element method is used to solve the equation▽u(x,y)+u(x,y)=f(x,y),domain is [0,1]×[0,1],the exact solution is uex=exp(x+y),using uex to get the f(x,y) and boundary condition(Dirichlet) .
When I use the second order, the error is greater than the first order. Can you tell me where I made a mistake? Are there any examples? Thank you very much.
my code :
// C++ include files that we need
#include
#include
#include <math.h>
#include "libmesh/tecplot_io.h"
// Basic include files needed for the mesh functionality.
#include "libmesh/libmesh.h"
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/vtk_io.h"
#include "libmesh/exodusII_io.h"
#include "libmesh/linear_implicit_system.h"
#include "libmesh/equation_systems.h"
// Define the Finite Element object.
#include "libmesh/fe.h"
// Define Gauss quadrature rules.
#include "libmesh/quadrature_gauss.h"
#include "libmesh/quadrature_trap.h"
// Define useful datatypes for finite element
// matrix and vector components.
#include "libmesh/sparse_matrix.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/dense_vector.h"
#include "libmesh/elem.h"
#include "libmesh/enum_solver_package.h"
// Define the DofMap, which handles degree of freedom
// indexing.
#include "libmesh/dof_map.h"
#include "libmesh/exact_error_estimator.h"
#include "libmesh/exact_solution.h"
#include "libmesh/exodusII_io.h"
// Bring in everything from the libMesh namespace
using namespace libMesh;
void assemble(EquationSystems& es, const std::string& system_name);
Number Exact_solution(const Point & p,
const Parameters &, // EquationSystem parameters, not needed
const std::string &, // sys_name, not needed
const std::string &); // unk_name, not needed);
Gradient Exact_derivative(const Point & p,
const Parameters &, // EquationSystems parameters, not needed
const std::string &, // sys_name, not needed
const std::string &); // unk_name, not needed);
// Function prototype for the exact solution.
Real exact_solution(const Real x, const Real y, const Real z = 0.);
int main(int argc, char** argv)
{
// Initialize libraries, like in example 2.
LibMeshInit init(argc, argv);
ExactSolution exact_sol(equation_systems);
exact_sol.attach_exact_value(Exact_solution);
exact_sol.attach_exact_deriv(Exact_derivative);
// Compute the error.
exact_sol.compute_error("TEST", "u");
// Print out the error values
libMesh::out << "L2-Error is: "
<< exact_sol.l2_error("TEST", "u")
<< std::endl;
libMesh::out << "Initial H1-Error is: "
<< exact_sol.h1_error("TEST", "u")
<< std::endl;
equation_systems.get_system("TEST").solution;
#ifdef LIBMESH_HAVE_EXODUS_API
if (!FIRST) /* CONSTANT, MONOMIAL_VEC */
{
// Warn about trivial solution for CONSTANT approximations (Poisson must be at least C1)
libMesh::out << "\nWARNING: The elemental vector order is CONSTANT. The solution" << std::endl
<< "written out to 'out_constant.e' will be trivial." << std::endl;
#endif
// All done.
return 0;
}
void assemble(EquationSystems& es,
const std::string& libmesh_dbg_var(system_name))
{
libmesh_assert_equal_to(system_name, "TEST");
}
Number Exact_solution(const Point & p,
const Parameters &, // parameters, not needed
const std::string &, // sys_name, not needed
const std::string &) // unk_name, not needed
{
return exp(p(0)+p(1));
}
Gradient Exact_derivative(const Point & p,
const Parameters &, // parameters, not needed
const std::string &, // sys_name, not needed
const std::string &) // unk_name, not needed
{
// Gradient value to be returned.
Gradient gradu;
return gradu;
}
Real exact_solution(const Real x,
const Real y,
const Real z = 0.)
{
return exp(x + y);
}
Beta Was this translation helpful? Give feedback.
All reactions