]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use dVector instead of vector<...> and fix a serious bug.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 7 Apr 1998 14:42:00 +0000 (14:42 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 7 Apr 1998 14:42:00 +0000 (14:42 +0000)
git-svn-id: https://svn.dealii.org/trunk@154 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/source/dofs/dof_handler.cc
deal.II/deal.II/source/numerics/base.cc

index fa51a436f33a0b3438c090a3e3f19c213b6a4e3c..bbf555b741d8c9fe4ac105756e9e611406ff603f 100644 (file)
@@ -1466,10 +1466,10 @@ unsigned int DoFHandler<2>::max_transfer_entries (const unsigned int max_level_d
 
 
 template <int dim>
-void DoFHandler<dim>::distribute_cell_to_dof_vector (const vector<double> &cell_data,
-                                                    dVector              &dof_data) const {
-  Assert (cell_data.size()==tria->n_active_cells(),
-         ExcWrongSize (cell_data.size(), tria->n_active_cells()));
+void DoFHandler<dim>::distribute_cell_to_dof_vector (const dVector &cell_data,
+                                                    dVector       &dof_data) const {
+  Assert (cell_data.n()==tria->n_active_cells(),
+         ExcWrongSize (cell_data.n(), tria->n_active_cells()));
 
                                   // assign the right size to the output vector
   dof_data.reinit (n_dofs());
@@ -1491,7 +1491,7 @@ void DoFHandler<dim>::distribute_cell_to_dof_vector (const vector<double> &cell_
        {
                                           // sum up contribution of the
                                           // present_cell to this dof
-         dof_data(dof_indices[i]) += cell_data[present_cell];
+         dof_data(dof_indices[i]) += cell_data(present_cell);
                                           // note that we added another
                                           // summand
          ++touch_count[dof_indices[i]];
index 69e88f9bfbc68c98d1369a72fd31d99a7bc230c5..993c0e6156df11fdc5032c17d41b25f6d5ed73f3 100644 (file)
@@ -6,6 +6,7 @@
 #include <grid/tria_iterator.h>
 #include <basic/data_io.h>
 #include <basic/function.h>
+#include <fe/quadrature.h>
 
 #include "../../../mia/control.h"
 #include "../../../mia/vectormemory.h"
@@ -148,16 +149,15 @@ void ProblemBase<dim>::solve () {
 
 template <int dim>
 void ProblemBase<dim>::integrate_difference (const Function<dim>      &exact_solution,
-                                            vector<double>           &difference,
+                                            dVector                  &difference,
                                             const Quadrature<dim>    &q,
                                             const FiniteElement<dim> &fe,
                                             const NormType           &norm) const {
   Assert ((tria!=0) && (dof_handler!=0), ExcNoTriaSelected());  
   Assert (fe == dof_handler->get_selected_fe(), ExcInvalidFE());
 
-  difference.erase (difference.begin(), difference.end());
-  difference.reserve (tria->n_cells());
-
+  difference.reinit (tria->n_active_cells());
+  
   FEValues<dim>::UpdateStruct update_flags;
   update_flags.q_points   = true;
   update_flags.jacobians  = true;
@@ -176,7 +176,7 @@ void ProblemBase<dim>::integrate_difference (const Function<dim>      &exact_sol
                                        endc = dof_handler->end();
   Triangulation<dim>::active_cell_iterator tria_cell = dof_handler->get_tria().begin_active();
   
-  for (; cell != endc; ++cell, ++tria_cell) 
+  for (unsigned int index=0; cell != endc; ++cell, ++tria_cell, ++index)
     {
       double diff=0;
                                       // initialize for this cell
@@ -209,6 +209,7 @@ void ProblemBase<dim>::integrate_difference (const Function<dim>      &exact_sol
                                             // and assign that to the vector
                                             // \psi.
            const unsigned int n_dofs = fe.total_dofs;
+           const unsigned int n_q_points = q.n_quadrature_points;
            const dFMatrix & shape_values = fe_values.get_shape_values();
            vector<double>   dof_values;
            cell->get_dof_values (solution, dof_values);
@@ -221,7 +222,7 @@ void ProblemBase<dim>::integrate_difference (const Function<dim>      &exact_sol
                                       psi);
                                             // then subtract finite element
                                             // solution
-           for (unsigned int j=0; j<n_dofs; ++j) 
+           for (unsigned int j=0; j<n_q_points; ++j) 
              for (unsigned int i=0; i<n_dofs; ++i)
                psi[j] -= dof_values[i]*shape_values(i,j);
 
@@ -277,7 +278,7 @@ void ProblemBase<dim>::integrate_difference (const Function<dim>      &exact_sol
       
                                       // append result of this cell
                                       // to the end of the vector
-      difference.push_back (diff);
+      difference(index) = diff;
     };
 };
 

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.