Mathematicaを用いてf(x,y)=cをみたす座標(x,y)のリストを手に入れたいです。

はじめまして。物理学が専門の大学院生です。Mathematica初心者ですがアドバイスをもらえると嬉しいです。

現在、ある滑らかな2変数関数f(x,y)=cをみたす座標(x,y)の集合を手に入れる必要があります。もちろんContourplotを使えばすぐに等高線図は書けます。しかし、欲しいのは等高線プロットではなく、f(x,y)=cを満たす座標(x,y)の集合です。そこで私はMathematicaで次のようにコードをかきました。

Table[{y, FindRoot[x^2 + y^2 == 4, {x, 1}]}, {y, -1, 1, 0.1}]

いま、話をわかりやすくするために、実定数cを4として,2変数関数f(x,y)をf(x,y)=x^2 + y^2としていますが、実際はもっと複雑な方程式です（微分方程式ではありません）。そのため、NsolveではなくFindRootを用いてます。このコードを実行すると、x^2 + y^2 == 4を満たす数値解(p,q)が{p,x→q}という形で20個求まります（実際は20個じゃ足りないのでもっと増やします）。

質問１．求まった数値解(p,q)の集合をリストとしてエクスポートしたいのですが、どうしたらよいでしょうか。このままでは、リストは(p,x→q)という風に変なところに矢印がついています。この矢印を消して、txtファイルやdatファイルに数値解の集合をエクスポートしたいのですが、その方法を教えてくれませんか。

質問２．上の方法では、明らかに等高線を引くにあたって求め欠けている解があります。具体的にはyが負の解が一切リストアップされません。この問題を解決するには解の探索の初期値を変える必要があると思いますが、他に方法はないのでしょうか。正直、上のやり方はなんだかとてもどんくさい気がして、f(x,y)=cをみたす(x,y)の集合をリストとして手に入れるもっと簡単な方法はないのでしょうか。
（できればMathematicaを用いて求めたいのですが、C言語やFortranで2分法などのアルゴリズムを用いて(x,y)の座標の集合を手に入れるのもありかなぁと思っています）。

よろしくお願いします。

No.1 ベストアンサー 回答者： rabbit_cat 回答日時： 2016/09/01 23:35

FindInstanceという関数があります。

http://reference.wolfram.com/language/ref/FindIn …

FindInstance[x^2 + y^2 == 4, {x,y}, Reals, 100]
とか

x->1 の -> は置換規則を表しているので
x ./ x->1
って置換(ReplaceAll）してやると、1 になります。

この回答へのお礼

ありがとうございます！！
ものすごく便利な関数ですね！

とりあえず、
Export["x.dat",
List[x /. FindInstance[x^2 + y^2 == 4, {x, y}, Reals, 1000] //
N ], "Data"]
Export["y.dat",
List[y /. FindInstance[x^2 + y^2 == 4, {x, y}, Reals, 1000] //
N ], "Data"]
とすることでf(x,y)=4をみたす解の集合(x,y)のxをx.datにyをy.datにエクスポートすることができました。

ただ、理想を言えば、解の集合を一つのファイル内で
x1 y1
x2 y2
x3 y3
というように、並べた状態で手に入れたいです（１段ずつ改行もさせておきたい）。その方が何かと便利で…。一応、エクセルにはりつけたりして手動でこの形にすることはできるのですが、これから大きな数値計算をするので毎回手動でやるのは大変です。

なんとか良い方法ないでしょうか。

お礼日時：2016/09/02 11:05