第5回　ベクトル・複素数の演算(2)

複素数の計算

今回は複素数の扱いを通して，データ構造とアルゴリズムを組み合わせ，計算処理を進める基本を学ぶ．
また，次回以降のフーリエ変換に対する準備でもあるので，しっかりと理解してほしい．

まず，構造体による複素数の表現を見てみよう． 


struct Complex
{
    double re;      /*  実部  */
    double im;      /*  虚部  */
};

この構造体を元にして，複素数の四則演算や指数・対数などを関数として準備しておけば，後の計算に便利なツールとなる．
次回，フーリエ変換の計算を行うが，ここでは複素数の四則演算に加えて，複素数の指数関数が必須となるので，その準備も行う．
（複素数の指数関数は周期関数となる．）

それでは，以下の複素数の基礎を思い出してみよう．

ここでは，複素数を z = x + j y とする． j は虚数単位で，j2 = -1 である．

複素数の演算について
［和］ z1 + z2 = (x1 + j y1) + (x2 + j y2 ) = (x1 + x2) + j (y1 + y2)
［差］ z1 - z2 = (x1 + j y1) - (x2 + j y2 ) = (x1 - x2) + j (y1 - y2)
［積］ z1 z2 = (x1 x2 - y1 y2) + j (x1 y2 + x2 y1)
［商］ z1 / z2 = (x1 x2 + y1 y2) / (x2 x2 + y2 y2) + j ( x2 y1 - x1 y2) / (x2 x2 + y2 y2)

［共役複素数(Conjugate)］ z = x - j y

［絶対値］ |z| = (z z)1/2 = (x2 + y2 )1/2

［偏角］ θ = arg(z) = tan-1 (y/x)　= atan (y/x)

［指数関数］(重要)
オイラーの式 e jθ = cosθ + j sinθ を用いて
ez = ex + j y = ex ( cos y + j sin y )

［累乗］
整数 n に対して，
zn = rn ( cos nθ + j sin nθ )
ただし，r = |z| = ( x2 + y2 ) 1/2, θ = arg(z)

［n乗根］
整数 n に対して，
z1/n = r1/n{ cos( (θ+2kπ)/n ) + j sin( (θ+2kπ)/n ) }
一般に，n乗根はn個あるため，k = 0, 1,・・・, n-1 となる．
例えば，8乗根ならn = 8で，k = 0から7まで．
関数に複素数 z と，n および k を渡す．

n乗根のうち，平方根（n = 2 の場合）であれば，以下のような公式がある．


上の式は初期化その他の問題で，微小でも平方根内が負の値を持つ可能性がある．
そのような場合，（math.hで利用可能な標準ライブラリ関数の）sqrt関数で計算するとエラーが発生するので，引数が負にならないように調整をすると安全なプログラムとなる．

本日の課題

実際の計算には上記の全てがあると便利だが，今回の実習ではフーリエ変換に必要な複素数の

  1. 実数2つを，複素数に変換
  2. 四則演算（和，差，積，商）
  3. 共役
  4. 絶対値，偏角
  5. 指数関数

を求める関数をつくろう．
特に，指数関数に使われるオイラーの式は制御工学や伝熱工学，電磁気学などにおいて頻繁に現れる．

複素数計算関数群が出来たら，以下のような動作確認用のmain 関数を作り，作成した各関数を呼び出して，計算が正しく行われるか確認せよ．

関数を含めたソースファイルを提出せよ．ファイル名は 153R000000-06-1.cpp とする．000000の部分は学生番号．

動作確認が済んだら，main 関数以外は別ファイル（ヘッダファイル complex.h および 実装ファイル complex.cpp）にすると良い．
後半のヒントを参照すること．
分割ができたら，全てのソースファイルを提出せよ． ファイル名は 153R000000-06-1.cpp, complex.h, complex.cpp とする．Makefile も作成したら，あわせて提出．


#include <stdio.h>
#include <math.h>
#include "complex.h"   /* 分割コンパイルを行う場合 */

/* 複素数型構造体の定義 */
/* ヘッダファイルに書くと良い． */
struct Complex
{
    double re;
    double im;
};

/* ここに関数を記述．確認ができたら別ファイル（.h および .cpp）に実装すると良い． */


int main()
{
    Complex z1 = {1.0, -1.0};
    Complex z2 = {-3.0, 2.0};
    Complex z3;
    /*  int ???;         */  /*  必要に応じて */
    /*  double ???;      */  /*  必要に応じて */

    /*
    このようにキーボード入力させるのも良い．

    printf("z1を入力：");
    z1 = Cinp();
    
    printf("z2を入力：");
    z2 = Cinp();
    */

    /* いろいろ計算 */
    z3 = Cadd(z1, z2);
    
    /* 表示関数を使って，結果を確認 */
    Cdisp(z3);

    return 0;
}


ヒント：ヘッダファイル，実装ファイルの例

一般的には，汎用性のある関数群は別ファイル（ヘッダファイルおよび実装ファイル）に分割する．
例えば

  1. 複素数の「構造体定義」および，複素数の演算関数一式の「宣言」が記載されているヘッダファイル（例えば complex.h）
  2. 複素数の演算関数一式の「実際の処理」を記述した実装ファイル（例えば complex.cpp）
  3. それらを利用するmain関数の入ったソースファイル（例えば153R000000-06-main.cpp）

の3つのファイルを生成し，コンパイルは

z:¥.windows2000> bcc32 153R000000-06-main.cpp complex.cpp

などとする．ファイル名は各自の学生番号に応じて適宜変更すること．
さらに可能であれば，Makefileを使用して自動化すると良い．（必須では無いが，できると効率Up!）

ヘッダファイル complex.h の例


struct Complex
{
    double re;
    double im;
};

Complex ToComplex(double x, double y);    /* 実数2つから複素数を生成 */

void Cdisp(Complex z);                    /* 複素数を画面に表示 */
Complex Cinp(void);                       /* キーボードから実部，虚部を入力して複素数を返す */

Complex Cadd(Complex a, Complex b);       /* 複素数 a, b の和を返す */
Complex Csub(Complex a, Complex b);       /* 差 */
Complex Cmul(Complex a, Complex b);       /* 積 */
Complex Cdiv(Complex a, Complex b);       /* 商 */

Complex Conj(Complex z);                  /* 共役な複素数を返す */

double Cabs(Complex z);                   /* 絶対値 */
double Carg(Complex z);                   /* 偏角 */

Complex Cexp(Complex z);             　　 /* 指数関数 */

実装ファイル complex.cpp の例


#include <stdio.h>
#include <math.h>
#include "complex.h"   /* ファイル名は各自で変更 */

Complex ToComplex(double x, double y)   //  実数x,yから複素数を作り，戻り値として返す関数
{
    Complex z;

    z.re = x;
    z.im = y;
    return z;
}

void Cdisp(Complex z)   //  実部と虚部を画面に表示
{
    //  ここに処理を記述
}

Complex Cinp(void)      // キーボードから実部，虚部を入力して返す
{
    //  ここに処理を記述
}


Complex Cadd(Complex a, Complex b)      //      複素数の和の例．
{
    Complex z;

    z.re = a.re + b.re;
    z.im = a.im + b.im;
    return z;
}

Complex Csub(Complex a, Complex b)
{
    //  ここに処理を記述
}

Complex Cmul(Complex a, Complex b)
{
    //  ここに処理を記述
}

Complex Cdiv(Complex a, Complex b)
{
    //  ここに処理を記述
}


Complex Conj(Complex z)
{
　　　　//  ここに処理を記述
}



double Cabs(Complex z)
{
    //  ここに処理を記述
}
    
double Carg(Complex z)
{
    //  ここに処理を記述．関数はatan2を使うと良い
}

Complex Cexp(Complex z)   // 指数関数
{
    //  ここに処理を記述
}