問題タブ [numerical-methods]

C ++のスパースSVDソルバーを知っている人はいますか? 私の問題には、列/行がゼロになっている可能性のある条件の悪い行列が含まれています。私のデータは、Harwell-Boeing スパース形式である uBLAS マトリックスに格納されています。

見つけるのに苦労しています：

SVD ソルバー

1. スパース行列を操作できる SVD ソルバー。Lapackはこれを行うことができないようですか？関数に疎行列を渡し、疎行列を出力したいと考えています。
2. 結果を再結合する方法... x=b(A^-1) から xs を読み取ることができるように。これは x=(b)(v.(d^-1).(u^t)) になると思います

GSLから次の2つのステップを再現したいと考えています

また、C++ で FORTRAN ライブラリをラップする方法もわかりません。PROPACK c/c++ バインディングはどこにありますか?

編集 1: PROPACK に問題があります。PROPACK はスパース行列を出力しますか? 「V(LDV,KMAX): DOUBLE PRECISION array」として V を出力するようです。そうでないことを意味するのはどれですか？

ほとんどのアプリケーション（特に統計アプリケーション）では、薄いSVDで十分です。ただし、行列（およびその共役）の零空間のオルソベーシスを取得するために、完全なSVDが必要になる場合があります。Rのsvd（）はシンバージョンのみを返すようです。フルバージョンを作成することは可能ですか？代替案はありますか？

uBlas の圧縮されたマトリックスから行を削除する最良の方法は何だろうと思っていましたか? ベクトル内の圧縮された列データを取得する方法があるかどうか疑問に思っていましたか? このページhttp://netlib.org/linalg/html_templates/node92.htmlを見ていますが、これらの要素にアクセスする方法が見つかりません。

私の行列は次のように宣言されています：

私の最終的な目標は、このメソッドを使用してすべてのゼロ行を削除することです。

編集: 解決策には、ublas::compressed_matrix 構造を内部的に変更することが含まれると思われます。ゼロ以外のエントリの数は 1k 程度になると想定しても安全だと思います。行を削除するメソッドをおそらく 20 回呼び出すことを望んでいます。

私はthisの存在を認識しており、このトピックについては this を使用しています。しかし、今回はPythonでの実際の実装についてまとめたいと思います。

私の唯一の問題は、肘のポイントがコードのさまざまなインスタンス化から変化しているように見えることです。この投稿に示されている 2 つのプロットを観察してください。視覚的には似ているように見えますが、エルボー ポイントの値が大幅に変更されています。両方の曲線は、20 回の異なる実行の平均から生成されました。それでも、エルボーポイントの値には大きな変化があります。値が特定の範囲内に収まるようにするには、どのような予防措置を講じることができますか?

私の試みを以下に示します。

私は本当に大きな数でいくつかの数学をやっています（私はPythonを使用していますが、この質問はPython固有ではありません）。1つの値について、私は私に与える式を持っていますf(t) = Pr(X < t)。この式を使用してを取得したいと思いますPr(X >= t) = 1 - f(t)。f(t)戻り値はゼロに非常に近いため、のlog( f(t) )代わりにログ変換と保存を使用していますf(t)。私log( f(t) )は-1e5かそこらのオーダーです。

乗算の場合、これは非常にうまく機能します。log( f(t) * g ) = log( f(t) ) + log(g)。

しかし、 ;log( 1 - f(t) )だけを使用して計算することは非常に困難です。log( f(t) )もちろん、保存して計算した値を一時的にべき乗することもできますが、がゼロに非常に近いため、値log( 1 - exp( log( f(t) ) )が返されます。log( 1 - 0.0 ) = 0.0log( f(t) )

「なぜ気にするのですか？ゼロに近い場合は、1から1を引いた値が1に非常に近くなります」と質問するかもしれません。さて、それはあなたが作った良い点です。あなたはスマートクッキーです。

問題は、これを使用して値をランク付けしたいので、一方がであるか、もう一方がlog(0.999)であるかを本当に気にしlog(0.9999)ます。log( f(t) )また、「では、ランク付けしてから、順序を逆にして、のランク付けを取得してみませんか」と質問することもできますlog( 1 - f(t) )。繰り返しになりますが、私はあなたの質問がどれほど素晴らしいかを指摘せずにはいられません。本当にお話できて光栄です。

しかし、ここに問題があります。私は単にランク付けしたくはありません1 - f(t)。私は実際にに基づいてランク付けしたいと思いPr(X >= t) * g(t) = (1 - f(t)) g(t)ます。ログを取った後、私は取得しlog( 1 - f(t) ) + log( g(t) )ます; 単独でランク付けしf(t)ても、正しい答えは得られません。

過去に私は計算するための小さなPython関数を書きlog(a + b)ましlog(a)たlog(b)：

最初にそれらを正規化して互いに接近させ、次にそれらが接近したときに指数化することで役立ちます。

残念ながら、引き算で同じトリックを実行することはできませんでした。これは、それらが非常に離れているために、互いに近づきlog(1)、近づく正規化係数がないためです。log( f(t) )

誰かがこれを解決する方法を知っていますか？それはそのような古典的な種類の問題のようです。log(1-x)から私に与えることができるビットレベルで動作する巧妙な関数があることを本当に期待/期待/祈っていますlog(x)。また、それがどのように機能するかを知っているなら、私は本当に、本当に知りたいです。

乾杯！オリバー

Monniauxの優れた記事「浮動小数点演算を検証する際の落とし穴」には、予期しない数値動作の例が示されています。ただし、ほとんどの例は、拡張精度FPUまたは拡張精度演算（積和演算など）に依存しています。符号付きゼロとNaNを含む比較に加えて、拡張高精度ハードウェアが無効になっている場合の驚くべき数値動作の良い例はありますか？

私は、パウエルの方法を使用してf(x)=0; の解を見つけるEigenhybrd1のminpackのポートを使用しています。この場合、ヤコビアンは数値的に計算されます。

ソルバーはxtol、パラメーター、相対許容誤差を受け入れます。これは、デフォルトsqrt(eps)で、特定の数値型 (double私の場合) に対応しています。

ほとんどの場合、 のような残差ベクトル ノルムで満足のいく解が得られ1.446575e-9ます。ただし、523.543 などの規範を与える解が受け入れられることもありますが、これはxtol実際には何を意味するのか疑問に思います。現在の残差ノルムが比較される基準値は何ですか?

ActiveState Recipesサイトには、Pythonで内部収益率を実装する機能があります。

このコードは正しい値を返します（少なくともExcelに対してチェックしたいくつかの例では）が、その理由を知りたいと思います。

• ニュートン法（導関数なし）または割線法（1回の反復のみを追跡する）の実装ではないようです。
• 特に、最初のキャッシュフロー要素としての投資変数の定義（およびその後の使用）は私を混乱させます。

何か案は？

GotoBLASという優れた線形ソルバーを次に示します。ダウンロード可能で、ほとんどのコンピューティング プラットフォームで実行できます。私の質問は、このソルバーを Mathematica カーネルとリンクして、LinearSolveのように呼び出すことができる簡単な方法はありますか? ほとんどの人が確実に同意することの 1 つは、非常に大規模な線形システムがある場合は、業界標準の線形ソルバーで解決する方がよいということです。組み込みのソルバーは、非常に大きな問題には適していません。

Mathematica 8 ではより優れたコンパイル機能とライブラリ リンク機能が提供されているため、これらのソルバーの一部を Mathematica 内で使用することが期待できます。問題は、ソース コードの調整をほとんど必要としないか、それとも高度なウィザードである必要があるかです。このフォーラムでは、GotoBLAS のようないくつかの優れたオープン ソース プログラムと Mathematica のリンクを開始し、意見を交換することができます。経験の浅い人でもプロのユーザーから洞察を得ることができ、最終的にははるかに強力な Mathematica を手に入れることができます. これは、増え続ける Mathematica コミュニティのためのオープンなプロジェクトであり、Mathematica 8 のこれらの新しく導入された機能を将来のユーザーのために透過的に文書化できるプラットフォームです。

Mathematica 内から GotoBLAS を実行する方法について、確かなアイデアを提供してくださる方がいらっしゃることを願っています。新しいコンパイルおよびライブラリ リンク機能は、通常、十分に文書化されていないため、一般的なユーザーによって頻繁に使用されることはありません。この質問は、Mathematica のこれらの新しい機能を文書化するためのおもちゃの例として機能します。経験豊富なフォーラム メンバーによるこの方向への支援は、私のような新規ユーザのモチベーションを本当に高めてくれるだけでなく、Mathematica の数値計算機能を拡張するための非常に役立つことを教えてくれます。

私は一般的にMatlab/プログラミングの初心者です。切り捨て誤差が確実に $< 0.5 \times 10 ^{ -5}$ を実行し、反復回数と根の推定値を出力します。

これが私のコンピューターを壊したように見える私の試みです...

私は反復回数のカウントビットを追加していません（これを行う方法がよくわかりません）、すでに行き詰まっています。

助けてくれてありがとう。