IMSL_FACTOR, int **p_pvt, float **p_factor 。入力/出力) int **p_pvt 。入力/出力)
分解 ッ 列 含 長さn 配列へ インタ ア 戻 必要 領域 imsl_f_lin_sol_gen_band によ 割 当 ます 一般的に int *p_pvt が宣言さ
&p_pvt が引数 し 使用さ ます float **p_factor 。入力/出力)
列 ッ 交換 A LU 分解 含 大 さ (2nlca + nuca +1) xn 配列へ
インタ ア 戻 必要 領域 imsl_f_lin_sol_gen _band によ 割
当 ます 一般的に float *p_factor が宣言さ &p_factorが引数 し 使用さ
ます
IMSL_FACTOR_USER, int pvt[], float factor[] 。入力/出力) int pvt[] 。入力/出力)
分解 ッ 列 含 大 さ n ー 割 当 配列 float factor[] 。入力/出力)
A LU 分解 含 大 さ (2nlca+nuca+1) ×n ー 割 当 配列 こ 配列
厳密 下 角部分 L 構成す に 角部分 U 構成す に必要 情報が 含ま ます A が必要 け factor a 最初 (2nlca + nuca + 1) ×n
領域 共有します
IMSL_SOLVE が指定さ こ パ ータ 入力 す そ 他 場合 出力に ま
す
IMSL_CONDITION, float *condition 。出力)
行列 A L1 条件数 推定値 含 へ インタ こ プ ン プ
ン IMSL_SOLVE_ONLY 共に使用す こ ませ
IMSL_FACTOR_ONLY
部分 ッ 交換法 A LU 分解 計算します IMSL_FACTOR_ONLY 用い 場合 IMSL_FACTOR 又 IMSL_FACTOR_USER い が必要に ます 引数 b 無視 さ imsl_f_lin_sol_gen_band 戻 値 NULL に ます
IMSL_SOLVE_ONLY
imsl_f_lin_sol_gen_band によ 事前に計算さ た LU 分解が与え Ax = b
解 ます フ Ax = b 解 imsl_f_lin_sol_gen_band によ 指定さ
います IMSL_SOLVE_ONLY 使用す 引数 IMSL_FACTOR_USER が必要に 引数 a 無視さ ます
IMSL_BLOCKING_FACTOR, int block_factor 入力
こ ロッキング因子 block_factor 32 よ 大 く設定し いけませ フ ：block_factor = 1
説明
関数 imsl_f_lin_sol_gen_band 実帯係数行列 A 持 線形方程式 解 ます こ 最初に Du
Crozそ 他 1990年 ロック化 LU 分解ア に基 い A LU 分解 計算します ベ
3 BLAS 呼び出し イン イン ープ 置 換えます ロッキング因子 block_factor
フ １ すが 32 越え い任意 正 値に ッ す こ が可能 す
よ 簡単 2 方程式 y = L-1b x = U-1y 解くこ によ 線形方程式 解 得ます 線形方程式
又 逆行列 解 求 A L1 条件数 推定値が Higham(1988年)にあ Hager 方法に対す
Higham 修正 使用し 計算さ ます 推定した条件数が 1/ε ここ ε ン精度 よ 大
い場合 警告 ッ ージが出ます こ A 非常に小さい変化が解 x 大 い変化 生 こ 示し います 分解 角部分 U が ロ 対角要素 持 場合 関数 imsl_f_lin_sol_gen_band
ー ます
例題 1
ここ 次 4元方程式 解くこ こ 関数 最も簡単 使用法 説明します
2x
1– x
2= 3
-3x
1+ x
2+ 2x
3= 1 -x
3+ 2x
4= 11 2x
3+ x
4= -2
#include <imsl.h>
int main () {
int n = 4;
int nuca = 1;
int nlca = 1;
float *x;
/* Note that a is in band storage mode */
float a[] = {0.0, -1.0, -2.0, 2.0, 2.0, 1.0, -1.0, 1.0,
-3.0, 0.0, 2.0, 0.0};
float b[] = {3.0, 1.0, 11.0, -2.0};
x = imsl_f_lin_sol_gen_band (n, a, nlca, nuca, b, 0);
imsl_f_write_matrix ("Solution x, of Ax = b", 1, n, x, 0);
}
出力
Solution x, of Ax = b
1 2 3 4
2 1 -3 4
例題 2
こ 例題 最初 例題 ータ 使用し Ax = b 解 ます 今回 分解 行 い LU
再計算せ にATx = b 解 ます
#include <imsl.h>
int main() {
int n = 4;
int nlca = 1;
int nuca = 1;
int *pivot;
f_complex *x;
f_complex *factor;
/* Note that a is in band storage mode */
f_complex a[] =
{{0.0, 0.0}, {4.0, 0.0}, {-2.0, 2.0}, {-4.0, -1.0}, {-2.0, -3.0}, {-0.5, 3.0}, {3.0, -3.0}, {1.0, -1.0}, {6.0, 1.0}, {1.0, 1.0}, {0.0, 2.0}, {0.0, 0.0}};
f_complex b[] =
{{-10.0, -5.0}, {9.5, 5.5}, {12.0, -12.0}, {0.0, 8.0}};
/* Solve Ax = b and return LU */
x = imsl_c_lin_sol_gen_band (n, a, nlca, nuca, b, IMSL_FACTOR, &pivot, &factor,
0);
imsl_c_write_matrix ("solution of Ax = b", n, 1, x, 0);
imsl_free (x);
/* Use precomputed LU to solve ctrans(A)x = b */
x = imsl_c_lin_sol_gen_band (n, a, nlca, nuca, b, IMSL_FACTOR_USER, pivot, factor,
IMSL_TRANSPOSE, 0);
imsl_c_write_matrix ("solution of ctrans(A)x = b", n, 1, x, 0);
}
出力
Solution of Ax = b
1 2 3 4
2 1 -3 4 Solution of trans(A)x = b
1 2 3 4
-6 -5 -1 -0
警告 ー
IMSL_ILL_CONDITIONED 入力行列 非常に悪条件 す L1 条件数 逆数 推定値
"rcond" = # す 解 恐 く正確 あ ませ
重大 ー
IMSL_SINGULAR_MATRIX こ 入力行列 特異 す
lin_sol_gen_band ( 複素数 )
複素一般線形帯方程式 Ax = b 解 ます プ ン 引数 使用し 関連す 計算 行うこ が
ます こ に 部分 ッ 選択 使用す A LU 分解 計算 AH x = b 解 又 A LU 分
解が与え た Ax = b 解 計算 があ ます
梗概
#include <imsl.h>
f_complex *imsl_c_lin_sol_gen_band (int n, f_complex a[], int nlca, int nuca, f_complex b[], …, 0)
d_complex型 関数 imsl_z_lin_sol_gen_band
必須 引数
int n ( 入力)
行列 行数 列数
f_complex a [ ] ( 入力)
帯格納 ー n × n 帯係数行列 含 大 さ (nlca + nuca + 1) × n 配列
int nlca ( 入力) a 下共対角項 数 int nuca ( 入力)
a 共対角項 数
f_complex b [ ] ( 入力)
右辺 含 大 さ n 配列
戻 値
線形方程式 Ax = b 解 x へ インタ こ 領域 解放に imsl_free 使用します 解が得 い場合 NULL が返さ ます
プ ン引数 梗概
#include <imsl.h>
f_complex *imsl_c_lin_sol_gen_band (int n, f_complex a[], int nlca, int nuca, f_complex b[],
IMSL_TRANSPOSE,
IMSL_RETURN_USER, f_complex x[],
IMSL_FACTOR, int **p_pvt, f_complex **p_factor, IMSL_FACTOR_USER, int pvt[], f_complex factor[], IMSL_CONDITION, float *condition,
IMSL_FACTOR_ONLY, IMSL_SOLVE_ONLY, 0)
プ ン引数
IMSL_TRANSPOSE
AH x = b 解く
フ ： Ax = b 解く
IMSL_RETURN_USER, f_complex x[] ( 出力)
解 x 含 長さ n ー 割 当 配列
IMSL_FACTOR, int **p_pvt, f_complex **p_factor (入力/出力) int **p_pvt (入力/出力)
分解 ッ 列 含 長さ n 配列へ インタ ア 戻 必要 領域 が imsl_c_lin_sol_gen_band によ 割 当 ます 一般的に int *p_pvt が 宣言さ &p_pvt が引数 し 使用さ ます
f_complex **p_factor (入力/出力)
列 ッ 交換によ A LU 分解 含 大 さ (2nlca + nuca+ 1) × n 配列へ
インタ ア 戻 必要 領域 imsl_c_lin_sol_gen_band によ
割 当 ます 一般的にf_complex *p_factor が宣言さ &p_factor が引数
し 使用さ ます
IMSL_FACTOR_USER, int pvt[], f_complex factor[] (入力/出力) int pvt[] (入力/出力)
分解 ッ 列 含 大 さ n ー 割 当 配列 f_complex factor[] (入力/出力)
A LU 分解 含 大 さ (2nlca + nuca + 1) × n ー 割 当 配列 A が 必要 け factor a 最初 (2nlca + nuca + 1) ×n 領域 共有し ます
IMSL_SOLVE_ONLY が指定さ こ パ ータ 入力 す そ 他 場合
出力に ます
IMSL_CONDITION, float *condition (出力)
行列 A L1 条件数 推定値 含 へ インタ こ プ ン
プ ンIMSL_SOLVE_ONLY 共に使用す こ ませ
IMSL_FACTOR_ONLY
部分 ッ 交換法 A LU 分解 計算します IMSL_FACTOR_ONLY が使用さ IMSL_FACTOR 又 IMSL_FACTOR_USER い が必要に ます 引数 b そ 時 無視さ imsl_c_lin_sol_gen_band 戻 値 NULL に ます IMSL_SOLVE_ONLY
imsl_c_lin_sol_gen_band によ 事前に計算さ た LU 分解が与え Ax = b
解 ます フ Ax = b 解 imsl_c_lin_sol_gen_band によ
指定さ います IMSL_SOLVE_ONLY が使用さ 引数 IMSL_FACTOR_USER が 必要に 引数 a 無視さ ます
説明
関数 imsl_c_lin_sol_gen_band 複素帯係数行列 A 持 線形方程式 解 ます 最初に ー ングさ た部分 ッ 法 利用し A LU 分解 計算します ー ングさ た部分 ッ 選択 各列が同 L∝ 持 し ー ングす いう点 部分 ッ 選択 異 い ます U が ロ 対角要素 持 場合 因子分解 ーに ます こ A が特異もしく 特異 行列に近い場合に起こ ます
よ 簡単 2 方程式 y = L-1 b x = U-1 y 解くこ によ 線形方程式 解 得 こ が ます 線形方程式 又 逆行列 解 求 A L1条件数 推定値 Higham(1988年) によ Hager
方法に対す Higham 修正 使用し 計算さ ます 推定した条件数が 1/ε ここ ε ン精度 よ 大 い場合 警告 ッ ージが出さ ます こ A 非常に小さい変化が解 x 大 い変化 生 こ 示し います
分解 角部分 U が ロ 対角要素 持 場合 関数imsl_c_lin_sol_gen_band ーに ます imsl_c_lin_sol_gen_band LINPACK ーチンCGBFA 元にし います Dongarra
1979年 参照し く さい CGBFA ー ングし い部分 ッ 選択 使用し います
例題 1
次 線形方程式 解 ます
i i
i i
x x x x
i i
i i
i i
i i
8 12 12
5 . 5 5 . 9
5 10
1 2 0
0
1 4 3 3 1
0
0 2 2 3 5 . 0 6
0 0 4
3 2
3 2 1 0
#include <imsl.h>
int main() {
int n = 4;
int nlca = 1;
int nuca = 1;
f_complex *x;
/* Note that a is in band storage mode */
f_complex a[] =
{{0.0, 0.0}, {4.0, 0.0}, {-2.0, 2.0}, {-4.0, -1.0}, {-2.0, -3.0}, {-0.5, 3.0}, {3.0, -3.0}, {1.0, -1.0}, {6.0, 1.0}, {1.0, 1.0}, {0.0, 2.0}, {0.0, 0.0}};
f_complex b[] =
{{-10.0, -5.0}, {9.5, 5.5}, {12.0, -12.0}, {0.0, 8.0}};
x = imsl_c_lin_sol_gen_band (n, a, nlca, nuca, b, 0);
imsl_c_write_matrix ("Solution, x, of Ax = b", n, 1, x, 0);
}
出力
Solution, x, of Ax = b 1 ( 3, -0) 2 ( -1, 1) 3 ( 3, 0) 4 ( -1, 1)
例題 2
最初 例題 ータ 使 Ax = b 解 ます こ 例題 因子分解 返し LU 再計算
せ に AH x = b 解 ます
#include <imsl.h>
int main() {
int n = 4;
int nlca = 1;
int nuca = 1;
int *pivot;
f_complex *x;
f_complex *factor;
/* Note that a is in band storage mode */
f_complex a[] =
{{0.0, 0.0}, {4.0, 0.0}, {-2.0, 2.0}, {-4.0, -1.0}, {-2.0, -3.0}, {-0.5, 3.0}, {3.0, -3.0}, {1.0, -1.0}, {6.0, 1.0}, {1.0, 1.0}, {0.0, 2.0}, {0.0, 0.0}};
f_complex b[] =
{{-10.0, -5.0}, {9.5, 5.5}, {12.0, -12.0}, {0.0, 8.0}};
/* Solve Ax = b and return LU */
x = imsl_c_lin_sol_gen_band (n, a, nlca, nuca, b, IMSL_FACTOR, &pivot, &factor,
0);
imsl_c_write_matrix ("solution of Ax = b", n, 1, x, 0);
imsl_free (x);
/* Use precomputed LU to solve ctrans(A)x = b */
x = imsl_c_lin_sol_gen_band (n, a, nlca, nuca, b, IMSL_FACTOR_USER, pivot, factor,
IMSL_TRANSPOSE, 0);
imsl_c_write_matrix ("solution of ctrans(A)x = b", n, 1, x, 0);
}
出力
solution of Ax = b 1 ( 3, -0) 2 ( -1, 1) 3 ( 3, 0) 4 ( -1, 1)
solution of ctrans(A)x = b 1 ( 5.58, -2.91) 2 ( -0.48, -4.67) 3 ( -6.19, 7.15) 4 ( 12.60, 30.20)
警告 ー
IMSL_ILL_CONDITIONED こ 入力行列 非常に悪条件す L1 条件数 逆数 推定値
"rcond"=# す 解 恐 く正確 あ ませ
重大 ー
IMSL_SINGULAR_MATRIX 入力行列 特異 す