From: Guido Kanschat Date: Fri, 19 Aug 2011 13:45:18 +0000 (+0000) Subject: add trace operators X-Git-Tag: v8.0.0~3651 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d592f4c719b0bc680f82ece4b1fc39e577aa560c;p=dealii.git add trace operators git-svn-id: https://svn.dealii.org/trunk@24117 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/integrators/maxwell.h b/deal.II/include/deal.II/integrators/maxwell.h index 54a3f16999..b7a6c334ba 100644 --- a/deal.II/include/deal.II/integrators/maxwell.h +++ b/deal.II/include/deal.II/integrators/maxwell.h @@ -37,7 +37,7 @@ namespace LocalIntegrators * \partial_3 u_2 - \partial_2 u_3 \\ * \partial_1 u_3 - \partial_3 u_1 \\ * \partial_2 u_1 - \partial_1 u_2 - * \end{pmatrix} + * \end{pmatrix}. * @f] * * In two space dimensions, the curl is obtained by extending a vector @@ -46,7 +46,7 @@ namespace LocalIntegrators * curl of a vector function and the vector curl of a scalar * function. The current implementation exchanges the sign and we have: * @f[ - * \nabla \times \mathbf u = \partial_1 u_2 - \partial 2 u_1 + * \nabla \times \mathbf u = \partial_1 u_2 - \partial 2 u_1, * \qquad * \nabla \times p = \begin{pmatrix} * \partial_2 p \\ -\partial_1 p @@ -156,8 +156,8 @@ namespace LocalIntegrators * \int_Z \nabla\!\times\! u \cdot * \nabla\!\times\! v \,dx * @f] + * in weak form. * - * @ingroup Integrators * @author Guido Kanschat * @date 2011 */ @@ -203,14 +203,15 @@ namespace LocalIntegrators } /** - * The curl operator + * The matrix for the curl operator * @f[ * \int_Z \nabla\!\times\! u \cdot v \,dx. * @f] * - * This is the standard curl operator in 3D and the scalar curl in 2D. - * - * @ingroup Integrators + * This is the standard curl operator in 3D and the scalar curl in + * 2D. The vector curl operator can be obtained by exchanging test and + * trial functions. + * * @author Guido Kanschat * @date 2011 */ @@ -247,6 +248,212 @@ namespace LocalIntegrators } } } + + /** + * The matrix for weak boundary + * condition of Nitsche type for + * the tangential component in + * Maxwell systems. + * + * @f[ + * \int_F \biggl( 2\gamma + * (u\times n) (v\times n) - + * (u\times n)(\nu \nabla\times + * v) - (v\times + * n)(\nu \nabla\times u) + * \biggr) + * @f] + * + * @author Guido Kanschat + * @date 2011 + */ + template + void nitsche_curl_matrix ( + FullMatrix& M, + const FEValuesBase& fe, + double penalty, + double factor = 1.) + { + const unsigned int n_dofs = fe.dofs_per_cell; + + AssertDimension(fe.get_fe().n_components(), dim); + AssertDimension(M.m(), n_dofs); + AssertDimension(M.n(), n_dofs); + + // Depending on the + // dimension, the cross + // product is either a scalar + // (2d) or a vector + // (3d). Accordingly, in the + // latter case we have to sum + // up three bilinear forms, + // but in 2d, we don't. Thus, + // we need to adapt the loop + // over all dimensions + const unsigned int d_max = (dim==2) ? 1 : dim; + + for (unsigned k=0;k& n = fe.normal_vector(k); + for (unsigned i=0;i + void tangential_trace_matrix ( + FullMatrix& M, + const FEValuesBase& fe, + double factor = 1.) + { + const unsigned int n_dofs = fe.dofs_per_cell; + + AssertDimension(fe.get_fe().n_components(), dim); + AssertDimension(M.m(), n_dofs); + AssertDimension(M.n(), n_dofs); + + // Depending on the + // dimension, the cross + // product is either a scalar + // (2d) or a vector + // (3d). Accordingly, in the + // latter case we have to sum + // up three bilinear forms, + // but in 2d, we don't. Thus, + // we need to adapt the loop + // over all dimensions + const unsigned int d_max = (dim==2) ? 1 : dim; + + for (unsigned k=0;k& n = fe.normal_vector(k); + for (unsigned i=0;i + inline void ip_curl_matrix ( + FullMatrix& M11, + FullMatrix& M12, + FullMatrix& M21, + FullMatrix& M22, + const FEValuesBase& fe1, + const FEValuesBase& fe2, + const double pen, + const double factor1 = 1., + const double factor2 = -1.) + { + const unsigned int n_dofs = fe1.dofs_per_cell; + + AssertDimension(fe1.get_fe().n_components(), dim); + AssertDimension(fe2.get_fe().n_components(), dim); + AssertDimension(M11.m(), n_dofs); + AssertDimension(M11.n(), n_dofs); + AssertDimension(M12.m(), n_dofs); + AssertDimension(M12.n(), n_dofs); + AssertDimension(M21.m(), n_dofs); + AssertDimension(M21.n(), n_dofs); + AssertDimension(M22.m(), n_dofs); + AssertDimension(M22.n(), n_dofs); + + const double nu1 = factor1; + const double nu2 = (factor2 < 0) ? factor1 : factor2; + const double penalty = .5 * pen * (nu1 + nu2); + + // Depending on the + // dimension, the cross + // product is either a scalar + // (2d) or a vector + // (3d). Accordingly, in the + // latter case we have to sum + // up three bilinear forms, + // but in 2d, we don't. Thus, + // we need to adapt the loop + // over all dimensions + const unsigned int d_max = (dim==2) ? 1 : dim; + + for (unsigned k=0;k& n = fe1.normal_vector(k); + for (unsigned i=0;i