Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
derivative_approximation.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2000 - 2021 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_derivative_approximation_h
17#define dealii_derivative_approximation_h
18
19#include <deal.II/base/config.h>
20
23
25
27#include <deal.II/fe/mapping.h>
28
30
31#include <deal.II/lac/vector.h>
32#ifdef _MSC_VER
34#endif
35#include <utility>
36
38
162{
178 template <int dim, class InputVector, int spacedim>
179 void
181 const DoFHandler<dim, spacedim> &dof,
182 const InputVector & solution,
184 const unsigned int component = 0);
185
189 template <int dim, class InputVector, int spacedim>
190 void
192 const InputVector & solution,
194 const unsigned int component = 0);
195
213 template <int dim, class InputVector, int spacedim>
214 void
216 const DoFHandler<dim, spacedim> &dof,
217 const InputVector & solution,
219 const unsigned int component = 0);
220
224 template <int dim, class InputVector, int spacedim>
225 void
227 const InputVector & solution,
229 const unsigned int component = 0);
230
244 template <int dim, int spacedim, class InputVector, int order>
245 void
247 const Mapping<dim, spacedim> & mapping,
248 const DoFHandler<dim, spacedim> &dof,
249 const InputVector & solution,
250#ifndef _MSC_VER
252#else
254 &cell,
255#endif
256 Tensor<order, dim> &derivative,
257 const unsigned int component = 0);
258
262 template <int dim, int spacedim, class InputVector, int order>
263 void
265 const DoFHandler<dim, spacedim> &dof,
266 const InputVector & solution,
267#ifndef _MSC_VER
269#else
271 &cell,
272#endif
273 Tensor<order, dim> &derivative,
274 const unsigned int component = 0);
275
279 template <int dim, int order>
280 double
281 derivative_norm(const Tensor<order, dim> &derivative);
282
287 int,
288 int,
289 << "The output vector needs to have a size equal "
290 "to the number of active cells of your triangulation "
291 "but has length "
292 << arg1 << "There are " << arg2
293 << " active cells in your triangulation.");
298 "While computing a finite difference approximation to "
299 "derivatives, the algorithm encountered a cell on which "
300 "the number of linearly "
301 "independent directions that span the matrix Y (discussed "
302 "in the documentation of the DerivativeApproximation "
303 "class) is not equal to dim. The matrix Y then is "
304 "rank deficient and can not be inverted. A common reason "
305 "why this might be happening is if a cell has neither "
306 "left/right (or up/down, or front/back) neighbors, for "
307 "example because the mesh is too coarse.");
308} // namespace DerivativeApproximation
309
310
311
313
314#endif
Abstract base class for mapping classes.
Definition mapping.h:317
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcInsufficientDirections()
static ::ExceptionBase & ExcVectorLengthVsNActiveCells(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
typename ActiveSelector::active_cell_iterator active_cell_iterator
void approximate_gradient(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InputVector &solution, Vector< float > &derivative_norm, const unsigned int component=0)
double derivative_norm(const Tensor< order, dim > &derivative)
void approximate_derivative_tensor(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InputVector &solution, const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell, Tensor< order, dim > &derivative, const unsigned int component=0)
void approximate_second_derivative(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InputVector &solution, Vector< float > &derivative_norm, const unsigned int component=0)