15#ifndef dealii_integrators_maxwell_h
16#define dealii_integrators_maxwell_h
142 normal[2] * (
g2[1] -
g0[2]) + normal[1] * (
g1[0] -
g0[1]);
144 normal[0] * (
g0[2] -
g1[0]) + normal[2] * (
g2[1] -
g1[2]);
146 normal[1] * (
g1[0] -
g2[1]) + normal[0] * (
g0[2] -
g2[0]);
166 const double factor = 1.)
168 const unsigned int n_dofs = fe.dofs_per_cell;
183 const unsigned int d_max = (dim == 2) ? 1 : dim;
185 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
187 const double dx = factor * fe.JxW(
k);
188 for (
unsigned int i = 0; i < n_dofs; ++i)
189 for (
unsigned int j = 0;
j < n_dofs; ++
j)
190 for (
unsigned int d = 0; d <
d_max; ++d)
192 const unsigned int d1 = (d + 1) % dim;
193 const unsigned int d2 = (d + 2) % dim;
195 const double cv = fe.shape_grad_component(i,
k,
d2)[
d1] -
196 fe.shape_grad_component(i,
k,
d1)[
d2];
197 const double cu = fe.shape_grad_component(
j,
k,
d2)[
d1] -
198 fe.shape_grad_component(
j,
k,
d1)[
d2];
222 const unsigned int n_dofs = fe.dofs_per_cell;
231 const unsigned int d_max = (dim == 2) ? 1 : dim;
233 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
235 const double dx = fe.JxW(
k) * factor;
236 for (
unsigned int i = 0; i <
t_dofs; ++i)
237 for (
unsigned int j = 0;
j < n_dofs; ++
j)
238 for (
unsigned int d = 0; d <
d_max; ++d)
240 const unsigned int d1 = (d + 1) % dim;
241 const unsigned int d2 = (d + 2) % dim;
243 const double vv =
fetest.shape_value_component(i,
k, d);
244 const double cu = fe.shape_grad_component(
j,
k,
d2)[
d1] -
245 fe.shape_grad_component(
j,
k,
d1)[
d2];
272 const unsigned int n_dofs = fe.dofs_per_cell;
288 const unsigned int d_max = (dim == 2) ? 1 : dim;
290 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
292 const double dx = factor * fe.JxW(
k);
294 for (
unsigned int i = 0; i < n_dofs; ++i)
295 for (
unsigned int j = 0;
j < n_dofs; ++
j)
296 if (fe.get_fe().has_support_on_face(i,
face_no) &&
297 fe.get_fe().has_support_on_face(
j,
face_no))
299 for (
unsigned int d = 0; d <
d_max; ++d)
301 const unsigned int d1 = (d + 1) % dim;
302 const unsigned int d2 = (d + 2) % dim;
304 const double cv = fe.shape_grad_component(i,
k,
d2)[
d1] -
305 fe.shape_grad_component(i,
k,
d1)[
d2];
306 const double cu = fe.shape_grad_component(
j,
k,
d2)[
d1] -
307 fe.shape_grad_component(
j,
k,
d1)[
d2];
309 fe.shape_value_component(i,
k,
d1) * n[
d2] -
310 fe.shape_value_component(i,
k,
d2) * n[
d1];
312 fe.shape_value_component(
j,
k,
d1) * n[
d2] -
313 fe.shape_value_component(
j,
k,
d2) * n[
d1];
333 const unsigned int n_dofs = fe.dofs_per_cell;
349 const unsigned int d_max = (dim == 2) ? 1 : dim;
351 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
353 const double dx = factor * fe.JxW(
k);
355 for (
unsigned int i = 0; i < n_dofs; ++i)
356 for (
unsigned int j = 0;
j < n_dofs; ++
j)
357 for (
unsigned int d = 0; d <
d_max; ++d)
359 const unsigned int d1 = (d + 1) % dim;
360 const unsigned int d2 = (d + 2) % dim;
362 const double v = fe.shape_value_component(i,
k,
d1) * n(
d2) -
363 fe.shape_value_component(i,
k,
d2) * n(
d1);
364 const double u = fe.shape_value_component(
j,
k,
d1) * n(
d2) -
365 fe.shape_value_component(
j,
k,
d2) * n(
d1);
367 M(i,
j) += dx *
u * v;
396 const unsigned int n_dofs =
fe1.n_dofs_per_cell();
423 const unsigned int d_max = (dim == 2) ? 1 : dim;
425 for (
unsigned int k = 0;
k <
fe1.n_quadrature_points; ++
k)
427 const double dx =
fe1.JxW(
k);
429 for (
unsigned int i = 0; i < n_dofs; ++i)
430 for (
unsigned int j = 0;
j < n_dofs; ++
j)
431 for (
unsigned int d = 0; d <
d_max; ++d)
433 const unsigned int d1 = (d + 1) % dim;
434 const unsigned int d2 = (d + 2) % dim;
438 fe1.shape_grad_component(i,
k,
d1)[
d2];
441 fe2.shape_grad_component(i,
k,
d1)[
d2];
451 fe1.shape_value_component(
j,
k,
d1) * n(
d2) -
452 fe1.shape_value_component(
j,
k,
d2) * n(
d1);
454 -
fe2.shape_value_component(
j,
k,
d1) * n(
d2) +
455 fe2.shape_value_component(
j,
k,
d2) * n(
d1);
457 fe1.shape_value_component(i,
k,
d1) * n(
d2) -
458 fe1.shape_value_component(i,
k,
d2) * n(
d1);
460 -
fe2.shape_value_component(i,
k,
d1) * n(
d2) +
461 fe2.shape_value_component(i,
k,
d2) * n(
d1);
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define AssertDimension(dim1, dim2)
#define DEAL_II_NOT_IMPLEMENTED()
void nitsche_curl_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const unsigned int face_no, double penalty, double factor=1.)
void curl_curl_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
void tangential_trace_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double factor=1.)
void curl_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, double factor=1.)
Tensor< 1, dim > tangential_curl(const Tensor< 1, dim > &g0, const Tensor< 1, dim > &g1, const Tensor< 1, dim > &g2, const Tensor< 1, dim > &normal)
void ip_curl_matrix(FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const double pen, const double factor1=1., const double factor2=-1.)
Tensor< 1, dim > curl_curl(const Tensor< 2, dim > &h0, const Tensor< 2, dim > &h1, const Tensor< 2, dim > &h2)
Library of integrals over cells and faces.