/usr/include/deal.II/numerics/error_estimator.h is in libdeal.ii-dev 6.3.1-1.1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 | //---------------------------------------------------------------------------
// $Id: error_estimator.h 18849 2009-05-15 17:59:46Z bangerth $
// Version: $Name$
//
// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 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.
//
//---------------------------------------------------------------------------
#ifndef __deal2__error_estimator_h
#define __deal2__error_estimator_h
#include <base/config.h>
#include <base/exceptions.h>
#include <base/function.h>
#include <base/multithread_info.h>
#include <dofs/function_map.h>
#include <map>
DEAL_II_NAMESPACE_OPEN
template <int, int> class DoFHandler;
template <int, int> class Mapping;
template <int> class Quadrature;
namespace hp
{
template <int, int> class DoFHandler;
template <int> class QCollection;
}
/**
* Implementation of the error estimator by Kelly, Gago, Zienkiewicz
* and Babuska. This error estimator tries to approximate the error
* per cell by integration of the jump of the gradient of the
* solution along the faces of each cell. It can be understood as a
* gradient recovery estimator; see the survey of Ainsworth for a
* complete discussion.
*
* It seem as if this error estimator should only be valid for linear trial
* spaces, and there are indications that for higher order trial spaces the
* integrals computed here show superconvergence properties, i.e. they tend
* to zero faster than the error itself, thus ruling out the values as error
* indicators.
*
* The error estimator really only estimates the error for the generalized
* Poisson equation $-\nabla\cdot a(x) \nabla u = f$ with either Dirichlet
* boundary conditions or generalized Neumann boundary conditions involving
* the conormal derivative $a\frac{du}{dn} = g$.
*
* The error estimator returns a vector of estimated errors per cell which
* can be used to feed the <tt>Triangulation<dim>::refine_*</tt> functions. This
* vector contains elements of data type @p float, rather than @p double,
* since accuracy is not so important here, and since this can save rather
* a lot of memory, when using many cells.
*
*
* <h3>Implementation</h3>
*
* In principle, the implementation of the error estimation is simple: let
* \f[
* \eta_K^2 = \frac h{24} \int_{\partial K} \left[a \frac{\partial u_h}{\partial n}\right]^2 do
* \f]
* be the error estimator for cell $K$. $[\cdot]$ denotes the jump of the
* argument at the face. In the paper of Ainsworth, $h$ is divided by $24$,
* but this factor is a bit esoteric, stemming from interpolation estimates
* and stability constants which may hold for the Poisson problem, but may
* not hold for more general situations. In the implementation, this factor
* is considered, but may lead to wrong results. You may scale the vector
* appropriately afterwards.
*
* To perform the integration, use is made of the FEFaceValues and
* FESubfaceValues classes. The integration is performed by looping
* over all cells and integrating over faces that are not yet treated.
* This way we avoid integration on faces twice, once for each time we
* visit one of the adjacent cells. In a second loop over all cells, we
* sum up the contributions of the faces (which are the integrated
* square of the jumps) of each cell and take the square root.
*
* The integration is done using a quadrature formula on the face.
* For linear trial functions (FEQ1), the QGauss2 or even the
* QMidpoint rule will suffice. For higher order elements, it is
* necessary to utilize higher order quadrature formulae as well.
*
* We store the contribution of each face in a @p map, as provided by the
* C++ standard library, with the iterator pointing to that face being the
* key into the map. In fact, we do not store the indicator per face, but
* only the integral listed above. When looping the second time over all
* cells, we have to sum up the contributions of the faces, multiply them
* with $\frac h{24}$ and take the square root. By doing the multiplication
* with $h$ in the second loop, we avoid problems to decide with which $h$
* to multiply, that of the cell on the one or that of the cell on the other
* side of the face.
*
* $h$ is taken to be the greatest length of the diagonals of the cell. For
* more or less uniform cells without deformed angles, this coincides with
* the diameter of the cell.
*
*
* <h3>Vector-valued functions</h3>
*
* If the finite element field for which the error is to be estimated
* is vector-valued, i.e. the finite element has more than one
* component, then all the above can be applied to all or only some
* components at the same time. The main function of this class takes
* a list of flags denoting the components for which components the
* error estimator is to be applied; by default, it is a list of only
* @p trues, meaning that all variables shall be treated.
*
* In case the different components of a field have different
* physical meaning (for example velocity and pressure in the Stokes
* equations), it would be senseless to use the same coefficient for
* all components. In that case, you can pass a function with as many
* components as there are components in the finite element field and
* each component of the error estimator will then be weighted by the
* respective component in this coefficient function. In the other
* case, when all components have the same meaning (for example the
* displacements in Lame's equations of elasticity), you can specify
* a scalar coefficient which will then be used for all components.
*
*
* <h3>%Boundary values</h3>
*
* If the face is at the boundary, i.e. there is no neighboring cell to which
* the jump in the gradiend could be computed, there are two possibilities:
* <ul>
* <li> The face belongs to a Dirichlet boundary. Then the face is not
* considered, which can be justified looking at a dual problem technique and
* should hold exactly if the boundary can be approximated exactly by the
* finite element used (i.e. it is a linear boundary for linear finite elements,
* quadratic for isoparametric quadratic elements, etc). For boundaries which
* can not be exactly approximated, one should consider the difference
* $z-z_h$ on the face, $z$ being a dual problem's solution which is zero at
* the true boundary and $z_h$ being an approximation, which in most cases
* will be zero on the numerical boundary. Since on the numerical boundary
* $z$ will not be zero in general, we would get another term here, but this
* one is neglected for practical reasons, in the hope that the error made
* here will tend to zero faster than the energy error we wish to estimate.
*
* Though no integration is necessary, in the list of face contributions we
* store a zero for this face, which makes summing up the contributions of
* the different faces to the cells easier.
*
* <li> The face belongs to a Neumann boundary. In this case, the
* contribution of the face $F\in\partial K$ looks like
* \f[ \int_F \left|g-a\frac{\partial u_h}{\partial n}\right|^2 ds \f]
* where $g$ is the Neumann boundary function. If the finite element is
* vector-valued, then obviously the function denoting the Neumann boundary
* conditions needs to be vector-valued as well.
*
* <li> No other boundary conditions are considered.
* </ul>
* The object describing the boundary conditions is obtained from the
* triangulation.
*
* Thanks go to Franz-Theo Suttmeier for clarifications about boundary
* conditions.
*
*
* <h3>Handling of hanging nodes</h3>
*
* The integration along faces with hanging nodes is quite tricky, since one
* of the elements has to be shifted one level up or down. See the
* documentation for the FESubFaceValues class for more information about
* technical issues regarding this topic.
*
* In praxi, since we integrate over each face only once, we do this when we
* are on the coarser one of the two cells adjacent to a subface (a subface
* is defined to be the child of a face; seen from the coarse cell, it is a
* subface, while seen from the refined cell it is one of its faces). The
* reason is that finding neighborship information is a bit easier then, but
* that's all practical reasoning, nothing fundamental.
*
* Since we integrate from the coarse side of the face, we have the mother
* face readily at hand and store the result of the integration over that
* mother face (being the sum of the integrals along the subfaces) in the
* abovementioned map of integrals as well. This consumes some memory more
* than needed, but makes the summing up of the face contributions to the
* cells easier, since then we have the information from all faces of all
* cells at hand and need not think about explicitly determining whether
* a face was refined or not. The same applies for boundary faces, see
* above.
*
*
* <h3>Multiple solution vectors</h3>
*
* In some cases, for example in time-dependent problems, one would
* like to compute the error estimates for several solution vectors
* on the same grid at once, with the same coefficients, boundary
* condition object, etc, e.g. for the solutions on several
* successive time steps. One could then call the functions of this
* class several times for each solution. However, the largest factor
* in the computation of the error estimates (in terms of computing
* time) is initialization of FEFaceValues and FESubFaceValues
* objects, and iterating through all faces and subfaces. If the
* solution vectors live on the same grid, this effort can be reduced
* significantly by treating all solution vectors at the same time,
* initializing the FEFaceValues objects only once per cell and for
* all solution vectors at once, and also only looping through the
* triangulation only once. For this reason, besides the @p estimate
* function in this class that takes a single input vector and
* returns a single output vector, there is also a function that
* accepts several in- and output vectors at the same time.
*
* @ingroup numerics
* @author Wolfgang Bangerth, 1998, 1999, 2000, 2004, 2006; parallelization by Thomas Richter, 2000
*/
template <int dim, int spacedim=dim>
class KellyErrorEstimator
{
public:
/**
* Implementation of the error
* estimator described above. You
* may give a coefficient, but
* there is a default value which
* denotes the constant
* coefficient with value
* one. The coefficient function
* may either be a scalar one, in
* which case it is used for all
* components of the finite
* element, or a vector-valued
* one with as many components as
* there are in the finite
* element; in the latter case,
* each component is weighted by
* the respective component in
* the coefficient.
*
* You might give a list of
* components you want to
* evaluate, in case the finite
* element used by the
* DoFHandler object is
* vector-valued. You then have
* to set those entries to true
* in the bit-vector
* @p component_mask for which the
* respective component is to be
* used in the error
* estimator. The default is to
* use all components, which is
* done by either providing a
* bit-vector with all-set
* entries, or an empty
* bit-vector.
*
* The @p subdomain_id parameter
* indicates whether we shall compute
* indicators for all cells (in case its
* value is the default,
* <tt>numbers::invalid_unsigned_int</tt>),
* or only for the cells belonging to a
* certain subdomain with the given
* indicator. The latter case is used for
* parallel computations where all
* processor nodes have the global grid
* stored, and could well compute all the
* indicators for all cells themselves,
* but where it is more efficient to have
* each process compute only indicators
* for the cells it owns, and have them
* exchange the resulting information
* afterwards. This is in particular true
* for the case where meshes are very
* large and computing indicators for @em
* every cells is too expensive, while
* computing indicators for only local
* cells is acceptable. Note that if you
* only ask for the indicators of a
* certain subdomain to be computed, you
* must nevertheless make sure that this
* function has access to the correct
* node values of @em all degrees of
* freedom. This is since the function
* needs to access DoF values on
* neighboring cells as well, even if
* they belong to a different subdomain.
*
* The @p material_id parameter has a
* similar meaning: if not set to its
* default value, it means that
* indicators will only be computed for
* cells with this particular material
* id.
*
* The @p n_threads parameter used to
* indicate the number of threads to be
* used to compute the error
* estimator. This parameter is now
* ignored, with the number of threads
* determined automatically. The
* parameter is retained for
* compatibility with old versions of the
* library.
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<dim, spacedim> &mapping,
const DH &dof,
const Quadrature<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Calls the @p estimate
* function, see above, with
* <tt>mapping=MappingQ1@<dim@>()</tt>.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const Quadrature<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Same function as above, but
* accepts more than one solution
* vector and returns one error
* vector for each solution
* vector. For the reason of
* existence of this function,
* see the general documentation
* of this class.
*
* Since we do not want to force
* the user of this function to
* copy around their solution
* vectors, the vector of
* solution vectors takes
* pointers to the solutions,
* rather than being a vector of
* vectors. This makes it simpler
* to have the solution vectors
* somewhere in memory, rather
* than to have them collected
* somewhere special. (Note that
* it is not possible to
* construct of vector of
* references, so we had to use a
* vector of pointers.)
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<dim, spacedim> &mapping,
const DH &dof,
const Quadrature<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Calls the @p estimate
* function, see above, with
* <tt>mapping=MappingQ1@<dim@>()</tt>.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const Quadrature<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<dim, spacedim> &mapping,
const DH &dof,
const hp::QCollection<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const hp::QCollection<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<dim, spacedim> &mapping,
const DH &dof,
const hp::QCollection<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const hp::QCollection<dim-1> &quadrature,
const typename FunctionMap<dim>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<dim> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Exception
*/
DeclException0 (ExcInvalidBoundaryIndicator);
/**
* Exception
*/
DeclException0 (ExcInvalidComponentMask);
/**
* Exception
*/
DeclException0 (ExcInvalidCoefficient);
/**
* Exception
*/
DeclException0 (ExcInvalidBoundaryFunction);
/**
* Exception
*/
DeclException2 (ExcIncompatibleNumberOfElements,
int, int,
<< "The number of elements " << arg1 << " and " << arg2
<< " of the vectors do not match!");
/**
* Exception
*/
DeclException0 (ExcInvalidSolutionVector);
/**
* Exception
*/
DeclException0 (ExcNoSolutions);
};
/**
* This is a specialization of the general template for 1d. The implementation
* is sufficiently different for 1d to justify this specialization. The basic
* difference between 1d and all other space dimensions is that in 1d, there
* are no faces of cells, just the vertices between line segments, so we have
* to compute the jump terms differently. However, this class offers exactly
* the same public functions as the general template, so that a user will not
* see any difference.
*
* @author Wolfgang Bangerth, 1998, 2004.
*/
template <int spacedim>
class KellyErrorEstimator<1,spacedim>
{
public:
/**
* Implementation of the error
* estimator described above. You
* may give a coefficient, but
* there is a default value which
* denotes the constant
* coefficient with value
* one. The coefficient function
* may either be a scalar one, in
* which case it is used for all
* components of the finite
* element, or a vector-valued
* one with as many components as
* there are in the finite
* element; in the latter case,
* each component is weighted by
* the respective component in
* the coefficient.
*
* You might give a list of components
* you want to evaluate, in case the
* finite element used by the DoFHandler
* object is vector-valued. You then have
* to set those entries to true in the
* bit-vector @p component_mask for which
* the respective component is to be used
* in the error estimator. The default is
* to use all components, which is done
* by either providing a bit-vector with
* all-set entries, or an empty
* bit-vector. All the other parameters
* are as in the general case used for 2d
* and higher.
*
* The estimator supports multithreading
* and splits the cells to
* <tt>multithread_info.n_default_threads</tt>
* (default) threads. The number of
* threads to be used in multithreaded
* mode can be set with the last
* parameter of the error estimator.
* Multithreading is not presently
* implemented for 1d, but we retain the
* respective parameter for compatibility
* with the function signature in the
* general case.
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<1,spacedim> &mapping,
const DH &dof,
const Quadrature<0> &quadrature,
const FunctionMap<1>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Calls the @p estimate
* function, see above, with
* <tt>mapping=MappingQ1<1>()</tt>.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const Quadrature<0> &quadrature,
const FunctionMap<1>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Same function as above, but
* accepts more than one solution
* vectors and returns one error
* vector for each solution
* vector. For the reason of
* existence of this function,
* see the general documentation
* of this class.
*
* Since we do not want to force
* the user of this function to
* copy around their solution
* vectors, the vector of
* solution vectors takes
* pointers to the solutions,
* rather than being a vector of
* vectors. This makes it simpler
* to have the solution vectors
* somewhere in memory, rather
* than to have them collected
* somewhere special. (Note that
* it is not possible to
* construct of vector of
* references, so we had to use a
* vector of pointers.)
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<1,spacedim> &mapping,
const DH &dof,
const Quadrature<0> &quadrature,
const FunctionMap<1>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Calls the @p estimate
* function, see above, with
* <tt>mapping=MappingQ1<1>()</tt>.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const Quadrature<0> &quadrature,
const FunctionMap<1>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<1,spacedim> &mapping,
const DH &dof,
const hp::QCollection<0> &quadrature,
const typename FunctionMap<1>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const hp::QCollection<0> &quadrature,
const typename FunctionMap<1>::type &neumann_bc,
const InputVector &solution,
Vector<float> &error,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const Mapping<1,spacedim> &mapping,
const DH &dof,
const hp::QCollection<0> &quadrature,
const typename FunctionMap<1>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Equivalent to the set of functions
* above, except that this one takes a
* quadrature collection for hp finite
* element dof handlers.
*/
template <typename InputVector, class DH>
static void estimate (const DH &dof,
const hp::QCollection<0> &quadrature,
const typename FunctionMap<1>::type &neumann_bc,
const std::vector<const InputVector *> &solutions,
std::vector<Vector<float>*> &errors,
const std::vector<bool> &component_mask = std::vector<bool>(),
const Function<1> *coefficients = 0,
const unsigned int n_threads = multithread_info.n_default_threads,
const unsigned int subdomain_id = numbers::invalid_unsigned_int,
const unsigned int material_id = numbers::invalid_unsigned_int);
/**
* Exception
*/
DeclException0 (ExcInvalidBoundaryIndicator);
/**
* Exception
*/
DeclException0 (ExcInvalidComponentMask);
/**
* Exception
*/
DeclException0 (ExcInvalidCoefficient);
/**
* Exception
*/
DeclException0 (ExcInvalidBoundaryFunction);
/**
* Exception
*/
DeclException2 (ExcIncompatibleNumberOfElements,
int, int,
<< "The number of elements " << arg1 << " and " << arg2
<< " of the vectors do not match!");
/**
* Exception
*/
DeclException0 (ExcInvalidSolutionVector);
/**
* Exception
*/
DeclException0 (ExcNoSolutions);
};
DEAL_II_NAMESPACE_CLOSE
#endif
|