Epanechnikov acceleration added to Krylov method
This commit is contained in:
@@ -9,8 +9,7 @@ librule(
|
||||
"multi_index_util.h",
|
||||
"quick_prune_lpr.h",
|
||||
"relative_prune_lpr.h"],
|
||||
deplibs = ["mlpack/series_expansion:series_expansion",
|
||||
"fastlib:fastlib_int"]
|
||||
deplibs = ["fastlib:fastlib_int"]
|
||||
)
|
||||
|
||||
# The library build rule for Krylov-subspace based local polynomial
|
||||
@@ -26,8 +25,7 @@ librule(
|
||||
"krylov_lpr_test.h",
|
||||
"lpr_util.h",
|
||||
"naive_lpr.h"],
|
||||
deplibs = ["mlpack/series_expansion:series_expansion",
|
||||
"fastlib:fastlib_int"] # dependency
|
||||
deplibs = ["fastlib:fastlib_int"] # dependency
|
||||
)
|
||||
|
||||
# The binary executable rule for Krylov-subspace based local
|
||||
@@ -37,7 +35,6 @@ binrule(
|
||||
sources = ["krylov_lpr_main.cc"],
|
||||
headers = [],
|
||||
deplibs = [":krylov_lpr",
|
||||
"mlpack/series_expansion:series_expansion",
|
||||
"fastlib:fastlib_int"]
|
||||
)
|
||||
|
||||
@@ -46,7 +43,6 @@ binrule(
|
||||
sources = ["dense_lpr_main.cc"],
|
||||
headers = [],
|
||||
deplibs = [":dense_lpr",
|
||||
"mlpack/series_expansion:series_expansion",
|
||||
"fastlib:fastlib_int"]
|
||||
)
|
||||
|
||||
|
||||
@@ -15,11 +15,6 @@
|
||||
#include "matrix_util.h"
|
||||
#include "multi_index_util.h"
|
||||
#include "fastlib/fastlib.h"
|
||||
#include "mlpack/series_expansion/farfield_expansion.h"
|
||||
#include "mlpack/series_expansion/local_expansion.h"
|
||||
#include "mlpack/series_expansion/mult_farfield_expansion.h"
|
||||
#include "mlpack/series_expansion/mult_local_expansion.h"
|
||||
#include "mlpack/series_expansion/kernel_aux.h"
|
||||
|
||||
/** @brief A computation class for dual-tree based local polynomial
|
||||
* regression.
|
||||
|
||||
@@ -36,7 +36,7 @@ class KrylovLpr {
|
||||
////////// Private Type Declarations //////////
|
||||
|
||||
/** @brief The internal query tree type used for the computation. */
|
||||
typedef BinarySpaceTree< DHrectBound<2>, Matrix, KrylovLprQStat >
|
||||
typedef BinarySpaceTree< DHrectBound<2>, Matrix, KrylovLprQStat<TKernel> >
|
||||
QueryTree;
|
||||
|
||||
/** @brief The internal reference tree type used for the
|
||||
@@ -273,7 +273,9 @@ class KrylovLpr {
|
||||
* @param qnode The query node.
|
||||
*/
|
||||
void FinalizeQueryTreeLanczosMultiplier_
|
||||
(QueryTree *qnode, const ArrayList<bool> &exclude_query_flag,
|
||||
(QueryTree *qnode, const Matrix &qset,
|
||||
const ArrayList<bool> &exclude_query_flag,
|
||||
const Matrix ¤t_lanczos_vectors,
|
||||
Matrix &lanczos_prod_l, Matrix &lanczos_prod_e,
|
||||
Vector &lanczos_prod_used_error, Vector &lanczos_prod_n_pruned,
|
||||
Matrix &neg_lanczos_prod_e, Matrix &neg_lanczos_prod_u,
|
||||
@@ -357,11 +359,6 @@ class KrylovLpr {
|
||||
// query point expansion to get the regression estimate.
|
||||
regression_estimates[i] = la::Dot(row_length_, query_pt_solution,
|
||||
query_point_expansion.ptr());
|
||||
|
||||
Vector query_pt_solution_vector;
|
||||
solution_vectors_e.MakeColumnVector(i, &query_pt_solution_vector);
|
||||
query_pt_solution_vector.PrintDebug();
|
||||
|
||||
}
|
||||
}
|
||||
|
||||
@@ -484,10 +481,35 @@ class KrylovLpr {
|
||||
// The second phase solves the least squares problem: (B^T W(q) B)
|
||||
// z(q) = B^T W(q) Y for each query point q.
|
||||
printf("Starting Phase 2...\n");
|
||||
SolveLeastSquaresByKrylov_(qroot, qset, right_hand_sides_e,
|
||||
/*
|
||||
for(index_t q = 0; q < qset.n_cols(); q++) {
|
||||
Matrix qset_single_alias, right_hand_sides_e_single_alias,
|
||||
solution_vectors_e_single_alias;
|
||||
qset_single_alias.Alias(qset.GetColumnPtr(q), qset.n_rows(), 1);
|
||||
right_hand_sides_e_single_alias.Alias(right_hand_sides_e.GetColumnPtr(q),
|
||||
right_hand_sides_e.n_rows(), 1);
|
||||
solution_vectors_e_single_alias.Alias(solution_vectors_e.GetColumnPtr(q),
|
||||
solution_vectors_e.n_rows(), 1);
|
||||
|
||||
// This is hack - construct a query tree out of only the current
|
||||
// query point.
|
||||
QueryTree *qroot_single = tree::MakeKdTreeMidpoint<QueryTree>
|
||||
(qset_single_alias, leaflen, NULL, NULL);
|
||||
|
||||
SolveLeastSquaresByKrylov_
|
||||
(qroot_single, qset_single_alias, right_hand_sides_e_single_alias,
|
||||
solution_vectors_e_single_alias);
|
||||
delete qroot_single;
|
||||
}
|
||||
*/
|
||||
SolveLeastSquaresByKrylov_(qroot, qset, right_hand_sides_e,
|
||||
solution_vectors_e);
|
||||
|
||||
printf("Phase 2 completed...\n");
|
||||
|
||||
// Delete the query tree.
|
||||
delete qroot;
|
||||
|
||||
// Proceed with the third phase of the computation to output the
|
||||
// final regression value.
|
||||
printf("Starting Phase 3...\n");
|
||||
|
||||
@@ -308,9 +308,9 @@ void KrylovLpr<TKernel, TPruneRule>::DualtreeRightHandSidesCanonical_
|
||||
else {
|
||||
|
||||
// Declare references to the query stats.
|
||||
KrylovLprQStat &q_stat = qnode->stat();
|
||||
KrylovLprQStat &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat &q_right_stat = qnode->right()->stat();
|
||||
KrylovLprQStat<TKernel> &q_stat = qnode->stat();
|
||||
KrylovLprQStat<TKernel> &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat<TKernel> &q_right_stat = qnode->right()->stat();
|
||||
|
||||
// Push down postponed bound changes owned by the current query
|
||||
// node to the children of the query node
|
||||
@@ -403,7 +403,7 @@ void KrylovLpr<TKernel, TPruneRule>::FinalizeQueryTreeRightHandSides_
|
||||
Matrix &right_hand_sides_l, Matrix &right_hand_sides_e,
|
||||
Vector &right_hand_sides_used_error, Vector &right_hand_sides_n_pruned) {
|
||||
|
||||
KrylovLprQStat &q_stat = qnode->stat();
|
||||
KrylovLprQStat<TKernel> &q_stat = qnode->stat();
|
||||
|
||||
if(qnode->is_leaf()) {
|
||||
for(index_t q = qnode->begin(); q < qnode->end(); q++) {
|
||||
@@ -436,8 +436,8 @@ void KrylovLpr<TKernel, TPruneRule>::FinalizeQueryTreeRightHandSides_
|
||||
}
|
||||
else {
|
||||
|
||||
KrylovLprQStat &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat &q_right_stat = qnode->right()->stat();
|
||||
KrylovLprQStat<TKernel> &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat<TKernel> &q_right_stat = qnode->right()->stat();
|
||||
|
||||
// Push down approximations
|
||||
la::AddTo(q_stat.postponed_ll_vector_l_,
|
||||
|
||||
@@ -153,8 +153,8 @@ void KrylovLpr<TKernel, TPruneRule>::DualtreeSolverCanonical_
|
||||
Vector &neg_lanczos_prod_used_error, Vector &neg_lanczos_prod_n_pruned) {
|
||||
|
||||
// Variables for storing changes due to a prune.
|
||||
double delta_used_error, delta_n_pruned, delta_neg_used_error,
|
||||
delta_neg_n_pruned;
|
||||
double delta_used_error = 0, delta_n_pruned = 0, delta_neg_used_error = 0,
|
||||
delta_neg_n_pruned = 0;
|
||||
Vector delta_l, delta_e, delta_neg_u, delta_neg_e;
|
||||
delta_l.Init(row_length_);
|
||||
delta_e.Init(row_length_);
|
||||
@@ -179,7 +179,7 @@ void KrylovLpr<TKernel, TPruneRule>::DualtreeSolverCanonical_
|
||||
|
||||
// try finite difference pruning first
|
||||
if(TPruneRule::PrunableKrylovSolver
|
||||
(internal_relative_error_,
|
||||
(internal_relative_error_,
|
||||
rnode->stat().sum_reference_point_expansion_norm_,
|
||||
qnode, rnode, dsqd_range, kernel_value_range,
|
||||
negative_dot_product_range, positive_dot_product_range,
|
||||
@@ -200,7 +200,28 @@ void KrylovLpr<TKernel, TPruneRule>::DualtreeSolverCanonical_
|
||||
|
||||
return;
|
||||
}
|
||||
|
||||
// For the Epanechnikov kernel, we can prune using the far field
|
||||
// moments if the maximum distance between the two nodes is within
|
||||
// the bandwidth! This if-statement does not apply to the Gaussian
|
||||
// kernel, so I need to fix in the future!
|
||||
if(rnode->stat().min_bandwidth_kernel.bandwidth_sq() >= dsqd_range.hi &&
|
||||
rnode->count() > dimension_ * dimension_) {
|
||||
|
||||
la::AddTo(delta_l, &(qnode->stat().postponed_ll_vector_l_));
|
||||
qnode->stat().postponed_ll_vector_n_pruned_ += delta_n_pruned;
|
||||
|
||||
la::AddTo(delta_neg_u, &(qnode->stat().postponed_neg_ll_vector_u_));
|
||||
qnode->stat().postponed_neg_ll_vector_n_pruned_ += delta_neg_n_pruned;
|
||||
|
||||
// Add the pruned reference node to the node list
|
||||
qnode->stat().epanechnikov_pruned_reference_nodes_.AddBackItem(rnode);
|
||||
|
||||
// Keep track of the far-field prunes.
|
||||
num_epanechnikov_prunes_++;
|
||||
return;
|
||||
}
|
||||
|
||||
// for leaf query node
|
||||
if(qnode->is_leaf()) {
|
||||
|
||||
@@ -241,9 +262,9 @@ void KrylovLpr<TKernel, TPruneRule>::DualtreeSolverCanonical_
|
||||
else {
|
||||
|
||||
// Declare references to the query stats.
|
||||
KrylovLprQStat &q_stat = qnode->stat();
|
||||
KrylovLprQStat &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat &q_right_stat = qnode->right()->stat();
|
||||
KrylovLprQStat<TKernel> &q_stat = qnode->stat();
|
||||
KrylovLprQStat<TKernel> &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat<TKernel> &q_right_stat = qnode->right()->stat();
|
||||
|
||||
// Push down postponed bound changes owned by the current query
|
||||
// node to the children of the query node.
|
||||
@@ -415,25 +436,6 @@ void KrylovLpr<TKernel, TPruneRule>::DotProductBetweenTwoBounds_
|
||||
reference_node_directional_bound.lo;
|
||||
}
|
||||
} // End of looping over each component...
|
||||
|
||||
/*
|
||||
for(index_t d = 0; d <= dimension_; d++) {
|
||||
printf("Lanczos vector: [%g %g]\n", lanczos_vectors_bound.get(d).lo,
|
||||
lanczos_vectors_bound.get(d).hi);
|
||||
if(d > 0) {
|
||||
printf("Reference: [%g %g]\n", rnode->bound().get(d - 1).lo,
|
||||
rnode->bound().get(d - 1).hi);
|
||||
}
|
||||
else {
|
||||
printf("Reference: [1 1]\n");
|
||||
}
|
||||
}
|
||||
|
||||
printf("Bounds: %g %g %g %g\n\n", negative_dot_product_range.lo,
|
||||
negative_dot_product_range.hi, positive_dot_product_range.lo,
|
||||
positive_dot_product_range.hi);
|
||||
exit(0);
|
||||
*/
|
||||
}
|
||||
|
||||
template<typename TKernel, typename TPruneRule>
|
||||
@@ -493,21 +495,85 @@ void KrylovLpr<TKernel, TPruneRule>::InitializeQueryTreeLanczosVectorBound_
|
||||
|
||||
template<typename TKernel, typename TPruneRule>
|
||||
void KrylovLpr<TKernel, TPruneRule>::FinalizeQueryTreeLanczosMultiplier_
|
||||
(QueryTree *qnode, const ArrayList<bool> &exclude_query_flag,
|
||||
(QueryTree *qnode, const Matrix &qset,
|
||||
const ArrayList<bool> &exclude_query_flag,
|
||||
const Matrix ¤t_lanczos_vectors,
|
||||
Matrix &lanczos_prod_l, Matrix &lanczos_prod_e,
|
||||
Vector &lanczos_prod_used_error, Vector &lanczos_prod_n_pruned,
|
||||
Matrix &neg_lanczos_prod_e, Matrix &neg_lanczos_prod_u,
|
||||
Vector &neg_lanczos_prod_used_error, Vector &neg_lanczos_prod_n_pruned) {
|
||||
|
||||
KrylovLprQStat &q_stat = qnode->stat();
|
||||
KrylovLprQStat<TKernel> &q_stat = qnode->stat();
|
||||
|
||||
if(qnode->is_leaf()) {
|
||||
for(index_t q = qnode->begin(); q < qnode->end(); q++) {
|
||||
|
||||
// Form the Epanechnikov moments on the fly here and evaluate the
|
||||
// expansions.
|
||||
ArrayList< ArrayList < EpanKernelMomentInfo > > moments;
|
||||
moments.Init(row_length_);
|
||||
for(index_t i = 0; i < row_length_; i++) {
|
||||
moments[i].Init(row_length_);
|
||||
for(index_t j = 0; j < row_length_; j++) {
|
||||
moments[i][j].Init(dimension_);
|
||||
}
|
||||
}
|
||||
|
||||
// Temporary variable for storing the multiindex expansion of a
|
||||
// reference point.
|
||||
Vector reference_point_expansion;
|
||||
reference_point_expansion.Init(row_length_);
|
||||
|
||||
for(index_t n = 0; n < q_stat.
|
||||
epanechnikov_pruned_reference_nodes_.size(); n++) {
|
||||
|
||||
// The current reference node in the list.
|
||||
ReferenceTree *rnode = q_stat.epanechnikov_pruned_reference_nodes_[n];
|
||||
|
||||
for(index_t r = rnode->begin(); r < rnode->end(); r++) {
|
||||
|
||||
// Get the pointer to the reference point.
|
||||
Vector r_col;
|
||||
rset_.MakeColumnVector(r, &r_col);
|
||||
|
||||
// Compute the reference point expansion.
|
||||
MultiIndexUtil::ComputePointMultivariatePolynomial
|
||||
(dimension_, lpr_order_, r_col.ptr(),
|
||||
reference_point_expansion.ptr());
|
||||
|
||||
for(index_t j = 0; j < row_length_; j++) {
|
||||
for(index_t i = 0; i < row_length_; i++) {
|
||||
moments[j][i].Add(reference_point_expansion[j] *
|
||||
reference_point_expansion[i],
|
||||
kernels_[r].bandwidth_sq(), r_col);
|
||||
}
|
||||
}
|
||||
|
||||
} // end of iterating over each reference point.
|
||||
|
||||
} // end of iterating over each pruned reference node.
|
||||
|
||||
// The matrix to store the evaluated moments at each query point.
|
||||
Matrix evaluated_moments;
|
||||
evaluated_moments.Init(row_length_, row_length_);
|
||||
Vector evaluated_moments_times_lanczos_vector;
|
||||
evaluated_moments_times_lanczos_vector.Init(row_length_);
|
||||
|
||||
// Iterate over each query point.
|
||||
for(index_t q = qnode->begin(); q < qnode->end(); q++) {
|
||||
|
||||
if(exclude_query_flag[q]) {
|
||||
continue;
|
||||
}
|
||||
|
||||
// Get the current query point.
|
||||
Vector q_col;
|
||||
qset.MakeColumnVector(q, &q_col);
|
||||
|
||||
// Get the pointer to the current lanczos vector owned by the
|
||||
// current query point.
|
||||
Vector q_current_lanczos_vector;
|
||||
current_lanczos_vectors.MakeColumnVector(q, &q_current_lanczos_vector);
|
||||
|
||||
// Get the column vectors accumulating the sums to update.
|
||||
double *q_lanczos_prod_l = lanczos_prod_l.GetColumnPtr(q);
|
||||
double *q_lanczos_prod_e = lanczos_prod_e.GetColumnPtr(q);
|
||||
@@ -523,12 +589,36 @@ void KrylovLpr<TKernel, TPruneRule>::FinalizeQueryTreeLanczosMultiplier_
|
||||
q_neg_lanczos_prod_e);
|
||||
la::AddTo(row_length_, (q_stat.postponed_neg_ll_vector_u_).ptr(),
|
||||
q_neg_lanczos_prod_u);
|
||||
}
|
||||
|
||||
// Evaluate the Epanechnikov moments.
|
||||
for(index_t i = 0; i < row_length_; i++) {
|
||||
for(index_t j = 0; j < row_length_; j++) {
|
||||
evaluated_moments.set(j, i, moments[j][i].ComputeKernelSum(q_col));
|
||||
}
|
||||
}
|
||||
|
||||
// Now compute the product between the evaluated moments and the
|
||||
// Lanczos vector owned by this query point.
|
||||
la::MulOverwrite(evaluated_moments, q_current_lanczos_vector,
|
||||
&evaluated_moments_times_lanczos_vector);
|
||||
|
||||
// Now accumulate the sum depending on the negativity or the
|
||||
// positivity of each component.
|
||||
for(index_t i = 0; i < row_length_; i++) {
|
||||
if(evaluated_moments_times_lanczos_vector[i] > 0) {
|
||||
q_lanczos_prod_e[i] += evaluated_moments_times_lanczos_vector[i];
|
||||
}
|
||||
else {
|
||||
q_neg_lanczos_prod_e[i] += evaluated_moments_times_lanczos_vector[i];
|
||||
}
|
||||
}
|
||||
|
||||
} // end of iterating over each query point.
|
||||
}
|
||||
else {
|
||||
|
||||
KrylovLprQStat &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat &q_right_stat = qnode->right()->stat();
|
||||
KrylovLprQStat<TKernel> &q_left_stat = qnode->left()->stat();
|
||||
KrylovLprQStat<TKernel> &q_right_stat = qnode->right()->stat();
|
||||
|
||||
// Push down approximations
|
||||
la::AddTo(q_stat.postponed_ll_vector_l_,
|
||||
@@ -549,13 +639,25 @@ void KrylovLpr<TKernel, TPruneRule>::FinalizeQueryTreeLanczosMultiplier_
|
||||
la::AddTo(q_stat.postponed_neg_ll_vector_u_,
|
||||
&(q_right_stat.postponed_neg_ll_vector_u_));
|
||||
|
||||
// Push down Epanechnikov pruned reference nodes.
|
||||
for(index_t i = 0; i < q_stat.
|
||||
epanechnikov_pruned_reference_nodes_.size(); i++) {
|
||||
|
||||
q_left_stat.epanechnikov_pruned_reference_nodes_.
|
||||
AddBackItem(q_stat.epanechnikov_pruned_reference_nodes_[i]);
|
||||
q_right_stat.epanechnikov_pruned_reference_nodes_.
|
||||
AddBackItem(q_stat.epanechnikov_pruned_reference_nodes_[i]);
|
||||
}
|
||||
q_stat.epanechnikov_pruned_reference_nodes_.Resize(0);
|
||||
|
||||
// Recurse both branches of the query node.
|
||||
FinalizeQueryTreeLanczosMultiplier_
|
||||
(qnode->left(), exclude_query_flag,
|
||||
(qnode->left(), qset, exclude_query_flag, current_lanczos_vectors,
|
||||
lanczos_prod_l, lanczos_prod_e, lanczos_prod_used_error,
|
||||
lanczos_prod_n_pruned, neg_lanczos_prod_e, neg_lanczos_prod_u,
|
||||
neg_lanczos_prod_used_error, neg_lanczos_prod_n_pruned);
|
||||
FinalizeQueryTreeLanczosMultiplier_
|
||||
(qnode->right(), exclude_query_flag,
|
||||
(qnode->right(), qset, exclude_query_flag, current_lanczos_vectors,
|
||||
lanczos_prod_l, lanczos_prod_e, lanczos_prod_used_error,
|
||||
lanczos_prod_n_pruned, neg_lanczos_prod_e, neg_lanczos_prod_u,
|
||||
neg_lanczos_prod_used_error, neg_lanczos_prod_n_pruned);
|
||||
@@ -631,7 +733,7 @@ void KrylovLpr<TKernel, TPruneRule>::SolveLeastSquaresByKrylov_
|
||||
|
||||
// Main iteration of the SYMMLQ algorithm - repeat until
|
||||
// "convergence"...
|
||||
for(index_t num_iter = 0; num_iter < row_length_; num_iter++) {
|
||||
for(index_t num_iter = 0; num_iter < sqrt(row_length_); num_iter++) {
|
||||
|
||||
// Determine how many queries are in the Krylov loop.
|
||||
int num_queries_in_krylov_loop = 0;
|
||||
@@ -640,7 +742,6 @@ void KrylovLpr<TKernel, TPruneRule>::SolveLeastSquaresByKrylov_
|
||||
num_queries_in_krylov_loop++;
|
||||
}
|
||||
}
|
||||
printf("%d queries are alive...\n", num_queries_in_krylov_loop);
|
||||
if(num_queries_in_krylov_loop == 0) {
|
||||
break;
|
||||
}
|
||||
@@ -667,16 +768,15 @@ void KrylovLpr<TKernel, TPruneRule>::SolveLeastSquaresByKrylov_
|
||||
neg_lanczos_prod_u, neg_lanczos_prod_used_error,
|
||||
neg_lanczos_prod_n_pruned);
|
||||
FinalizeQueryTreeLanczosMultiplier_
|
||||
(qroot, query_should_exit_the_loop,
|
||||
(qroot, qset, query_should_exit_the_loop, current_lanczos_vectors,
|
||||
lanczos_prod_l, lanczos_prod_e, lanczos_prod_used_error,
|
||||
lanczos_prod_n_pruned, neg_lanczos_prod_e, neg_lanczos_prod_u,
|
||||
neg_lanczos_prod_used_error, neg_lanczos_prod_n_pruned);
|
||||
printf("Finished multiplying Lanczos...\n");
|
||||
|
||||
// Compute v_tilde_mat (the residue after applying the linear
|
||||
// operator the current Lanczos vector).
|
||||
la::AddOverwrite(lanczos_prod_e, neg_lanczos_prod_e, &v_tilde_mat);
|
||||
|
||||
|
||||
/*
|
||||
printf("Positive matrix: %g\n",
|
||||
MatrixUtil::EntrywiseLpNorm(lanczos_prod_e, 1));
|
||||
@@ -706,7 +806,7 @@ void KrylovLpr<TKernel, TPruneRule>::SolveLeastSquaresByKrylov_
|
||||
// vector and v_tilde vector).
|
||||
double alpha = la::Dot(row_length_, current_lanczos_vector,
|
||||
v_tilde_mat_column);
|
||||
|
||||
|
||||
// Subtract the component of the current Lanczos vector (a form
|
||||
// of Gram-Schmidt orthogonalization.)
|
||||
la::AddExpert(row_length_, -alpha, current_lanczos_vector,
|
||||
@@ -720,7 +820,7 @@ void KrylovLpr<TKernel, TPruneRule>::SolveLeastSquaresByKrylov_
|
||||
for(index_t i = 0; i < row_length_; i++) {
|
||||
previous_lanczos_vector[i] = current_lanczos_vector[i];
|
||||
}
|
||||
|
||||
|
||||
// Set a new current Lanczos vector based on v_tilde_mat_column.
|
||||
// A potential place to watch out for division by zero!!
|
||||
if(beta_vec[q] > 0) {
|
||||
|
||||
@@ -4,170 +4,6 @@
|
||||
#error "This file is not a public header file!"
|
||||
#endif
|
||||
|
||||
/** @brief The node statistics used for the query tree.
|
||||
*/
|
||||
class KrylovLprQStat {
|
||||
|
||||
public:
|
||||
|
||||
////////// Member Variables //////////
|
||||
|
||||
/** @brief The lower bound on the norm of the vector computation.
|
||||
*/
|
||||
double ll_vector_norm_l_;
|
||||
|
||||
/** @brief The upper bound on the used error for approximating the
|
||||
* positive components of the vector computation.
|
||||
*/
|
||||
double ll_vector_used_error_;
|
||||
|
||||
/** @brief The lower bound on the portion of the reference set
|
||||
* pruned for the query points owned by this node.
|
||||
*/
|
||||
double ll_vector_n_pruned_;
|
||||
|
||||
/** @brief The lower bound on the norm of the negative components
|
||||
* of the vector computation.
|
||||
*/
|
||||
double neg_ll_vector_norm_l_;
|
||||
|
||||
/** @brief The upper bound on the used error for approximating the
|
||||
* negative components of the vector computation.
|
||||
*/
|
||||
double neg_ll_vector_used_error_;
|
||||
|
||||
/** @brief The lower bound on the portion of the reference set
|
||||
* pruned for the query points owned by this node for the
|
||||
* negative components.
|
||||
*/
|
||||
double neg_ll_vector_n_pruned_;
|
||||
|
||||
/** @brief The lower bound vector offset passed from the above on
|
||||
* each sum component of the vector owned by this node.
|
||||
*/
|
||||
Vector postponed_ll_vector_l_;
|
||||
|
||||
/** @brief This stores the portion pruned by finite difference for
|
||||
* each sum component.
|
||||
*/
|
||||
Vector postponed_ll_vector_e_;
|
||||
|
||||
ArrayList<EpanKernelMomentInfo> postponed_moment_ll_vector_e_;
|
||||
|
||||
/** @brief The amount of used error passed down from above for
|
||||
* approximating the positive components of the vector sum.
|
||||
*/
|
||||
double postponed_ll_vector_used_error_;
|
||||
|
||||
/** @brief The portion of the reference set pruned for approximating
|
||||
* the positive components of the vector sum passed down
|
||||
* from above.
|
||||
*/
|
||||
double postponed_ll_vector_n_pruned_;
|
||||
|
||||
/** @brief This stores the portion pruned by finite difference for
|
||||
* each negative sum component of the vector owned by this
|
||||
* node.
|
||||
*/
|
||||
Vector postponed_neg_ll_vector_e_;
|
||||
|
||||
/** @brief The upper bound vector offset passed from above on each
|
||||
* negative sum component of the right hand sides owned by
|
||||
* this node.
|
||||
*/
|
||||
Vector postponed_neg_ll_vector_u_;
|
||||
|
||||
/** @brief The amount of used error passed down from above for
|
||||
* approximating the negative components of the vector sum.
|
||||
*/
|
||||
double postponed_neg_ll_vector_used_error_;
|
||||
|
||||
/** @brief The portion of the reference set pruned for approximating
|
||||
* the negative components of the vector sum passed down
|
||||
* from above.
|
||||
*/
|
||||
double postponed_neg_ll_vector_n_pruned_;
|
||||
|
||||
/** @brief The bounding box for the Lanczos vectors. */
|
||||
DHrectBound<2> lanczos_vectors_bound_;
|
||||
|
||||
////////// Constructor/Destructor //////////
|
||||
|
||||
/** @brief The constructor which does not do anything. */
|
||||
KrylovLprQStat() {}
|
||||
|
||||
/** @brief The destructor which does not do anything. */
|
||||
~KrylovLprQStat() {}
|
||||
|
||||
////////// Functions during the tree construction //////////
|
||||
|
||||
/** @brief Resets all bounds to zero.
|
||||
*/
|
||||
void Reset() {
|
||||
ll_vector_norm_l_ = 0;
|
||||
ll_vector_used_error_ = 0;
|
||||
ll_vector_n_pruned_ = 0;
|
||||
neg_ll_vector_norm_l_ = 0;
|
||||
neg_ll_vector_used_error_ = 0;
|
||||
neg_ll_vector_n_pruned_ = 0;
|
||||
postponed_ll_vector_l_.SetZero();
|
||||
postponed_ll_vector_e_.SetZero();
|
||||
postponed_ll_vector_used_error_ = 0;
|
||||
postponed_ll_vector_n_pruned_ = 0;
|
||||
postponed_neg_ll_vector_e_.SetZero();
|
||||
postponed_neg_ll_vector_u_.SetZero();
|
||||
postponed_neg_ll_vector_used_error_ = 0;
|
||||
postponed_neg_ll_vector_n_pruned_ = 0;
|
||||
lanczos_vectors_bound_.Reset();
|
||||
|
||||
for(index_t i = 0; i < postponed_moment_ll_vector_e_.size(); i++) {
|
||||
postponed_moment_ll_vector_e_[i].Reset();
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Allocate and initialize memory for the given dimension.
|
||||
*
|
||||
* @param dimension The dimensionality.
|
||||
*/
|
||||
void AllocateMemory(int dimension) {
|
||||
|
||||
// For local polynomial regression order p, each vector contains
|
||||
// (D + p) choose D numbers.
|
||||
int lpr_order = fx_param_int_req(NULL, "lpr_order");
|
||||
int matrix_dimension =
|
||||
(int) math::BinomialCoefficient(dimension + lpr_order, dimension);
|
||||
|
||||
postponed_ll_vector_l_.Init(matrix_dimension);
|
||||
postponed_ll_vector_e_.Init(matrix_dimension);
|
||||
postponed_moment_ll_vector_e_.Init(matrix_dimension);
|
||||
for(index_t i = 0; i < postponed_moment_ll_vector_e_.size(); i++) {
|
||||
postponed_moment_ll_vector_e_[i].Init(dimension);
|
||||
}
|
||||
postponed_neg_ll_vector_e_.Init(matrix_dimension);
|
||||
postponed_neg_ll_vector_u_.Init(matrix_dimension);
|
||||
|
||||
lanczos_vectors_bound_.Init(matrix_dimension);
|
||||
}
|
||||
|
||||
/** @brief Computing the statistics for a leaf node involves
|
||||
* explicitly running over the points owned by the node.
|
||||
*/
|
||||
void Init(const Matrix &dataset, index_t start, index_t count) {
|
||||
|
||||
// Allocate all memory required for the statistics.
|
||||
AllocateMemory(dataset.n_rows());
|
||||
}
|
||||
|
||||
void Init(const Matrix &dataset, index_t start, index_t count,
|
||||
const KrylovLprQStat &left_stat,
|
||||
const KrylovLprQStat &right_stat) {
|
||||
|
||||
// Allocate all memory required for the statatistics.
|
||||
AllocateMemory(dataset.n_rows());
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
/** @brief The node statistics used for the reference tree.
|
||||
*/
|
||||
template<typename TKernel>
|
||||
@@ -305,3 +141,176 @@ class KrylovLprRStat {
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
/** @brief The node statistics used for the query tree.
|
||||
*/
|
||||
template<typename TKernel>
|
||||
class KrylovLprQStat {
|
||||
|
||||
public:
|
||||
|
||||
////////// Member Variables //////////
|
||||
|
||||
/** @brief The lower bound on the norm of the vector computation.
|
||||
*/
|
||||
double ll_vector_norm_l_;
|
||||
|
||||
/** @brief The upper bound on the used error for approximating the
|
||||
* positive components of the vector computation.
|
||||
*/
|
||||
double ll_vector_used_error_;
|
||||
|
||||
/** @brief The lower bound on the portion of the reference set
|
||||
* pruned for the query points owned by this node.
|
||||
*/
|
||||
double ll_vector_n_pruned_;
|
||||
|
||||
/** @brief The lower bound on the norm of the negative components
|
||||
* of the vector computation.
|
||||
*/
|
||||
double neg_ll_vector_norm_l_;
|
||||
|
||||
/** @brief The upper bound on the used error for approximating the
|
||||
* negative components of the vector computation.
|
||||
*/
|
||||
double neg_ll_vector_used_error_;
|
||||
|
||||
/** @brief The lower bound on the portion of the reference set
|
||||
* pruned for the query points owned by this node for the
|
||||
* negative components.
|
||||
*/
|
||||
double neg_ll_vector_n_pruned_;
|
||||
|
||||
/** @brief The lower bound vector offset passed from the above on
|
||||
* each sum component of the vector owned by this node.
|
||||
*/
|
||||
Vector postponed_ll_vector_l_;
|
||||
|
||||
/** @brief This stores the portion pruned by finite difference for
|
||||
* each sum component.
|
||||
*/
|
||||
Vector postponed_ll_vector_e_;
|
||||
|
||||
ArrayList<EpanKernelMomentInfo> postponed_moment_ll_vector_e_;
|
||||
|
||||
/** @brief The amount of used error passed down from above for
|
||||
* approximating the positive components of the vector sum.
|
||||
*/
|
||||
double postponed_ll_vector_used_error_;
|
||||
|
||||
/** @brief The portion of the reference set pruned for approximating
|
||||
* the positive components of the vector sum passed down
|
||||
* from above.
|
||||
*/
|
||||
double postponed_ll_vector_n_pruned_;
|
||||
|
||||
/** @brief This stores the portion pruned by finite difference for
|
||||
* each negative sum component of the vector owned by this
|
||||
* node.
|
||||
*/
|
||||
Vector postponed_neg_ll_vector_e_;
|
||||
|
||||
/** @brief The upper bound vector offset passed from above on each
|
||||
* negative sum component of the right hand sides owned by
|
||||
* this node.
|
||||
*/
|
||||
Vector postponed_neg_ll_vector_u_;
|
||||
|
||||
/** @brief The amount of used error passed down from above for
|
||||
* approximating the negative components of the vector sum.
|
||||
*/
|
||||
double postponed_neg_ll_vector_used_error_;
|
||||
|
||||
/** @brief The portion of the reference set pruned for approximating
|
||||
* the negative components of the vector sum passed down
|
||||
* from above.
|
||||
*/
|
||||
double postponed_neg_ll_vector_n_pruned_;
|
||||
|
||||
/** @brief The bounding box for the Lanczos vectors. */
|
||||
DHrectBound<2> lanczos_vectors_bound_;
|
||||
|
||||
/** @brief The list of reference nodes that were pruned using the
|
||||
* Epanechnikov series expansion.
|
||||
*/
|
||||
ArrayList<BinarySpaceTree< DHrectBound<2>, Matrix, KrylovLprRStat<TKernel> > *> epanechnikov_pruned_reference_nodes_;
|
||||
|
||||
////////// Constructor/Destructor //////////
|
||||
|
||||
/** @brief The constructor which does not do anything. */
|
||||
KrylovLprQStat() {}
|
||||
|
||||
/** @brief The destructor which does not do anything. */
|
||||
~KrylovLprQStat() {}
|
||||
|
||||
////////// Functions during the tree construction //////////
|
||||
|
||||
/** @brief Resets all bounds to zero.
|
||||
*/
|
||||
void Reset() {
|
||||
ll_vector_norm_l_ = 0;
|
||||
ll_vector_used_error_ = 0;
|
||||
ll_vector_n_pruned_ = 0;
|
||||
neg_ll_vector_norm_l_ = 0;
|
||||
neg_ll_vector_used_error_ = 0;
|
||||
neg_ll_vector_n_pruned_ = 0;
|
||||
postponed_ll_vector_l_.SetZero();
|
||||
postponed_ll_vector_e_.SetZero();
|
||||
postponed_ll_vector_used_error_ = 0;
|
||||
postponed_ll_vector_n_pruned_ = 0;
|
||||
postponed_neg_ll_vector_e_.SetZero();
|
||||
postponed_neg_ll_vector_u_.SetZero();
|
||||
postponed_neg_ll_vector_used_error_ = 0;
|
||||
postponed_neg_ll_vector_n_pruned_ = 0;
|
||||
lanczos_vectors_bound_.Reset();
|
||||
|
||||
for(index_t i = 0; i < postponed_moment_ll_vector_e_.size(); i++) {
|
||||
postponed_moment_ll_vector_e_[i].Reset();
|
||||
}
|
||||
epanechnikov_pruned_reference_nodes_.Resize(0);
|
||||
}
|
||||
|
||||
/** @brief Allocate and initialize memory for the given dimension.
|
||||
*
|
||||
* @param dimension The dimensionality.
|
||||
*/
|
||||
void AllocateMemory(int dimension) {
|
||||
|
||||
// For local polynomial regression order p, each vector contains
|
||||
// (D + p) choose D numbers.
|
||||
int lpr_order = fx_param_int_req(NULL, "lpr_order");
|
||||
int matrix_dimension =
|
||||
(int) math::BinomialCoefficient(dimension + lpr_order, dimension);
|
||||
|
||||
postponed_ll_vector_l_.Init(matrix_dimension);
|
||||
postponed_ll_vector_e_.Init(matrix_dimension);
|
||||
postponed_moment_ll_vector_e_.Init(matrix_dimension);
|
||||
for(index_t i = 0; i < postponed_moment_ll_vector_e_.size(); i++) {
|
||||
postponed_moment_ll_vector_e_[i].Init(dimension);
|
||||
}
|
||||
postponed_neg_ll_vector_e_.Init(matrix_dimension);
|
||||
postponed_neg_ll_vector_u_.Init(matrix_dimension);
|
||||
|
||||
lanczos_vectors_bound_.Init(matrix_dimension);
|
||||
|
||||
epanechnikov_pruned_reference_nodes_.Init();
|
||||
}
|
||||
|
||||
/** @brief Computing the statistics for a leaf node involves
|
||||
* explicitly running over the points owned by the node.
|
||||
*/
|
||||
void Init(const Matrix &dataset, index_t start, index_t count) {
|
||||
|
||||
// Allocate all memory required for the statistics.
|
||||
AllocateMemory(dataset.n_rows());
|
||||
}
|
||||
|
||||
void Init(const Matrix &dataset, index_t start, index_t count,
|
||||
const KrylovLprQStat &left_stat,
|
||||
const KrylovLprQStat &right_stat) {
|
||||
|
||||
// Allocate all memory required for the statatistics.
|
||||
AllocateMemory(dataset.n_rows());
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
@@ -274,7 +274,7 @@ class RelativePruneLpr {
|
||||
delta_neg_n_pruned =
|
||||
rnode->stat().sum_reference_point_expansion_norm_;
|
||||
|
||||
// check pruning condition
|
||||
// check pruning condition
|
||||
return (delta_used_error <= allowed_err &&
|
||||
delta_neg_used_error <= neg_allowed_err);
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user