連立一次方程式の数値解法

準備

1. プログラミングIII用の「ソリューション」を作成する。
2. その下に連立方程式の数値解法のプログラミング用「プロジェクト」を作成する。以降の実習においても章ごとにプロジェクトを作成するようにしよう。ソースファイルが増えると管理が大変なので。
3. テキストのサンプルプログラムをダウンロード、各自のフォルダに解凍し、利用できるようにしておく。
4. 以下の「事前学習」〜「課題３．１のプログラムでテスト」までについて，行ったことを報告しなさい．

1., 2. については、前準備の節に書いてあります。

作業報告は，Word等で作成したファイル(画面のスナップショットを含める）をPDF形式に出力し，e-Learningのページ（下記URL）ヘアップロードしなさい．

https://moodle.yamanashi.ac.jp/2016/course/view.php?id=5288

レポートサンプル

【締切：10/5の午後６時】

10/12

1. program21.c，program22.cのdvector, dmatrix関数でのメモリーの利用に関する説明を理解する．補足説明 (ここまで30分程度)
2. 課題1-1, 1-2を行い，レポートを提出する．

テキスト

独自に作成したテキスト(PDF版）（学内からのアクセスのみ） （印刷して授業で配布しますが、オンラインで見たい場合はここから取得してください。）

ガウス消去法のプログラム課題の実行

事前学習 (テキスト§2.1, 2.2のプログラム理解)

３章のプログラムでは、１１行目から２２行目にかけて定義してある関数（行列やベクトル変数の領域確保など）は、多くのプログラムに共通して用いられる関数であり、テキストには書いてない．（筆者のページからダウンロードしたプログラムには付属している。）

したがって， 独自にはじめからプログラムを書く場合には、これらの関数の中身を書き加える必要がある。（プログラムは２．１〜２．３にある。）

この部分（プログラム２．１〜２．）について，プログラムをよく読み，行なっていることについて理解しておくこと．（C言語を使いこなす上での最重要ポイントの一つであるポインタについての復習にもなる．）

(もし，2.1, 2.2のプログラムにあるdvector(), dmatrix()の定義について説明が必要であれば，質問してほしい．メモリの確保のされ方を図示できるようになっておくことが望ましい．）

プログラムを読むだけではなく，実行してみることをすすめる．

テキストの課題3.1のプログラムでテスト

• テキストのプログラム3.1を作成しビルドする。ビルドが成功すると「出力」欄に

1> prog3ex.vcxproj -> U:\prog3ee\prog3ex\Debug\prog3ex.exe
========== ビルド: 1 正常終了、0 失敗、0 更新不要、0 スキップ ==========

のようにメッセージがでる。exeファイルが置かれたフォルダを確認し、そのフォルダにデータファイル(input31.dat)を置き、実行してみる。

正常に動作すれば、output31.datができるはずである。VisualStudioで「ファイル→開く→ファイル」を選択し、output31.datを開いて内容を確認する。

ファイルに出力するfprintfの代わりに，コンソールに出力するprintfに書き換えて，結果を確認しても良い．

課題１-1

プリントテキスト（資料１）の課題1-1を行いなさい。（締め切り：10月12日）

• (追加指示) 時間があれば，R2とI0の数値をファイルにセーブし，Excel等で図にしてレポートに含めなさい．
• さらに時間があれば，（線形代数の復習であるが），正則ではない行列の例を自分で考えてデータファイルを作成し，結果が得られないことを確かめなさい．
• レポートは、テキストかPDF形式で提出してください。実習中に、Wordやテキストエディタ（Visual Studioでテキストファイルを開いても良い）を開いて、行ったことをその都度記入していくことをすすめる。
  • Wordであれば、「保存と送信」メニューにてPDFファイルへ変換できる。テキストファイルは、Windowsの印刷機能を用いて、PDFファイルへの出力を行う。

E-Learningのページへ提出してください．

課題1-2

プリントテキストの課題1-2を考えてみなさい．（締め切り： 10月12日）

（注）締切に時刻が書いてない場合は、原則として、その日の２３時５０分です。

補足

実際には

授業では、もっとも基本的な数値解法のみを実習したが、行列の特徴、たとえば対称行列とかほとんどの要素がゼロの「疎な行列」とかの場合ごとに効率的な方法が考案されてる。それを見極め、適切な数値計算ライブラリ関数を選択することが実際には大事な作業になる。

広く用いられているオープンソースの数値計算ライブラリとしてGSLがある。下記のWebページを参照してみてほしい。

https://www.researchgate.net/publication/288933699_How_to_use_GSL_with_Visual_C_Visual_C_deGSLwoshiufangfa

VC++特有なこと

コンパイルするときに、fopenをより安全な関数にしたとのことで、fopen_sの使用を推奨する警告がでる。今回は無視してもよいが、警告がでないようにするには次のようにする。

fopen_sでは、

  if ((fopen_s(&fin, "input31.dat", "r")) == 0)

同様に,fprintfはfprintf_sに、fscanfはfscan_sに変える。

Scilabで解く

Scilab (matlabでも同様）では次のようにlinsolve(A,b)という関数を用いるだけで解を得ることができます．

// 線型方程式の解を求める
A = [1, 2, 1, 1;...
     4, 5, -2, 4;...
    4, 3, -3, 1;...
    2, 1, 1, 3]

b = [1; 7; 12; -2];

[x0, kerA] = linsolve(A,b)  // Ax +b = 0 の解を求める。（符号に注意）

disp(x0)

matlabの場合

matlabでは、コメントは%で始まります。また、linsolveはAx=bを解く仕様になっているので、Scilabとは符号が逆になります。 継続行指定、コマンドの最後のセミコロンの要不要なども若干異なります。

% 線型方程式の解を求める
A = [1, 2, 1, 1;...
     4, 5, -2, 4;...
    4, 3, -3, 1;...
    2, 1, 1, 3];

b = [-1; -7; -12; 2];
% Ax=bを解く
[x0, kerA] = linsolve(A,b);

disp(x0)