風と翼 Windandwing は、c++の数値計算から風と翼の計算へ を目的とした数値計算を行っています。

第1章 c++ STL　および　関数の種類

数値計算で、c++ を利用した 計算法を記載する。
この章は、備忘録的な内容が主であるので、具体的な構成や区分けは、記載しながら、考えていく。

---

STLの適用例

1、曲線のゼロとなる点や交点の算出法
　曲線のゼロとなる点や交点(両曲線の差がゼロとなる点)の算出は、ゼロとなる点の曲線の前後を掛ると、負になる事より算出できる。この性質を利用して、算出する。

例　第3章　方程式求根

//区間を細分した時、ある区間内にゼロ点あるとその前後で、f(x)の値が変化する。
//区分が一つ増加部との積、v[n]*v[n+1]は必ず負値になる。その区間以外は全て正である。
//この考え方で、zero点区間をSTLを利用し求めて、解を算出する。
これを STL/transform で行う。

区間を [ a , b ] とする。

double h=(b-a)/n;

vector<double>v1(n+1,0.0);//　n+1の数は位置bでのf(b)値を含めるため。
for(int i=0;i<=n;i++) v1[i]=f_bisection()(a+i*h);

vector<double>v2(n+1,0.0); //f(a+(n-1)*h)*f(b)の積を求めるため。
transform(v1.begin(),v1.end()-1,v1.begin()+1,v2.begin(),multiplies<double>());//積をv2に代入。

// 負値になる区間を算出する。
vector<double>::iterator p1;// find_if()のitaratorである。
p1=find_if(v2.begin(),v2.end(),minusvalue);//p1で負値の個所を得る。
int near0=p1-v2.begin();


2、一つ飛びの繰り返し文 for 文の手法 (偶数番号あるいは奇数番号だけの for 文）

例　第5章 直交関数 多項式級数展開(Taylor 展開) 　sin() と cos() の計算

for(int L=k_Le;L>=0;L=L-2) // L=L-2 で、一つ飛びで減少させる。
{

　for(int i=L*(L+1)/2;i<m+L*(L+1)/2+1;i++)
　{
　fx_Le.push_back(hx_Le[i]);
　fx_Ch.push_back(hx_Ch[i]);
　}
　fx_Le.resize(fx_Le.size());
　fx_Ch.resize(fx_Ch.size());
}


3、n 行 n 列の行列を STL/vector で表現している時
　　1)単位行列 u の対角要素の 1.0 は、以下の方法で得る事が出来る。
　　　　vector<double>u;
　　　　u.assign(n*n,0.0);
　　　　for(int j0=0;j0<n*n;j0=j0+n+1) u[j0]=1.0;

　　2)正則行列が STL/vector ans とした場合の転置行列は、以下で得る事が出来る。
　　　int n=(int)sqrt((double)ans.size());//サイズ n*nのｎの算出法
　　　vector<double>transposed_determinant(n*n,0.0);
　　　for(int j=0;j<n;j++) for(int i=0;i<n;i++) transposed_determinant[i+j*n]=ans[j+n*i];

　　3)partial_sum による vector の全ての要素の積　これは行数が多い。従来法 ( for( i ) x*=y[i] ) が良い。
　　　vector<double>det_value(n,0.0);
　　　vector<double>temp_vector(n,0.0);
　　　for(int i=0;i<n;i++) temp_vector[i]=data0[i+i*n];
　　　partial_sum(temp_vector.begin(),temp_vector.end(),det_value.begin(), multiplies<double>());
　　　double det_val_d_partial_sum=det_value[n-1]; // det_val_d_partial_sum で全ての要素の積

　　　例、c++ 数値計算 8 行列式(余因子、逆行列)、行列の 2、行列式 7、連立一次方程式の解法による逆行列


4、m 行 m 列の行列 u と v の内積　STL / inner_product 利用による
　　STL / inner_product で、行列 u と v の内積を求める場合は、行列ｖを転置行列のように位置を変える必要
　　がある。この位置変換を行い、内積を求める計算を示す。 位置変換の目的が転置行列 transposed matrix )
　　生成ではないので、この位置変換の行列を chnged matrix 　と表示した。

　　int m=3;
　　vector<double>u(m*m,0.0);vector<double>v(m*m,0.0);vector<double>w(m*m,0.0);

　　u[0]=2.0;u[1]=3.0;u[2]=4.0;u[3]=5.0;u[4]=6.0;u[5]=7.0;u[6]=8.0;u[7]=9.0;u[8]=1.0;
　　v[0]=2.0;v[1]=5.0;v[2]=8.0;v[3]=3.0;v[4]=1.0;v[5]=6.0;v[6]=4.0;v[7]=7.0;v[8]=9.0;

　　vector<double>changed_v(m*m,0.0);//STL/inner product 利用のために aij を aji に位置変換
　　for(int j=0;j<m;j++) for(int i=0;i<m;i++) changed_v[i+j*m]=v[m*i+j];

　　for(int j=0;j<m;j++) for(int i=0;i<m;i++)
　　 w[i+j*m]=inner_product(u.begin()+j*m,u.begin()+j*m+m,changed_v.begin()+i*m,0.0);

　　ofstream out_w("w_value.csv",ios::out |ios::binary);//29,41,70,56,80,139,47,56,127
　　out_w<<"wの要素 "<<endl;
　　for(int i=0; i<m*m; i++) out_w<<"w["<<i<<"]=,"<<w[i]<<endl;
　　out_w.close();


5、収束条件　ノルム法
　　1）ノルム法
　　　ノルム法が良いと感じた。
　　
　　 double gg=0.0;//ノルム法
　　 for(int i=0; i<n;i++) gg=gg+fabs(a[(n+1)*i]-b[(n+1)*i]);
　　 if(gg<eps) break;

　　2）条件の論理積
　　　収束の条件を論理積で行ったが、コードが適切でないのか、最後まで、条件を満たす事が出来なかった。

　　for(int i=0; i<n-1; i++) // 論理積
　　{
　　if(fabs(a[(n+1)*(0+i)]-b[(n+1)*(0+i)])<eps && fabs(a[(n+1)*(1+i)]-b[(n+1)*(1+i)])<eps ))
　　break;
　　}


6、vector<complex<double>> u の (0.0,0.0) の 初期化
　　u(n) の 初期化の方法として、以下が適切である。
　　
　　予め、以下を宣言する。

　　const complex<double> real_unit(1.0,0.0),imag_unit(0.0,1.0),zero_unit(0.0,0.0);
　　
　　文中では、以下で、定義する。

　　vector<complex<double>>u(n,zero_unit);
　　vector<complex<double>>v(n,real_unit);
　　vector<complex<double>>w(n,imag_unit);


7、多次元配列の vector<complex<double>> v の (0.0,0.0) の 初期化

　7-1) complex<double> v2[m][n] について
　　　　vector<vector<complex<double>>> v2(m,vector<complex<double>>(n,zero_unit));

　7-2) complex<double> v3[m][n][p] について(以下 2 行は、一つの文である。）
　vector<vector<vector<complex<double>>>> v3(m, vector<vector<complex<double>>>(n,
　vector<complex<double>>(p,zero_unit)));

　7-3) 数値の代入
　for(int k=0;k<p;k++) for(int j=0;j<m;j++) for(int i=0;i<n;i++) v3[i][j][k]=u*v*w;


8、関数定義文とClass文に関して

　8-1)、関数定義文の例 本ホームページ c+数値計算 1 の三次方程式で使用
　complex<double> csum1(complex<double> d1,complex<double> d2){;return ( pow(d1,1/3.0) + pow (d2,1/3.0));}

　8-2)、Class 文による関数定義　本ホームページ c++ 数値計算 1 の三次方程式で使用
　class binary_csum2
　{
　　public: complex<double> operator()(complex<double> x,complex<double> y)
　　{
　　　return pow(x,1/3.0)+pow(y,1/3.0);
　　}
　};


9、文字列の数値変換方法の制約
　　入力を文字列として、それを数値に変換する方法を本ホームページ c++ 数値計算 1 の三次方程式に示し
　ている。その部分を下記に示す。しかし、これは、実数の場合のみに使用できて、 複素数 (1.0,2.0) のよ
　うな ( が ある場合には適用できない。(,,,,) の記号が変換出来ないためである。

　　文字列の数値変換　indata=strtod((const char*)string_data0.c_str(),NULL);//文字列を数値に変更。


　tring string_data0;//データを文字列で読み、strtodで数値に変換し、vectorに格納する。
　vector<double>a;//a[0],a[1],a[2]となる。
　double indata;

　ifstream input0("in3eq.csv", ios::in | ios::binary);//
　if(!input0.is_open())
　{
　　cout<<"Cannot open file　\n";exit(1);
　}
　while(!input0.eof()) //eof()関数で最後までデータを読込む。
　{
　　input0>>string_data0;
　　if(!input0.eof())
　　{
　　　indata=strtod((const char*)string_data0.c_str(),NULL);//文字列を数値に変更。
　　　a.push_back(indata);
　　}
　}
　input0.close();
　a.resize(a.size());

---

STLの種類

準備中

関数の種類

準備中



























　　　　　

2012年

この章の全体履歴
2012年1月　 　

2014年

この章の変更追加履歴
5)行列の積の記述追加と収束条件の方法を追記　2015/5/11
4)行列の積の方法を追記　2015/5/10
3)行列式のpartial_sum部に 1 行の式を追加　2014/5/14
2)行列式での経験を追記　2014/5/8
1)行列式での経験を記載　2014/5/4