解 x 含 長さ n ー 割 当 配列 IMSL_FACTOR, int **p_pvt, f_complex **p_factor (出力)
int **p_pvt (出力)
分解 ッ 列 含 長さ n 配列へ インタ ア 戻 必要 領域 imsl_c_lin_sol_gen によ 割 当 ます 一般にint *p_pvt が宣言さ &p_pvt が引数 し 使用さ ます
f_complex **p_factor (出力)
列 ッ 交換 A LU 分解 含 大 さn × n 配列へ インタ ア
戻 必要 領域 imsl_c_lin_sol_gen によ 割 当 ます こ 配列 下
角部分 L 構成す に また 角部分 U 構成す に必要 情報 含 います 一般にf_complex *p_factor が宣言さ &p_factor が引数 し 使用さ ま す
IMSL_FACTOR_USER, int pvt[], f_complex factor[] (入力/出力) int pvt[] (入力/出力)
分解 ッ 列 含 大 さ n ー 割 当 配列 f_complex factor[] (入力/出力)
A LU 分解 含 大 さ n × n ー 割 当 配列 こ 配列 下 角部分
L 構成す た に 角部分 U 構成す た に必要 情報 含 います
IMSL_SOLVE が指定さ こ パ ータ 入力 そ 他 場合 出力に
ます A が必要 け factor a 同 領域 共有します
IMSL_FAC_COL_DIM, int fac_col_dim (入力)
A LU 分解 含 列次元
フ ： fac_col_dim = n IMSL_INVERSE, f_complex **p_inva (出力)
行列 A 逆行列 含 大 さ n × n 配列へ インタ ア 戻 必要 領域 imsl_c_lin_sol_gen によ 割 当 ます 一般に f_complex *p_inva が 宣言さ &p_inva が引数 し 使用さ ます
IMSL_INVERSE_USER, f_complex inva[] (出力)
A 逆行列 含 大 さ n× n ー 割 当 配列
IMSL_INV_COL_DIM, int inva_col_dim (入力)
A 逆行列 含 配列 列次元
フ ： inva_col_dim = n IMSL_CONDITION, float *cond (出力)
行列 A L1 条件数 推定値 含 へ インタ こ プ ン プ
ン IMSL_SOLVE_ONLY 共に使用す こ ませ
IMSL_FACTOR_ONLY
部分 ッ 交換法 A LU 分解 計算します IMSL_FACTOR_ONLY 使用す 場 合 IMSL_FACTOR 又 IMSL_FACTOR_USER い が必要に ます 引数 b そ 時 無視さ imsl_c_lin_sol_gen 戻 値 NULL に ます
IMSL_SOLVE_ONLY
imsl_c_lin_sol_gen によ 前も 計算さ た LU 分解が与え た Ax = b 解
ます フ Ax = b 解 imsl_c_lin_sol_gen によ 指定さ ます
IMSL_SOLVE_ONLY 使用す 場合 引数IMSL_FACTOR_USER が必要に 引数 a 無視さ ます
IMSL_INVERSE_ONLY
行列 A 逆行列 計算します IMSL_INVERSE_ONLY 使用す 場合 IMSL_INVERSE 又 IMSL_INVERSE_USER が必要に ます 引数 b そ 時無視さ
imsl_c_lin_sol_gen 戻 値 NULL に ます
説明
関数 imsl_c_lin_sol_gen 複素係数行列 A 持 線形方程式 解 ます こ 最初に L-1A = U 部分 ッ 選択 A LU 分解 計算します F プ ン引数 IMSL_FACTOR 計算さ
行列 p_factor します 角行列 U F 角に格納さ ます F 厳密 下 角に
1 n1 n1 1 1
L LPL P 用い L-1 再構成す に必要 情報が含ま います
因子 P i L i 部分 ッ 交換によ 決定さ ます P i 行 i p_pvt[i-1] 交換さ た単 位行列 す Li 列 i 対角 下に挿入さ た Fji (j = i + 1, …, ni) 単位行列 す
線形方程式 解 2 よ 簡単 方程式 y = L-1 b x = U-1 y 解くこ によ 得 ます 線
形方程式 また 逆行列 解が得 A L1 条件数 推定値が Dongarra そ 他。1979年) 同
ア 使用し 計算さ ます 推定さ た条件数が 1/ε ここにε ン精度 よ 大 い 場合 警告 ッ ージが出ます こ A 非常に小さい変化が解 x 大 い変化 生 こ 示し
います 関数 imsl_c_lin_sol_gen 分解 角部分 U が ロ 対角要素 持 ー
ます
例題 1
ここ 次 3 元連立方程式 解 ます
(1 + i) x
1+ (2 + 3i) x
2+ (3 − 3i) x
3= 3 + 5i (2 + i) x
1+ (5 + 3i) x
2+ (7 − 5i) x
3= 22 + 10i (−2 + i) x
1+ (−4 + 4i) x
2+ (5 + 3i) x
3= −10 + 4i
#include <imsl.h>
f_complex a[] = {{1.0, 1.0}, {2.0, 3.0}, {3.0, -3.0}, {2.0, 1.0}, {5.0, 3.0}, {7.0, -5.0}, {-2.0, 1.0}, {-4.0, 4.0}, {5.0, 3.0}};
f_complex b[] = {{3.0, 5.0}, {22.0, 10.0}, {-10.0, 4.0}};
int main() {
int n = 3;
f_complex *x;
/* Solve Ax = b for x */
x = imsl_c_lin_sol_gen (n, a, b, 0);
/* Print x */
imsl_c_write_matrix ("Solution, x, of Ax = b", 1, n, x, 0);
}
出力
Solution, x, of Ax = b
1 2 3 ( 1, -1) ( 2, 4) ( 3, -0)
例題 2
こ 例題 転置問題 AH x = b 解い 部分 ッ 交換 A LU 分解 返します こ 例題 最初 例題 異 解 配列 インプログ に割 当 ます
#include <imsl.h>
f_complex a[] = {{1.0, 1.0}, {2.0, 3.0}, {3.0, -3.0}, {2.0, 1.0}, {5.0, 3.0}, {7.0, -5.0}, {-2.0, 1.0}, {-4.0, 4.0}, {5.0, 3.0}};
f_complex b[] = {{3.0, 5.0}, {22.0, 10.0}, {-10.0, 4.0}};
int main() {
int n = 3, pvt[3];
f_complex factor[9];
f_complex x[3];
/* Solve ctrans(A)*x = b for x */
imsl_c_lin_sol_gen (n, a, b, IMSL_TRANSPOSE, IMSL_RETURN_USER, x,
IMSL_FACTOR_USER, pvt, factor, 0);
/* Print x */
imsl_c_write_matrix ("Solution, x, of ctrans(A)x = b", 1, n, x, 0);
/* Print factors and pivot sequence */
imsl_c_write_matrix ("LU factors of A", n, n, factor, 0);
imsl_i_write_matrix ("Pivot sequence", 1, n, pvt, 0);
}
出力
Solution, x, of ctrans(A)x = b
1 2 3 ( -9.79, 11.23) ( 2.96, -3.13) ( 1.85, 2.47)
LU factors of A
1 2 3 1 ( -2.000, 1.000) ( -4.000, 4.000) ( 5.000, 3.000) 2 ( 0.600, 0.800) ( -1.200, 1.400) ( 2.200, 0.600) 3 ( 0.200, 0.600) ( -1.118, 0.529) ( 4.824, 1.294) Pivot sequence
1 2 3 3 3 3
例題 3
こ 例題 最初 例題 3 × 3 行列 A 逆行列 計算し 同 方程式 解 ます 行列積 C =
A-1A がチ ック し 計算さ ます 近似結果 C = I が得 ます
#include <imsl.h>
f_complex a[] = {{1.0, 1.0}, {2.0, 3.0}, {3.0, -3.0}, {2.0, 1.0}, {5.0, 3.0}, {7.0, -5.0}, {-2.0, 1.0}, {-4.0, 4.0}, {5.0, 3.0}};
f_complex b[] = {{3.0, 5.0}, {22.0, 10.0}, {-10.0, 4.0}};
int main() {
int n = 3;
f_complex *x;
f_complex *p_inva;
f_complex *C;
/* Solve Ax = b for x */
x = imsl_c_lin_sol_gen (n, a, b,
IMSL_INVERSE, &p_inva, 0);
/* Print solution */
imsl_c_write_matrix ("Solution, x, of Ax = b", 1, n, x, 0);
/* Print input and inverse matrices */
imsl_c_write_matrix ("Input A", n, n, a, 0);
imsl_c_write_matrix ("Inverse of A", n, n, p_inva, 0);
/* Check and print result */
C = imsl_c_mat_mul_rect ("A*B",
IMSL_A_MATRIX, n,n, p_inva, IMSL_B_MATRIX, n,n, a, 0);
imsl_c_write_matrix ("Product, inv(A)*A", n, n, C, 0);
}
出力
Solution, x, of Ax = b
1 2 3 ( 1, -1) ( 2, 4) ( 3, -0)
Input A
1 2 3 1 ( 1, 1) ( 2, 3) ( 3, -3) 2 ( 2, 1) ( 5, 3) ( 7, -5) 3 ( -2, 1) ( -4, 4) ( 5, 3) Inverse of A
1 2 3 1 ( 1.330, 0.594) ( -0.151, 0.028) ( -0.604, 0.613) 2 ( -0.632, -0.538) ( 0.160, 0.189) ( 0.142, -0.245) 3 ( -0.189, 0.160) ( 0.193, -0.052) ( 0.024, 0.042)
Product, inv(A)*A
1 2 3 1 ( 1, -0) ( -0, -0) ( -0, 0) 2 ( 0, 0) ( 1, 0) ( 0, -0) 3 ( -0, -0) ( -0, 0) ( 1, 0)
警告 ー
IMSL_ILL_CONDITIONED 入力行列 非常に悪条件 す L1 条件数 逆数
推定値 “rcond” = # す 解 恐 く正確 あ
ませ
重大 ー
IMSL_SINGULAR_MATRIX 入力行列 特異 す