changeset 44:4134a0f2423d

Finalise VTK implementation, biggest bugs squashed, linear wave equation works.
author Jordi Gutiérrez Hermoso <jordigh@gmail.com>
date Wed, 17 Mar 2010 01:25:28 -0600
parents 8aea477c7cf8
children 330d2813788b
files src/include/rbf.hpp src/include/vtkplot.hpp src/main.cpp src/utils.cpp src/vtkplot.cpp
diffstat 5 files changed, 115 insertions(+), 97 deletions(-) [+]
line wrap: on
line diff
--- a/src/include/rbf.hpp
+++ b/src/include/rbf.hpp
@@ -12,8 +12,8 @@
 #include "error.hpp"
 #include "func.hpp"
 
+namespace kwantix{
 
-namespace kwantix{
 ///An exception struct for RBFs when attempting to use the wrong dimension.
 struct badDimension : public badArgument{
   badDimension() {};
@@ -21,9 +21,6 @@
   size_t actual_dimension;
   size_t given_dimension;
 };
-}
-
-namespace kwantix{
 
 ///Base abstract class.
 class radial_basis_function : public kwantix::realfunc{
@@ -79,10 +76,8 @@
   static size_t dimension; 
 };
    
-}
+//Two important subclasses
 
-//Two important subclasses
-namespace kwantix{
 /// Piecewise smooth RBFs
 class piecewise_smooth_rbf : public radial_basis_function{
 public:
@@ -104,11 +99,10 @@
   ///Shape parameter
   static double eps; 
 };
-}
+
 
 
 //Specific rbf's follow...
-namespace kwantix{
 
 /// \f$r^n\f$ with \f$n\f$ odd
 class piecewise_polynomial : public piecewise_smooth_rbf{
@@ -135,9 +129,8 @@
 
 ///a common synonym.
 typedef piecewise_polynomial conical; 
-}
 
-namespace kwantix{
+
 /// \f$r^n \log(r)\f$ with \f$n\f$ even
 class thin_plate_spline : public piecewise_smooth_rbf{
 public:
@@ -158,10 +151,8 @@
   double d(double r) const;
   double d2(double r) const;
 };
-}
 
 
-namespace kwantix{
 /// \f$\sqrt{1+(\epsilon r)^2}\f$ with \f$\epsilon > 0\f$
 class multiquadric : public c_infty_rbf{
 public:
@@ -177,9 +168,8 @@
   double d(double r) const;
   double d2(double r) const;
 };
-}
 
-namespace kwantix{
+
 /// \f$1/\sqrt{1 + (\epsilon r)^2}\f$ with \f$\epsilon > 0\f$
 class inverse_multiquadric : public c_infty_rbf{
 public:
@@ -195,9 +185,8 @@
   double d(double r) const;
   double d2(double r) const;
 };
-}
 
-namespace kwantix{
+
 /// \f$1/(1+(\epsilon r)^2)\f$ with \f$\epsilon > 0 \f$
 class inverse_quadratic : public c_infty_rbf{
 public:
@@ -213,9 +202,8 @@
   double d(double r) const;
   double d2(double r) const;
 };
-}
 
-namespace kwantix{
+
 /// \f$ e^{- (\epsilon r)^2}\f$ with \f$\epsilon > 0\f$.
 class gaussian : public c_infty_rbf{
 public:
@@ -231,7 +219,7 @@
   double d(double r) const;
   double d2(double r) const;
 };
-}
 
+}//namespace kwantix
 
 #endif // __RBF_HPP__ 
--- a/src/include/vtkplot.hpp
+++ b/src/include/vtkplot.hpp
@@ -22,11 +22,11 @@
    */
   template<typename RBF>
   vtkplot(const interpolator<RBF>& u, 
+          double zmin_in, double zmax_in,
           bool offscreen_in = false);
   ///Defines the view direction when plotting.
   void set_view_direction(double x, double y, double z);
   ///Without changing the triangulation, update the values.
-  template<typename RBF>
   void update_values(const interp_values& vals_in);
   ///Enable offscreen plotting.
   void set_offscreen(bool offscreen_in = true);
@@ -35,6 +35,9 @@
   ///Actually do the plot, whether offscreen or onscreen.
   void plot() const;
 
+  ///Enable window interaction
+  void begin_interaction() const;
+
   ///The VTK callback for updating the surface points
   void AdjustPoints();
   ///Trampoline function for VTK callback
@@ -52,12 +55,18 @@
   ///Setup the pipeline to plot offscreen or not
   void InitOffscreen();
 
+  ///View direction
+  double xview, yview, zview;
+
+  ///Min and max ranges (for colour values)
+  double zmin, zmax;
+
   ///Values with which to update the surface data
-  interp_values& vals;
+  interp_values vals;
 
   /*! \name VTK pipeline endpoints
    *
-   * Endpoint
+   * Endpoint of the VTK pipeline
    */
   //@{
   vtkSmartPointer<vtkRenderer> ren1;
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -15,6 +15,7 @@
 #include "include/diff_op.hpp"
 #include "include/interpolator.hpp"
 #include "include/func.hpp"
+#include "include/vtkplot.hpp"
 
 using namespace kwantix;
 
@@ -71,16 +72,23 @@
     interpolator<RBF_TYPE> u(the_bvp);
     u.precompute_ev();
     interpolator<RBF_TYPE> u0 = u, u1 = u;
+    
+    u = 5*u;
+
+    vtkplot plotter(u,-0.1, 0.1);
+    plotter.set_view_direction(0,1,1);
+    plotter.begin_interaction();
 
     //Main loop
     size_t maxiter = 500;
     for(size_t n = 1; n <= maxiter; n++)
     {
-      save_iteration(u,n);
+      plotter.plot();
       u0 = u1;
       u1 = u;
       cout << "Iteration: " << n << endl;
       u.set_f( (u0 - 2*u1).at() );
+      plotter.update_values(u.at());
     }
 
   }
@@ -91,17 +99,3 @@
   }
 
 }
-
-void 
-save_iteration(const interpolator<RBF_TYPE>& u,
-               size_t n)
-{
-  using namespace std;
-  stringstream ss;
-  string filename = "results/u";
-  ss << filename << setw(5) << setfill('0') << n << ".map";
-  ss >> filename;
-  ofstream ofs(filename.c_str());
-  ofs << u;
-  ofs.close();
-}
--- a/src/utils.cpp
+++ b/src/utils.cpp
@@ -100,9 +100,8 @@
 void show_exception(kwantix::error exc){
   using namespace std;
 
-  cerr << exc.file << ": ""Caught an exception!" << endl;
+  cerr << "Caught an exception!" << endl;
   cerr << exc.file << ":" << exc.line << " error: "<< exc.reason << endl;
-  cerr << "From file " << exc.file << endl;
 }
 
 }
--- a/src/vtkplot.cpp
+++ b/src/vtkplot.cpp
@@ -1,4 +1,5 @@
 #include "include/vtkplot.hpp"
+#include "include/rbf.hpp"
 
 #include <vtkPoints.h>
 #include <vtkPointData.h>
@@ -50,8 +51,13 @@
 
 template<typename RBF>
 vtkplot::vtkplot(const interpolator<RBF>& u, 
+                 double zmin_in, 
+                 double zmax_in,
                  bool offscreen_in)
-  :hash(u.rbfs_hash), offscreen(offscreen)
+  : xview(1), yview(1), zview(1), 
+    zmin(zmin_in), zmax(zmax_in),
+    vals(u.at()), 
+    hash(u.rbfs_hash), basename("surf"), offscreen(offscreen_in)
 {
   if(u.thebvp -> get_domain() -> get_dimension() != 2)
   {
@@ -64,22 +70,33 @@
     throw exc;
   }
   
-  SetupPipeline(u, offscreen);
+  SetupPipeline(u);
 }
 
+//Instantiations
+template vtkplot::vtkplot(const interpolator<conical>&, double, double, bool);
+
 void vtkplot::set_offscreen(bool offscreen_in)
 {
   offscreen = offscreen_in;
   InitOffscreen();
 }
 
-template<typename RBF>
 void vtkplot::update_values(const interp_values& vals_in)
 {
   vals = vals_in;
   programmableFilter -> Modified();
 }
 
+void vtkplot::set_view_direction(double x,
+                                 double y,
+                                 double z)
+{
+  xview = x; yview = y; zview = z;
+  ren1->GetActiveCamera()->SetPosition(xview, yview, zview);
+  ren1->ResetCamera();
+}
+
 void vtkplot::trampoline(void* this_ptr)
 {
   vtkplot* instance = static_cast<vtkplot*>(this_ptr);
@@ -93,19 +110,22 @@
 	vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
   vtkSmartPointer<vtkDoubleArray> scalars = vtkSmartPointer<vtkDoubleArray>::New();
 
+  size_t idx = 1;
   for(auto i = u.thebvp -> get_domain() -> get_interior().begin();
       i != u.thebvp -> get_domain() -> get_interior().end(); i++)
   {
-    double z = 0;
-    vtkIdType idx = points -> InsertNextPoint( *i(1), *i(2), 0);
+    double z = vals(idx);
+    vtkIdType idx = points -> InsertNextPoint( (*i)(1), (*i)(2), 0);
     scalars -> InsertTuple(idx, &z);
+    idx++;
   }
   for(auto i = u.thebvp -> get_domain() -> get_boundary().begin();
       i != u.thebvp -> get_domain() -> get_boundary().end(); i++)
   {
-    double z = 0;
-    vtkIdType idx = points -> InsertNextPoint( *i(1), *i(2), 0);
+    double z = vals(idx);
+    vtkIdType idx = points -> InsertNextPoint( (*i)(1), (*i)(2), 0);
     scalars -> InsertTuple(idx, &z);
+    idx++;
   }
  
 	//add the grid points to a polydata object
@@ -144,7 +164,7 @@
 	vtkSmartPointer<vtkPolyDataMapper> contMapper 
     = vtkSmartPointer<vtkPolyDataMapper>::New();
   contMapper->SetInput(normals -> GetOutput() );
-  contMapper->SetScalarRange(0,1);
+  contMapper->SetScalarRange(zmin,zmax);
   contMapper->SetLookupTable(lut);
   contMapper->ImmediateModeRenderingOn();
 
@@ -163,6 +183,55 @@
   }
 	renWin->AddRenderer(ren1);
 
+	// add the actors to the scene
+	ren1->SetBackground(0.0,0.0,0.0);
+  ren1->GetActiveCamera()->SetViewUp(0,0,1);
+  ren1->GetActiveCamera()->SetPosition(xview, yview, zview);
+  ren1->GetActiveCamera()->SetParallelProjection(1);
+	ren1->AddActor(contActor);
+  ren1->ResetCamera();
+ 
+}
+
+//Instantiations
+template void vtkplot::SetupPipeline(const interpolator<conical>&);
+
+void vtkplot::plot() const
+{
+  if(offscreen)
+  {
+    static size_t i = 0;
+
+    using namespace std;
+    string filename = basename;
+    stringstream ss; string i_str;
+    ss << setw(5) << setfill('0')  << i;
+    ss >> i_str;
+    filename += i_str + ".png";
+   
+    // Write to file
+    vtkSmartPointer<vtkWindowToImageFilter> windowToImageFilter = 
+      vtkSmartPointer<vtkWindowToImageFilter>::New();
+    windowToImageFilter->SetInput(renWin);
+    windowToImageFilter->Update();
+ 
+    vtkSmartPointer<vtkPNGWriter> writer 
+      = vtkSmartPointer<vtkPNGWriter>::New();
+    writer->SetFileName(filename.c_str() );
+    writer->SetInput(windowToImageFilter->GetOutput());
+    writer->Write();
+   
+    i++;
+  }
+  else
+  {
+    // render an image (lights and cameras are created automatically)
+    renWin->Render();
+  }
+}
+
+void vtkplot::begin_interaction() const
+{
   // an interactor
   vtkSmartPointer<vtkRenderWindowInteractor> iren 
     = vtkSmartPointer<vtkRenderWindowInteractor>::New();
@@ -170,57 +239,16 @@
 
   // Initialize must be called prior to creating timer events.
   iren->Initialize();
-  iren->CreateRepeatingTimer(10);
- 
-	// add the actors to the scene
-	ren1->SetBackground(0.0,0.0,0.0);
-  ren1->GetActiveCamera()->SetViewUp(0,0,1);
-  ren1->GetActiveCamera()->SetPosition(1, 1, 1);
-  ren1->GetActiveCamera()->SetParallelProjection(1);
-	ren1->AddActor(contActor);
-  ren1->ResetCamera();
- 
-  if(offscreen)
-  {
-    for(size_t i = 0; i < 1000; i++)
-    {
-      using namespace std;
-      string filename = "foo";
-      stringstream ss; string i_str;
-      ss << setw(4) << setfill('0')  << i;
-      ss >> i_str;
-      filename += i_str + ".png";
-   
-      // Write to file
-      vtkSmartPointer<vtkWindowToImageFilter> windowToImageFilter = 
-        vtkSmartPointer<vtkWindowToImageFilter>::New();
-      windowToImageFilter->SetInput(renWin);
-      windowToImageFilter->Update();
- 
-      vtkSmartPointer<vtkPNGWriter> writer 
-        = vtkSmartPointer<vtkPNGWriter>::New();
-      writer->SetFileName(filename.c_str() );
-      writer->SetInput(windowToImageFilter->GetOutput());
-      writer->Write();
 
-      programmableFilter -> Modified();
-    }
-  }
-  else
-  {
+  //The style
+  vtkSmartPointer<vtkInteractorStyleTerrain> terrain_style 
+    = vtkSmartPointer<vtkInteractorStyleTerrain>::New();
+   iren -> SetInteractorStyle(terrain_style);
 
-    // render an image (lights and cameras are created automatically)
-    renWin->Render();
-    //The style
-    vtkSmartPointer<vtkInteractorStyleTerrain> terrain_style 
-      = vtkSmartPointer<vtkInteractorStyleTerrain>::New();
- 
-    // Set the style
-    iren -> SetInteractorStyle(terrain_style);
-  
-    // begin mouse interaction
-    iren->Start(); 
-  }
+
+  renWin->Render();  
+  // begin mouse interaction
+  iren->Start(); 
 }
 
 void vtkplot::InitOffscreen()