The source file meshless_interpolation_function.h with comments:

        #ifndef LIBMESH_MESHLESS_INTERPOLATION_FUNCTION_H
        #define LIBMESH_MESHLESS_INTERPOLATION_FUNCTION_H
        
Local Includes
        #include "libmesh/function_base.h"
        #include "libmesh/meshfree_interpolation.h"
        #include "libmesh/threads.h"
        
        
C++ includes
        #include <cstddef>
        
        namespace libMesh
        {
        
        
        
Forward Declarations
        template <typename T>
        class DenseVector;
        
        
------------------------------------------------------------ MeshlessInterpolationFunction class definition
        class MeshlessInterpolationFunction : public FunctionBase<Number>
        {
        private:
          const MeshfreeInterpolation &_mfi;
          mutable std::vector<Point> _pts;
          mutable std::vector<Number> _vals;
          Threads::spin_mutex &_mutex;
        
        public:
        
          /**
           * Constructor.  Requires a \p \pMeshlessInterpolation object.
           */
          MeshlessInterpolationFunction (const MeshfreeInterpolation &mfi,
        				 Threads::spin_mutex &mutex) :
          _mfi (mfi),
          _mutex(mutex)
          {}
        
        
          /**
           * The actual initialization process.
           */
          void init ();
        
          /**
           * Clears the function.
           */
          void clear ();
        
          /**
           * Returns a new deep copy of the function.
           */
          virtual AutoPtr<FunctionBase<Number> > clone () const;
        
          /**
           * @returns the value at point \p p and time
           * \p time, which defaults to zero.
           */
          Number operator() (const Point& p,
        		     const Real time=0.);
        
          /**
           * Like before, but returns the values in a
           * writable reference.
           */
          void operator() (const Point& p,
        		   const Real time,
        		   DenseVector<Number>& output);
        
        };
        
        
        
------------------------------------------------------------ MeshlessInterpolationFunction inline methods
        inline
        Number MeshlessInterpolationFunction::operator() (const Point& p,
        						  const Real /* time */)
        {
          _pts.clear();
          _pts.push_back(p);
          _vals.resize(1);
        
          Threads::spin_mutex::scoped_lock lock(_mutex);
        
          _mfi.interpolate_field_data (_mfi.field_variables(),
        			       _pts, _vals);
        
          return _vals.front();
        }
        
        
        
        inline
        void MeshlessInterpolationFunction::operator() (const Point& p,
        						const Real time,
        						DenseVector<Number>& output)
        {
          output.resize(1);
          output(0) = (*this)(p,time);
          return;
        }
        
        
        
        inline
        void MeshlessInterpolationFunction::init ()
        {
        }
        
        
        
        inline
        void MeshlessInterpolationFunction::clear ()
        {
        }
        
        
        
        inline
        AutoPtr<FunctionBase<Number> >
        MeshlessInterpolationFunction::clone () const
        {
          return AutoPtr<FunctionBase<Number> > (new MeshlessInterpolationFunction (_mfi, _mutex) );
        }
        
        
        } // namespace libMesh
        
        
        #endif // LIBMESH_MESHLESS_INTERPOLATION_FUNCTION_H



The source file miscellaneous_ex8.C with comments:

        #include "libmesh/libmesh.h"
        #include "libmesh/meshfree_interpolation.h"
        #include "libmesh/radial_basis_interpolation.h"
        #include "libmesh/mesh.h"
        #include "libmesh/equation_systems.h"
        #include "libmesh/numeric_vector.h"
        #include "libmesh/tecplot_io.h"
        #include "libmesh/threads.h"
        #include "meshless_interpolation_function.h"
        
C++ includes
        #include <cstdlib>
        
        
Bring in everything from the libMesh namespace
        using namespace libMesh;
        
        
        void create_random_point_cloud (const unsigned int Npts,
        				std::vector<Point> &pts,
        				const Real max_range = 10)
        {
          std::cout << "Generating "<< Npts << " point cloud...";
          pts.resize(Npts);
        
          for (size_t i=0;i<Npts;i++)
            {
              pts[i](0) = max_range * (std::rand() % 1000) / Real(1000);
              pts[i](1) = max_range * (std::rand() % 1000) / Real(1000);
              pts[i](2) = max_range * (std::rand() % 1000) / Real(1000);
            }
          std::cout << "done\n";
        }
        
        
        
        Real exact_solution_u (const Point &p)
        {
          const Real
            x = p(0),
            y = p(1),
            z = p(2);
        
          return (x*x*x   +
        	  y*y*y*y +
        	  z*z*z*z*z);
        }
        
        
        
        Real exact_solution_v (const Point &p)
        {
          const Real
            x = p(0),
            y = p(1),
            z = p(2);
        
          return (x*x   +
        	  y*y +
        	  z*z*z);
        }
        
        Number exact_value (const Point& p,
                            const Parameters&,
                            const std::string&,
                            const std::string&)
        {
          return exact_solution_v(p);
        }
        
We now define the function which provides the initialization routines for the "Convection-Diffusion" system. This handles things like setting initial conditions and boundary conditions.
        void init_sys(EquationSystems& es,
                      const std::string& system_name)
        {
Get a reference to the Convection-Diffusion system object.
          System & system =
            es.get_system<System>(system_name);
        
          system.project_solution(exact_value, NULL, es.parameters);
        }
        
        
        
        
        int main(int argc, char** argv)
        {
Skip this example if we do not meet certain requirements
          libmesh_example_assert(3 <= LIBMESH_DIM, "3D support");
        #ifndef LIBMESH_HAVE_EIGEN
          libmesh_example_assert(false, "--enable-eigen");
        #endif
        #ifndef LIBMESH_HAVE_ZLIB_H
          libmesh_example_assert(false, "--enable-zlib");
        #endif
        
Initialize libMesh.
          LibMeshInit init (argc, argv);
          {
Demonstration case 1
            {
              std::vector<Point>       tgt_pts;
              std::vector<Number>      tgt_data_idi, tgt_data_rbi;
              std::vector<std::string> field_vars;
        
              field_vars.push_back("u");
              field_vars.push_back("v");
        
              InverseDistanceInterpolation<3> idi (init.comm(),
        					   /* n_interp_pts = */ 8,
        					   /* power =        */ 2);
        
              RadialBasisInterpolation<3> rbi (init.comm());
        
              idi.set_field_variables (field_vars);
              rbi.set_field_variables (field_vars);
        
              create_random_point_cloud (100,
        				 idi.get_source_points());
        
        
Explicitly set the data values we will interpolate from
              {
        	const std::vector<Point> &src_pts  (idi.get_source_points());
        	std::vector<Number>      &src_vals (idi.get_source_vals());
        
        	src_vals.clear(); src_vals.reserve(2*src_pts.size());
        
        	for (std::vector<Point>::const_iterator pt_it=src_pts.begin();
        	     pt_it != src_pts.end(); ++pt_it)
        	  {
        	    src_vals.push_back (exact_solution_u (*pt_it));
        	    src_vals.push_back (exact_solution_v (*pt_it));
        	  }
              }
        
give rbi the same info as idi
              rbi.get_source_points() = idi.get_source_points();
              rbi.get_source_vals()   = idi.get_source_vals();
        
              idi.prepare_for_use();
              rbi.prepare_for_use();
        
              std::cout << idi;
        
Interpolate to some other random points, and evaluate the result
              {
        	create_random_point_cloud (10,
        				   tgt_pts);
        
tgt_pts = rbi.get_source_points();

                idi.interpolate_field_data (field_vars,
        				    tgt_pts,
        				    tgt_data_idi);
        
        	rbi.interpolate_field_data (field_vars,
        				    tgt_pts,
        				    tgt_data_rbi);
        
              	std::vector<Number>::const_iterator
        	  v_idi=tgt_data_idi.begin(),
        	  v_rbi=tgt_data_rbi.begin();
        
        	for (std::vector<Point>::const_iterator  p_it=tgt_pts.begin();
        	     p_it!=tgt_pts.end(); ++p_it)
        	  {
        	    std::cout << "\nAt target point " << *p_it
        		      << "\n u_interp_idi="   << *v_idi
        		      << ", u_interp_rbi="    << *v_rbi
        		      << ", u_exact="         << exact_solution_u(*p_it);
        	    ++v_idi;
        	    ++v_rbi;
        	    std::cout << "\n v_interp_idi=" << *v_idi
        		      << ", v_interp_rbi="  << *v_rbi
        		      << ", v_exact="       << exact_solution_v(*p_it)
        		      << std::endl;
        	    ++v_idi;
        	    ++v_rbi;
        	  }
              }
            }
        
        
Demonstration case 2
            {
              Mesh mesh_a(init.comm()), mesh_b(init.comm());
        
              mesh_a.read("struct.ucd.gz"); mesh_b.read("unstruct.ucd.gz");
        
Create equation systems objects.
              EquationSystems
        	es_a(mesh_a), es_b(mesh_b);
        
              System
        	&sys_a = es_a.add_system<System>("src_system"),
        	&sys_b = es_b.add_system<System>("dest_system");
        
              sys_a.add_variable ("Cp", FIRST);
              sys_b.add_variable ("Cp", FIRST);
        
              sys_a.attach_init_function (init_sys);
              es_a.init();
        
Write out the initial conditions.
              TecplotIO(mesh_a).write_equation_systems ("src.dat",
        						es_a);
        
              InverseDistanceInterpolation<3> idi (init.comm(),
        					   /* n_interp_pts = */ 4,
        					   /* power =        */ 2);
              RadialBasisInterpolation<3> rbi (init.comm());
        
              std::vector<Point>  &src_pts  (idi.get_source_points());
              std::vector<Number> &src_vals (idi.get_source_vals());
              std::vector<std::string> field_vars;
              field_vars.push_back("Cp");
              idi.set_field_variables(field_vars);
        
We now will loop over every node in the source mesh and add it to a source point list, along with the solution
              {
        	MeshBase::const_node_iterator nd  = mesh_a.local_nodes_begin();
        	MeshBase::const_node_iterator end = mesh_a.local_nodes_end();
        
        	for (; nd!=end; ++nd)
        	  {
        	    const Node *node(*nd);
        	    src_pts.push_back(*node);
        	    src_vals.push_back(sys_a.current_solution(node->dof_number(0,0,0)));
        	  }
        
        	rbi.set_field_variables(field_vars);
        	rbi.get_source_points() = idi.get_source_points();
        	rbi.get_source_vals()   = idi.get_source_vals();
              }
        
We have only set local values - prepare for use by gathering remote gata
              idi.prepare_for_use();
              rbi.prepare_for_use();
        
Create a MeshlessInterpolationFunction that uses our InverseDistanceInterpolation object. Since each MeshlessInterpolationFunction shares the same InverseDistanceInterpolation object in a threaded environment we must also provide a locking mechanism.
              {
        	Threads::spin_mutex mutex;
        	MeshlessInterpolationFunction mif(idi, mutex);
        
project the solution onto system b
                es_b.init();
        	sys_b.project_solution (&mif);
        
Write the result
                TecplotIO(mesh_b).write_equation_systems ("dest_idi.dat",
        						  es_b);
              }
        
Create a MeshlessInterpolationFunction that uses our RadialBasisInterpolation object. Since each MeshlessInterpolationFunction shares the same RadialBasisInterpolation object in a threaded environment we must also provide a locking mechanism.
              {
        	Threads::spin_mutex mutex;
        	MeshlessInterpolationFunction mif(rbi, mutex);
        
project the solution onto system b
                sys_b.project_solution (&mif);
        
Write the result
                TecplotIO(mesh_b).write_equation_systems ("dest_rbi.dat",
        						  es_b);
              }
            }
        
        
        
          }
          return 0;
        }



The source file meshless_interpolation_function.h without comments:

 
  #ifndef LIBMESH_MESHLESS_INTERPOLATION_FUNCTION_H
  #define LIBMESH_MESHLESS_INTERPOLATION_FUNCTION_H
  
  #include "libmesh/function_base.h"
  #include "libmesh/meshfree_interpolation.h"
  #include "libmesh/threads.h"
  
  
  #include <cstddef>
  
  namespace libMesh
  {
  
  
  
  template <typename T>
  class DenseVector;
  
  
  class MeshlessInterpolationFunction : public FunctionBase<Number>
  {
  private:
    const MeshfreeInterpolation &_mfi;
    mutable std::vector<Point> _pts;
    mutable std::vector<Number> _vals;
    Threads::spin_mutex &_mutex;
  
  public:
  
    /**
     * Constructor.  Requires a \p \pMeshlessInterpolation object.
     */
    MeshlessInterpolationFunction (const MeshfreeInterpolation &mfi,
  				 Threads::spin_mutex &mutex) :
    _mfi (mfi),
    _mutex(mutex)
    {}
  
  
    /**
     * The actual initialization process.
     */
    void init ();
  
    /**
     * Clears the function.
     */
    void clear ();
  
    /**
     * Returns a new deep copy of the function.
     */
    virtual AutoPtr<FunctionBase<Number> > clone () const;
  
    /**
     * @returns the value at point \p p and time
     * \p time, which defaults to zero.
     */
    Number operator() (const Point& p,
  		     const Real time=0.);
  
    /**
     * Like before, but returns the values in a
     * writable reference.
     */
    void operator() (const Point& p,
  		   const Real time,
  		   DenseVector<Number>& output);
  
  };
  
  
  
  inline
  Number MeshlessInterpolationFunction::operator() (const Point& p,
  						  const Real /* time */)
  {
    _pts.clear();
    _pts.push_back(p);
    _vals.resize(1);
  
    Threads::spin_mutex::scoped_lock lock(_mutex);
  
    _mfi.interpolate_field_data (_mfi.field_variables(),
  			       _pts, _vals);
  
    return _vals.front();
  }
  
  
  
  inline
  void MeshlessInterpolationFunction::operator() (const Point& p,
  						const Real time,
  						DenseVector<Number>& output)
  {
    output.resize(1);
    output(0) = (*this)(p,time);
    return;
  }
  
  
  
  inline
  void MeshlessInterpolationFunction::init ()
  {
  }
  
  
  
  inline
  void MeshlessInterpolationFunction::clear ()
  {
  }
  
  
  
  inline
  AutoPtr<FunctionBase<Number> >
  MeshlessInterpolationFunction::clone () const
  {
    return AutoPtr<FunctionBase<Number> > (new MeshlessInterpolationFunction (_mfi, _mutex) );
  }
  
  
  } // namespace libMesh
  
  
  #endif // LIBMESH_MESHLESS_INTERPOLATION_FUNCTION_H



The source file miscellaneous_ex8.C without comments:

 
  #include "libmesh/libmesh.h"
  #include "libmesh/meshfree_interpolation.h"
  #include "libmesh/radial_basis_interpolation.h"
  #include "libmesh/mesh.h"
  #include "libmesh/equation_systems.h"
  #include "libmesh/numeric_vector.h"
  #include "libmesh/tecplot_io.h"
  #include "libmesh/threads.h"
  #include "meshless_interpolation_function.h"
  
  #include <cstdlib>
  
  
  using namespace libMesh;
  
  
  void create_random_point_cloud (const unsigned int Npts,
  				std::vector<Point> &pts,
  				const Real max_range = 10)
  {
    std::cout << "Generating "<< Npts << " point cloud...";
    pts.resize(Npts);
  
    for (size_t i=0;i<Npts;i++)
      {
        pts[i](0) = max_range * (std::rand() % 1000) / Real(1000);
        pts[i](1) = max_range * (std::rand() % 1000) / Real(1000);
        pts[i](2) = max_range * (std::rand() % 1000) / Real(1000);
      }
    std::cout << "done\n";
  }
  
  
  
  Real exact_solution_u (const Point &p)
  {
    const Real
      x = p(0),
      y = p(1),
      z = p(2);
  
    return (x*x*x   +
  	  y*y*y*y +
  	  z*z*z*z*z);
  }
  
  
  
  Real exact_solution_v (const Point &p)
  {
    const Real
      x = p(0),
      y = p(1),
      z = p(2);
  
    return (x*x   +
  	  y*y +
  	  z*z*z);
  }
  
  Number exact_value (const Point& p,
                      const Parameters&,
                      const std::string&,
                      const std::string&)
  {
    return exact_solution_v(p);
  }
  
  void init_sys(EquationSystems& es,
                const std::string& system_name)
  {
    System & system =
      es.get_system<System>(system_name);
  
    system.project_solution(exact_value, NULL, es.parameters);
  }
  
  
  
  
  int main(int argc, char** argv)
  {
    libmesh_example_assert(3 <= LIBMESH_DIM, "3D support");
  #ifndef LIBMESH_HAVE_EIGEN
    libmesh_example_assert(false, "--enable-eigen");
  #endif
  #ifndef LIBMESH_HAVE_ZLIB_H
    libmesh_example_assert(false, "--enable-zlib");
  #endif
  
    LibMeshInit init (argc, argv);
    {
      {
        std::vector<Point>       tgt_pts;
        std::vector<Number>      tgt_data_idi, tgt_data_rbi;
        std::vector<std::string> field_vars;
  
        field_vars.push_back("u");
        field_vars.push_back("v");
  
        InverseDistanceInterpolation<3> idi (init.comm(),
  					   /* n_interp_pts = */ 8,
  					   /* power =        */ 2);
  
        RadialBasisInterpolation<3> rbi (init.comm());
  
        idi.set_field_variables (field_vars);
        rbi.set_field_variables (field_vars);
  
        create_random_point_cloud (100,
  				 idi.get_source_points());
  
  
        {
  	const std::vector<Point> &src_pts  (idi.get_source_points());
  	std::vector<Number>      &src_vals (idi.get_source_vals());
  
  	src_vals.clear(); src_vals.reserve(2*src_pts.size());
  
  	for (std::vector<Point>::const_iterator pt_it=src_pts.begin();
  	     pt_it != src_pts.end(); ++pt_it)
  	  {
  	    src_vals.push_back (exact_solution_u (*pt_it));
  	    src_vals.push_back (exact_solution_v (*pt_it));
  	  }
        }
  
        rbi.get_source_points() = idi.get_source_points();
        rbi.get_source_vals()   = idi.get_source_vals();
  
        idi.prepare_for_use();
        rbi.prepare_for_use();
  
        std::cout << idi;
  
        {
  	create_random_point_cloud (10,
  				   tgt_pts);
  
  
  	idi.interpolate_field_data (field_vars,
  				    tgt_pts,
  				    tgt_data_idi);
  
  	rbi.interpolate_field_data (field_vars,
  				    tgt_pts,
  				    tgt_data_rbi);
  
        	std::vector<Number>::const_iterator
  	  v_idi=tgt_data_idi.begin(),
  	  v_rbi=tgt_data_rbi.begin();
  
  	for (std::vector<Point>::const_iterator  p_it=tgt_pts.begin();
  	     p_it!=tgt_pts.end(); ++p_it)
  	  {
  	    std::cout << "\nAt target point " << *p_it
  		      << "\n u_interp_idi="   << *v_idi
  		      << ", u_interp_rbi="    << *v_rbi
  		      << ", u_exact="         << exact_solution_u(*p_it);
  	    ++v_idi;
  	    ++v_rbi;
  	    std::cout << "\n v_interp_idi=" << *v_idi
  		      << ", v_interp_rbi="  << *v_rbi
  		      << ", v_exact="       << exact_solution_v(*p_it)
  		      << std::endl;
  	    ++v_idi;
  	    ++v_rbi;
  	  }
        }
      }
  
  
      {
        Mesh mesh_a(init.comm()), mesh_b(init.comm());
  
        mesh_a.read("struct.ucd.gz"); mesh_b.read("unstruct.ucd.gz");
  
        EquationSystems
  	es_a(mesh_a), es_b(mesh_b);
  
        System
  	&sys_a = es_a.add_system<System>("src_system"),
  	&sys_b = es_b.add_system<System>("dest_system");
  
        sys_a.add_variable ("Cp", FIRST);
        sys_b.add_variable ("Cp", FIRST);
  
        sys_a.attach_init_function (init_sys);
        es_a.init();
  
        TecplotIO(mesh_a).write_equation_systems ("src.dat",
  						es_a);
  
        InverseDistanceInterpolation<3> idi (init.comm(),
  					   /* n_interp_pts = */ 4,
  					   /* power =        */ 2);
        RadialBasisInterpolation<3> rbi (init.comm());
  
        std::vector<Point>  &src_pts  (idi.get_source_points());
        std::vector<Number> &src_vals (idi.get_source_vals());
        std::vector<std::string> field_vars;
        field_vars.push_back("Cp");
        idi.set_field_variables(field_vars);
  
        {
  	MeshBase::const_node_iterator nd  = mesh_a.local_nodes_begin();
  	MeshBase::const_node_iterator end = mesh_a.local_nodes_end();
  
  	for (; nd!=end; ++nd)
  	  {
  	    const Node *node(*nd);
  	    src_pts.push_back(*node);
  	    src_vals.push_back(sys_a.current_solution(node->dof_number(0,0,0)));
  	  }
  
  	rbi.set_field_variables(field_vars);
  	rbi.get_source_points() = idi.get_source_points();
  	rbi.get_source_vals()   = idi.get_source_vals();
        }
  
        idi.prepare_for_use();
        rbi.prepare_for_use();
  
        {
  	Threads::spin_mutex mutex;
  	MeshlessInterpolationFunction mif(idi, mutex);
  
  	es_b.init();
  	sys_b.project_solution (&mif);
  
  	TecplotIO(mesh_b).write_equation_systems ("dest_idi.dat",
  						  es_b);
        }
  
        {
  	Threads::spin_mutex mutex;
  	MeshlessInterpolationFunction mif(rbi, mutex);
  
  	sys_b.project_solution (&mif);
  
  	TecplotIO(mesh_b).write_equation_systems ("dest_rbi.dat",
  						  es_b);
        }
      }
  
  
  
    }
    return 0;
  }



The console output of the program:

make[4]: Entering directory `/net/spark/workspace/roystgnr/libmesh/git/devel/examples/miscellaneous/miscellaneous_ex8'
***************************************************************
* Running Example miscellaneous_ex8:
*  mpirun -np 4 example-devel  -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc
***************************************************************
 
*** Warning, This code is untested, experimental, or likely to see future API changes: ../../../include/libmesh/radial_basis_interpolation.h, line 79, compiled Apr 19 2013 at 11:51:51 ***
Generating 100 point cloud...done
bounding box is 
(x,y,z)=(    0.34,     0.19,     0.11)
(x,y,z)=(    9.65,     9.96,     9.96)
*** Warning, This code is untested, experimental, or likely to see future API changes: ./include/libmesh/radial_basis_functions.h, line 83, compiled Apr 19 2013 at 11:42:53 ***
bounding box is 
(x,y,z)=(    0.34,     0.19,     0.11)
(x,y,z)=(    9.65,     9.96,     9.96)
r_bbox = 16.7078
rbf(r_bbox/2) = 0.1875
MeshfreeInterpolation
 n_source_points()=400
 n_field_variables()=2
  variables = u v 
Generating 10 point cloud...done
*** Warning, This code is untested, experimental, or likely to see future API changes: ../src/solution_transfer/meshfree_interpolation.C, line 283, compiled Apr 19 2013 at 11:42:51 ***
*** Warning, This code is untested, experimental, or likely to see future API changes: ../src/solution_transfer/radial_basis_interpolation.C, line 158, compiled Apr 19 2013 at 11:42:53 ***

At target point (x,y,z)=(     0.9,     6.84,     3.76)
 u_interp_idi=2931.24, u_interp_rbi=3006.88, u_exact=2941.14
 v_interp_idi=100.789, v_interp_rbi=101.434, v_exact=100.753

At target point (x,y,z)=(    5.42,     9.36,     1.07)
 u_interp_idi=6246.44, u_interp_rbi=6877.55, u_exact=7836.06
 v_interp_idi=115.044, v_interp_rbi=113.128, v_exact=118.211

At target point (x,y,z)=(    4.45,     7.56,     1.79)
 u_interp_idi=2958.51, u_interp_rbi=3427.22, u_exact=3373.03
 v_interp_idi=70.7865, v_interp_rbi=81.8618, v_exact=82.6914

At target point (x,y,z)=(    4.18,     8.87,     4.12)
 u_interp_idi=4632.05, u_interp_rbi=6651.05, u_exact=7450.19
 v_interp_idi=156.888, v_interp_rbi=162.735, v_exact=166.084

At target point (x,y,z)=(    3.48,     1.72,     6.59)
 u_interp_idi=19625.3, u_interp_rbi=10933.5, u_exact=12479.6
 v_interp_idi=356.174, v_interp_rbi=293.021, v_exact=301.26

At target point (x,y,z)=(    0.09,     3.36,      2.1)
 u_interp_idi=498.469, u_interp_rbi=806.722, u_exact=168.297
 v_interp_idi=39.6305, v_interp_rbi=21.7629, v_exact=20.5587

At target point (x,y,z)=(    3.42,     5.87,     2.06)
 u_interp_idi=971.312, u_interp_rbi=1072.91, u_exact=1264.38
 v_interp_idi=48.0353, v_interp_rbi=53.0049, v_exact=54.8951

At target point (x,y,z)=(    3.01,     7.13,     3.72)
 u_interp_idi=4955.74, u_interp_rbi=3311.65, u_exact=3324.05
 v_interp_idi=162.394, v_interp_rbi=110.872, v_exact=111.376

At target point (x,y,z)=(    3.21,     2.55,     8.19)
 u_interp_idi=26700.4, u_interp_rbi=37098, u_exact=36923.8
 v_interp_idi=466.577, v_interp_rbi=565.731, v_exact=566.16

At target point (x,y,z)=(    5.99,     7.21,     9.04)
 u_interp_idi=53370.1, u_interp_rbi=62518.5, u_exact=63290.2
 v_interp_idi=755.935, v_interp_rbi=822.984, v_exact=826.627
bounding box is 
(x,y,z)=(      -1,        0,        0)
(x,y,z)=(       0,        1,        1)
bounding box is 
(x,y,z)=(      -1,        0,        0)
(x,y,z)=(       0,        1,        1)
r_bbox = 1.73205
rbf(r_bbox/2) = 0.1875

 -------------------------------------------------------------------------------------------------------------------
| Processor id:   0                                                                                                 |
| Num Processors: 4                                                                                                 |
| Time:           Fri Apr 19 11:52:17 2013                                                                          |
| OS:             Linux                                                                                             |
| HostName:       spark.ices.utexas.edu                                                                             |
| OS Release:     2.6.32-279.22.1.el6.x86_64                                                                        |
| OS Version:     #1 SMP Tue Feb 5 14:33:39 CST 2013                                                                |
| Machine:        x86_64                                                                                            |
| Username:       roystgnr                                                                                          |
| Configuration:  ../configure  '--enable-everything'                                                               |
|  'METHODS=devel'                                                                                                  |
|  '--prefix=/h2/roystgnr/libmesh-test'                                                                             |
|  'CXX=distcc /usr/bin/g++'                                                                                        |
|  'CC=distcc /usr/bin/gcc'                                                                                         |
|  'FC=distcc /usr/bin/gfortran'                                                                                    |
|  'F77=distcc /usr/bin/gfortran'                                                                                   |
|  'PETSC_DIR=/opt/apps/ossw/libraries/petsc/petsc-3.3-p2'                                                          |
|  'PETSC_ARCH=gcc-system-mkl-gf-10.3.12.361-mpich2-1.4.1p1-cxx-opt'                                                |
|  'SLEPC_DIR=/opt/apps/ossw/libraries/slepc/slepc-3.3-p2-petsc-3.3-p2-cxx-opt'                                     |
|  'TRILINOS_DIR=/opt/apps/ossw/libraries/trilinos/trilinos-10.12.2/sl6/gcc-system/mpich2-1.4.1p1/mkl-gf-10.3.12.361'|
|  'VTK_DIR=/opt/apps/ossw/libraries/vtk/vtk-5.10.0/sl6/gcc-system'                                                 |
|  'HDF5_DIR=/opt/apps/ossw/libraries/hdf5/hdf5-1.8.9/sl6/gcc-system'                                               |
 -------------------------------------------------------------------------------------------------------------------
 ------------------------------------------------------------------------------------------------------------
| libMesh Performance: Alive time=17.3246, Active time=17.2494                                               |
 ------------------------------------------------------------------------------------------------------------
| Event                          nCalls    Total Time  Avg Time    Total Time  Avg Time    % of Active Time  |
|                                          w/o Sub     w/o Sub     With Sub    With Sub    w/o S    With S   |
|------------------------------------------------------------------------------------------------------------|
|                                                                                                            |
|                                                                                                            |
| DofMap                                                                                                     |
|   add_neighbors_to_send_list() 2         0.0080      0.004017    0.0091      0.004570    0.05     0.05     |
|   create_dof_constraints()     2         0.0033      0.001668    0.0033      0.001668    0.02     0.02     |
|   distribute_dofs()            2         0.0157      0.007827    0.0465      0.023253    0.09     0.27     |
|   dof_indices()                7528      0.0272      0.000004    0.0272      0.000004    0.16     0.16     |
|   prepare_send_list()          2         0.0000      0.000013    0.0000      0.000013    0.00     0.00     |
|   reinit()                     2         0.0269      0.013447    0.0269      0.013447    0.16     0.16     |
|                                                                                                            |
| EquationSystems                                                                                            |
|   build_solution_vector()      3         0.0086      0.002855    0.0250      0.008318    0.05     0.14     |
|                                                                                                            |
| InverseDistanceInterpolation<>                                                                             |
|   construct_kd_tree()          4         0.0064      0.001610    0.0064      0.001610    0.04     0.04     |
|   interpolate_field_data()     4219      0.0128      0.000003    0.0128      0.000003    0.07     0.07     |
|                                                                                                            |
| Mesh                                                                                                       |
|   find_neighbors()             2         0.0209      0.010456    0.0345      0.017260    0.12     0.20     |
|   read()                       2         0.0407      0.020368    0.0407      0.020368    0.24     0.24     |
|   renumber_nodes_and_elem()    4         0.0020      0.000507    0.0020      0.000507    0.01     0.01     |
|                                                                                                            |
| MeshCommunication                                                                                          |
|   broadcast()                  2         0.0106      0.005299    0.0252      0.012591    0.06     0.15     |
|   compute_hilbert_indices()    4         0.0723      0.018084    0.0723      0.018084    0.42     0.42     |
|   find_global_indices()        4         0.0092      0.002311    0.1095      0.027365    0.05     0.63     |
|   parallel_sort()              4         0.0017      0.000415    0.0267      0.006686    0.01     0.16     |
|                                                                                                            |
| MeshOutput                                                                                                 |
|   write_equation_systems()     3         0.0001      0.000040    0.0997      0.033224    0.00     0.58     |
|                                                                                                            |
| MeshfreeInterpolation                                                                                      |
|   gather_remote_data()         4         0.0003      0.000085    0.0009      0.000237    0.00     0.01     |
|                                                                                                            |
| MetisPartitioner                                                                                           |
|   partition()                  2         0.0937      0.046840    0.1419      0.070953    0.54     0.82     |
|                                                                                                            |
| Parallel                                                                                                   |
|   allgather()                  20        0.0030      0.000150    0.0034      0.000170    0.02     0.02     |
|   broadcast()                  48        0.0044      0.000091    0.0043      0.000090    0.03     0.03     |
|   max(scalar)                  243       0.0045      0.000018    0.0045      0.000018    0.03     0.03     |
|   max(vector)                  56        0.0009      0.000016    0.0032      0.000058    0.01     0.02     |
|   min(bool)                    289       0.0037      0.000013    0.0037      0.000013    0.02     0.02     |
|   min(scalar)                  237       0.9160      0.003865    0.9160      0.003865    5.31     5.31     |
|   min(vector)                  56        0.0010      0.000018    0.0058      0.000104    0.01     0.03     |
|   probe()                      72        0.0024      0.000033    0.0024      0.000033    0.01     0.01     |
|   receive()                    72        0.0004      0.000005    0.0028      0.000038    0.00     0.02     |
|   send()                       72        0.0003      0.000004    0.0003      0.000004    0.00     0.00     |
|   send_receive()               68        0.0003      0.000004    0.0028      0.000041    0.00     0.02     |
|   sum()                        24        0.0254      0.001059    0.0284      0.001181    0.15     0.16     |
|                                                                                                            |
| Parallel::Request                                                                                          |
|   wait()                       72        0.0001      0.000001    0.0001      0.000001    0.00     0.00     |
|                                                                                                            |
| Partitioner                                                                                                |
|   set_node_processor_ids()     2         0.0054      0.002697    0.0215      0.010748    0.03     0.12     |
|   set_parent_processor_ids()   2         0.0022      0.001107    0.0022      0.001107    0.01     0.01     |
|                                                                                                            |
| RadialBasisInterpolation<>                                                                                 |
|   interpolate_field_data()     4219      0.2936      0.000070    0.2936      0.000070    1.70     1.70     |
|   prepare_for_use()            2         15.3048     7.652407    15.3048     7.652407    88.73    88.73    |
|                                                                                                            |
| System                                                                                                     |
|   project_vector()             3         0.2461      0.082041    0.5710      0.190338    1.43     3.31     |
|                                                                                                            |
| TecplotIO                                                                                                  |
|   write_nodal_data()           3         0.0744      0.024791    0.0744      0.024791    0.43     0.43     |
 ------------------------------------------------------------------------------------------------------------
| Totals:                        17355     17.2494                                         100.00            |
 ------------------------------------------------------------------------------------------------------------

 
***************************************************************
* Done Running Example miscellaneous_ex8:
*  mpirun -np 4 example-devel  -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc
***************************************************************
make[4]: Leaving directory `/net/spark/workspace/roystgnr/libmesh/git/devel/examples/miscellaneous/miscellaneous_ex8'

Site Created By: libMesh Developers
Last modified: April 23 2013 04:19:30 UTC

Hosted By:
SourceForge.net Logo