C言語のマクローリン展開ローラン展開のコードについて

1/(z^2-1) のマクローリン展開と、z = 1 の周りを 0<|z-1|<2 でローラン展開したときの、近似値をチェックするために以下のようにCでコーディングしたのですが、ローラン展開の実部が怪しいです。関数の使い方に誤りがあるのでしょうか？　コンパイラは gcc です。
　ちなみに、十進BASICでは問題ありませんでした。

z = 0.5 + 0.02i で確認

1/(z^2-1)
-1.33168 -0.0354924i

-1 - z^2 - z^4 - z^6 - z^8
-1.33048 -0.03495i

1/(2*(z-1)) - 1/4 + (z-1)/8 - (z-1)^2/16 + (z-1)^3/32
-1.07262 -0.0366546i
↑がおかしい。
----------------------------------------------------------------
#include <stdio.h>
#include <complex.h>

// 複素数の表示
void cprint(double complex z)
{
　　//printf("%g% + gi\n", creal(z), cimag(z));
　　if (cimag(z) < 0)
　　　　printf("%g　%gi\n", creal(z), cimag(z));
　　else
　　　　printf("%g + %gi\n", creal(z), cimag(z));
}

//FuncO(z) = 1/(z^2-1)
double complex FuncO(double complex z)
{
　　double complex dmy = 1/(z*z-1);
　　return dmy;
}

//cpow(x, y) = x^y
//マクローリン展開
//FuncM(z) = -1 - z^2 - z^4 - z^6 - z^8
double complex FuncM(double complex z)
{
　　double complex dmy = -1 - cpow(z,2) - cpow(z,4) - cpow(z,6) - cpow(z,8);
　　return dmy;
}

//z=1 (0<|z-1|<2)でローラン展開
//FuncR(z) = 1/(2*(z-1)) - 1/4 + (z-1)/8 - (z-1)^2/16 + (z-1)^3/32
double complex FuncR(double complex z)
{
　　double complex dmy = 1/(2*(z-1)) - 1/4 + (z-1)/8 - cpow(z-1,2)/16 - cpow(z-1,3)/32;
　　//double complex dmy = 1/(2*(z-1)) - 1/4 + (z-1)/8 - (z-1)*(z-1)/16 - (z-1)*(z-1)*(z-1)/32;
　　//としてもいっしょ
　　return dmy;
}

int main(void)
{
　　double complex a = 0.5 + 0.02 * I;
　　printf("1/(z^2-1)\n");
　　cprint(FuncO(a));
　　printf("-1 - z^2 - z^4 - z^6 - z^8\n");
　　cprint(FuncM(a));
　　printf("1/(2*(z-1)) - 1/4 + (z-1)/8 - (z-1)^2/16 + (z-1)^3/32\n");
　　cprint(FuncR(a));
　　
　　return 0;
}

質問者からの補足コメント

• > double complex dmy = 1/(2*(z-1)) - 1/4 + (z-1)/8 - cpow(z-1,2)/16 - cpow(z-1,3)/32;
  　第5項が負になっていました。
  
  　　double complex dmy = 1/(2*(z-1)) - 1/4 + (z-1)/8 - cpow(z-1,2)/16 + cpow(z-1,3)/32;
  
  が正しいです。結果も少し違いますが、実部はやはりおかしいです。
  1/(2*(z-1)) - 1/4 + (z-1)/8 - (z-1)^2/16 + (z-1)^3/32
  
  -1.08039 -0.0357176i

    補足日時：2022/12/15 16:01

No.2 ベストアンサー 回答者： tatsumaru77 回答日時： 2022/12/15 17:31

ローラン展開についてはわかりませんが、

double complex dmy = 1/(2*(z-1)) - 1/4 + (z-1)/8 - cpow(z-1,2)/16 - cpow(z-1,3)/32;
を
double complex dmy = 1.0/(2.0*(z-1.0)) - 1.0/4.0 + (z-1.0)/8.0 - cpow(z-1.0,2.0)/16.0 - cpow(z-1.0,3.0)/32.0;

にしてみてください。
結果が
-1.07262 -0.0366546i
から
-1.32262 -0.0366546i
になります。

この回答へのお礼

丁寧な回答まことにありがとうございました。

double complex

ですもんね。

お礼日時：2022/12/15 18:23

No.3 回答者： Tacosan 回答日時： 2022/12/15 17:34

試してみようと手で打っていて気付いた. 1/4 って 0 じゃん....

なお「ローラン展開」の「ローラン」は Laurent だね.

この回答へのお礼

> なお「ローラン展開」の「ローラン」は Laurent だね
ああ、そうでした。ありがとう。

お礼日時：2022/12/15 18:22

No.1 回答者： osamuy 回答日時： 2022/12/15 15:25

^はべき乗じゃなくて、排他的論理和演算子かと。

この回答へのお礼

回答ありがとうございます。
> ^はべき乗じゃなくて、排他的論理和演算子かと。
　いや、それは表示用です。プログラムでのべき乗は複素数専用の cpow を使っています。
　ローラン展開の数式は

　　1/(2(z-1)) - 1/4 + (z-1)/8 - (z-1)^2/16 + (z-1)^3/32

なので、プログラムでは

double complex dmy
= 1/(2*(z-1)) - 1/4 + (z-1)/8 - cpow(z-1,2)/16 - cpow(z-1,3)/32;

したのですけれど、期待された結果が出ないので、この使い方に問題はないのかということをお聞きしたいのです。

お礼日時：2022/12/15 15:54