scipy.odrを使用して、測定したデータ ポイントをモデルに適合させようとしています。ヤコビアンを指定せずにフィッティングを行うとうまくいきますが、一部のモデルでは悪い結果が生じますが、一部のモデルのみで独自のヤコビアンを指定すると、セグメンテーション エラーが発生します。
#!/usr/bin/python2
from __future__ import print_function
import numpy as np
import scipy.odr as odr
tspace = np.linspace(0, 10, 100)
def fun(p, x):
a, x0 = p
r = a * (x - x0)**2
return r
def param_jacobian(p, x, *args):
a, x0 = p
r = np.array([
(x.T - x0)**2,
-2*a*x0*(x.T - x0),
])
return r
def value_jacobian(p, x, *args):
a, x0 = p
r = np.array([
2*a*(x.T - x0),
])
return r
correct = fun((0.5, 3), tspace)
p0 = (0.1, 1)
model = odr.Model(fun, fjacb=param_jacobian, fjacd=value_jacobian)
data = odr.Data(tspace, correct)
fitter = odr.ODR(data, model, beta0=p0)
fitter.set_job(0, 2)
fitter.run()
fitter.output.pprint()
この例は、python 2.7、64 ビット、scipy 0.9.0 で segfault でクラッシュします。fitter.set_job(0, 0)
またはを呼び出してカスタム ヤコビアンの使用を切り替えると、機能します (このfitter.set_job(0, 1)
例では当てはめの結果は良好ですが、他の例ではそうではありません)。
私の間違いはどこですか?それは私の間違いですか?
更新:もう一度その問題に遭遇しました。fit スクリプトで gdb と valgrind を実行しました。gdb はmalloc_consolidate
、モジュールのインポート試行によって深く呼び出されたへのスタックトレースを示していsyslog
ます... 私には奇妙に思えますが、valgrindparam_jacobian
が最初に呼び出された直後に次の 2 倍前に表示されるため、それは無関係だと思います (これは後です)。value_jacobian
)。結果を取得しようとすると、セグメンテーション違反が発生します。セグメンテーション違反はヒープの破損によるものであり、これはodrpack.so
.
==32764== Source and destination overlap in memcpy(0x14239ee0, 0x14239ee0, 32)
==32764== at 0x4A09306: memcpy@@GLIBC_2.14 (mc_replace_strmem.c:653)
==32764== by 0x1B3663EB: fcn_callback (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32764== by 0x1B3A1821: doddrv_ (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32764== by 0x1B3A282A: dodcnt_ (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32764== by 0x1B3A355A: dodrc_ (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32764== by 0x1B36A452: odr (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32764== by 0x3A41249382: PyObject_Call (abstract.c:2529)
==32764== by 0x3A412DA8F6: PyEval_CallObjectWithKeywords (ceval.c:3967)
==32764== by 0x3A412D89E9: builtin_apply (bltinmodule.c:195)
==32764== by 0x3A412E03DA: PyEval_EvalFrameEx (ceval.c:4098)
==32764== by 0x3A412E075D: PyEval_EvalFrameEx (ceval.c:4184)
==32764== by 0x3A412E19A4: PyEval_EvalCodeEx (ceval.c:3330)
さらに後で、*_jacobian
関数が 2 回目に呼び出された後:
==32100== Invalid read of size 8
==32100== at 0x18189BFB: gen_output (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32100== by 0x1818C4BA: odr (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32100== by 0x3A41249382: PyObject_Call (abstract.c:2529)
==32100== by 0x3A412DA8F6: PyEval_CallObjectWithKeywords (ceval.c:3967)
==32100== by 0x3A412D89E9: builtin_apply (bltinmodule.c:195)
==32100== by 0x3A412E03DA: PyEval_EvalFrameEx (ceval.c:4098)
==32100== by 0x3A412E075D: PyEval_EvalFrameEx (ceval.c:4184)
==32100== by 0x3A412E19A4: PyEval_EvalCodeEx (ceval.c:3330)
==32100== by 0x3A412DFF02: PyEval_EvalFrameEx (ceval.c:4194)
==32100== by 0x3A412E19A4: PyEval_EvalCodeEx (ceval.c:3330)
==32100== by 0x3A412E1AD1: PyEval_EvalCode (ceval.c:689)
==32100== by 0x3A412FBD5B: run_mod (pythonrun.c:1361)
==32100== Address 0x10207b40 is 0 bytes inside a block of size 80 free'd
==32100== at 0x4A0662E: free (vg_replace_malloc.c:366)
==32100== by 0xC1050FF: ??? (in /usr/lib64/python2.7/site-packages/numpy/core/multiarray.so)
==32100== by 0xC110199: ??? (in /usr/lib64/python2.7/site-packages/numpy/core/multiarray.so)
==32100== by 0x18189A7E: gen_output (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32100== by 0x1818C4BA: odr (in /usr/lib64/python2.7/site-packages/scipy/odr/__odrpack.so)
==32100== by 0x3A41249382: PyObject_Call (abstract.c:2529)
==32100== by 0x3A412DA8F6: PyEval_CallObjectWithKeywords (ceval.c:3967)
==32100== by 0x3A412D89E9: builtin_apply (bltinmodule.c:195)
==32100== by 0x3A412E03DA: PyEval_EvalFrameEx (ceval.c:4098)
==32100== by 0x3A412E075D: PyEval_EvalFrameEx (ceval.c:4184)
==32100== by 0x3A412E19A4: PyEval_EvalCodeEx (ceval.c:3330)
==32100== by 0x3A412DFF02: PyEval_EvalFrameEx (ceval.c:4194)
同様のエラー(Invalid read of size 8
同じ関数内から)が数回発生します。
この (ほとんど無関係な) bugreportLD_PRELOAD
に投稿された小さな memcpy の代替品を試してみましたが、役に立ちませんでした。