From 7892e8cd06c1f938f49a8bf45b33544568e3d9f8 Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 24 Jan 2008 17:17:19 +0000 Subject: [PATCH] Add no_flux_04. git-svn-id: https://svn.dealii.org/trunk@15670 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/deal.II/no_flux_04.cc | 144 +++++++++++++++++++++++++++ tests/deal.II/no_flux_04/cmp/generic | 6 ++ 2 files changed, 150 insertions(+) create mode 100644 tests/deal.II/no_flux_04.cc create mode 100644 tests/deal.II/no_flux_04/cmp/generic diff --git a/tests/deal.II/no_flux_04.cc b/tests/deal.II/no_flux_04.cc new file mode 100644 index 0000000000..361a3d09d5 --- /dev/null +++ b/tests/deal.II/no_flux_04.cc @@ -0,0 +1,144 @@ +//---------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2007, 2008 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------------------------------------------------- + + +// check the creation of no-flux boundary conditions for a finite +// element that consists of only a single set of vector components +// (i.e. it has dim components) +// +// like no_flux_04 but apply the constraints to a vector field to see +// whether the result looks alright + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + + +// a function that shows something useful on the surface of a sphere +template +class RadialFunction : public Function +{ + public: + RadialFunction() : Function (dim) {} + + virtual void vector_value (const Point &p, + Vector &v) const + { + Assert (v.size() == dim, ExcInternalError()); + + switch (dim) + { + case 2: + v(0) = p[0] + p[1]; + v(1) = p[1] - p[0]; + break; + case 3: + v(0) = p[0] + p[1]; + v(1) = p[1] - p[0]; + v(2) = p[2]; + break; + default: + Assert (false, ExcNotImplemented()); + } + } +}; + + + +template +void test (const Triangulation& tr, + const FiniteElement& fe) +{ + DoFHandler dof(tr); + dof.distribute_dofs(fe); + + deallog << "FE=" << fe.get_name() + << std::endl; + + std::set boundary_ids; + boundary_ids.insert (0); + + ConstraintMatrix cm; + VectorTools::compute_no_normal_flux_constraints (dof, 0, boundary_ids, cm); + cm.close (); + + DoFHandler dh (tr); + dh.distribute_dofs (fe); + + Vector v(dh.n_dofs()); + VectorTools::interpolate (dh, RadialFunction(), v); + + cm.distribute (v); + + std::ofstream o("x.vtk"); + DataOut data_out; + data_out.attach_dof_handler (dh); + + std::vector + data_component_interpretation + (dim, DataComponentInterpretation::component_is_part_of_vector); + + data_out.add_data_vector (v, "x", + DataOut::type_dof_data, + data_component_interpretation); + data_out.build_patches (fe.degree); + + data_out.write_vtk (o); +} + + + +template +void test_hyper_sphere() +{ + Triangulation tr; + GridGenerator::hyper_ball(tr); + + static const HyperBallBoundary boundary; + tr.set_boundary (0, boundary); + + tr.refine_global(2); + + for (unsigned int degree=1; degree<6-dim; ++degree) + { + FESystem fe (FE_Q(degree), dim); + test(tr, fe); + } +} + + +int main() +{ + std::ofstream logfile ("no_flux_04/output"); + deallog << std::setprecision (2); + deallog << std::fixed; + deallog.attach(logfile); + deallog.depth_console (0); + deallog.threshold_double(1.e-12); + + test_hyper_sphere<2>(); + test_hyper_sphere<3>(); +} diff --git a/tests/deal.II/no_flux_04/cmp/generic b/tests/deal.II/no_flux_04/cmp/generic new file mode 100644 index 0000000000..22fa52f745 --- /dev/null +++ b/tests/deal.II/no_flux_04/cmp/generic @@ -0,0 +1,6 @@ + +DEAL::FE=FESystem<2>[FE_Q<2>(1)^2] +DEAL::FE=FESystem<2>[FE_Q<2>(2)^2] +DEAL::FE=FESystem<2>[FE_Q<2>(3)^2] +DEAL::FE=FESystem<3>[FE_Q<3>(1)^3] +DEAL::FE=FESystem<3>[FE_Q<3>(2)^3] -- 2.39.5