WIP: Remove assumption that state(0) is start_state.

Signed-off-by: Matthew Emmett <[email protected]>

Author: memmett

include/pfasst/encap/encap_sweeper.hpp

@@ -108,7 +108,7 @@ namespace pfasst
           }
 
           if (coarse) {
-            size_t num_fas = this->quadrature->left_is_node() ? num_nodes -1 : num_nodes;
+            size_t num_fas = this->quadrature->left_is_node() ? num_nodes - 1 : num_nodes;
             for (size_t m = 0; m < num_fas; m++) {
               this->fas_corrections.push_back(this->get_factory()->create(pfasst::encap::solution));
             }
@@ -190,18 +190,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 +276,6 @@ 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);
         }
 
         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);
         }
         //! @}
     };

include/pfasst/encap/poly_interp.hpp

@@ -46,14 +46,14 @@ namespace pfasst
           auto fine_factory = fine.get_factory();
 
           auto crse_delta = crse_factory->create(solution);
-          this->restrict(crse_delta, fine.get_state(0));
-          crse_delta->saxpy(-1.0, crse.get_state(0));
+          this->restrict(crse_delta, fine.get_start_state());
+          crse_delta->saxpy(-1.0, crse.get_start_state());
 
           auto fine_delta = fine_factory->create(solution);
           this->interpolate(fine_delta, crse_delta);
-          fine.get_state(0)->saxpy(-1.0, fine_delta);
+          fine.get_start_state()->saxpy(-1.0, fine_delta);
 
-          fine.evaluate(0);
+          fine.reevaluate(true);
         }
 
 
@@ -68,6 +68,11 @@ namespace pfasst
             tmat = pfasst::quadrature::compute_interp<time>(fine.get_nodes(), crse.get_nodes());
           }
 
+          if (interp_initial) {
+            this->interpolate_initial(dst, src);
+            return;
+          }
+
           size_t nfine = fine.get_nodes().size();
           size_t ncrse = crse.get_nodes().size();
 
@@ -80,20 +85,15 @@ namespace pfasst
           for (size_t m = 0; m < ncrse; m++) { fine_delta[m] = fine_factory->create(solution); }
 
           auto crse_delta = crse_factory->create(solution);
-          size_t m0 = interp_initial ? 0 : 1;
-          for (size_t m = m0; m < ncrse; m++) {
+          for (size_t m = 0; m < ncrse; m++) {
             crse_delta->copy(crse.get_state(m));
             crse_delta->saxpy(-1.0, crse.get_saved_state(m));
             interpolate(fine_delta[m], crse_delta);
           }
 
-          if (!interp_initial) {
-            fine_delta[0]->zero();
-          }
-
           fine.get_state(0)->mat_apply(fine_state, 1.0, tmat, fine_delta, false);
 
-          for (size_t m = m0; m < nfine; m++) { fine.evaluate(m); }
+          fine.reevaluate();
         }
 
 
@@ -112,7 +112,8 @@ namespace pfasst
         {
           auto& crse = as_encap_sweeper(dst);
           auto& fine = as_encap_sweeper(src);
-          this->restrict(crse.get_state(0), fine.get_state(0));
+          this->restrict(crse.get_start_state(), fine.get_start_state());
+          crse.reevaluate(true);
         }
 
 
@@ -128,17 +129,21 @@ namespace pfasst
           auto const num_crse = crse_nodes.size();
           auto const num_fine = fine_nodes.size();
 
+          if (restrict_initial) {
+            this->restrict_initial(dst, src);
+            return;
+          }
+
           int trat = (int(num_fine) - 1) / (int(num_crse) - 1);
 
-          int m0 = restrict_initial ? 0 : 1;
-          for (size_t m = m0; m < num_crse; m++) {
+          for (size_t m = 0; m < num_crse; m++) {
             if (crse_nodes[m] != fine_nodes[m * trat]) {
               throw NotImplementedYet("coarse nodes must be nested");
             }
             this->restrict(crse.get_state(m), fine.get_state(m * trat));
           }
 
-          for (size_t m = m0; m < num_crse; m++) { crse.evaluate(m); }
+          crse.reevaluate();
         }