今週のプログラムその4

今週はif文の簡単な使い方

ポイント

  1. if文は条件によってプログラムの実行範囲を変えるための文である。 FORTRANでは
          if(条件1)then
            実行部分1
          else if(条件2)then
            実行部分2
          else
            実行部分3
          endif
    
    という形式で書く。
    この場合、条件1が満たされれば実行部分1が実行される。条件1が満たされず条件2が満たされれば実行部分2が実行される。条件1も2も満たされなければ実行部分3が実行される。
    実行部分は一文だけプログラムでも長いプログラムでもかまわない。また
    else if(条件)then 実行部分
    はひとつとは限らず、異なる条件をいくつ並べてもよい。その場合、条件判断は上から順に行なわれ、最初に満たされた条件に対応する実行部分が実行される。else以下はどの条件も満たされなかった場合に実行される文である。
    なお、if(条件)thenendifだけが必須であり、else ifelseに関する部分はなくても構わない。
    条件を記述する関係演算子には、値が等しい(.eq.)、等しくない(.ne.)、左辺が右辺より大きい(.gt.)、小さい(.lt.)、左辺が右辺以上(.ge.)、 以下(.le.)がある。複数の条件の論理和(.or.)・論理積(.and.)等も使える。詳細は参考書にあたること。 同じプログラムをCで書くと
          if(条件1){
            実行部分1
          }
          else if(条件2){
            実行部分2
          }
          else{
            実行部分3
          }
    
    となる。
    条件を記述する関係演算子には、値が等しい(==)、等しくない(!=)、左辺が右辺より大きい(>)、小さい(<)、左辺が右辺以上(>=)、 以下(<=)がある。複数の条件の論理和(||)・論理積(&&)等も使える。詳細は参考書にあたること。
  2. プログラム9中では絶対値を計算する関数を使っている。FORTRANではabs(x)、Cではfabs(x)という関数である。同じようにsin(x),cos(x)などの関数も用意されている。どういう関数があるかは解説書を参照のこと。
    なお、Cではこういった数学関数を使うために#include<math.h>というおまじないが最初に必要。また、コンパイル時には
    #cc -o Prog prog.c -lm
    
    と数学関数ライブラリを使うための宣言-lmが必要になる場合がある (OSや環境の設定によっては、なくてもいい場合もある)。
  3. program 9で使う入力データはひとつの2*2行列(自明なので例は省略)
なお、プログラム自体は

まずはC版

program 9
..................
#include<stdio.h>
#include<math.h>                数学関数ライブラリを使うためのおまじない
/* 2*2行列の逆行列 */
main()
{
  float a[2][2],ainv[2][2],check[2][2],ad,bc,det;
  int i,j,k;

/*  行列に値を読み込む */
  for( i=0; i<2; ++i){
    for( j=0; j<2; ++j){
      scanf("%f",&a[i][j]);
    }
  }

/*  行列式の計算 */
  ad = a[0][0]*a[1][1];
  bc = a[0][1]*a[1][0];
  det = ad-bc;

/*  行列式が0に近いと桁落ちで計算精度が悪くなるので計算中止 */
  if( fabs((ad-bc)/ad) < 0.00001){          (ad-bc)/adの絶対値が 0.00001より小さかったら{}の中を実行 
    printf("matrix is near singular\n");
    for( i=0; i<2; ++i){
      for( j=0; j<2; ++j){
        printf("  %f",a[i][j]);
      }
      printf("\n");
    }
    printf("with determinant = %f\n", det );
    return;                                 プログラムの実行を終了して処理をOSに戻す
  }

/*  精度が大丈夫なら逆行列を計算 */
  else{                                (ad-bc)/adの絶対値が 0.00001以上なら{}の中を実行
    ainv[0][0]=a[1][1]/det;
    ainv[0][1]=-a[0][1]/det;
    ainv[1][0]=-a[1][0]/det;
    ainv[1][1]=a[0][0]/det;
  }

/* 出力 */
    printf("matrix a:\n");
    for( i=0; i<2; ++i){
      for( j=0; j<2; ++j){
        printf("  %f",a[i][j]);
      }
      printf("\n");
    }
    printf("inverse of a:\n");
    for( i=0; i<2; ++i){
      for( j=0; j<2; ++j){
        printf("  %f",ainv[i][j]);
      }
      printf("\n");
    }

/* 念のため、元の行列と求めた逆行列の積を計算 */
  for( i=0; i<2; ++i){
    for( j=0; j<2; ++j){
    check[i][j]=0;                             実は初期化はここでやってもよい
      for( k=0; k<2; ++k){
        check[i][j] += ainv[i][k]*a[k][j];
      }
    }
  }

/* ちゃんと単位行列になっているかどうか出力 */      
    printf("check :\n");
    for( i=0; i<2; ++i){
      for( j=0; j<2; ++j){
        printf("  %f",check[i][j]);
      }
      printf("\n");
    }
}

FORTRAN版

program 9
..................
C 2*2行列の逆行列
      program main
      real a(2,2),ainv(2,2),check(2,2)

C  ファイルから行列に値を読み込む
      do i=1,2
        read(*,*)(a(i,j),j=1,2)
      enddo

C  行列式の計算
      ad = a(1,1)*a(2,2)
      bc = a(1,2)*a(2,1)
      det = ad-bc

C  行列式が0に近いと桁落ちで計算精度が悪くなるので計算中止
      if( abs((ad-bc)/ad) .lt. 0.00001)then   (ad-bc)/adの絶対値が 0.00001より小さかったら以下を実行
        write(*,*)"matrix is near singular"
        do i=1,2
          write(*,*)(a(i,j),j=1,2)
        enddo
        write(*,*)"with determinant = ", det 
        stop                                  プログラムの実行を終了
   
C  精度が大丈夫なら逆行列を計算
      else                                   (ad-bc)/adの絶対値が 0.00001より大きかったら以下を実行
        ainv(1,1)=a(2,2)/det
        ainv(1,2)=-a(1,2)/det
        ainv(2,1)=-a(2,1)/det
        ainv(2,2)=a(1,1)/det
      endif

C 出力
      write(*,*)"matrix a:"
      do i=1,2
        write(*,*)(a(i,j),j=1,2)
      enddo
      write(*,*)"inverse of a:"
      do i=1,2
        write(*,*)(ainv(i,j),j=1,2)
      enddo

C 念のため、元の行列と求めた逆行列の積を計算
      do i=1,2
        do j=1,2
          check(i,j)=0
          do k=1,2
            check(i,j)=check(i,j)+a(i,k)*ainv(k,j)
          enddo
        enddo
      enddo

C ちゃんと単位行列になっているかどうか出力       
      write(*,*)"check :"
      do i=1,2
        write(*,*)(check(i,j),j=1,2)
      enddo

      end