octaveのquadの評価関数に、パラメタ変数入りの関数を使って、パラメタを外から与える。

[2011-02-24追記:
octave 3にて無名関数が追加され

ret = @(c) quad(@(c,x) c * x^2, a, b)

などと簡潔に書けるようになりました.このretをさらにquadに渡して使うなどの多段呼び出しはやはりうまくいかないみたいですが.]

octaveの積分用関数quadは、

quad("f", a, b);

や、

quad("f(x)", a, b);

などとして呼び出すことができます。

 

しかし、この関数が

function y=f(c,x)
  f(c,x) = c * x**2;
endfunction

のようにパラメータを含む場合は、

a=1;b=2;
quad("f(2,x)", a, b);

ということはできても、

a=1;b=2;c=2;
quad("f(c, x)", a, b);

などという呼び出しはできません。
次のように、「cが定義されていません。」とエラーがでます。

error: `c' undefined near line 1 column 38
error: evaluating argument list element number 1
error: evaluating assignment expression near line 1, column 34
error: called from `quad_fcn_N'
error: quad: evaluation of user-supplied function failed

これを回避する方法の一つとして、sprintfでcを評価後の文字列を与えるという手があります。

quad(sprintf("f(%f, x)",c), a, b);

しかし、これの問題はcがベクトルのときに働かないということです。
cがベクトルの時には、自分で要素抽出してやって、関数に与え、返り値をひとつひとつベクトルに納め直すという方法が考えられますが、めんどくさいです。

c = 1:10;
function y = f(c)
  for i=1:length(c)
    y(i) = quad(sprintf("f(%f, x)", c(i)), a, b);
  end
endfunction

めんどくさい上、演算に時間がかかるので、これでフィットなどしようものならどんだけーです。

quadの第一引数は、ヘッダを見ればどうなっているかわかるわけですが、見ずに推測するならば、feval的なことをしていると思われます。なので、fevalにうまくベクトルを渡すことを考えればいいと思うのですが・・。
現段階で、他の方法は思いつきません。

  • このエントリーをはてなブックマークに追加
  • Pocket

コメントを残す