From: guido Date: Tue, 28 Jun 2005 22:39:02 +0000 (+0000) Subject: projection compiles, but results are wrong X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=31e555e476af1d05e1b4fc9c234ccd2cf35d0d1e;p=dealii-svn.git projection compiles, but results are wrong git-svn-id: https://svn.dealii.org/trunk@10960 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 3bba9f9d9e..090cf49d20 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -500,23 +500,21 @@ FETools::compute_embedding_matrices(const FiniteElement& fe, const unsigned int nq = q_fine.n_quadrature_points; FEValues fine (mapping, fe, q_fine, - update_q_points | update_JxW_values | update_values); - - /** - * We search for the polynomial on - * the small cell, being equal to - * the coarse polynomial in all - * quadrature points. - * - * First build the matrix for this - * least squares problem. This - * contains the values of the fine - * cell polynomials in the fine - * cell grid points. - * - * This matrix is the same for all - * children. - */ + update_q_points | update_JxW_values | update_values); + + // We search for the polynomial on + // the small cell, being equal to + // the coarse polynomial in all + // quadrature points. + + // First build the matrix for this + // least squares problem. This + // contains the values of the fine + // cell polynomials in the fine + // cell grid points. + + // This matrix is the same for all + // children. fine.reinit(tria.begin_active()); FullMatrix A(nq*nd, n); for (unsigned int d=0;d& fe, -//TODO[GK]: this function is presently not instantiated, and probably doesn't even compile. +//TODO[GK]: this function does not work yet. template void FETools::compute_projection_matrices(const FiniteElement& fe, FullMatrix* matrices) { -//TODO[GK]: clean up this function and make it work - Assert (false, ExcNotImplemented()); - + Assert(false, ExcNotImplemented()); const unsigned int nc = GeometryInfo::children_per_cell; const unsigned int n = fe.dofs_per_cell; const unsigned int nd = fe.n_components(); @@ -608,23 +604,10 @@ FETools::compute_projection_matrices(const FiniteElement& fe, Assert(matrices[i].m() == n, ExcDimensionMismatch(matrices[i].m(),n)); } - /* - * Set up two meshes, one with a - * single reference cell and the - * other refined, together with - * DoFHandler and finite elements. - */ -//TODO[GK]: This should be changed just like the function above, i.e. not use a DoFHandler at all (and use the Triangulation alone), and to use only a single triangulation. this would make it significantly more efficient - Triangulation tr_coarse; - Triangulation tr_fine; - GridGenerator::hyper_cube (tr_coarse, 0, 1); - GridGenerator::hyper_cube (tr_fine, 0, 1); - tr_fine.refine_global(1); - DoFHandler dof_coarse(tr_coarse); - dof_coarse.distribute_dofs(fe); - DoFHandler dof_fine(tr_fine); - dof_fine.distribute_dofs(fe); - + Triangulation tr; + GridGenerator::hyper_cube (tr, 0, 1); + tr.refine_global(1); + MappingCartesian mapping; QGauss q_fine(degree+1); const unsigned int nq = q_fine.n_quadrature_points; @@ -634,43 +617,88 @@ FETools::compute_projection_matrices(const FiniteElement& fe, FEValues fine (mapping, fe, q_fine, update_q_points | update_JxW_values | update_values); - typename DoFHandler::active_cell_iterator coarse_cell - = dof_coarse.begin_active(); - typename DoFHandler::active_cell_iterator fine_cell; + typename Triangulation::cell_iterator coarse_cell + = tr.begin(0); + typename Triangulation::cell_iterator fine_cell; // Compute the coarse level mass // matrix - coarse.reinit(dof_coarse); + coarse.reinit(coarse_cell); FullMatrix A(n, n); for (unsigned int k=0;k H(A); Vector v_coarse(n); Vector v_fine(n); - unsigned int cell_number = 0; - for (fine_cell = dof_fine.begin_active(); - fine_cell != dof_fine.end(); - ++fine_cell, ++cell_number) + for (unsigned int cell_number=0;cell_number::children_per_cell;++cell_number) { + FullMatrix &this_matrix = matrices[cell_number]; + // Compute right hand side, // which is a fine level basis // function tested with the // coarse level functions. - Assert(false, ExcNotImplemented()); + fine.reinit(coarse_cell->child(cell_number)); + Quadrature q_coarse (fine.get_quadrature_points(), + fine.get_JxW_values()); + FEValues coarse (mapping, fe, q_coarse, update_values); + coarse.reinit(coarse_cell); + + // Build RHS + + // Outer loop over all fine + // grid shape functions phi_j + for (unsigned int j=0;j (const FiniteElement &, FullMatrix*); +template +void FETools::compute_projection_matrices +(const FiniteElement &, FullMatrix*); + template void FETools::interpolate (const DoFHandler &, const Vector &,