質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

新規登録して質問してみよう
ただいま回答率
85.50%
C

C言語は、1972年にAT&Tベル研究所の、デニス・リッチーが主体となって作成したプログラミング言語です。 B言語の後継言語として開発されたことからC言語と命名。そのため、表記法などはB言語やALGOLに近いとされています。 Cの拡張版であるC++言語とともに、現在世界中でもっとも普及されているプログラミング言語です。

Q&A

解決済

1回答

1841閲覧

オイラーの等式について

lack_un

総合スコア58

C

C言語は、1972年にAT&Tベル研究所の、デニス・リッチーが主体となって作成したプログラミング言語です。 B言語の後継言語として開発されたことからC言語と命名。そのため、表記法などはB言語やALGOLに近いとされています。 Cの拡張版であるC++言語とともに、現在世界中でもっとも普及されているプログラミング言語です。

0グッド

0クリップ

投稿2016/07/12 08:47

###前提・実現したいこと
オイラーの等式で、e^iπ+1=0を証明するプログラムです。
具体的には、実部と虚部にわけて計算し、それらを最後に足して答えを出そうというものでしたが、答え自体がでません。

まず、<complex.h>の使い方が合っているのか、と思ったのですが

###発生している問題・エラーメッセージ

エラーはないが、計算ができない。

C

1#include <stdio.h> 2#include <math.h> 3#include <complex.h> 4 5typedef struct{ 6 double re; 7 double im; 8}COMPLEX; 9 10COMPLEX compbeki(COMPLEX a); 11int kaijo(int n); 12 13int main(){ 14 15 COMPLEX z1; 16 COMPLEX z2; 17 z2 = compbeki(z1); 18 19 printf("方程式の答えは%fです。",z2.re+z2.im); 20 21 return 0; 22} 23 24COMPLEX compbeki(COMPLEX a){ 25 COMPLEX b; 26 int i; 27 double x,y; 28 x=1; 29 y=0; 30 31 for(i=1; i<10000; i++){ 32 33 if(i%2 == 0){ 34 x= x+ ((-1.0)/kaijo(i))*pow(M_PI,i*2); 35 y= y+ ((-1.0)/kaijo(2*i+1))*pow(M_PI,2*i-1); 36 } 37 else{ 38 x= x+ ((1.0)/kaijo(i))*pow(M_PI,i*2); 39 y= y+ ((1.0)/kaijo(2*i+1))*pow(M_PI,2*i-1); 40 } 41 } 42 43 b.re = creal(x); 44 b.im = cimag(y); 45 46 return b; 47} 48 49int kaijo(int n){ 50 int i; 51 int seki; 52 seki=0; 53 54 for(i=0; i<n; i++){ 55 seki = n*kaijo(n-1); 56 } 57 58 return seki; 59}

###試したこと
今見てみると、for文章の中がかなりおかしい、というのがわかるのですが、どうも頭をひねっても答えが思いつきませんでした..。
こういう風に、実部と虚部にわけて足して演算する形でよいのでしょうか。

###補足情報(言語/FW/ツール等のバージョンなど)
より詳細な情報

気になる質問をクリップする

クリップした質問は、後からいつでもMYページで確認できます。

またクリップした質問に回答があった際、通知やメールを受け取ることができます。

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

guest

回答1

0

ベストアンサー

一番おかしな所は階乗の所です。for文の使い方というよりも、intで大きな階乗を扱うのが不味いです。階乗は大きくなるのがはやく、13!=6,227,020,800 であり、2**32=4,294,967,296 を越えます。現在のほとんどのパソコンではintは32bitですので、すぐに限界になってしまいます。

そのほか<complex.h>をincludeしているのに使ってないとか、ちょっと色々怪しいです。参考までに作ってみたので見比べてみてください。(C11オンリー、gccに"-std=c11"をつけないと無理かも知れない。)

C

1/** 2 * オイラーの等式の計算 3 * $e^{i\pi} = \sum_{n=0}^{\infty}\frac{\pi^n}{n!}i^n = -1$ 4*/ 5 6#include <complex.h> 7#include <limits.h> 8#include <stdio.h> 9#include <tgmath.h> 10 11void print_complex(long double complex c); 12 13int main(void) 14{ 15 long double complex eipi = 1; 16 long double complex prev_eipi = 1; 17 long double pi = acos((long double)-1); 18 // $\frac{\pi^n}{n!}$ 19 // 最初の \pi^0 が 0 になるように pi で割っておく 20 long double pin_n = 1 / pi; 21 for (long long n = 0; n < LLONG_MAX; n++) { 22 pin_n *= pi; 23 if (n != 0) { 24 pin_n /= n; 25 } 26 // 計算できないほど小さいなら終わりにする。 27 if (pin_n == 0) { 28 break; 29 } 30 prev_eipi = eipi; 31 switch (n % 4) { 32 case 0: 33 // $n^(4k+0) = (n^4)^k = 1$ 34 eipi += pin_n; 35 break; 36 case 1: 37 // $n^(4k+1) = (n^4)^ki = i$ 38 eipi += pin_n * I; 39 break; 40 case 2: 41 // $n^(4k+2) = (n^4)^ki^2 = -1$ 42 eipi += -pin_n; 43 break; 44 case 3: 45 // $n^(4k+3) = (n^4)^ki^3 = -i$ 46 eipi += -pin_n * I; 47 break; 48 } 49 printf("%lld: ", n); 50 print_complex(eipi); 51 printf("\n"); 52 // long doubleの限界を超えて同じになるなら終わりにする。 53 if (prev_eipi == eipi) { 54 break; 55 } 56 } 57 printf("%s = ", "e^(i*pi)+1"); 58 print_complex(eipi); 59 printf("\n"); 60 return 0; 61} 62 63void print_complex(long double complex c) 64{ 65 long double c_r = creal(c); 66 long double c_i = cimag(c); 67 if (c_i == 0) { 68 printf("%Lf", c_r); 69 } else if (c_r == 0) { 70 printf("%Lfi", c_i); 71 } else if (c_i > 0) { 72 printf("%Lg+%Lgi", c_r, c_i); 73 } else { 74 printf("%Lg%Lgi", c_r, c_i); 75 } 76}

Q/A

  • $...$は? -> mathjaxを使うときのtexの記法。Qiitaあたりで使うと数式になるよ。https://stackedit.io/editorとかで確かめてみると良いかも。
  • <tgmath.h>って何? -> _Genericを使ったmath。floatでもdoubleでlong doubleでもcomplexでも同じ関数名でできるようになる。
  • long doubleって意味あるの? -> 4倍精度浮動小数点数対応のCPUとコンパイラならきっと…たぶん。
  • 途中で止めてる? -> long doubledoubleと同じぐらいの環境では、34回目あたりで表示の限界になり、49回目で足しても値が変わらなくなる。これ以上は無意味なので、止めてます。
  • 階乗するところは別の関数にしないの? -> 前に足した部分へ π/n かければ次になるのでそのままforで一緒に回している。関数にする場合はメモ化とかしないと遅くなる。
  • M_PIは使わないの? -> M_PIはCの標準じゃ無いから、acos(-1)でわざわざ計算している。tgmath.hのせいであらかじめlong doubelにキャストしないとlong dobuleで計算してくれない。
  • 最終的にはどうなるの? -> 手元では e^(i*pi)+1 = -2.27573e-19-1.43047e-19i となった。限りなく0に近くなるので間違ってはいないと思う。

投稿2016/07/12 12:03

raccy

総合スコア21733

バッドをするには、ログインかつ

こちらの条件を満たす必要があります。

lack_un

2016/07/12 21:29

回答ありがとうございます。 なかなかむずかしい...、ちょっと理解に時間がかかりそうですが、しっかり理解したのち自分でもソースを書いてみたいと思います。
guest

あなたの回答

tips

太字

斜体

打ち消し線

見出し

引用テキストの挿入

コードの挿入

リンクの挿入

リストの挿入

番号リストの挿入

表の挿入

水平線の挿入

プレビュー

15分調べてもわからないことは
teratailで質問しよう!

ただいまの回答率
85.50%

質問をまとめることで
思考を整理して素早く解決

テンプレート機能で
簡単に質問をまとめる

質問する

関連した質問