fourth_error_estimators.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 #include "libmesh/libmesh_config.h"
20 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
21 
22 
23 // C++ includes
24 #include <algorithm> // for std::fill
25 #include <cstdlib> // *must* precede <cmath> for proper std:abs() on PGI, Sun Studio CC
26 #include <cmath> // for sqrt
27 
28 
29 // Local Includes
30 #include "libmesh/libmesh_common.h"
32 #include "libmesh/error_vector.h"
33 #include "libmesh/fe_base.h"
35 #include "libmesh/elem.h"
36 #include "libmesh/system.h"
37 
38 #include "libmesh/dense_vector.h"
39 #include "libmesh/tensor_tools.h"
40 
41 
42 namespace libMesh
43 {
44 
45 
46 void
48  ErrorVector&,
49  bool)
50 {
51  // We'll need second derivatives for Laplacian jump computation
52  fe_fine->get_d2phi();
53  fe_coarse->get_d2phi();
54 }
55 
56 
57 
58 void
60 {
61  Real error = 1.e-30;
62  unsigned int n_qp = fe_fine->n_quadrature_points();
63  unsigned int n_fine_dofs = Ufine.size();
64  unsigned int n_coarse_dofs = Ucoarse.size();
65 
66  unsigned int dim = fine_elem->dim();
67 
68  std::vector<std::vector<RealTensor> > d2phi_coarse = fe_coarse->get_d2phi();
69  std::vector<std::vector<RealTensor> > d2phi_fine = fe_fine->get_d2phi();
70  std::vector<Real> JxW_face = fe_fine->get_JxW();
71 
72  for (unsigned int qp=0; qp != n_qp; ++qp)
73  {
74  // Calculate solution gradients on fine and coarse elements
75  // at this quadrature point
76  Number laplacian_fine = 0., laplacian_coarse = 0.;
77 
78  for (unsigned int i=0; i != n_coarse_dofs; ++i)
79  {
80  laplacian_coarse += d2phi_coarse[i][qp](0,0) * Ucoarse(i);
81  if (dim > 1)
82  laplacian_coarse += d2phi_coarse[i][qp](1,1) * Ucoarse(i);
83  if (dim > 2)
84  laplacian_coarse += d2phi_coarse[i][qp](2,2) * Ucoarse(i);
85  }
86 
87  for (unsigned int i=0; i != n_fine_dofs; ++i)
88  {
89  laplacian_fine += d2phi_fine[i][qp](0,0) * Ufine(i);
90  if (dim > 1)
91  laplacian_fine += d2phi_fine[i][qp](1,1) * Ufine(i);
92  if (dim > 2)
93  laplacian_fine += d2phi_fine[i][qp](2,2) * Ufine(i);
94  }
95 
96 
97  // Find the jump in the Laplacian
98  // at this quadrature point
99  const Number jump = laplacian_fine - laplacian_coarse;
100  const Real jump2 = TensorTools::norm_sq(jump);
101 
102  // Accumulate the jump integral
103  error += JxW_face[qp] * jump2;
104  }
105 
106  // Add the h-weighted jump integral to each error term
107  fine_error =
108  error * fine_elem->hmax() * error_norm.weight(var);
109  coarse_error =
110  error * coarse_elem->hmax() * error_norm.weight(var);
111 }
112 
113 } // namespace libMesh
114 
115 #else // defined (LIBMESH_ENABLE_SECOND_DERIVATIVES)
116 
118 
119 namespace libMesh
120 {
121 
122 void
124  ErrorVector&,
125  bool)
126 {
127  libMesh::err << "Error: LaplacianErrorEstimator requires second "
128  << "derivative support; try configuring libmesh with "
129  << "--enable-second" << std::endl;
130  libmesh_error();
131 }
132 
133 
134 void
136 {
137  libMesh::err << "Error: LaplacianErrorEstimator requires second "
138  << "derivative support; try configuring libmesh with "
139  << "--enable-second" << std::endl;
140  libmesh_error();
141 }
142 
143 } // namespace libMesh
144 
145 #endif // defined (LIBMESH_ENABLE_SECOND_DERIVATIVES)

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:05 UTC

Hosted By:
SourceForge.net Logo