P1673: `vector_sum_of_squares`: Rescaling algorithm may not be correct
まだ誰も着手していません。
評価
調査の方向性
include/experimental/__p1673_bits/blas1_vector_sum_of_squares.hpp と vector_sum_of_squares の再スケーリングループから始めます。quotient と scale の更新を、P1673 のセクション 16.9.2.7 およびリンクされている LAPACK dlassq のリファレンスと比較します。実装が、初期値を含め、記載されたスケーリング係数の不変条件を満たせば完了です。
索引モデルが issue の本文から書いたものです。
説明
@mhoemmen
I was looking at this and it seems to me this is wrong?
Scalar scale = init.scaling_factor;
Scalar ssq = init.scaled_sum_of_squares;
for (extents<>::size_type i = 0; i < x.extent(0); ++i) {
if (abs(x(i)) != 0.0) {
const auto absxi = abs(x(i));
const auto quotient = scale / absxi;
if (scale < absxi) {
ssq = Scalar(1.0) + ssq * quotient * quotient;
scale = absxi;
}
else {
ssq = ssq + quotient * quotient;
}
}
}
-
why does scale change along the way?
-
also, shouldn't
scale / absxi;instead beabsxi/scale? Looking at the spec, section 16.9.2.7:Constraints: For all i in the domain of v, and for f and ssq of type T, the expression ssq = ssq + (abs(x[i]) / f)*(abs(x[i]) / f) is well formed.It says
abs(x[i])/fshould be well-formed, but notf/abs(x[i]). -
also, if things are correct, in theory one should have that:
scaling_factor * scaling_factor * scaled_sum_of_squares = sum of squares of abs(x[i]) plus init.scaling_factor * init.scaling_factor * init.scaled_sum_of_squares.see dlassq. But that does not work for the impl above.
So i tried something on my end:
{
using std::abs;
using std::max;
T scale = init.scaling_factor;
for (std::size_t i = 0; i < x.extent(0); ++i) {
scale = max(scale, abs(x(i)));
}
T ssq = (init.scaling_factor*init.scaling_factor*init.scaled_sum_of_squares)/(scale*scale);
T s= {};
for (std::size_t i = 0; i < x.extent(0); ++i) {
const auto absxi = abs(x(i));
const auto quotient = absxi/scale;
ssq = ssq + quotient * quotient;
s += absxi*absxi;
}
std::experimental::linalg::sum_of_squares_result<T> result;
result.scaled_sum_of_squares = ssq;
result.scaling_factor = scale;
// verify that things are consistent according to definition
const auto lhs = scale*scale*ssq;
const auto rhs = s+init.scaling_factor*init.scaling_factor*init.scaled_sum_of_squares;
std::cout << lhs << " " << rhs << '\n';
EXPECT_NEAR(lhs, rhs, 1e-9);
}
and this works, in the sense that it satisfies:
scaling_factor * scaling_factor * scaled_sum_of_squares = sum of squares of abs(x[i]) plus init.scaling_factor *
init.scaling_factor * init.scaled_sum_of_squares.
- 主要言語
- C++
- スター
- 167
- フォーク
- 31
- PR マージ指標
- 30日以内にマージされた PR はありません
環境構築
このプロジェクトには開発コンテナ、Dockerfile、コントリビューションガイドがありません。まず README を読み、一般的な手順ははじめてのコントリビューションガイドを参照してください。
はじめの一歩
- issue を最後まで読み、次にプロジェクトのコントリビューションガイドを読みます。
- 着手することを issue にコメントします — 二人が同じ作業をするのを防げます。
- リポジトリをフォークし、ブランチを切って変更します。
- issue 番号を参照したプルリクエストを送ります。
kokkos/stdBLAS のほかの issue
-
難易度 3/5 1〜2日 初心者へのやさしさ 48/100
-
難易度 3/5 1〜2日 初心者へのやさしさ 48/100
-
難易度 5/5 1週間以上 初心者へのやさしさ 35/100
-
[linalg] NB commentsオープン
難易度 5/5 1週間以上 初心者へのやさしさ 15/100
-
難易度 4/5 3〜5日 初心者へのやさしさ 35/100
似ている issue
-
Status: Needs confirmation Status: Needs triage
難易度 2/5 1〜3時間 初心者へのやさしさ 64/100
メンテナーはふだん 1 日以内に返信
-
bug-unconfirmed
難易度 2/5 1〜3時間 初心者へのやさしさ 82/100
メンテナーはふだん 1 日以内に返信
-
難易度 2/5 1〜3時間 初心者へのやさしさ 76/100
メンテナーはふだん 1 日以内に返信
-
A-combat P-low S-ready-to-start task
難易度 2/5 1〜3時間 初心者へのやさしさ 77/100
メンテナーはふだん 1 日以内に返信
-
bug intrinsic medium priority rt_lib
難易度 2/5 1〜3時間 初心者へのやさしさ 75/100
メンテナーはふだん 1 日以内に返信