本文整理汇总了C++中orthotope类的典型用法代码示例。如果您正苦于以下问题:C++ orthotope类的具体用法?C++ orthotope怎么用?C++ orthotope使用的例子?那么恭喜您, 这里精选的类代码示例或许可以为您提供帮助。
在下文中一共展示了orthotope类的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于我们的系统推荐出更棒的C++代码示例。
示例1: write_matrix_to_octave_file
void write_matrix_to_octave_file(
orthotope<T> const& A
, std::string const& name
)
{
BOOST_ASSERT(2 == A.order());
std::size_t const rows = A.extent(0);
std::size_t const cols = A.extent(1);
std::ofstream file(name + ".mat");
BOOST_ASSERT(file.is_open());
file << "# name: " << name << "\n"
<< "# type: matrix\n"
<< "# rows: " << rows << "\n"
<< "# columns: " << cols << "\n";
for (std::size_t i = 0; i < rows; ++i)
{
for (std::size_t j = 0; j < cols; ++j)
{
T const v = (compare_floating(0.0, A(i, j), 1e-6) ? 0.0 : A(i, j));
file << " " << v;
}
file << "\n";
}
file.close();
}
开发者ID:STEllAR-GROUP,项目名称:hpx_historic,代码行数:32,代码来源:explicit_qr_algorithm.cpp
示例2: random_symmetric_matrix
inline void random_symmetric_matrix(
orthotope<T>& A
, std::size_t seed
)
{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(A.hypercube());
std::size_t const n = A.extent(0);
boost::random::mt19937_64 rng(seed);
boost::random::uniform_int_distribution<> dist(-100, 100);
for (std::size_t l = 0; l < n; ++l)
A(l, l) = T(dist(rng));
for (std::size_t l = 0; l < (n - 1); ++l)
{
for (std::size_t i = l + 1; i < n; ++i)
{
A(i, l) = T(dist(rng));
A(l, i) = A(i, l);
}
}
}
开发者ID:STEllAR-GROUP,项目名称:hpx_historic,代码行数:25,代码来源:explicit_qr_algorithm.cpp
示例3: euclidean_norm
inline T euclidean_norm(
orthotope<T> const& w
)
{
BOOST_ASSERT(1 == w.order());
T sum = T();
for (std::size_t i = 0; i < w.extent(0); ++i)
sum += (w(i) * w(i));
return std::sqrt(sum);
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:13,代码来源:householder_qr_factorization.cpp
示例4: matrix_add
inline void matrix_add(
orthotope<T> const& A
, orthotope<T> const& B
)
{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(2 == B.order());
BOOST_ASSERT(A.hypercube());
BOOST_ASSERT(B.hypercube());
std::size_t const n = A.extent(0);
for (std::size_t i = 0; i < n; ++i)
for (std::size_t j = 0; j < n; ++j)
A(i, j) += B(i, j);
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:16,代码来源:householder_qr_factorization.cpp
示例5: random_matrix
inline void random_matrix(
orthotope<T>& A
, std::size_t seed
)
{
BOOST_ASSERT(2 == A.order());
std::size_t const n = A.extent(0);
std::size_t const m = A.extent(1);
boost::random::mt19937_64 rng(std::time(0));
boost::random::uniform_int_distribution<> dist(-100, 100);
for (std::size_t i = 0; i < n; ++i)
for (std::size_t j = 0; j < m; ++j)
A(i, j) = T(dist(rng));
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:17,代码来源:householder_qr_factorization.cpp
示例6: householders_tri_factor
void householders_tri_factor(
orthotope<T>& A
, std::size_t block_size
, T eps = 1e-8
)
{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(A.hypercube());
std::size_t const n = A.extent(0);
for (std::size_t l = 0; l < (n - 2); ++l)
{
boost::int16_t sign = -compute_sign(A(l + 1, l));
T alpha = 0.0;
for (std::size_t j = (l + 1); j < n; ++j)
alpha += (A(j, l) * A(j, l));
if (alpha < eps)
continue;
alpha = sign * std::sqrt(alpha);
T const r = std::sqrt(0.5 * ((alpha * alpha) - (alpha * A(l + 1, l))));
orthotope<T> w({n});
w(l + 1) = (A(l + 1, l) - alpha) / (2.0 * r);
for (std::size_t j = (l + 2); j < n; ++j)
w(j) = A(j, l) / (2.0 * r);
orthotope<T> H = compute_H(w);
/// A_l = H * A_l_minus_1 * H
orthotope<T> H_A_l_minus_1 = blocked_matrix_multiply(H, A, block_size);
A = blocked_matrix_multiply(H_A_l_minus_1, H, block_size);
}
}
开发者ID:STEllAR-GROUP,项目名称:hpx_historic,代码行数:42,代码来源:explicit_qr_algorithm.cpp
示例7: matrix_multiply
inline orthotope<T> matrix_multiply(
orthotope<T> const& A
, orthotope<T> const& B
)
{
BOOST_ASSERT(A.order() == 2);
BOOST_ASSERT(B.order() == 2);
BOOST_ASSERT(A.extent(1) == B.extent(0));
// (n x m) * (m x p) = (n x p)
std::size_t const n = A.extent(0);
std::size_t const m = A.extent(1);
std::size_t const p = B.extent(1);
orthotope<T> C({n, p});
for (std::size_t i = 0; i < n; ++i)
{
for (std::size_t j = 0; j < p; ++j)
{
T sum = T();
for (std::size_t l = 0; l < m; ++l)
sum += A(i, l) * B(l, j);
C(i, j) = sum;
}
}
return C;
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:31,代码来源:householder_qr_factorization.cpp
示例8: print
inline void print(
orthotope<T> const& A
, std::string const& name
)
{
BOOST_ASSERT(2 == A.order());
for (std::size_t i = 0; i < A.extent(0); ++i)
{
if (i == 0)
std::cout << (boost::format("%- 8s = [ ") % name);
else
std::cout << (boost::format("%|11T |[ "));
for (std::size_t j = 0; j < A.extent(1); ++j)
output_float(A(i, j));
std::cout << "]\n";
}
std::cout << "\n";
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:22,代码来源:householder_qr_factorization.cpp
示例9: compute_H
inline orthotope<T> compute_H(
orthotope<T> const& w
)
{
BOOST_ASSERT(1 == w.order());
std::size_t const n = w.extent(0);
orthotope<T> H({n, n});
for (std::size_t i = 0; i < n; ++i)
{
for (std::size_t j = 0; j < n; ++j)
{
if (i == j)
H(i, j) = 1 - 2 * (w(i) * w(j));
else
H(i, j) = 0 - 2 * (w(i) * w(j));
}
}
return H;
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:23,代码来源:householder_qr_factorization.cpp
示例10: compute_sigma
inline T compute_sigma(
orthotope<T> const& R
, std::size_t n
, std::size_t l
)
{
BOOST_ASSERT(2 == R.order());
T sum = T();
for (std::size_t i = l; i < n; ++i)
sum += (R(i, l) * R(i, l));
return std::sqrt(sum);
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:15,代码来源:householder_qr_factorization.cpp
示例11: blocked_matrix_multiply
inline orthotope<T> blocked_matrix_multiply(
orthotope<T> const& A
, orthotope<T> const& B
, std::size_t block_size
)
{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(2 == B.order());
BOOST_ASSERT(A.hypercube());
BOOST_ASSERT(B.hypercube());
BOOST_ASSERT(0 != block_size);
BOOST_ASSERT(0 == (A.extent(0) % block_size));
std::size_t const n = A.extent(0);
orthotope<T> C({n, n});
// TODO: Figure out how large this will be and do a reserve.
std::vector<hpx::lcos::future<void> > stop_list;
for (std::size_t i = 0; i < n; i += block_size)
{
for (std::size_t j = 0; j < n; j += block_size)
{
orthotope<T> C_sub(C, {block_size, block_size}, {i, j});
matrix_mutex mtx;
for (std::size_t l = 0; l < n; l += block_size)
{
orthotope<T> A_sub(A, {block_size, block_size}, {i, l})
, B_sub(B, {block_size, block_size}, {l, j});
stop_list.push_back(
hpx::async<multiply_and_add_action>(
hpx::find_here(), C_sub, A_sub, B_sub, mtx
).get_future());
}
}
}
hpx::lcos::wait(stop_list);
return C;
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:44,代码来源:householder_qr_factorization.cpp
示例12: matrix_equal
inline bool matrix_equal(
orthotope<T> const& A
, orthotope<T> const& B
)
{
// Are A and B both 2 dimensional?
if ((2 != A.order()) || (2 != B.order()))
return false;
// Do A and B have the same dimensions?
if ((A.extent(0) != B.extent(0)) || (A.extent(1) != B.extent(1)))
return false;
std::size_t const n = A.extent(0);
std::size_t const m = A.extent(1);
for (std::size_t i = 0; i < n; ++i)
for (std::size_t j = 0; j < m; ++j)
if (!compare_floating(A(i, j), B(i, j)))
return false;
return true;
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:23,代码来源:householder_qr_factorization.cpp
示例13: householders
void householders(
orthotope<T> const& A
)
{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(A.hypercube());
std::size_t const n = A.extent(0);
orthotope<T> R = A.copy();
orthotope<T> Q({n, n});
for (std::size_t l = 0; l < n; ++l)
Q(l, l) = 1.0;
for (std::size_t l = 0; l < (n - 1); ++l)
{
T const sigma = compute_sigma(R, n, l);
boost::int16_t const sign = compute_sign(R(l, l));
#if defined(HPXLA_DEBUG_HOUSEHOLDERS)
std::cout << (std::string(80, '#') + "\n")
<< "ROUND " << l << "\n\n";
print(sigma, "sigma");
print(sign, "sign");
#endif
orthotope<T> w({n});
w(l) = R(l, l) + sign * sigma;
for (std::size_t i = (l + 1); i < w.extent(0); ++i)
w(i) = R(i, l);
#if defined(HPXLA_DEBUG_HOUSEHOLDERS)
print(w, "u");
#endif
T const w_norm = euclidean_norm(w);
for (std::size_t i = l; i < n; ++i)
w(i) /= w_norm;
#if defined(HPXLA_DEBUG_HOUSEHOLDERS)
print(w, "v");
#endif
orthotope<T> H = compute_H(w);
#if defined(HPXLA_DEBUG_HOUSEHOLDERS)
print(H, "H");
#endif
R = matrix_multiply(H, R);
Q = matrix_multiply(Q, H);
for (std::size_t i = l + 1; i < n; ++i)
R(i, l) = 0;
}
#if defined(HPXLA_DEBUG_HOUSEHOLDERS)
std::cout << std::string(80, '#') << "\n";
#endif
print(A, "A");
print(Q, "Q");
print(R, "R");
#if defined(HPXLA_DEBUG_HOUSEHOLDERS)
check_QR(A, Q, R);
#endif
}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:74,代码来源:householder_qr_factorization.cpp
示例14: check_QR
void check_QR(
orthotope<T> const& A
, orthotope<T> const& Q
, orthotope<T> const& R
)
{ // {{{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(2 == Q.order());
BOOST_ASSERT(2 == R.order());
BOOST_ASSERT(A.hypercube());
BOOST_ASSERT(Q.hypercube());
BOOST_ASSERT(R.hypercube());
std::size_t const n = A.extent(0);
BOOST_ASSERT(n == Q.extent(0));
BOOST_ASSERT(n == R.extent(0));
///////////////////////////////////////////////////////////////////////////
/// Make sure Q * R equals A.
orthotope<T> QR = matrix_multiply(Q, R);
for (std::size_t l = 0; l < n; ++l)
{
for (std::size_t i = 0; i < n; ++i)
{
if (!compare_floating(A(l, i), QR(l, i), 1e-6))
std::cout << "WARNING: QR[" << l << "][" << i << "] (value "
<< QR(l, i) << ") is not equal to A[" << l << "]["
<< i << "] (value " << A(l, i) << ")\n";
}
}
///////////////////////////////////////////////////////////////////////////
/// Make sure R is an upper triangular matrix.
for (std::size_t l = 0; l < (n - 1); ++l)
{
for (std::size_t i = l + 1; i < n; ++i)
{
if (!compare_floating(0.0, R(i, l), 1e-6))
std::cout << "WARNING: R[" << i << "][" << l << "] is not 0 "
"(value is " << R(i, l) << "), R is not an upper "
"triangular matrix\n";
}
}
///////////////////////////////////////////////////////////////////////////
/// Make sure Q is orthogonal. A matrix is orthogonal if its transpose is
/// equal to its inverse:
///
/// Q^T = Q^-1
///
/// This implies that:
///
/// Q^T * Q = Q * Q^T = I
///
/// We use the above formula to verify Q's orthogonality.
orthotope<T> QT = Q.copy();
// Transpose QT.
for (std::size_t l = 0; l < (n - 1); ++l)
for (std::size_t i = l + 1; i < n; ++i)
std::swap(QT(l, i), QT(i, l));
// Compute Q^T * Q and store the result in QT.
QT = matrix_multiply(Q, QT);
for (std::size_t l = 0; l < n; ++l)
{
for (std::size_t i = 0; i < n; ++i)
{
// Diagonals should be 1.
if (l == i)
{
if (!compare_floating(1.0, QT(l, i), 1e-6))
std::cout << "WARNING: (Q^T * Q)[" << l << "][" << i << "] "
"is not 1 (value is " << QT(l, i) << "), Q is "
"not an orthogonal matrix\n";
}
// All other entries should be 0.
else
{
if (!compare_floating(0.0, QT(l, i), 1e-6))
std::cout << "WARNING: (Q^T * Q)[" << l << "][" << i << "] "
"is not 0 (value is " << QT(l, i) << "), Q is "
"not an orthogonal matrix\n";
}
}
}
} // }}}
开发者ID:fpelliccioni,项目名称:hpx,代码行数:91,代码来源:householder_qr_factorization.cpp
示例15: qr_eigenvalue
std::vector<std::complex<T> > qr_eigenvalue(
orthotope<T> const& A
, std::size_t max_iterations
, std::size_t block_size
, T const& tolerance = 1.0
)
{
BOOST_ASSERT(2 == A.order());
BOOST_ASSERT(A.hypercube());
std::size_t const n = A.extent(0);
/*
std::vector<std::complex<T> > evs;
evs.reserve(n);
*/
std::complex<T> const nan_(std::numeric_limits<T>::quiet_NaN()
, std::numeric_limits<T>::quiet_NaN());
std::vector<std::complex<T> > evs(n, nan_), old(n, nan_);
orthotope<T> Ak = A.copy(), R, Q;
householders_tri_factor(Ak, block_size);
write_matrix_to_octave_file(Ak, "hess_A0");
std::size_t iterations = 0;
while (true)
{
/*
T const mu = Ak(n, n);
if (0 != iterations)
{
for (std::size_t i = 0; i < (n - 1); ++i)
Ak(i, i) -= mu;
Ak(n, n) = 0.0;
}
*/
householders_qr_factor(Ak, Q, R, block_size);
Ak = blocked_matrix_multiply(R, Q, block_size);
/*
if (0 != iterations)
for (std::size_t i = 0; i < n; ++i)
Ak(i, i) += mu;
*/
/*
bool pseudo_upper_triangular = true;
for (std::size_t j = 0; j < (n - 1); ++j)
{
// Make sure we're in Hessenberg form.
for (std::size_t i = j + 2; i < n; ++i)
{
if (!compare_floating(0.0, Ak(i, j), 1e-6))
pseudo_upper_triangular = false;
}
/// Check for convergence. Either we converge to 2x2 complex
/// conjugates eigenvalues, which take the form:
///
/// [ a b ]
/// [ c a ]
///
/// Where b * c < 0. Or, we converge to real eigenvalues which take
/// the form:
///
/// [ a ]
/// [ 0 ]
///
// Determine if we've failed to converge to a real eigenvalue.
if (!compare_floating(0.0, Ak(j, j + 1), 1e-6))
{
// Determine if we've failed to converge to a pair of complex
// eigenvalues.
if (!compare_floating(Ak(j, j), Ak(j + 1, j + 1), 1e-6))
pseudo_upper_triangular = false;
}
}
*/
bool converged = true;
// std::cout << "ITERATION " << iterations << "\n";
for (std::size_t j = 0; j < n; ++j)
{
if (j != n)
{
// Check for complex eigenvalues.
if (!compare_floating(0.0, Ak(j + 1, j), 1e-6))
//.........这里部分代码省略.........
开发者ID:STEllAR-GROUP,项目名称:hpx_historic,代码行数:101,代码来源:explicit_qr_algorithm.cpp
注:本文中的orthotope类示例由纯净天空整理自Github/MSDocs等源码及文档管理平台,相关代码片段筛选自各路编程大神贡献的开源项目,源码版权归原作者所有,传播和使用请参考对应项目的License;未经允许,请勿转载。 |
请发表评论