From: Guido Kanschat Date: Thu, 29 Sep 2005 12:45:08 +0000 (+0000) Subject: new constructor using auxiliary functions X-Git-Tag: v8.0.0~13051 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=835c27889b594fbb073013d9915cda62dd6df23d;p=dealii.git new constructor using auxiliary functions git-svn-id: https://svn.dealii.org/trunk@11552 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc b/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc index ca491a306b..34d5650077 100644 --- a/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc +++ b/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc @@ -45,18 +45,54 @@ FE_RaviartThomasNodal::FE_RaviartThomasNodal (const unsigned int deg) std::vector(dim,true))) { Assert (dim >= 2, ExcImpossibleInDim(dim)); + const unsigned int n_dofs = this->dofs_per_cell; this->mapping_type = this->independent_on_cartesian; // These must be done first, since // they change the evaluation of // basis functions + + // Set up the generalized support + // points initialize_unit_support_points (deg); - initialize_node_matrix(); + + //Now compute the inverse node + //matrix, generating the correct + //basis functions from the raw + //ones. + + // We use an auxiliary matrix in + // this function. Therefore, + // inverse_node_matrix is still + // empty and shape_value_component + // returns the 'raw' shape values. + FullMatrix M(n_dofs, n_dofs); + FETools::compute_node_matrix(M, *this); + this->inverse_node_matrix.reinit(n_dofs, n_dofs); + this->inverse_node_matrix.invert(M); + // From now on, the shape functions + // will be the correct ones, not + // the raw shape functions anymore. for (unsigned int i=0; i::children_per_cell; ++i) this->prolongation[i].reinit (this->dofs_per_cell, this->dofs_per_cell); FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + + std::vector > + face_embeddings(1<<(dim-1), FullMatrix(this->dofs_per_face, + this->dofs_per_face)); + FETools::compute_face_embedding_matrices(*this, &face_embeddings[0], 0, 0); + interface_constraints.reinit((1<<(dim-1)) * this->dofs_per_face, + this->dofs_per_face); + unsigned int target_row=0; + for (unsigned int d=0;d void FE_RaviartThomasNodal::initialize_node_matrix () { - const unsigned int n_dofs = this->dofs_per_cell; - // We use an auxiliary matrix in - // this function. Therefore, - // inverse_node_matrix is still - // empty and shape_value_component - // returns the 'raw' shape values. - FullMatrix M(n_dofs, n_dofs); - FETools::compute_node_matrix(M, *this); - this->inverse_node_matrix.reinit(n_dofs, n_dofs); - this->inverse_node_matrix.invert(M); }