x^2+y^2
数学の計算に使用したい式がいくつかあります。私がやりたいことの 1 つは、式の偏微分を取ることです。
したがって、に対するf(x,y) = x^2 + y^2
の部分がになる場合、 に対する の部分は になります。f
x
2x
y
2y
有限差分法を使用して厄介な関数を作成しましたが、浮動小数点の精度に関する多くの問題に直面しています。たとえば、1.99234
代わりに2
. 記号微分をサポートするライブラリはありますか? 他の提案はありますか?
x^2+y^2
数学の計算に使用したい式がいくつかあります。私がやりたいことの 1 つは、式の偏微分を取ることです。
したがって、に対するf(x,y) = x^2 + y^2
の部分がになる場合、 に対する の部分は になります。f
x
2x
y
2y
有限差分法を使用して厄介な関数を作成しましたが、浮動小数点の精度に関する多くの問題に直面しています。たとえば、1.99234
代わりに2
. 記号微分をサポートするライブラリはありますか? 他の提案はありますか?
私はそのようなライブラリをいくつかの異なる言語で実装しましたが、残念ながら C では実装していません。多項式 (和、積、変数、定数、累乗) だけを扱っている場合は、かなり簡単です。三角関数も悪くない。より複雑なものは、時間をかけて他の人のライブラリをマスターする方がよいでしょう。
自分で巻くことにした場合、生活を簡素化するためのいくつかの提案があります。
不変のデータ構造 (純粋に機能的なデータ構造) を使用して式を表現します。
Hans Boehm のガベージ コレクターを使用して、メモリを管理します。
線形和を表すには、有限マップ (バイナリ検索ツリーなど) を使用して、各変数をその係数にマップします。
Luaを C コードに埋め込んで計算を行う場合は、Lua コードをhttp://www.cs.tufts.edu/~nr/drop/luaに置きます。優れた機能の 1 つは、シンボリック式を取得して微分し、結果を Lua にコンパイルできることです。もちろん、ドキュメントは何も見つかりません:-(
(エラーを最小限に抑えるという意味で)数値微分を「正しく」することは非常に難しい場合があります。開始するには、数値微分に関する数値レシピのセクションを確認することをお勧めします。
無料の記号数学パッケージについては、GiNaCを参照してください。また、自己完結型の純粋なPython記号数学パッケージであるSymPyをチェックすることもできます。SymPyは、Pythonコマンドラインからインタラクティブに使用できるため、探索がはるかに簡単であることがわかります。
商用の面では、MathematicaとMapleの両方にCAPIがあります。ライブラリを使用するには、インストール/ライセンスされたバージョンのプログラムが必要ですが、どちらも、目的のタイプの記号の区別を簡単に行うことができます。
2 つの簡単な方法で、数値微分の精度を向上させることができます。
より小さなデルタを使用してください。前後の値を使用したようです1e-2
。から始めて1e-8
、小さな痛みや助けになるかどうかをテストします。明らかに、マシンの精度に近づきすぎることはできません-約1e-16
2倍です。
前方 (または後方) 差異ではなく、中央差異を使用します。すなわちdf_dx =(f(x+delta) - f(x-delta)) / (2.0*delta)
、より高い切り捨て項のキャンセルを行うため、中心差分推定値の誤差は、delta^2
前方差分のデルタではなく次数です。http://en.wikipedia.org/wiki/Finite_differenceを参照してください
独自のパッケージを作成するよりも既存のパッケージを活用する方が確かに簡単ですが、独自のパッケージを作成することに決めており、C ++の暗い部分について学ぶ準備ができている場合は、Boost.Protoを使用できます。 Boostから独自のライブラリを設計するまで。
基本的に、Boost.Protoを使用するx * x + y * y
と、式テンプレート(基本的にはネストされたstruct
sを使用したその式の解析ツリーの表現)などの有効なC ++式を変換し、後でその解析ツリーに対して任意の計算を実行できます。それを呼び出すproto::eval()
ことによって時間。デフォルトでproto::eval()
は、ツリーを直接実行されたかのように評価するために使用されますが、代わりにシンボリック導関数を取得するように各関数または演算子の動作を変更できない理由はありません。
これは問題に対する非常に複雑な解決策になりますが、それでも、C++テンプレートメタプログラミング手法を使用して独自の式テンプレートを作成するよりもはるかに簡単です。
これが実際に使用したい種類の関数である場合、クラス ライブラリを作成するのは簡単です。係数と指数を含む 1 つの項から始めます。項のコレクションで構成される多項式があります。
対象の数学的メソッド (たとえば、add/sub/mul/div/differentiate/integrate) のインターフェイスを定義すると、GoF Composite パターンが表示されます。Term と Polynomial の両方がそのインターフェイスを実装します。多項式は、そのコレクション内の各項を単純に反復します。
これは、C / C ++ではなくLispに適用されるため、ちょっと脇に置いてありますが、他の人が同様のタスクを探している場合や、C /C++で同様のものを自分で実装するためのアイデアを得ることができます。SICPには、lispの主題に関するいくつかの講義があります。
Lispでは、それは非常に簡単です(そして、強力なパターンマッチングとポリモーフィック型を備えた他の関数型言語では)。Cでは、同じパワーを得るには、おそらく列挙型と構造体を多用する必要があります(割り当て/割り当て解除は言うまでもありません)。1時間以内にocamlで必要なものを確実にコーディングすることができます-タイピング速度が制限要因だと思います。Cが必要な場合は、実際にCからocamlを呼び出すことができます(その逆も可能です)。
導関数の 1 次のみを計算することは、実装するのが非常に簡単です。しかし、それを速くするのは芸術です。を含むクラスが必要です
次に、加算や減算などの演算子と、この演算の基本的でよく知られているルールを実装する sin() などの関数を記述します。
高次導関数を計算するには、切り捨てられたテイラー級数を使用して行う必要があります。上記のクラスをそれ自体に適用することもできます。値と導関数の値の型はテンプレート引数である必要があります。しかし、これは導関数の計算と保存を複数回行うことを意味します。
切り詰められたテイラー級数 -- これには 2 つのライブラリが利用可能です: