output

frame/pvfmm/lapack.h:/*!  DGETRI computes the inverse of a matrix using the LU factorization

extern void dgetrf_(int *M, int *N, double *A, int *LDA, int *IPIV, int *INFO);
/*!  DGETRI computes the inverse of a matrix using the LU factorization
      *  computed by DGETRF.
      *
      *  This method inverts U and then computes inv(A) by solving the system
      *  inv(A)*L = inv(U) for inv(A).
      *
      *  See http://www.netlib.org/lapack/double/dgetri.f for more information
      */
extern void dgetri_(int *N, double *A, int *LDA, int *IPIV, double *WORK, int *LWORK, int *INFO);
}


template<typename T, typename TREE>
hmlpError_t Solve( TREE &tree, Data<T> &input )
{
  using NODE = typename TREE::NODE;

  const bool AUTO_DEPENDENCY = true;
  const bool USE_RUNTIME     = true;

  /** copy input to output */
  auto *output = new Data<T>( input.row(), input.col() );

  SolverTreeViewTask<NODE>             treeviewtask;
  MatrixPermuteTask<true,  NODE> forwardpermutetask;
  MatrixPermuteTask<false, NODE> inversepermutetask;
  /** Sherman-Morrison-Woodbury */
  SolveTask<NODE, T>      solvetask1;
  /** ULV */
  ULVForwardSolveTask<NODE, T>   ulvforwardsolvetask;
  ULVBackwardSolveTask<NODE, T>  ulvbackwardsolvetask;

  /** attach the pointer to the tree structure */
  tree.setup.input  = &input;
  tree.setup.output = output;

  if ( tree.setup.do_ulv_factorization )
  {
    /** clean up all dependencies on tree nodes */
    RETURN_IF_ERROR( tree.dependencyClean() );
    tree.traverseDown( treeviewtask );
    tree.traverseLeafs( forwardpermutetask );
    tree.traverseUp( ulvforwardsolvetask );
    tree.traverseDown( ulvbackwardsolvetask );
    if ( USE_RUNTIME ) hmlp_run();

    /** clean up all dependencies on tree nodes */
    RETURN_IF_ERROR( tree.dependencyClean() );
    tree.traverseLeafs( inversepermutetask );
    if ( USE_RUNTIME ) hmlp_run();
  }
  else
  {
    /** clean up all dependencies on tree nodes */
    RETURN_IF_ERROR( tree.dependencyClean() );
    tree.traverseDown( treeviewtask );
    tree.traverseLeafs( forwardpermutetask );
    tree.traverseUp( solvetask1 );
    if ( USE_RUNTIME ) hmlp_run();
    /** clean up all dependencies on tree nodes */
    RETURN_IF_ERROR( tree.dependencyClean() );
    tree.traverseLeafs( inversepermutetask );
    if ( USE_RUNTIME ) hmlp_run();
  }

  /** delete buffer space */
  delete output;

  return HMLP_ERROR_SUCCESS;
}; /** end Solve() */


  tree.setup.output = output;

  /** w = inv( K + lambda * I ) * u where w is weight and u is potential.
   and rhs is u, the potential*/
  Solve( tree, rhs );


    void Factorize( Data<T> &Kaa ) 


  auto *tree_ptr = gofmm::Compress( K, NN, splitter, rkdtsplitter, config );
  auto &tree = *tree_ptr;
  /** Examine accuracies. */
  auto error = gofmm::SelfTesting( tree, 100, cmd.nrhs );


/** */
template<class ARGUMENT, class NODEDATA>
class Tree


template<typename T, typename SPDMATRIX>
class SimpleGOFMM
{
  public:

    SimpleGOFMM( SPDMATRIX &K, T stol, T budget )
    {
      tree_ptr = Compress( K, stol, budget );
    };

    ~SimpleGOFMM()
    {
      if ( tree_ptr ) delete tree_ptr;
    };

    void Multiply( Data<T> &y, Data<T> &x )
    {
      //hmlp::Data<T> weights( x.col(), x.row() );

      //for ( size_t j = 0; j < x.col(); j ++ )
      //  for ( size_t i = 0; i < x.row(); i ++ )
      //    weights( j, i ) = x( i, j );


      y = gofmm::Evaluate( *tree_ptr, x );
      //auto potentials = hmlp::gofmm::Evaluate( *tree_ptr, weights );

      //for ( size_t j = 0; j < y.col(); j ++ )
      //  for ( size_t i = 0; i < y.row(); i ++ )
      //    y( i, j ) = potentials( j, i );

    };

  private:

    /** GOFMM tree */
    tree::Tree<
      gofmm::Argument<SPDMATRIX, centersplit<SPDMATRIX, 2, T>, T>, 
      gofmm::NodeData<T>> *tree_ptr = NULL; 

}; /** end class SimpleGOFMM */



// return type: ptr to tree::Tree (tree: namespace, Tree: class name)
../gofmm/gofmm.hpp:*Compress( SPDMATRIX &K, T stol, T budget, size_t m, size_t k, size_t s )
../gofmm/gofmm.hpp:  return Compress<SPLITTER, RKDTSPLITTER>

template<class ARGUMENT, class NODEDATA>
class Tree
{
  public:                
    typedef Node<ARGUMENT, NODEDATA> NODE;
};


template<typename ARGUMENT, typename NODEDATA>
class Node : public ReadWrite
{ };


Data<T> proj;

Data
T getvalue( size_t i ) { return (*this)[ i ]; };
T getvalue( size_t i, size_t j ) { return (*this)( i, j ); };