From 3313a43a9968d7110739a6b44a010c150d40a2db Mon Sep 17 00:00:00 2001 From: kanschat Date: Mon, 15 Aug 2011 12:28:26 +0000 Subject: [PATCH] implement output of face elements git-svn-id: https://svn.dealii.org/trunk@24073 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/integrators/patches.h | 62 +++++++++++++ .../deal.II/numerics/mesh_worker_output.h | 91 +++++++++++-------- 2 files changed, 113 insertions(+), 40 deletions(-) create mode 100644 deal.II/include/deal.II/integrators/patches.h diff --git a/deal.II/include/deal.II/integrators/patches.h b/deal.II/include/deal.II/integrators/patches.h new file mode 100644 index 0000000000..68d3c9aeb6 --- /dev/null +++ b/deal.II/include/deal.II/integrators/patches.h @@ -0,0 +1,62 @@ +//--------------------------------------------------------------------------- +// $Id: laplace.h 23983 2011-07-29 22:45:15Z kanschat $ +// +// Copyright (C) 2010, 2011 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. +// +//--------------------------------------------------------------------------- +#ifndef __deal2__integrators_patches_h +#define __deal2__integrators_patches_h + + +#include +#include +#include +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +namespace LocalIntegrators +{ +/** + * @brief Integrators writing patches with values in quadrature points + * + * @ingroup Integrators + * @author Guido Kanschat + * @date 2011 + */ + namespace Patches + { + template + inline + void + points_and_values( + Table<2, double>& result, + const FEValuesBase& fe, + const VectorSlice > >& input) + { + const unsigned int n_comp = fe.get_fe().n_components(); + AssertVectorVectorDimension(input, n_comp, fe.n_quadrature_points); + AssertDimension(result.n_rows(), fe.n_quadrature_points); + AssertDimension(result.n_cols(), n_comp+dim); + + for (unsigned int k=0;k #include #include +#include #include #include @@ -46,15 +47,22 @@ namespace MeshWorker /** * Initialize for writing * n data vectors. The + * number of points is the + * number of quadrature + * points in a single + * direction in a tensor + * product formula. It must + * match the number in the + * actual Quadrature used to + * create the patches. The * total number of data * vectors produced is n+dim * and the first dim will be * the space coordinates of * the points. */ - void initialize (const unsigned int dim, - const unsigned int n_vectors, - const unsigned int n_points); + void initialize (const unsigned int n_points, + const unsigned int n_vectors); /** * Set the stream #os to @@ -72,19 +80,14 @@ namespace MeshWorker * object used later for * assembling. * - * The second parameter is - * used to distinguish - * between the data used on - * cells and boundary faces - * on the one hand and - * interior faces on the - * other. Interior faces may - * require additional data - * being initialized. + * The info object refers to + * a cell if + * !face, or + * else to an interior or + * boundary face. */ template - void initialize_info(DoFInfo& info, - bool interior_face); + void initialize_info(DoFInfo& info, bool face); /** * Assemble the local values @@ -121,7 +124,6 @@ namespace MeshWorker */ void write_endl () const; - unsigned int patch_dimension; unsigned int n_vectors; unsigned int n_points; @@ -162,11 +164,9 @@ namespace MeshWorker inline void - GnuplotPatch::initialize (const unsigned int dim, - const unsigned int np, + GnuplotPatch::initialize (const unsigned int np, const unsigned int nv) { - patch_dimension = dim; n_vectors = nv; n_points = np; } @@ -181,10 +181,12 @@ namespace MeshWorker template inline void - GnuplotPatch::initialize_info(DoFInfo& info, bool) + GnuplotPatch::initialize_info(DoFInfo& info, bool face) { - AssertDimension(patch_dimension,dim); - info.initialize_quadrature(n_points, n_vectors); + if (face) + info.initialize_quadrature(Utilities::fixed_power(n_points), n_vectors+dim); + else + info.initialize_quadrature(Utilities::fixed_power(n_points), n_vectors+dim); } @@ -192,31 +194,40 @@ namespace MeshWorker inline void GnuplotPatch::assemble(const DoFInfo& info) { - const unsigned int n_q_points = info.n_quadrature_points(); + const unsigned int np = info.n_quadrature_points(); const unsigned int nv = info.n_quadrature_values(); - const unsigned int square_root = static_cast(std::pow(n_q_points, 1./dim)+.5); - for (unsigned int k1=0; k1(std::pow(np, 1./patch_dim)+.5); + + // If patches are 1D, end the + // patch after a row, else and + // it after a square + const unsigned int row_length2 = (patch_dim==1) ? row_length : (row_length*row_length); + + for (unsigned int k=0; k inline void - GnuplotPatch::assemble(const DoFInfo& /*info1*/, const DoFInfo& /*info2*/) + GnuplotPatch::assemble(const DoFInfo& info1, const DoFInfo& info2) { - Assert(false, ExcNotImplemented()); + assemble(info1); + assemble(info2); } } } -- 2.39.5