Merge pull request #132 from memmett/memmett/state0

LGTM. Fix assumption that `state(0)` is the same as `start_state`.

Author: torbjoernk

examples/advection_diffusion/advection_diffusion_sweeper.hpp

@@ -133,7 +133,6 @@ class AdvectionDiffusionSweeper
 
       auto n = this->get_controller()->get_step();
       auto k = this->get_controller()->get_iteration();
-      CLOG(INFO, "data_values") << n << k << qend;
       LOG(INFO) << "err:" << n << k << max << "(" << qend.size() << "," << predict << ")";
 
       this->errors.insert(pair<pair<size_t, size_t>, double>(pair<size_t, size_t>(n, k), max));

examples/advection_diffusion/vanilla_sdc.cpp

@@ -23,14 +23,13 @@ error_map run_vanilla_sdc(double abs_residual_tol)
 
   const size_t nsteps = pfasst::config::get_value<size_t>("num_steps", 4);
   const double dt     = pfasst::config::get_value<double>("delta_step", 0.01);
-  const size_t nnodes = pfasst::config::get_value<size_t>("num_nodes", 5);
+  const size_t nnodes = pfasst::config::get_value<size_t>("num_nodes", 3);
   const size_t ndofs  = pfasst::config::get_value<size_t>("spatial_dofs", 64);
   const size_t niters = pfasst::config::get_value<size_t>("num_iter", 4);
   const pfasst::quadrature::QuadratureType quad_type = \
-    pfasst::config::get_value<pfasst::quadrature::QuadratureType>("nodes_type", pfasst::quadrature::QuadratureType::GaussLobatto);
+    pfasst::config::get_value<pfasst::quadrature::QuadratureType>("nodes_type", pfasst::quadrature::QuadratureType::GaussLegendre);
 
-   // auto quad    = pfasst::quadrature::quadrature_factory(nnodes, pfasst::quadrature::QuadratureType::GaussLobatto);
-  auto quad    = pfasst::quadrature::quadrature_factory(nnodes-2, pfasst::quadrature::QuadratureType::GaussLegendre);
+  auto quad    = pfasst::quadrature::quadrature_factory(nnodes, quad_type);
   auto factory = make_shared<pfasst::encap::VectorFactory<double>>(ndofs);
   auto sweeper = make_shared<AdvectionDiffusionSweeper<>>(ndofs);
 
@@ -56,12 +55,8 @@ error_map run_vanilla_sdc(double abs_residual_tol)
 #ifndef PFASST_UNIT_TESTING
 int main(int argc, char** argv)
 {
-  pfasst::encap::enable_config_options();
-  // First we want to enable command line options for the Advection-Diffusion Sweeper ...
   AdvectionDiffusionSweeper<>::enable_config_options();
-  // ... then we initialize all options other default options and parse given parameters ...
   pfasst::init(argc, argv);
-  pfasst::encap::set_logfile_from_options();
 
   run_vanilla_sdc(0.0);
 }

examples/boris/boris_sweeper.hpp

@@ -380,7 +380,7 @@ namespace pfasst
 
             // building rules for Q_E and Q_I:
             //  Q_E is striclty lower diagonal matrix with delta nodes of column index
-            //  Q_I is lower diagonal matrix with first row and column all zero and delta nodes of 
+            //  Q_I is lower diagonal matrix with first row and column all zero and delta nodes of
             //      column index minus one
             Matrix<time> qe_mat = Matrix<time>(nnodes, nnodes);
             qe_mat.fill(time(0.0));
@@ -423,7 +423,7 @@ namespace pfasst
             this->energy_evals.front() = this->energy_evals.back();
           }
 
-          virtual void evaluate(size_t m) override
+          void evaluate(size_t m)
           {
       //       Vector3d<scalar> vel; vel.fill(scalar(0.0));
       //       Vector3d<scalar> B; B.fill(scalar(0.0));
@@ -448,6 +448,13 @@ namespace pfasst
             this->f_evals++;
           }
 
+          virtual void reevaluate(bool initial) override
+          {
+            for (size_t m = 0; m < this->get_quadrature()->get_num_nodes(); m++) {
+              this->evaluate(m);
+            }
+          }
+
           virtual void predict(bool initial) override
           {
             UNUSED(initial);

include/pfasst/encap/encap_sweeper.hpp

@@ -108,19 +108,16 @@ namespace pfasst
           }
 
           if (coarse) {
-            size_t num_fas = this->quadrature->left_is_node() ? num_nodes -1 : num_nodes;
-            for (size_t m = 0; m < num_fas; m++) {
+            for (size_t m = 0; m < num_nodes; m++) {
               this->fas_corrections.push_back(this->get_factory()->create(pfasst::encap::solution));
             }
           }
-
         }
 
         //! @{
         virtual void spread() override
         {
           for (size_t m = 1; m < this->quadrature->get_num_nodes(); m++) {
-            //            this->get_state(m)->copy(this->start_state);
             this->state[m]->copy(this->state[0]);
           }
         }
@@ -190,18 +187,12 @@ namespace pfasst
         }
 
         /**
-         * evaluates the right hand side at given time node
-         *
-         * This evaluates the right hand side at the given time node with index `m` as returned by
-         * pfasst::encap::EncapSweeper::get_nodes:
-         *
-         * @param[in] m index of the time node to evaluate at
+         * Re-evaluate function values.
          *
          * @note This method must be implemented in derived sweepers.
          */
-        virtual void evaluate(size_t m)
+        virtual void reevaluate(bool initial_only=false)
         {
-          UNUSED(m);
           throw NotImplementedYet("sweeper");
         }
 
@@ -282,8 +273,9 @@ namespace pfasst
         virtual void recv(ICommunicator* comm, int tag, bool blocking) override
         {
           this->start_state->recv(comm, tag, blocking);
-          // XXX
-          this->state.front()->copy(this->start_state);
+          if (this->quadrature->left_is_node()) {
+            this->state[0]->copy(this->start_state);
+          }
         }
 
         virtual void broadcast(ICommunicator* comm) override
@@ -292,8 +284,6 @@ namespace pfasst
             this->start_state->copy(this->end_state);
           }
           this->start_state->broadcast(comm);
-          // XXX
-          this->state.front()->copy(this->start_state);
         }
         //! @}
     };