３次方程式の求解プログラム（ニュートン法）について！

質問

初心者です！
ニュートン法での３次方程式の求解プログラムを作成しようとしています。
ですがよく分かりません(-_-;)
どなたか教えていただけないでしょうか。。。

回答No.6

funaho · ベストアンサー

プログラム書きました。
これで、どうでしょう？

/* 3次方程式の求解 */
#include<stdio.h>
#include<math.h>
#include<float.h>

#define PI 3.14159265358979323846

#define REAL 0 /* 実根 */
#define IMAGE 1 /* 虚数根 */

int main(void){

double a, b, c, d; /* 3次方程式の係数 */
double p, q;
double D;
double x1, x2, x3; /* 根の実部 */
double y1, y2, y3;
double Im; /* 根の虚部 */
double u, v;
double d_u, d_v;
double phi;

int kind; /* 根の種類: REAL, IMAGE */

// ここから
printf("3次方程式(ax^3+bx^2+cx+d=0)の各係数：\n");
printf("a>>>");
scanf("%lf", &a);
printf("b>>>");
scanf("%lf", &b);
printf("c>>>");
scanf("%lf", &c);
printf("d>>>");
scanf("%lf", &d);

// ここまで消して、a, b, c, d値を代入する式を書いてください。

printf("3次方程式(%gx^3%+gx^2%+gx%+g=0)の解は\n", a, b, c, d );
if( a == 0.0 ){
printf("a(%g)がゼロなので，データエラーとみなし，\n", a );
printf("解を求めずに終了します．\n");
return 1; /* 終了コードを1にして，プログラムを終了 */
}

p = ( 3*a*c - b*b ) / ( 9*a*a );
q = ( 2*b*b*b - 9*a*b*c + 27*a*a*d ) / ( 54*a*a*a );

D = q*q + p*p*p; /* Dには誤差が含まれてしまうことに注意 */

if( fabs(D) <= DBL_EPSILON ){ /* Dを数値的にゼロかどうかチェックする */
/* 3つの重根(2重根・3重根)を含む実根 */
y1 = 2*pow(-q, 1.0/3.0);
y2 = y3 = -y1/2;
kind = REAL;
}else if( D < 0 ){ /* 3つの異なる実根 */
phi = acos( -q / pow( -p*p*p, 1.0/3.0 ) );
y1 = 2*sqrt(-p)*cos(phi/3);
y2 = -2*sqrt(-p)*cos(phi/3+PI/3);
y3 = -2*sqrt(-p)*cos(phi/3-PI/3);
kind = REAL;
}else{ /* 1つの実根，2つの共役複素数の解 */
d_u = -q+sqrt(D);
d_v = -q-sqrt(D);

if( d_u >= 0 ) u = pow( d_u, 1.0/3.0 );
else u = -pow( -d_u, 1.0/3.0 );
if( d_v >= 0 ) v = pow( d_v, 1.0/3.0 );
else v = -pow( -d_v, 1.0/3.0 );
y1 = u+v;
y2 = y3 = -y1/2;
Im = sqrt(3)/2*(u-v);
kind = IMAGE;
}

x1 = y1 - b/(3*a);
x2 = y2 - b/(3*a);
x3 = y3 - b/(3*a);

if( kind == REAL ){
printf("x1=%g, x2=%g, x3=%g\n", x1, x2, x3 );
}else{
printf("x1=%g, x2=(%g + %g i), x3=(%g - %g i)\n",
x1, x2, Im, x3, Im );
}

return 0;
}

補足

プログラムまで作っていただいてどうもすみません！
プログラムが↓の部分辺りからよくわからないのですが
p = ( 3*a*c - b*b ) / ( 9*a*a );
q = ( 2*b*b*b - 9*a*b*c + 27*a*a*d ) / ( 54*a*a*a );

あと、初期値も表示したいのですが(^_^;)
まだ、始めたばかりの初心者なので理解悪くてすみません！

回答No.7

funaho

>プログラムが↓の部分辺りからよくわからないのですが
>p = ( 3*a*c - b*b ) / ( 9*a*a );
>q = ( 2*b*b*b - 9*a*b*c + 27*a*a*d ) / ( 54*a*a*a );
ここは、()の中の3*a*cとb*bを先に計算してから引いてます。下も同じような感じです。
これは、演算子が優先順位が高い方から計算されるからです。

演算子の優先順位関して
http://www9.plala.or.jp/sgwr-t/c/sec14.html#s14-5

>あと、初期値も表示したいのですが(^_^;)
プログラム内に記述済み。
これです。「printf("3次方程式(%gx^3%+gx^2%+gx%+g=0)の解は\n", a, b, c, d );」

お礼

どうもありがとうございました！
参考にして頑張ってみたいと思います。

回答No.5

funaho

>ax^3+bx^2+cx+d=0の場合についても解が求められ、初期値は入力しなくていいようにしたいのですが^_^;
入力をやめて、初期値を設定するようにすれば、問題ないと思います。

補足

わざわざありがとうございます。
a.b.c.dを入力しニュートン法で解を求める。初期値は入力しないという決まりで作成しなくてはいけないのです(>_<)

回答No.4

swimeye

ニュートン法を使わずに解の公式から解析的に求めたほうが早いのでは？

３次方程式に限定するとどうなるかかんがえてませんが、
ニュートン法は収束性を保障しないので、
任意の式に対して、かならず解を得られるように
初期値を自動で決めることはできないような気がします

回答No.3

funaho

これですか。

参考URL：

http://www.tonami-h.tym.ed.jp/kadai/H17suu2.pdf

補足

ax^3+bx^2+cx+d=0の場合についても解が求められ、初期値は入力しなくていいようにしたいのですが^_^;
無理ですかね・・・？
何度もすみません<(_ _)>

回答No.2

funaho

ここらへんを参考にしてみては。

参考URL：

http://www.ma.noda.tus.ac.jp/u/tg/comp/L11.pdf

お礼

わざわざありがとうございます<m(__)m>
現在、
１．f(x) = ax＾3 + bx＾2 + cx + dのa.b.c.dの値を入力する
２．判別式Dによって極値の数を求める
３．極値の数により場合分けをして初期値の位置を決める（初期値もプログラムが決定）
４．Newton法により解を求める
という考えでいます。
参考に載っていた数値計算プログラムの例とは少し違うように思うのですが、説明が難しくてなかなか理解できません(>_<)

回答No.1

επιστημη（@episteme）

なにが/どうわからないのかがわからないので答えようがありません。

補足

説明不足でした...すみません(>_<)
以下が私の考えです
１．f(x) = ax＾3 + bx＾2 + cx + dのa.b.c.dの値を入力する
２．判別式Dによって極値の数を求める
３．極値の数により場合分けをして初期値の位置を決める（初期値もプログラムが決定）
４．Newton法により解を求める
この３・４辺りが分からずに困っています。