]> https://gitweb.dealii.org/ - dealii.git/commitdiff
local mg works for Q1
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 20 Apr 2000 14:42:12 +0000 (14:42 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 20 Apr 2000 14:42:12 +0000 (14:42 +0000)
git-svn-id: https://svn.dealii.org/trunk@2761 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/mg_smoother.templates.h
deal.II/deal.II/include/multigrid/multigrid.h
deal.II/deal.II/include/multigrid/multigrid.templates.h
deal.II/deal.II/source/multigrid/multigrid.cc

index cc98583b82c2e14bce010a58eb7be2aba89f47d0..b786b24addc2877c8ce0c969d4279ed76695bca7 100644 (file)
@@ -27,6 +27,7 @@ MGSmootherRelaxation<number>::smooth (const unsigned int   level,
   for(unsigned i=0;i<get_steps();++i)
     {
       (*matrix)[level].residual(h1, u, rhs);
+      set_zero_interior_boundary(level,h1);
       ((*matrix)[level].*relaxation)(h2,h1,omega);
       set_zero_interior_boundary(level,h2);
       u.add(h2);
index cd98244072345059a0c609017c84f8b9c2033c2c..7cfa21e1aee799d10b4b2ca6886f9c2d116ad894 100644 (file)
@@ -221,7 +221,7 @@ class MGTransferPrebuilt : public MGTransferBase
                                          * while row indices belong to the
                                          * child cell, i.e. the fine level.
                                          */
-    vector<SparseMatrix<float> > prolongation_matrices;
+    vector<SparseMatrix<double> > prolongation_matrices;
 };
 
 
index 820d2bdc7d02369eca920ade19e8755e8c674cbc..77d7be4f8b49ed8f4a35b409bae39798181681a2 100644 (file)
@@ -38,15 +38,18 @@ Multigrid<dim>::copy_to_mg (const Vector<number>&)
 template <int dim>
 template <typename number>
 void
-Multigrid<dim>::copy_to_mg (const Vector<number>& src)
+Multigrid<dim>::copy_to_mg (const Vector<number>& osrc)
 {
+                                  // Make src a real finite element function
+  Vector<number> src = osrc;
+//  constraints->distribute(src);
+
   const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell;
   const unsigned int dofs_per_face = mg_dof_handler->get_fe().dofs_per_face;
 
                                   // set the elements of the vectors
                                   // on all levels to zero
   defect.clear();
-//  constraints->condense(src);
   
   vector<unsigned int> global_dof_indices (dofs_per_cell);
   vector<unsigned int> level_dof_indices  (dofs_per_cell);
@@ -102,14 +105,16 @@ Multigrid<dim>::copy_to_mg (const Vector<number>& src)
          for (unsigned int i=0; i<dofs_per_cell; ++i)
            defect[level](level_dof_indices[i]) = src(global_dof_indices[i]);
 
-                                          // Delete values on refinement edge,
-                                          // since restriction will add them again.
          for (unsigned int face_n=0; face_n<GeometryInfo<dim>::faces_per_cell; ++face_n)
            {
              const MGDoFHandler<dim>::face_iterator face = level_cell->face(face_n);
              if (face->has_children())
                {
                  face->get_mg_dof_indices(level_face_indices);
+
+
+                                                  // Delete values on refinement edge,
+                                                  // since restriction will add them again.
                  for (unsigned int i=0; i<dofs_per_face; ++i)
                    defect[level](level_face_indices[i]) = 0.;
                };
index 4a7c0442af7c84fb6cc9670da43322ece448794a..a6878404c0726356867bc1df3438abd5a9e3a205 100644 (file)
@@ -131,7 +131,7 @@ void MGTransferPrebuilt::build_matrices (const MGDoFHandler<dim> &mg_dof)
          };
       prolongation_sparsities[level].compress ();
 
-      prolongation_matrices.push_back (SparseMatrix<float>());
+      prolongation_matrices.push_back (SparseMatrix<double>());
       prolongation_matrices[level].reinit (prolongation_sparsities[level]);
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.