
以下のプログラム
*********************************************
* 連立方程式の解法 ( ガウスの消去法 )
*********************************************/
#include <iostream> // for cout
#include <stdio.h> // for printf()
// 元の数定義
#define N 4 // 3
using namespace std;
/*
* 計算クラス
*/
class Calc
{
double a[N][N + 1];
// 各種変数
double d; // ダミー
int i, j, k; // LOOP インデックス
public:
// 連立方程式を解く(ガウスの消去法)
void calcGaussElimination();
};
/*
* 連立方程式を解く(ガウスの消去法)
*/
void Calc::calcGaussElimination()
{
// 係数
static double a[N][N + 1] = {
//{ 2.0, -3.0, 1.0, 5.0},
//{ 1.0, 1.0, -1.0, 2.0},
//{ 3.0, 5.0, -7.0, 0.0}
{ 1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{ 5.0, 4.0, 7.0, 1.0, -1.0},
{ 9.0, 7.0, 3.0, 5.0, 4.0}
};
// 元の連立方程式をコンソール出力
for (i = 0; i < N; i++) {
for (j = 0; j < N; j++)
printf("%+fx%d ", a[i][j], j + 1);
printf("= %+f\n", a[i][N]);
}
// 前進消去
for (k = 0; k < N -1; k++) {
for (i = k + 1; i < N; i++) {
d = a[i][k] / a[k][k];
for (j = k + 1; j <= N; j++)
a[i][j] -= a[k][j] * d;
}
}
// 後退代入
for (i = N - 1; i >= 0; i--) {
d = a[i][N];
for (j = i + 1; j < N; j++)
d -= a[i][j] * a[j][N];
a[i][N] = d / a[i][i];
}
// 結果出力
for (k = 0; k < N; k++)
printf("x%d = %f\n", k + 1, a[k][N]);
}
/*
* メイン処理
*/
int main()
{
try
{
// 計算クラスインスタンス化
Calc objCalc;
// 連立方程式を解く(ガウスの消去法)
objCalc.calcGaussElimination();
}
catch (...) {
cout << "例外発生!" << endl;
return -1;
}
// 正常終了
return 0;
}
を載せさせて頂いたプログラムのように配列ではなく、方程式の形のまま、Nを使わず、for文なし、scanfで座標を入力し、n次方程式を解くプログラムを作れるでしょうか?
#include <stdio.h>
int main(void) {
a *x1*x1 + b * x1 + c == 32;//(1)
d *x2*x2 + e * x2 + f == 40;//(2)
g *x3*x3 + h * x3 + i == 35;//(3)
d*x2*x2 - (a*x1*x1*(d*x2*x2 /a*x1*x1)) + e*x2 - (b*x1*(a*x2*x2 /d* x1*x1)) - f - (c*(a*x2*x2 /d*x1*x1));
printf("%d", d*x2*x2 - (a*x1*x1*(d*x2*x2 /a*x1*x1)) + e*x2 - (b*x1*(a*x2*x2 /d*x1*x1)) - f - (c*(a*x2*x2 /d*x1*x1)));
return 0;
}
A 回答 (15件中1~10件)
- 最新から表示
- 回答順に表示
No.15
- 回答日時:
9. リファクタリング
さて、最後です。
一応、ガウスの消去法のプログラムは完成しました。
ところで、前にも見た通り、「前進消去」と「後退代入」の中心アルゴリズムは「殆ど同じだ」と言う事を見ました。
そういうわけで、関数backward_substitution_auxは関数forward_elimination_auxを「コピペして」改造する、って方針を取ったわけですが。
ちょっと二つを並べてみましょう。
vector<vector<double>> forward_elimination_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
if (i == a.size()) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return forward_elimination_aux(a, u, k, i+1);
}
}
vector<vector<double>> backward_substitution_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
if (i == k) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return backward_substitution_aux(a, u, k, i+1);
}
}
この二つは本当に、「終了条件が違う」だけで全くと言って良い程「やってる内容は同じ」なんですね。
つまり言い換えると「無駄に関数が二つある」と言う事に他なりません。
まあ、気分の問題なんですが、このテの「機能が殆ど同じ」関数を二つ作るよりは、「まとめちゃって」、終了条件を「切り替える」機能を付けた方がマシなんじゃねーの、と言う話があります。
このテの作業を「リファクタリング」と呼びます。
注: 「ファクター: factor」は「まとめる」と言う意味で、リファクタリングはまさしく「まとめ直し」を意味します。ちなみに、数学での「因数分解」も英語で「ファクタリング: factoring」と言います。
どうせプログラミングしてるウチに「書き直し」は出てくるんで、書き直しをまずは怖がらない、と言うのが重要です。
そして、プログラムを書いてる時、「関数を書いていく」と言うのは単一の機能をそれぞれに「まとめて」書いていく事に他なりません。当然やっていくウチに「機能が重複していく」事もあり得るのです。それら「重複した機能」をまとめて直し(リファクタリング)していけば、汎用の関数になって、そのうち「個人のライブラリ」として重宝するようになるやもしれません。
ここでは、あくまで一例として、forward_elimination_auxとbackward_substitution_auxを「まとめて」同じ関数にしてしまいましょう。
まず今、終了条件が「二種類」ある事が分かっています。どうせ二種類なんで、もう一つ引数を増やして、bool値で「終了条件」が切り替わるようにしてみましょう。
vector<vector<double>> aux_engine(vector<vector<double>> a, vector<double>u, int k, int i, bool toggle = true);
そうすれば、終了条件を&&と||を使って、変数「toggle」を見て、合致した条件で行列aを返すように改造すれば、わざわざ補助関数を二種類も持たなくて良い事が分かります。
vector<vector<double>> aux_engine(vector<vector<double>> a, vector<double>u, int k, int i, bool toggle = true) {
// 終了条件を変数toggleを見るように変更する
if ((toggle == true && i == a.size()) || (toggle == false && i == k)) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return aux_engine(a, u, k, i+1, toggle); // 再帰するべき関数名を変更するのを忘れないように
}
}
これで終了、です。あとは関数forward_eliminationとbackward_substitutionからこいつを呼び出すようにすれば良い、です。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
// 補助関数
vector<vector<double>> aux_engine(vector<vector<double>> a, vector<double>u, int k, int i, bool toggle = true) {
if ((toggle == true && i == a.size()) || (toggle == false && i == k)) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return aux_engine(a, u, k, i+1, toggle);
}
}
// 前進消去
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n = 0) {
if (n == a.size()){
return a;
} else {
double pivot = a[n][n];
vector<double> result;
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;});
a[n] = result;
a = aux_engine(a, result, n, n+1); // 補助関数を呼び出す
return forward_elimination(a, n+1);
}
}
// 後退代入
vector<vector<double>> backward_substitution(vector<vector<double>> a, int n) {
if (n == 0) {
return a;
} else {
a = aux_engine(a, a[n], n, 0, false); // 補助関数をfalseを付けて呼び出す
return backward_substitution(a, n-1);
}
}
int main(int argc, char const* argv[]){
// 拡大係数行列
vector<vector<double>> a {
{1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{5.0, 4.0, 7.0, 1.0, -1.0},
{9.0, 7.0, 3.0, 5.0, 4.0}
};
vector<vector<double>> result;
result = backward_substitution(forward_elimination(a), a.size()-1);
for (auto v1: result) {
for (auto v2: v1) {
cout << v2 << " ";
}
cout << endl;
}
return 0;
}
と言うわけで「余計な部分」を削って短くして、本当に終了です。
「forやwhileを使わない」と言うのが分かりやすいか分かりづらいか、ってのは人に依りますが、取り敢えずこれでガウスの消去法は完成ですね。
お疲れ様でした。
No.14
- 回答日時:
続きです。
ガウスの消去法の全コードは、
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
// 前進消去の補助関数
vector<vector<double>> forward_elimination_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
if (i == a.size()) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return forward_elimination_aux(a, u, k, i+1);
}
}
// 前進消去
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n = 0) {
if (n == a.size()){
return a;
} else {
double pivot = a[n][n];
vector<double> result;
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;});
a[n] = result;
a = forward_elimination_aux(a, result, n, n+1);
return forward_elimination(a, n+1);
}
}
// 後退代入の補助関数
vector<vector<double>> backward_substitution_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
if (i == k) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return backward_substitution_aux(a, u, k, i+1);
}
}
// 後退代入
vector<vector<double>> backward_substitution(vector<vector<double>> a, int n) {
if (n == 0) {
return a;
} else {
a = backward_substitution_aux(a, a[n], n, 0);
return backward_substitution(a, n-1);
}
}
int main(int argc, char const* argv[]){
// 拡大係数行列
vector<vector<double>> a {
{1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{5.0, 4.0, 7.0, 1.0, -1.0},
{9.0, 7.0, 3.0, 5.0, 4.0}
};
vector<vector<double>> result;
result = backward_substitution(forward_elimination(a), a.size()-1); // ガウスの消去法
for (auto v1: result) {
for (auto v2: v1) {
cout << v2 << " ";
}
cout << endl;
}
return 0;
}
となります。
コンパイルして出力を見てみましょう。
➜ ~ ./a.out
1 0 0 0 1
0 1 0 0 3
-0 -0 1 0 -2
-0 -0 -0 1 -4
これは
x[1] = 1
x[2] = 3
x[3] = -2
x[4] = -4
を表してるんで、無事、「ガウスの消去法」を実装し終えた、って事を意味してますね。
「どういう表示を出力に回すのか」ってのは個人の好みなんで、main関数での「出力」をどういじるのか、ってのは任せましょう。ぶっちゃけ、そのへんは「計算」じゃないんで、どーでも良いです(笑)。
No.13
- 回答日時:
続きです。
8. 後退代入
さて、前項で、「前進消去」も「後退代入」も中心となるアルゴリズムは「殆ど同じになる」と言う事を見ました。
そして、こういう際には「関数の組み立て方針」も「ほぼ同じになる」と考えた方が見通しが良くなるわけです。
すなわち、二つの関数
vector<vector<double>> backward_substitution(vector<vector<double>> a, int n);
vector<vector<double>> backward_substitution_aux(vector<vector<double>> a, vector<double> u, int k, int i);
を作る、と言うのを作戦とし、また、「必要とする引数」も一致させる、と言う方針を取ります。そして関数backward_substitution_auxはbackward_substitutionの補助関数、とします。
では大本の後退代入用の関数、backward_substitutionから実装していきましょう。基本的にある引数nを取り、「1づつ引きながら」更新していき、nが0に到達した時点で、行列aを返して計算を脱出する、ってのを今まで見てきました。
これだけ、なんで大枠では次のように書けば良い、ってのがすぐ分かりますね。
vector<vector<double>> backward_substitution(vector<vector<double>> a, int n) {
if (n == 0) {
return a;
} else {
// ここで補助関数を呼び出す
return backward_substitution(a, n-1);
}
}
前進消去と違って、ピボットで割る、なんつー余計な作業がないので、これが大枠なんだ、ってのがすぐ分かるでしょう。ついでに言うと、これは「カウンタである第2引数」の値を1づつ小さくしていってるだけ、なんで、実は「何もしてない」に等しいです。
取り敢えずiostreamでも使って、「nの値が一つづつ小さくなっていってる」って事だけ確認してみて下さい。
注: つまり関数backward_substitutionは、「ループの外側」だけを関数として「敢えて」実装しただけだ、と言う事になる。前に書いた通り、これも末尾再帰なんで、コンパイラによっては「末尾再帰最適化」でアセンブリレベルでは単なるジャンプ構文に変換される。これが意味する事は、当然「人力で」for文に書き換える事が可能だと言う事を示唆していて、for文での「多重ループ」の外側が、実はこの関数backward_substitutionの「やってる事」だと言う事を意味してる。
次に補助関数backward_substitution_auxを実装しますが、これが実は「殆ど前進消去の計算の中核」(具体的に言うと補助関数forward_elimination_aux)」のアルゴリズムと変わらない、と言う事を見ました。
と言う事はforward_elimination_auxを「コピペして」改造すれば一丁上がり、って事を意味しています。
終了条件を良く見てforward_elimination_auxを「改造」しましょう。
vector<vector<double>> backward_substitution_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
// 当然関数名は変える
if (i == k) {
// 終了条件を変更
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return backward_substitution_aux(a, u, k, i+1); // 再帰する対象関数の名前を変更する事を忘れないように
}
}
基本的にはこれで実装終了、です。あとは、減数vectorの指定(当然n番目の行数なんでa[n]です)、被減数の番号は0から始まる、って事だけ留意してbackward_substitutionからbackward_substitution_auxを呼び出せば良い。
vector<vector<double>> backward_substitution(vector<vector<double>> a, int n) {
if (n == 0) {
return a;
} else {
a = backward_substitution_aux(a, a[n], n, 0); // ここで補助関数を呼び出す
return backward_substitution(a, n-1);
}
}
と言うわけで全実装が終了しました。
あとは「ガウスの消去法」を実装すれば良い、んですが、これは簡単ですよね。ガウスの消去法とは、
「「行列aに前進消去を適用」したものに後退代入を適用」したもの
です。そして、後退代入の第二引数nには「行列の行数から1引いた数」を与えれば良い、って事が分かっています(つまり、vectorの「大きさ」が4だとすれば、最後の要素は「3番目」ですね)。
つまり、ガウスの消去法は
backward_substitution(forward_elimination(a), a.size()-1);
と表現出来ます。
続きます。
No.12
- 回答日時:
続きです。
// forward_elimination
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
// forward_eliminationの補助関数
vector<vector<double>> forward_elimination_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
if (i == a.size()) {
return a;
} else {
vector<double> v = a[i];
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];});
a[i] = v;
return forward_elimination_aux(a, u, k, i+1);
}
}
// 前進消去
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n = 0) {
if (n == a.size()){
return a;
} else {
double pivot = a[n][n];
vector<double> result;
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;});
a[n] = result;
a = forward_elimination_aux(a, result, n, n+1); // 補助関数を呼び出して、結果を行列aに代入
return forward_elimination(a, n+1);
}
}
int main(int argc, char const* argv[]){
// 拡大係数行列
vector<vector<double>> a {
{1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{5.0, 4.0, 7.0, 1.0, -1.0},
{9.0, 7.0, 3.0, 5.0, 4.0}
};
vector<vector<double>> result;
result = forward_elimination(a);
for (auto v1: result) {
for (auto v2: v1) {
cout << v2 << " ";
}
cout << endl;
}
return 0;
}
コンパイルして実行してみると次のようになりますね。
➜ ~ ./a.out
1 -2 3 -4 5
0 1 14 -11 19
-0 -0 1 -0.857843 1.43137
-0 -0 -0 1 -4
はい、お疲れ様です。無事に「前進消去」が実装されました。
・・・とは言っても、まだ「後方代入」が残ってます。そしてこれが実装し終わらないと「ガウスの消去法」は完成しません(ここで止めて「ガウス・ジョルダンの消去法でござい」って言う手もありますけどね・笑)。
さぁ、どうするか。
実はここまでが「実装としては大変」なんですが、この後の「後方代入」は「前進消去」より遥かに簡単で、前進消去を実装してしまえばあとは鼻クソ同然なんです。
では、取り敢えず、今度は件のページ( https://www.mk-mode.com/octopress/2013/09/24/cpp … )で解説されている「後方代入」のアルゴリズムを確認してみましょう。
【後方代入】
1. b'[n]からb'[1]に向かって以下の操作を繰り返す。
(a) b'[i]を初期値とし、b'[i+1]からb'[n]までa'[i][j]×b'[j]の値を減算する。
(b) (a)で求められた値をa'[i][i]で除算する。
ふ〜む・・・・・・。何書いてるのかサッパリ分かりませんね(笑)。間違ってんじゃねぇの、これ(苦笑)?
そもそも (b)の行程でまたピボットで除算する必要があるんか・・・・・・?
ちょっと整理してみましょうか。
今、お題の行列aに対して前進消去し終わった連立方程式があるとしましょう。つまり、先程のアレですね。
x[1] -2*x[2] +3*x[3] -4*x[4] = 5
x[2] +14*x[3] -11*x[4] = 19
x[3] -0.8*x[4] = 1.43137
x[4] = -4
1回目の後方代入のプロセスを考えます。現時点、拡大係数行列として考えると、対角線上の係数a[i][i]は全て1になっていますね。当然です、そのために前進消去したわけですよね。
1回目の後方代入のプロセスだと減数用の行は4行目のx[4] = -4、最初の被減数行は1行目のx[1] -2*x[2] +3*x[3] -4*x[4] = 5 となっていて、そのx[4]を「消去」する為に4行目×(-4)を計算して1行目から引くわけです。
x[1] -2*x[2] +3*x[3] = 9
x[2] +14*x[3] -11*x[4] = 19
x[3] -0.8*x[4] = 1.43137
x[4] = -4
次は2行目から4行目の「−11倍」を引く。
x[1] -2*x[2] +3*x[3] = 9
x[2] +14*x[3] = 63
x[3] -0.8*x[4] = 1.43137
x[4] = -4
そのまた次は3行目から4行目の「-0.8倍」を引く。
x[1] -2*x[2] +3*x[3] = 9
x[2] +14*x[3] = 63
x[3] = 2.289213
x[4] = -4
これで1回目のプロセスが終了して、次は減数行が3行目になって、被減数である1行目と2行目に対して演算を行っていきます・・・・。
これで十分じゃね(笑)?
つまり、中心作業としては、先程見た、(transformを使った)
transform(a[i].begin(), a[i].end(), a[k].begin(), a[i].begin(), [](double x, double y){return x - y*a[i][k];})
と全く同じ事をやってる、と言う事ですね。結局i行とk行の間で減算を行ってる、ってのはまるで変わりありません(要するにi行-k行×a[i][k]を計算してる)。
前進消去との違いが何なのか、と言うと、一回の計算プロセスに於いてプログラム上では、
i) 減数行が1行目からではなく、k行目から始まってる(「前進消去」だと0->1->...と進んで行って、kがa.size()に到達した際に「計算終了」だが、逆に「後退代入」だとk->k-1->...とすすんで行って、kが0の時に「計算終了」となる)。
ii) 被減数行であるi行目は1づつ増えていくが、終了地点が違う(「前進消去」だと行数(a.size())に到達して終了、だけど「後退代入」の場合、iがkに到達した時点で「計算終了」となる)
殆どアルゴリズムが変わらないんで、「クッソ簡単に実装出来る」と言う予感しか無いでしょう。それでいいのです、その通りなんです。
続きます。
No.11
- 回答日時:
続きです。
さて、アルゴリズムの(b)をもう一度見てみましょう。
(b) i行-k行×a[i][k]/ピボットを行う
これを「行列の要素一個一個に対してやっていく」と言う事を件のページ( https://www.mk-mode.com/octopress/2013/09/24/cpp … )で書いているわけですが。
一方、実はここまで実装したforward_eliminationでは、既に
「各行の全要素をそれぞれの行に対応したピボットで割る」
って事を終わらせているんですね。
実際、これは我々が「ガウスの消去法」を手計算でやる時と同じです。元のプログラムと違って、人間が「ガウスの消去法」を計算する際に、「係数一つ一つ」に対して「バラバラに」考える、って事はあり得ないわけですよ。必ず「一行を」まとめて考えます。少なくとも中学校で連立方程式を習う際に、「式毎に」考える、って事を訓練されてるわけですね。
transformでvectorを「まとめて扱う」にはそういう利点があります。「人が手計算で扱う」ように行列を扱えるんで、結局「見通し」が良くなるのです。
ちょっと計算過程を(途中からですが)確認してみます。
# (a)と(b)のプロセスが一回終わったトコロ
x[1]-2*x[2]+3*x[3]-4*x[4]=5
x[2]+14*x[3]-11*x[4]=19
14*x[2]-8*x[3]+21*x[4]=-26
25*x[2]-24*x[3]+41*x[4]=-41
# 3行目-2行目×14を実行 -> 3行目-2行目×係数a[3][2]
x[1]-2*x[2]+3*x[3]-4*x[4]=5
x[2]+14*x[3]-11*x[4]=19
-204*x[3]+175*x[4]=-292
25*x[2]-24*x[3]+41*x[4]=-41
# 4行目-2行目×25を実行 -> 4行目-2行目×係数a[4][2]
x[1]-2*x[2]+3*x[3]-4*x[4]=5
x[2]+14*x[3]-11*x[4]=19
-204*x[3]+175*x[4]=-292
-374*x[3]+316*x[4]=-516
つまり、2周目のプロセスを見ると、手計算では「ここでピボットで割る」必要は全くなくって、単に
(b) i行目-k行目×a[i][k]
を行えば良い、と言う事を示してるのです。
従って、transformを使う際には、基本的には
transform(a[i].begin(), a[i].end(), a[k].begin(), a[i].begin(), [](double x, double y){return x - y*a[i][k];})
か、それに類する式を使えば良い、と言う作戦が成り立ちます(実際はC++のtransformの仕様のせいで、上の式そのままでは使えないんですが)。
そして、前進消去の為の関数forward_elimination(本)とその補助関数(補)の間には次のプロセスが成り立ちます。
本->補:(1行目-0行目×a[1][0]->2行目-0行目×a[2][0]->3行目-0行目×a[3][0]: 終了)->本->補:(2行目-1行目×a[2][1]->3行目-1行目×a[3][1]: 終了)->本->補:(3行目-2行目×a[3][2]: 終了)->本(終了)
ここで分かる事はまずは次の二つです。
i) 補助関数は最低でも引数に行列aとiとkの3つを取らないとならない。
ii) 補助関数は「1回」(b)のプロセスをやれば良いので、「計算全部」をやろうと考えない(本関数forward_eliminationが全体の面倒を見てくれる)。
i) に関して言うと、補助関数のカウンタとしては、補助関数が一回呼び出される毎に「kは全く動かない」です。動くのはiだけ、ですね。そしてiが「行列aのサイズ」に到達すればその時点での「計算が終了」した、って事なんで、制御は本関数であるforward_eliminationに戻ります。
注: なお、このようにkが一回の呼び出し内中で「全く動かない」場合、一々引数として持つのはどーよ、ってのがあるんですが、関数型言語やPythonなんかの場合、この補助関数をローカル関数にしてしまって、kをローカル関数内部から外側の関数に置かれた「外部変数」として参照させる、なんて事が出来るんですが、C++だと色々試してみたけど上手く行きませんでした。何故なら、恐らく、C++は「レキシカルスコープ」を持ってても「レキシカル変数」を持ってない為だと思われます。従って、上手い具合に「固定された」変数を参照させるには、補助関数の引数を長くしてでもその中に入れておかないといけないみたいです(誰か、上手くやる方法を知ってたら教えて下さい・笑)。
kが動かない、と言う事は当然、減数側の行、つまり減数vector、a[k]も動かない、って事ですね。動くのは被減数vector、a[i]だと言う事です。しかも、ループ内でa[k]を再計算させるのも面倒ですし、a[k]が何なのか、と言うと単純にforward_elimination内にあるvector、resultの事です。と言うわけで、こいつも引数として貰っちまいましょう(kとa[k]が別々にある、ってのもみっともないんですが、前述の通り、「レキシカル変数」がC++には存在しないようなんで、ある程度の「みっともなさ」はどうしようもない模様です)。
さて、最後に実装上の注意点です。先程書いた通りtransformの仕様がちょっと力量不足のせいか、多元vectorの「内側を弄る」ってのが出来ない模様です(出来るのかもしれませんが、発見出来ませんでした・笑)。
じゃあどうするか、と言うとa[i]を代入するvectorを作っておいて、操作して、またa[i]に戻す、と言うようなやり方をするしか(少なくとも今の時点では)ないようです。関数型言語やPythonだと「一行で済む」処理がC++だとどうしても「増えちゃう」んで、困ったモンですね(笑)。
と言うわけで、補助関数は次のようになります。
// forward_eliminationの補助関数
vector<vector<double>> forward_elimination_aux(vector<vector<double>> a, vector<double> u, int k, int i) {
if (i == a.size()) {
return a; // iがaの大きさに到達したらaを返す
} else {
vector<double> v = a[i]; // a[i]を代入するvectorを作る
transform(v.begin(), v.end(), u.begin(), v.begin(), [&](double x, double y) {return x - y * a[i][k];}); // v - u*a[i][k]を計算
a[i] = v; // 結果vをa[i]に代入
return forward_elimination_aux(a, u, k, i+1); // iを1増やして再帰
}
}
そして、減数vectorが例えばa[0]の時、非減数vectorはa[1]からスタートしてますね。と言う事は、k = nの時、i = n+1から補助関数が呼び出される、ってのも「ガウスの消去法」の計算プロセスから自明です。
従って、ここまでのプログラムは次のようになるでしょう。
続きます。
No.10
- 回答日時:
ここからちょっと余談です。
ここまで見てきたように再帰を書くのは鼻クソなんですが、同時に「再帰を書くのは凄く簡単なのになんでforやwhileなんてあるんだろ?」って疑問が出てくると思います。
理由はたった一つ。「再帰はforやwhileに比べると一般的には実行効率が悪い」と思われてるからです。半分真実で半分嘘なんですが。
より正確に言うと、「プログラミング言語の種類によって、再帰が効率的な言語もあれば非効率な言語もある」と言う事です。で、CやC++なんかは後者のプログラミング言語なんですね。
元々C++の前身であるC言語はハードウェアの動きに密接に関係してて、「アセンブリ言語の難解な記述を抽象化する」為に生まれました。アセンブリ言語レベルだと「繰り返し」はかいつまんで言うと、
「条件を満たせば次の処理に移行し、そうじゃなければ前にある特定の場所へ"ジャンプ"しろ」
と言うようなカンジで書くので、プログラミングがクソメンド臭いわけです。
クソメンド臭い記述を簡単にする為「抽象化」としてfor文やwhile文が生まれたんですね。
一方、「再帰」は今まで見てきた通り、原理的には「関数呼び出し」を駆使します。それでその「関数呼び出し」と言うのは通常コストがかかるんです。連続して関数を呼び出しまくる「再帰」は低レベルだとあまり効率が良くないんですね。
それがC/C++等で「よっぽどの事が無ければ」再帰を使わないで済むなら使わないようにしたい、と言う理由です。
他方、関数型言語の場合、そもそも「関数型言語」自体が抽象化の塊のようなモノで、ハードウェアの動きとはあまり密接に関係していません。元々存在自体がC/C++から見ると「非効率の塊」のような存在で(笑)、そもそも再帰しようがしまいが「コストがかかってる」ってのが基本大前提なんです。そしてこのテの言語だとハードウェアに於ける「負担」よりも「記述が簡単で済む」方が重要視されてるのです。C/C++で100行書かなきゃならない処理が10行で済むのなら、それが「関数型言語の正義」なんですね。そしてC/C++で100日かかる作業が10日で書き上がるのなら、そっちの方がハードウェアの負担よりマシだろう、と言う考え方が背景にあります。
両者ともに設計思想が根本的に違うわけです。
と言うわけで、今回は「forやwhileを使わないで」と言う要求なんで敢えて再帰でやってみましたが、効率的に言うと、一般的には「サイテーです」って言われる部類でしょう(笑)。
なお、例示のコードのように、
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n = 0) {
.....
return forward_elimination(a, n+1); // 自分自身「だけ」returnされる
}
と「returnで返す関数が自分自身しかない」再帰を特別に「末尾再帰」と呼びます。
✗末尾再帰じゃない例
int foo(int n = 0) {
.....
return foo(n+1) + 1; // 自分自身だけじゃなくって他の操作が加わってるのでこれは末尾再帰じゃない
}
で、C++の言語仕様では要求されてない(筈)ですが、コンパイラによっては「末尾再帰で書かれたコード」をアセンブリレベルでは「ジャンプ構文にコンパイルして」くれて、要するに、for文やwhile文で書かれたコードと全く同じアセンブリを吐き出してくれるスグレモノがあったりします。
こういう機能を「末尾再帰最適化」と言って、多くの関数型言語ではこれを利用してるケースが多い為、「再帰は必ずしも効率が悪いとは限らない」論拠になっています。
C++のコンパイラだと、GCCや、恐らくCLANGは「最適化を施せ」と指令したら末尾再帰最適化を実行してくれますね(少なくとも、Cで書かれた末尾再帰は最適化出来ます)。一方、Microsoft のVisual C++だと・・・どうだろ、末尾再帰最適化、って機能は恐らく無いんじゃないでしょうか。
本当は、「あるコンパイラ特有の機能をアテにして」プログラムを書くべきではない、とは思うんですが、コンパイラによっては必ずしも「末尾再帰なら」非効率なプログラムにはならないよ、と言う例がまさしく今回のプログラムになっています。
なお、7.のコードで<iostream>を使って「forward_eliminationの第2引数nが増えていく様」を観察して、再帰を視覚的に体験してみて下さい。これも「実験」です。
と言う辺りで本題に戻りましょう。
前進消去のプログラムを完成させるには、アルゴリズムで紹介されてる(b)を実装しないとなりません。
現在、7.の段階で「与えられた拡大係数行列の係数行列の対角線上を全て1にする」事に成功しました。
問題は「対角線上を1にする」その途中で(b)の命令をぶっこまないとならないわけです。
いつそれを行うか?つまり例えば、1行目のピボットによる演算と2行目のピボットによる演算の「途中で」(b)の処理を行わないとならないんですね。一体どこに「処理」をぶっこむべきでしょうか?
そうですね、自分自身をreturnする「その直前」に(b)の作業をぶっこまないとならない。
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n = 0) {
if (n == a.size()){
return a;
} else {
double pivot = a[n][n];
vector<double> result;
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;});
a[n] = result;
// 再帰の直前の「ここ」で何らかの処理を行う
return forward_elimination(a, n+1);
}
}
上だと「ここ」って簡単に書いてますが、(b)の作業と言うのは、ぶっちゃけた話、「前進消去」で一番面倒くさい部分です。要するに「1行で処理が終わる」程簡単ではない、って事ですね。
そうなると、ここを仮に「1行で処理が終わる」ように書くとしたら、ここで別に組み立てた関数を呼び出すしかない、と言う意味になります。要するにforward_eliminationの為の「補助関数」をどっかに作って、ここでそれを「呼び出せば」前進消去の実装は「終了した」と言えるわけです。
あるいは、フツーのC++プログラミングで考えるなら、先の余談で説明した通り、「実行効率を重視して」ここでfor/while文を書いて「必要な処理を行う」と言うのがスタンダードなテクニックになるでしょう(練習としてはこれも良い題材です)。
いずれにせよ、(b)自体も「繰り返し」が必要そうなんで、「forやwhileを使わない」のなら、やはり再帰関数に頼るしかないでしょう。
まずはちょっとアルゴリズム(b)を精査しましょう。
続きます。
No.9
- 回答日時:
続きです。
7. 再帰プログラミング
さて、6. の実験を通して、関数forward_eliminationの第二引数の値を0 <= i <= 3の範囲で色々と変えてみたら任意の行(i行目)の全ての要素が対応するピボットa[i][i]で割られてる事が分かりました。
と言う事は、
「「「「aの0行目に関数forward_eliminationを適用した結果」の1行目に関数forward_eliminationを適用した結果」の2行目に関数forward_eliminationを適用した結果」の3行目に関数forward_eliminationを適用した結果」
が「全ての行の全要素が、それぞれの行に対応したピボットで割られてる」結果を得る行為になると分かるでしょう。
これ、文章で書けばメンド臭いんですが、実際はforward_eliminationの第二引数の値を「一つづつ増やしていけば良い」と言ってるだけなんで、クッソ簡単なんですね。for文やwhile文で考えるより簡単です。
つまり、方針としては
関数forward_eliminationの中で、関数forward_eliminationの第二引数を1増やした状態をreturnする
と言う事を行います。このテクニックは「ある関数の中で自分自身を使う」ので、(自分に再び帰ると言う意味を込めて)「再帰プログラミング」と呼ばれるのです。
ちなみに、「マトモなプログラミング言語」なら、大体どんな言語でもこの「再帰プログラミング」は行えるので、別に特殊なテクニックでも何でもありません(関数型言語界隈だと特に人気があるプログラミング方法ではありますが)。
ところで、永遠に第二引数を増やしていくと当然与えられたvectorの「大きさ」を超えてしまうので(例えば6. の「実験」で第二引数に4を与えるととんでもない結果になるのが分かるでしょう)、「適当なトコで止めないと」いけません。じゃあ、いつ止めるのか?
これも簡単で、forward_eliminationの第二引数は0から初めて「vector aのサイズ(全行数)になったトコ」で止めれば良いので鼻クソですね(ちなみに、あるvector、aのサイズ/大きさはa.size()で取得出来ます)。その時点で「操作」し終わった二次元vector、aを返せば済むのです。
と言うわけで再帰プログラミングを施したforward_eliminationは次のようになります。
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n) {
if (n == a.size()){
return a; // nがvector aの大きさに到達したらaを返す(終了条件)
} else {
double pivot = a[n][n];
vector<double> result;
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;});
a[n] = result;
return forward_elimination(a, n+1); // 自分自身の第二引数を1増やして返す
}
}
元のforward_eliminationを「ちょっとだけ改造したら」こうなりますね。これで「再帰プログラミングになる」んで何て再帰は簡単なんでしょう。宣言した通り、「forやwhileより簡単で」全く悩む必要がありません。鼻クソですね。
せっかくなんでもうちょっとだけ改造してみます。行列aの操作は「0行目」から始める事は猿にだって分かります。従って、forward_eliminationを使う際にわざわざnに0を与える、ってのは「メンド臭い」ですね。
こういう場合、C++では「デフォルト引数」と言う機能を使います。デフォルト引数は関数が呼び出された際に「自動的に引数に値を与える」指定で、次のようなフォーマットで記述します。
// 例
int foo(int n = 0) { ....
上の例だと関数fooが呼び出された際に、nに自動的に0が代入されて呼び出されます。便利な機能ですね!(実はC言語にはこの機能がありません)
これを使うとプログラムは次のように書けます。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n = 0) {
// デフォルト引数を使って n = 0 から計算をスタート
if (n == a.size()){
return a;
} else {
double pivot = a[n][n];
vector<double> result;
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;});
a[n] = result;
return forward_elimination(a, n+1); // ただし、再帰呼出しではnに「1を加えた」引数を渡さないとならない
}
}
int main(int argc, char const* argv[]){
// 拡大係数行列
vector<vector<double>> a {
{1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{5.0, 4.0, 7.0, 1.0, -1.0},
{9.0, 7.0, 3.0, 5.0, 4.0}
};
vector<vector<double>> result;
result = forward_elimination(a); // 関数forward_eliminationを呼び出すには行列aを与えるだけで、n=0から計算がスタートする
for (auto v1: result) {
for (auto v2: v1) {
cout << v2 << " ";
}
cout << endl;
}
return 0;
}
実行結果は次のようになりますね。
1 -2 3 -4 5
-0.4 1 1.6 -0.6 1.8
0.714286 0.571429 1 0.142857 -0.142857
1.8 1.4 0.6 1 0.8
まだ「前進消去」と言う段階ではありませんが(まだアルゴリズムのbを実装してないから当たり前)、取り敢えずi行目の「全要素」がa[i][i]で割られてるのが分かるでしょうか。検算してみて下さい。
そして係数行列の対角線上の要素が全部1になっていますね。これはa[i][i]がピボットの為、演算結果のresult[i][i]は当然a[i][i]/a[i][i]なので1になるのです。
続きます。
No.8
- 回答日時:
続きです。
6. 前進消去を書いてみよう
と言うわけでここからは参照サイト( https://www.mk-mode.com/octopress/2013/09/24/cpp … )で説明されているアルゴリズムを参考にして「ガウスの消去法」を実装していきます。
なお、「アルゴリズムは参考に」しますが、「プログラミングスタイルは参考に」しません。個人的には、オブジェクト指向にはそれほど明るくないんですが、件のサイトの「基底クラス(計算クラス)を作って」「計算クラスを継承して」って実装方針があまり意味があるとは思えないんですよね。ぶっちゃけ、「かなり無駄な実装をしてる」と思います。
また、ここでは「オブジェクト指向でプログラミング」していくより「関数型プログラミング(に近い形)」で実装して行った方がSTLの「能力」を活かせる、と言う判断です。行列(二次元vector)を受け取り適当な操作をして行列(二次元vector)を返す関数を二つ書いて(内情は「前進消去」と「後退代入」)、それらを組み合わせて最終的に「ガウスの消去法」にしてしまう方が遥かに見通しが良い、と思います(まあ、この辺も極論個人の好み、ですけどね)。
と言うわけで、まずは「前進消去」を実装していきましょう。かつ、ここが実は最大の山場で、「ガウスの消去法」のプログラミングの70%〜80%の「労力」はここに集中してますし、ここさえ出来れば後は鼻クソレベルだと言っても過言ではありません。
まずはアルゴリズムの確認から。
【前進消去】
ピボットを1行1列からn-1行n-1列に移しながら以下の操作を繰り返す。
(a) a[k][k]をピボットにしi行についてa[i][k]/a[k][k]を求め、i行-k行×a[i][k]/a[k][k]を行う
なるほど、です。
ちなみに、この(a)は「二つの作業」が一つに書かれてるので、次のようにして見てみればより分かりやすいですね。
【前進消去】
ピボットを1行1列からn-1行n-1列に移しながら以下の操作を繰り返す。
(a) a[k][k]をピボットにしi行についてa[i][k]/ピボットを求める
(b) i行-k行×a[i][k]/ピボットを行う
さて、まずは(a)の作業を考えましょう。と言うより(a)は実は次のように読み替える事が出来ます。
(a) a[k][k]をピボットにしてi行として表現されたvectorの全要素をピボットで割れ
はい、ここで前項目で書いた「vectorをマッピング操作しろ」と同じ事を指定されてるのか分かるでしょうか。
そうすると、今、二次元vectorであるaが次のように定義されているとして、
vector<vector<double>> a{
{1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{5.0, 4.0, 7.0, 1.0, -1.0},
{9.0, 7.0, 3.0, 5.0, 4.0}
};
aのi番目もvectorなんでそれを対象としたマッピングを考えれば良いわけです。
例: a[0] => {1.0, -2.0, 3.0, -4.0, 5.0} を意味してる。これが行列の「0行目」を表す。
それでピボットに関して言うと、例えば0行目のピボットはa[0][0]で1、1行目のピボットはa[1][1]で5、・・・ってなるのは分かるでしょう。
そうするとまず「書かなければならない」操作用のコードは、transformとラムダ式を利用して「以下のような部分解を築き上げる」ってのがぼんやりと分かってきます。
double pivot = a[n][n]; // ピボットを定義
vector<double> result; // aのi番目のvectorを操作した結果を代入するvectorを定義する
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x/pivot;}); // aのi番目の全要素をピボットで割る
a[n] = result; // 結果をaのi番目に代入して書き換える
そして二次元vector、aを受け取って、操作済みのvector aを返す関数として、次のような雛形がでっち上げられる、ってのが分かるでしょう。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
vector<vector<double>> forward_elimination(vector<vector<double>> a, int n) {
double pivot = a[n][n]; // ピボットを定義
vector<double> result; // aのi番目のvectorを操作した結果を代入するvectorを定義する
transform(a[n].begin(), a[n].end(), back_inserter(result), [&](double x){return x / pivot;}); // aのi番目の全要素をピボットで割る
a[n] = result; // 結果をaのi番目に代入して書き換える
return a; // 処理が終わったaを返す
}
int main(int argc, char const* argv[]){
// 拡大係数行列
vector<vector<double>> a {
{1.0, -2.0, 3.0, -4.0, 5.0},
{-2.0, 5.0, 8.0, -3.0, 9.0},
{5.0, 4.0, 7.0, 1.0, -1.0},
{9.0, 7.0, 3.0, 5.0, 4.0}
};
vector<vector<double>> result; // 結果代入用の二次元vector
result = forward_elimination(a, 1); // 取り敢えずaの1行目を操作した結果を見てみる
for (auto v1: result) {
for (auto v2: v1) {
cout << v2 << " ";
}
cout << endl;
}
return 0;
}
上のコードの実行結果は以下のようになりますね。
1 -2 3 -4 5
-0.4 1 1.6 -0.6 1.8
5 4 7 1 -1
9 7 3 5 4
行列の二行目(vectorを使ったプログラミングでは一行目)の全要素がピボットa[1][1] = 5.0で割られているのが分かるでしょうか?従って結果の行列の要素result[1][1]が1になっていますね。
ところで、main関数内の
result = forward_elimination(a, 1);
のforward_eliminationの第二引数を0 <= i <= 3の範囲で色々と変えてみると「行列aのi行目の全要素がピボットa[i][i]で割られてる」のが分かると思います。あるいは実験目的にscanfを使ってforward_eliminationの第二引数に代入してみても良いかもしれません。いずれにせよ、実験してみて下さい。
さて、これで拡大係数行列aの「任意の行」のすべての要素を対応するピボットで「割る」事が出来る、と言う事は分かりました。しかし実際は、最終的には「すべての行の要素をそれぞれの行に対応するピボットで割る」と言うのが目的です。任意の行に付いて実験は成功ですが、「すべての行のケース」を何とかしないといけないのです。
ここで、for文やwhileに慣れていれば「あ、繰り返し構文使えば良いな」って事になるんで、もちろんその実装方法で構わないでしょう(練習目的で自分で書いてみて下さい)。
でも、今回はあくまで「forやwhileを使わないで」と言う提案なんで、それに伴った書き方をしてみましょう。
4つ目の武器を「再帰プログラミング」と言います。
続きます。
No.7
- 回答日時:
続きです。
5. ラムダ式(無名関数)
さて、先程のLispやPythonのコードを見ると、mapの他にもう一つ共通してる「言葉」があります。
;; Lisp(Scheme)
(map (lambda (x) (+ x 1)) ls)
(map (lambda (x y) (+ x y)) lsa lsb)
# Python
list(map(lambda x: x + 1, ls))
list(map(lambda x, y: x + y, lsa, lsb))
そうですね、lambda(ラムダ、と読みます)と書かれた部分です。ちょっと抜書しますか。
;; Lisp(Scheme)
(lambda (x) (+ x 1))
(lambda (x y) (+ x y))
# Python
lambda x: x + 1
lambda x, y: x + y
これがマッピング関数に「何をさせたいか(どんな計算を配列に適用させたいか)」を記述する部分で、lambdaと言う「キーワード」が指し示すように「ラムダ式」と呼びます。
ラムダ式は機能的には関数そのもの、なんですが、「関数としての名前」は持ってません。そんなわけで、「無名関数」等と呼んだりもします。「使い捨て目的の関数」なんですね。
例えば、ここではLispでの例を見ますが、
(map (lambda (x) (+ x 1)) ls) ;; lsは '(1, 2, 3, 4, 5) とする
と記述された際、ラムダ式の引数xには1から始まって順番に5まで次々と値が与えられていきます。
当然作用は x + 1、つまり引数に1を足して返される。次は2が入って3が返される・・・となっていくわけですね。
つまりここでのラムダ式の引数xは「lsの中身が順次投入されていく」と対象になっているわけです。
もう一つの例の場合、
(map (lambda (x y) (+ x y)) lsa lsb) ;; lsaは '(1, 2, 3, 4, 5)、lsbは'(6, 7, 8, 9, 10) とする
ラムダ式の引数はxとyが二つ、です。引数xはlsaに対応、引数yはlsbに対応してます。当然lsaとlsbからは「同時に」次々とラムダ式内の引数x, yに要素がぶち込まれていって、計算結果、つまりxとyが足し合わされた結果が返されていきます。こうやってリスト操作が成されていくわけです。
昨今の言語だと大体ラムダ式が入ってる事が多くて、ラムダ式は今や大人気ですね(C#や最新のJavaなんかにも入ってます)。
そんなラムダ式をC++も持っています。先程の例示コード
// C++
transform(ls.begin(), ls.end(), back_inserter(result), [](int x){return x + 1;});
transform(lsa.begin(), lsa.end(), lsb.begin(), lsa.begin(), [](int x, int y){return x + y;});
のこの部分がラムダ式です(lambdaと言うキーワードがないですが)。
// C++
[](int x){return x + 1;}
[](int x, int y){return x + y;}
C++ではこういう形で「どんな計算をさせたいか」を記述します。そしてC++でのラムダ式のフォーマットは基本的には次のようになっています。
[](引数の型 引数, ...) { 色んな処理...; return 結果;}
C++の通常の関数の書き方とほぼ同じですが、
i) 関数の返り値の型と関数の名前を記述しない
ii) 代わりに[]と言う記号が使われてる
と言う二点だけ、が違いますね。
基本的には、transformとラムダ式の組み合わせで行列の演算(ガウスの消去法)の面倒な演算の半分くらいは軽減出来ます。
ところで、C++のラムダ式(とマッピング関数の組み合わせ)には注意点が一つだけあります。
例えばLispやPythonで行う次のような演算、
;; Lisp(Scheme)での例
(let ((i 2))
(map (lambda (x) (* x i)) '(1 2 3 4 5))) ;; リスト(配列)の各要素をi倍する
'(2 4 6 8 10) ;; 結果
# Pythonでの例
>>> i = 2
>>> list(map(lambda x: x * i, [1, 2, 3, 4, 5])) # リスト(配列)の各要素をi倍する
[2, 4, 6, 8, 10] # 結果
を見てみます。両者とも配列(リスト)の各要素をi倍する、と言う計算を行ってて、その「i」自体はマッピング関数の「外側」で定義されています。
このような計算をする場合、LispでもSchemeでも問題なく計算を終える(つまり外部にある変数iを参照可能)んですが、このコードをC++に移植しようとするとちょっとした障害に突き当たります。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
int main(int argc, char const* argv[]){
int i = 2;
vector<int> ls {1, 2, 3, 4, 5};
vector<int> result;
transform(ls.begin(), ls.end(), back_inserter(result), [](int x){return x * i;}); // これだとコンパイルエラーが起きる
for (auto v: result) {
cout << v << " ";
}
cout << endl;
return 0;
}
上のC++でのコードをコンパイルしようとすると、思いっきりエラーが吐き出されるでしょう(一回実験してみて下さい)。
原因が何か、と言うと、これはマッピングであるtransformの外にある変数iを内部にあるラムダ式が参照出来ない、と言う事に由来します。LispやPythonでは出来るのにC++だと「自動で外側にある変数を探す事が出来ない」(これを専門的に「スコープが違う」と表現します)。
実はC++の場合、外部の変数を参照する場合、その変数を「明示的に参照せよ」と指定する必要があります。
例えば、上のコードの場合、「transformの外側にあるiを参照せよ」と言うんで次のように書くか、
transform(ls.begin(), ls.end(), back_inserter(result), [i](int x){return x * i;}); // ラムダ式冒頭の[]内にラムダ式内から参照すべき変数名を記述する
あるいは、「とにかく何でも良いんでtransformの外側にある変数を全部参照せよ」と指定するか、です。
transform(ls.begin(), ls.end(), back_inserter(result), [&](int x){return x * i;}); // ラムダ式冒頭の[]内にアンパサンド(&)を記述するとtransformの外にある変数をラムダ式内で取り敢えず全部参照出来る
ザーッとWeb上検索すると、C++のラムダ式は後者のスタイルで書かれる事が多い模様です。あるいはC++でラムダ式を記述をする場合、とにかく[&]で書き始める、ってのも面倒が無くて良いかもしれません。
いずれにせよ、下がLisp/Pythonでの上記のコードのC++への「正しい」移植です。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
int main(int argc, char const* argv[]){
int i = 2;
vector<int> ls {1, 2, 3, 4, 5};
vector<int> result;
transform(ls.begin(), ls.end(), back_inserter(result), [&](int x){return x * i;}); // こうすればC++コンパイラを正しく通る
for (auto v: result) {
cout << v << " ";
}
cout << endl;
return 0;
}
と言うわけでガウスの消去法をやっつける武器が取り敢えず3つ揃いました。復習ですが、その3つとは
・ vector
・ transform
・ ラムダ式
です。この3つの基本の武器装備で「forやwhileを使わずにガウスの消去法を書く」事に挑戦しましょう。
あと、もう2つほど「武器」が必要ですが、それはおいおい書いていきます。
ではやっていきましょう。続きます。
No.6
- 回答日時:
続きです。
4. transform
さて、行列を扱う際に何がメンド臭いのか、と言うと、
「要素一つ一つを順繰りに操作していかないといけない」
から、ですね。これに尽きます。
結局forが二重ループになるとか三重ループになるとか言うのは、要素を一つ一つ手繰って行ってるからですよね。
ところが考えてみると、「プログラム上は」配列の配列だったりするわけで、内側の配列を「まとめて」操作出来れば随分と見通しが良くなります。要素一個一個で考えていかなくて良い。
ところで、「関数型プログラミング」が得意な「関数型言語」と呼ばれる種類のプログラミング言語は、このテの「配列を一気に処理する」機能に長けてます。
例えばLisp言語族と言われるプログラム言語のファミリーではリスト(配列のようなもの)を一気に処理する事なんて良くやります。
;; Scheme言語の例
(define ls '(1 2 3 4 5)) ;; リスト(配列)を定義
(map (lambda (x) (+ x 1)) ls) ;; リスト(配列)の各要素に1を足す
'(2 3 4 5 6) # 結果: lsの各要素に1が足される
同様に、最近流行りの女の子、もとい、プログラミング言語のPythonなんかでも同様の事が出来ます。
# Pythonの例
>>> ls = [1, 2, 3, 4, 5] # リスト(配列)を定義
>>> list(map(lambda x: x + 1, ls)) # リスト(配列)の各要素に1を足す
[2, 3, 4, 5, 6]
>>>
細かいトコは今は端折りますが、いずれにせよ「forやwhileを使わなくても」望んだ結果が出ているのが分かるでしょう。
これが「抽象化」の力です。
そして、上記の2つのコードで「同じ名前」があるのが分かるでしょうか。そうですね、mapです。
このテのリスト(配列)を一気に処理する関数を「マッピング関数」と呼びます。そしてこれが「イテレータ」の一種です。
注: ちなみに、上記のPythonの例はちょっと「古い」書き方で、今Python使ってる人たちは大方
[x + 1 for x in ls]
と書く方が好みです。が「マッピング」を説明する為に敢えて「古い」書き方をしています。
さて、この「マッピング関数」があれば、随分と行列を「扱う」のがラクになると思いませんか?
STLでこの「マッピング関数」にあたるのがtransformです。インクルードするヘッダはalgorithm。
#include <algorithm>
そして上記の例をC++で書くと、こんなカンジになるでしょう。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
int main(int argc, char const* argv[]){
vector<int> ls {1, 2, 3, 4, 5};
vector<int> result;
transform(ls.begin(), ls.end(), back_inserter(result), [](int x){return x + 1;});
for (auto v: result) {
cout << v << " ";
}
cout << endl;
}
Lisp(Scheme)やPythonに比べるといささか長くなりますが、まあ、そのへんはコンパイラ型言語であるC++ならでは、ですね。出力をキチンと指定せなアカンですから。
まだ説明してないトコもありますが、重要なのは次の部分です。
transform(ls.begin(), ls.end(), back_inserter(result), [&](int x){return x + 1;});
これがLispやPythonでmapがしてるのとほぼ同じ事をやってます。
フォーマットは次の通り。
transform(vectorの始点指定, vectorの終点指定, back_inserter(代入目的のvector指定), したい処理);
LispやPythonのmapと違うところもあります。差異は次のようなモノです。
i) vectorを直接引数に取らない。代わりにvectorの始点と終点を引数に取るが、beginとendと言うメソッド(メンバ関数?)(注: これをC++では「何故か」イテレータと呼ぶ)でそれらを簡単に指定できる。
ii) 特にLispはそのまま返り値としてリストを返すが、C++のtransformは結果を別のvectorに「代入」する(別に代入目的のresultと言うvectorを宣言してるのはそのため)。
iI)に関して言うと、resultに「代入」する為にback_inserterと言う(C++用語で言う)「挿入イテレータ」を用いています。こいつを使う為にiteratorと言うヘッダをインクルードしないといけないのです。
さて、「やりたい処理」をどう書くか、と言うのは取り敢えず後回しにして、transformのもう一つの使い方を見てみましょう。
LispやPythonのmapは実は可変長引数を取れて、2つ以上のリスト(配列)を取って次のような操作が可能です。
;; Lisp(Scheme)での例
(define lsa '(1 2 3 4 5)) ;; 配列その1を設定
(define lsb '(6 7 8 9 10)) ;; 配列その2を設定
(map (lambda (x y) (+ x y)) lsa lsb) ;; 2つ配列の要素を順番に足し合わせる
'(7 9 11 13 15) ;; 結果
# Pythonでの例
lsa = [1, 2, 3, 4, 5] # 配列その1を設定
lsb = [6, 7, 8, 9, 10] # 配列その2を設定
>>> list(map(lambda x, y: x + y, lsa, lsb)) ;; 2つの配列の要素を順番に足し合わせる
[7, 9, 11, 13, 15]
>>>
これに似たような事がC++で出来ないか、と言う事ですね。これも出来れば「ガウスの消去法」で行同士の要素を「順番に」引いたり出来るわけじゃないですか。
取り敢えず結論から書くと、C++でtransformを使うと次のように書けます。
#include <iostream>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
int main(int argc, char const* argv[]){
vector<int> lsa {1, 2, 3, 4, 5};
vector<int> lsb {6, 7, 8, 9, 10};
transform(lsa.begin(), lsa.end(), lsb.begin(), lsa.begin(), [](int x, int y){return x + y;});
for (auto v: lsa) {
cout << v << " ";
}
cout << endl;
}
ぶっちゃけ、何だか良くわからない状態になっていますが(笑)、取り敢えず
i) フォーマットは
transform(一つ目のvectorの始点, 一つ目のvectorの終点, 二つ目のvectorの始点、一つ目のvectorの始点, やりたい処理)
になってる
ii) LispやPythonの場合、(リストにせよオブジェクトにせよ)「新しい値を返す」が、transformの場合、一つ目のvectorが計算結果によって「書き換えられる」
と言う二つが注意事項です。
特にii)に関して言うと、最初のlsaと「計算後の」lsaとは中身が全く変わっています。
これは自分でiostreamを使ってでも調べてみて下さい。
次は「やりたい処理」の書き方、です。
続きます。
お探しのQ&Aが見つからない時は、教えて!gooで質問しましょう!
似たような質問が見つかりました
- C言語・C++・C# C言語の課題が出たのですが自力でやっても分かりませんでした。 要素数がnであるint型の配列v2の並 3 2022/11/19 17:41
- C言語・C++・C# C 言語の Gauss Jordan 法について 2 2022/12/28 11:16
- C言語・C++・C# 至急教えてください。プログラミングの問題です。 最初に正の整数nの入力を受け付け、次に分数の分子と分 1 2022/07/19 17:03
- C言語・C++・C# 至急お願いします。プログラミングの問題です。 最初に正の整数nの入力を受け付け、次に分数の分子と分母 3 2022/07/19 17:09
- 数学 線形代数の対称行列についての問題がわからないです。 2 2023/01/08 14:59
- C言語・C++・C# 至急教えてください。プログラミングの問題です。 malloc関数を使ってください!お願いします! 最 1 2022/07/21 09:28
- その他(プログラミング・Web制作) Pythonにおける物理のシミュレーションでの単位変換について 2 2023/06/02 17:11
- Visual Basic(VBA) VBAプログラミング 2 2022/11/27 12:07
- C言語・C++・C# LU分解法のピボッティングについて(C言語/gcc-9) 3 2022/07/11 23:10
- C言語・C++・C# LU分解法のピボット選択機能実装について(C言語・gcc-9) 1 2022/07/22 15:20
関連するカテゴリからQ&Aを探す
おすすめ情報
デイリーランキングこのカテゴリの人気デイリーQ&Aランキング
-
プログラムでの数字につく”f”の...
-
float型とdouble型の変数の違い...
-
C言語を実行すると-infが出てき...
-
c言語で、繰り返し文の中で、0....
-
至急です! マクロ定義で #defi...
-
C言語 関数プロトタイプ宣言の...
-
C言語の型による処理速度の違い
-
Cプログラミングの問題です。ニ...
-
C言語初心者 構造体 課題について
-
c言語のプログラミングについて...
-
C言語(プログラミング)関連の質...
-
関数におけるif文とreturn文に...
-
型変換のitoaのaって?
-
浮動小数点数が表示されないん...
-
doubleの変数にintとintの割り...
-
DWORDの警告
-
C 開放してるのにエラー(doubl...
-
浮動小数点の定数
-
以下のプログラム ************...
-
c言語でユーザ関数を利用して複...
マンスリーランキングこのカテゴリの人気マンスリーQ&Aランキング
-
プログラムでの数字につく”f”の...
-
float型とdouble型の変数の違い...
-
C言語を実行すると-infが出てき...
-
C 開放してるのにエラー(doubl...
-
c言語で、繰り返し文の中で、0....
-
doubleの変数にintとintの割り...
-
至急です! マクロ定義で #defi...
-
C言語の型による処理速度の違い
-
C言語 関数プロトタイプ宣言の...
-
2次方程式の解を求めるプログ...
-
関数におけるif文とreturn文に...
-
doubleは常に%lfとするべきなのか
-
int とdoubleの比較
-
C言語のプログラムで#include<m...
-
C言語で-23乗を取り扱うには
-
データ数の多い構造体配列
-
指数の表示
-
C言語のpow関数の不具合
-
c言語のプログラミングについて...
-
c言語のコンパイルエラー canno...
おすすめ情報
最初の条件は無かったことに配列を使わず、 a *x1*x1 + b * x1 + c == 32のように表す。かつ、繰り返しを使わずに表現できないでしょうか?
もし配列を使わなくてはならないならば配列を使うかつ、forやwhileなどの繰り返しを使わないで表現して頂けないでしょうか?