3

私はScipyを混合統合と差別化に使用して学習しようとしましたが、最初のステップで次の問題が発生しました。

数値微分の場合、呼び出し可能な関数で機能するScipy関数は、私が正しければscipy.derivative（）だけのようです！？しかし、私はそれを扱うことができませんでした：

1番目）微分が行われるポイントを指定しない場合。たとえば、微分が積分の下にあり、その被積分関数に数値を割り当てるのは積分である場合。私ではありません。簡単な例として、Sageのノートブックでこのコードを試しました。

import scipy as sp
from scipy import integrate, derivative
var('y')
f=lambda x: 10^10*sin(x)
g=lambda x,y: f(x+y^2)
I=integrate.quad( sp.derivative(f(y),y, dx=0.00001, n=1, order=7) , 0, pi)[0]; show(I)
show( integral(diff(f(y),y),y,0,1).n() )

また、「警告：丸め誤差の発生が検出されたため、要求された許容値を達成できません。誤差が過小評価されている可能性があります」という警告が表示されます。「dx」を増やして「order」を減らしてもこの警告が続くので、この警告が何を意味するのかわかりません。

2番目）上記の例でg（x、y）のような多変数関数の導関数と、sp.derivative（g（x、y）、（x、0.5）、dx = 0.01、n = 1、order = 3）は、簡単に予想されるように、エラーを出します。

数値微分に関する上記の問題を解決する方法について、あなたからの連絡を楽しみにしています。よろしくお願いします

4

1 に答える 1

0

コードには奇妙な問題がいくつかあり、Pythonをブラッシュアップする必要があることを示唆しています。これらは合法的な構文ではないため、Pythonでこれらの定義をどのように作成したのかわかりません。

まず、古いバージョンのscipyを使用していると思います。最近のバージョン（少なくとも0.12 +から）では、が必要from scipy.misc import derivativeです。derivativescipyグローバル名前空間にありません。

第二に、var定義されていませんが、とにかく必要ではありません（最初にsympyをインポートして使用するつもりだったと思いますsympy.var('y')）。sinまた、math（または必要に応じてnumpy）からインポートされていません。showsympyまたはscipyでは有効な関数ではありません。

^はPythonのパワー演算子ではありません。あなたは**を意味しました

ここでは、記号演算と数値微積分演算の概念を混同しているようです。scipyは、シンボリックオブジェクトを含む式を数値的に区別しません。導関数の2番目の引数は、導関数（つまり数値）を取得したいポイントであると想定されています。あなたが数値微分をしようとしているとあなたが言うように、私はその目的のために問題を解決します。

from scipy import integrate
from scipy.misc import derivative
from math import *

f = lambda x: 10**10*sin(x)
df = lambda x: derivative(f, x, dx=0.00001, n=1, order=7)
I = integrate.quad( df, 0, pi)[0]

さて、この最後の式はあなたが言及した警告を生成し、返される値は絶対値で-0.0731642869874073でゼロにあまり近くありませんが、それはのスケールに比べて悪くはありませんf。有限差分での丸め誤差の問題を理解する必要があります。あなたの関数fは0から10^10の間の間隔で変化します！おそらく逆説的に思えdxますが、微分の値を小さくしすぎると、実際には丸め誤差が大きくなり、数値が不安定になる可能性があります。説明については、ここの2番目のグラフ（「丸め誤差と式誤差の両方が原因でhを選択するのが難しい例」）を参照してください：http：//en.wikipedia.org/wiki/Numerical_differentiation

実際、この場合は、0.001に増やす必要があります。df = lambda x: derivative(f, x, dx=0.001, n=1, order=7)

そうすれば、ひどい丸めをすることなく、安全に統合できます。

I=integrate.quad( df, 0, pi)[0]

から2番目の戻り値を破棄することはお勧めしませんquad。これは「結果の絶対誤差の推定値」であるため、何が起こったかの重要な検証です。この場合、I == 0.0012846582250212652であり、absエラーは〜0.00022であり、これは悪くありません（これは、まだゼロが含まれていないことを意味する間隔です）。たぶん、クワッドのdxと絶対許容誤差をもう少しいじると、さらに良い解決策が得られるでしょうが、うまくいけば、あなたはその考えを理解するでしょう。

gx2番目の問題では、y = 0.5に沿ってg（x、y）を表す適切なスカラー関数（これを呼び出す）を作成する必要があります（これはコンピューターサイエンスではカリー化と呼ばれます）。

g = lambda x, y: f(x+y**2)
gx = lambda x: g(x, 0.5)

derivative(gx, 0.2, dx=0.01, n=1, order=3)

x=0.2での導関数の値を示します。当然のことながら、の規模を考えるとその価値は非常に大きいですf。上で示したように、クワッドを使用して統合できます。

g自体を微分できるようにしたい場合は、別の数値微分関数が必要です。scipyやnumpyがこれをサポートしているとは思いませんが、2Dファインメッシュ（サイズdx）を作成し、numpy.gradientを使用することで、中央の差の計算をハックすることができます。私が知らない他のライブラリソリューションがあるかもしれませんが、私のPyDSToolソフトウェアにはdiffそれを行う関数が含まれていることを知っています（代わりに1つの配列引数を取るようにgを書き直した場合）。これはRidderの方法を使用し、NumericalRecipes擬似コードから着想を得ています。

于 2015-04-10T21:03:43.723 に答える