マクローリン展開でeの値を求めたい

Question

マクローリン展開を用いてeの値であるe=2.718...を求めたいのですが以下のプログラムではＮの値をどれだけ変化させてもe=2.718と出てきません。
どのように修正する必要があるのでしょうか？
ご回答宜しくお願いします。

#include <stdio.h>
#include <math.h>
#define N 5

double k(int x)
{
  int i,k=1,l=1;

  for(i=0;i<x;i++){
    k=k*l;
    l=l-1;
  }
  return k;
}

double myexp(double x,int n)
{
  int i=1;
  double y;
  double e;

  e=1;
  y=pow(x,i);

  for(i=1;i<N;i++){
    e=e+(y/(k(i)));
  }
  return e;
}

int main(void)
{
  double e;

  e=myexp(1,N);
  printf("e=%.10f\n",e);

  return 0;
}

Comments

[tiitoi]

pow(x,i); や k=k*l; など初期化されていない変数を計算で使用していますが、その場合の計算結果は動作未定義です。(つまり、やってはいけない)
初期化されていない変数が使われているというのがおかしいと思うので、そこをまず直してみてはどうでしょうか

[forza]

ご指摘有難う御座います。
１で初期化してみましたところN =3においてe=2.0 、N=5においてinfとなりました。
１での初期化が誤っているのでしょうか。

Answer 1

初期化されていない変数があります。

double k(int x)

int i, k, l;

double myexp(double x,int n)

int i;

何を入れるのが正しいのでしょうか?

Comments

[forza]

ご回答及びご指摘有難う御座います。
それぞれ1を代入して見たのですが正しいでしょうか。

Answer 2

先程のニュートン法に関する質問でもそうですが、まず数式で確認してから、コードに落とし込むとよいかと思います。

数式を見ると、for 文の中で x^k / k! を足していけばよいとわかりますね。
(0! = 1, x^0 = 1 に注意)

修正案

#include <math.h>
#include <stdio.h>
#define N 20

/**
   @brief k! を計算する。
 */
double factorial(int k)
{
    // 1 * 2 * ... * k を計算する。
    int ret = 1;
    for (int i = 1; i <= k; ++i)
        ret *= i;

    return ret;
}

/**
   @brief e^x を n 項までマクローリン展開した値を返す。
 */
double myexp(double x, int n)
{
    double e = 0;
    for (int k = 0; k < N; ++k)
        e += pow(x, k) / factorial(k);  // x^k / k!

    return e;
}

int main(void)
{
    double e = myexp(1, N);
    printf("e=%.10lf\n", e);  // e=2.7182818351
}

Comments

[forza]

いつもご回答して頂きまして有難う御座います。
数式をよく理解していない状態で更に理解していないプログラミングを書く事は早計でした。
お陰様でこの度も勉強させて頂きました。なるべく１人で解決できるように努めます。
有難う御座いました。

Answer 3

いろいろ、未代入の変数を使っています。myexpの引数のnも使ってないし。
この変数を何に使うかが途中で忘れちゃうんでしょうかね。
よく考え直しましょう。変数の意味を紙に書いては？

Comments

[forza]

ご回答及びご指摘有難う御座います。
一度初期化も含めて初めから見直してみたく思います。

Answer 4

他の方々も仰っておられるように変数の初期化がなされていないのが原因ですね。

それはさておき、階乗関数kを次のように修正してはいかがでしょうか?(変数kと混同するので関数名をfactorialと改めました。)

double factorial(int x)
{
  int i,k;

  k = 1;
  for(i=1;i<=x;i++){
    k *= i;
  }
  return k;
}

それとe = 1;ではなくe = 0;とする必要もあります。

Comments

[forza]

ご回答有難う御座います。
参考にさせて頂きましてe=0として考え直してみます。