
C言語で連立一次方程式の解き方の解き方を教えて!
こういった悩みにお答えします.
本記事の信頼性
- リアルタイムシステムの研究歴12年.
- 東大教員の時に,英語でOSの授業.
- 2012年9月~2013年8月にアメリカのノースカロライナ大学チャペルヒル校コンピュータサイエンス学部(2021年の世界大学学術ランキングで20位)で客員研究員として勤務.C言語でリアルタイムLinuxの研究開発.
- プログラミング歴15年以上,習得している言語: C/C++,Solidity/Vyper,Java,Python,Ruby,HTML/CSS/JS/PHP,MATLAB,Assembler (x64,ARM).
- 東大教員の時に,C++言語で開発した「LLVMコンパイラの拡張」,C言語で開発した独自のリアルタイムOS「Mcube Kernel」をGitHubにオープンソースとして公開.
こういった私から学べます.
本記事は行列を習得していることを前提とします.
-
-
【C言語】2次元配列で行列の四則演算(足し算,引き算,掛け算,割り算)【逆行列】
こういった悩みにお答えします. こういった私から学べます. 目次1 行列の四則演算(足し算,引き算,掛け算,割り算)2 行列の足し算と引き算3 行列の掛け算4 行列の割り算(逆行列の掛け算)4.1 2 ...
続きを見る
C言語で連立一次方程式の解法
C言語で連立一次方程式の解法を紹介します.
具体的には,以下の解法を紹介します.
- ガウスの消去法(掃き出し法)
- クラメルの公式(クラメルの法則)
ガウスの消去法(掃き出し法)
ガウスの消去法(掃き出し法)は,連立一次方程式を解くための多項式時間アルゴリズムです.
特に,行基本変形を用いて行列を行階段形に変形することをガウス・ジョルダンの消去法と呼びます.
行列が行階段形であるとは,以下の2つが成立することを意味します.
- 0でない成分を持つ行(少なくとも1つの成分が0でない行)が,0しか成分に持たない行よりも上に位置すること(もし0成分だけからなる行が存在する場合,それらは行列の最下部に配置)
- 主成分(行の最も左にある0でない成分で,ピボットとも呼ばれる)が,その行の上にある行の主成分よりも真に右側に位置すること(主成分は必ず1でなければならない場合もある.)
ガウス・ジョルダンの消去法のコードは以下になります.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 |
/* * Author: Hiroyuki Chishiro * License: 2-Clause BSD */ #include <stdio.h> #define N 3 int main(void) { double a[N][N + 1] = { {1.0, -2.0, 3.0, 5.0}, {2.0, -1.0, 1.0, 6.0}, {1.0, 3.0, -5.0, 2.0} }; double p, d; int i, j, k; for (i = 0; i < N; i++) { p = a[i][i]; for (j = 0; j < (N + 1); j++) { a[i][j] = a[i][j] / p; } for (j = 0; j < N; j++) { if (i != j) { d = a[j][i]; for (k = i; k < (N + 1); k++) { a[j][k] = a[j][k] - d * a[i][k]; } } } } for (i = 0; i < N; i++) { printf("x%d = %f\n", i + 1, a[i][N]); } return 0; } |
実行結果は以下になります.
こちらと同じ結果になることを確認しましょう.
1 2 3 4 5 |
$ gcc gauss_jordan.c $ a.out x1 = 6.000000 x2 = 17.000000 x3 = 11.000000 |
クラメルの公式(クラメルの法則)
クラメルの公式(クラメルの法則)は,未知数の数と方程式の数が同じで,かつ一意的に解ける線型方程式系を明示的に解く行列式の公式です.
クラメルの公式は逆行列が存在しない場合(det(A)が0の場合)は利用できないことに注意して下さい.
クラメルの公式のコードは以下になります.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 |
/* * Author: Hiroyuki Chishiro * License: 2-Clause BSD */ #include <stdio.h> #include <stdlib.h> #include <math.h> #define N 3 #define DBL_EPSILON 2.2204460492503131e-16 double get_det(double a[N][N]) { double det = a[0][0] * a[1][1] * a[2][2] + a[0][1] * a[1][2] * a[2][0] + a[0][2] * a[1][0] * a[2][1] - a[0][2] * a[1][1] * a[2][0] - a[0][1] * a[1][0] * a[2][2] - a[0][0] * a[1][2] * a[2][1]; return det; } int main(void) { double a[N][N + 1] = { {1.0, -2.0, 3.0, 5.0}, {2.0, -1.0, 1.0, 6.0}, {1.0, 3.0, -5.0, 2.0} }; double det = a[0][0] * a[1][1] * a[2][2] + a[0][1] * a[1][2] * a[2][0] + a[0][2] * a[1][0] * a[2][1] - a[0][2] * a[1][1] * a[2][0] - a[0][1] * a[1][0] * a[2][2] - a[0][0] * a[1][2] * a[2][1]; double b[N][N]; double x[N]; int i, j, k; if (fabs(det) <= DBL_EPSILON) { fprintf(stderr, "Error: cannot find inverse matrix.\n"); exit(1); } for (i = 0; i < N; i++) { for (j = 0; j < N; j++) { for (k = 0; k < N; k++) { if (i == k) { b[j][k] = a[j][N]; } else { b[j][k] = a[j][k]; } } } x[i] = get_det(b) / det; } for (i = 0; i < N; i++) { printf("x%d = %f\n", i + 1, x[i]); } return 0; } |
実行結果は以下になります.
1 2 3 4 5 |
$ gcc cramer.c $ a.out x1 = 6.000000 x2 = 17.000000 x3 = 11.000000 |
まとめ
C言語で連立一次方程式の解法として,ガウスの消去法(掃き出し法)とクラメルの公式(クラメルの法則)を紹介しました.
どちらも良く使うので,きちんと理解しておきましょう!
方程式の解の公式を学びたいあなたはこちらからどうぞ.
-
-
【C言語】一次方程式,二次方程式,三次方程式,四次方程式の解の公式
こういった悩みにお答えします. こういった私から学べます. 目次1 C言語で一次方程式,二次方程式,三次方程式,四次方程式の解の公式2 一次方程式の解の公式3 ニ次方程式の解の公式3.1 ニ次方程式の ...
続きを見る
C言語を独学で習得することは難しいです.
私にC言語の無料相談をしたいあなたは,公式LINE「ChishiroのC言語」の友だち追加をお願い致します.
私のキャパシティもあり,一定数に達したら終了しますので,今すぐ追加しましょう!
独学が難しいあなたは,元東大教員がおすすめするC言語を学べるオンラインプログラミングスクール5社で自分に合うスクールを見つけましょう.後悔はさせません!