タ
Imsl_symbolic_factor 構造体 詳細 次 通
パ ータ ータ型 説明
nzsub int ** Cholesky因子 非 ロ 非対角要素に対す
縮行指標 含 配列へ インタ
xnzsub int ** *nzsub に対す イン ック 含 長さ
n + 1 配列へ インタ Cholesky因子 列 j にあ 非 ロ 要素に対す 行指標 (*nzsub)[(*xnzsub)[j]]) 始ま 連 続領域に格納さ
maxsub int 指標 し 使わ 配列 *nzsub 要素
数 *nzsub 大 さ maxsub よ 大 いこ もあ こ に注意
xlnz int ** 長さ n + 1 配列へ インタ こ 配
列 *alnz 非 ロ 非対角要素 抽出
す に使わ 始 終わ イン ック
含 alnz に い プ ン引
数 IMSL_RETURN_NUMERIC_FACTOR 説 明 参照 因子行列 列 j に対し *alnz
始 終わ イン ック そ
(*xlnz)[j] (*xlnz)[j + 1] に格納 さ
maxlnz int Cholesky 因子 非 ロ非対角要素 数
perm int ** 置換ベク 含 長さ n 配列へ イ
ンタ
invp int ** 逆置換ベク 含 長さ n 配列へ
インタ
multifrontal_space int フロンタ 行列 タック た に必要 作
業用領域 大 さ チフロンタ 因子分 解が用い い場合 こ 値 ロに ッ さ
こ 構造体内 アロ ー さ た 解放す に 関数 imsl_free_symbolic_factor 使うこ
IMSL_SUPPLY_SYMBOLIC_FACTOR, Imsl_symbolic_factor *sym_factor (入力)
Imsl_symbolic_factor 型 構造体へ インタ こ 構造体 入力行列 ン ック分
解 含 IMSL_RETURN_SYMBOLIC_FACTOR プ ン 用いた imsl_c_lin_sol_posdef_coordinate によ 計算さ ます IMSL_SYMBOLIC_FACTOR_ONLY
入力行列 ン ック分解 計算します 引数 b 無視さ ます IMSL_RETURN_NUMERIC_FACTOR, Imsl_c_numeric_factor *num_factor (出力)
戻 入力行列 数値分解 含 Imsl_c_numeric_factor 型 構造体へ インタ
Imsl_c_numeric_factor構造体 詳細 次 通 パ ータ ータ型 説明
nzsub int ** Cholesky因子 非 ロ 非対角要素に対す 行指標 含 配列へ インタ こ 配列 長さ nz 割 当
が 配列 全 要素 使うこ い xnzsub int ** nzsub に対す 指標 含 長さ n + 1 配列へ
インタ Cholesky因子 列 j にあ 非 ロ 要素に対 す 行指標 nzsub[xnzsub[j]] 始ま 連続領域 に格納さ
xlnz int ** 長さ n + 1 配列へ インタ こ 配列 alnz
非 ロ 非対角要素 抽出す に使わ 始 終わ 指標 含 因子行列 列 j に対し
alnz 始 終わ 指標 そ xlnz[j]
xlnz[j + 1] に格納さ
alnz f_complex ** 長さ nz - n 配列 Cholesky因子 非 ロ 非対 角要素 含
perm int ** 置換ベク 含 長さ n 配列へ インタ
diag f_complex ** Cholesky因子 対角要素 含 長さ n 配列へ イ ンタ
L a 数値的 因子 します に述べた構造体 L 対角要素 diag に格納さ ます L 非 ロ 非対角要素 alnzに格納さ ます 列jに対し alnz L 非 ロ要素 抽出
す に使わ 始 終わ 指標 そ xlnz[j] xlnz[j + 1] に格納さ ます L
非 ロ 要素 列指標 nzsub にあ ます xnzsub[i] 列iに対し L 行指標 抽出す
始 nzsub 指標 含 ます 次 コー 構造体 成分 因子行列L 下 角
再構築す や 方 示します：
Imsl_c_numeric_factor numfctr;
. .
for (i = 0; i < n; i++){
L[i][i] = (*numfctr.diag)[i];
if((*numfctr.xlnz)[i] > (nz-n)+1) break;
start = (*numfctr.xlnz)[i]-1;
stop = (*numfctr.xlnz)[i+1]-1;
k = (*numfctr.xnzsub)[i]-1;
for (j = start; j < stop; j++){
L[(*numfctr.nzsub)[k]-1][i] = (*numfctr.alnz)[j];
k++;
} }
こ 構造体内 アロ ー さ た 解放す に 関数 imsl_c_free_numeric_factor 使うこ
IMSL_SUPPLY_NUMERIC_FACTOR, Imsl_c_numeric_factor *num_factor 。入力)
Imsl_c_numeric_factor型 構造体 指す インタ こ 構造体
IMSL_RETURN_NUMERIC_FACTOR プ ン
imsl_c_lin_sol_posdef_coordinate によ 計算さ た入力行列 数値的分解 含 います こ 構造体 プ ン引数 IMSL_RETURN_NUMERIC_FACTOR に説明さ います
こ 構造体内 アロ ー さ た 解放す に 関数 imsl_c_free_numeric_factor 使うこ
IMSL_NUMERIC_FACTOR_ONLY,
入力行列 数値的分解 計算し 返します 引数 b 無視さ ます IMSL_SOLVE_ONLY,
A 数値的 又 ン ック分解が与え Ax = b 解 ます こ プ ン IMSL_SUPPLY_NUMERIC_FACTOR 又 IMSL_SUPPLY_SYMBOLIC_FACTOR 使用 必 要 します
IMSL_MULTIFRONTAL_FACTORIZATION,
チフロンタ 法 使用し 数値的分解 実行します フ 疎 縮格納方式に
基 い 標準 分解 計算します
IMSL_RETURN_USER, f_complex x[] 。出力)
解 x 含 長さ n ー 割 当 配列
IMSL_SMALLEST_DIAGONAL_ELEMENT, float *small_element 。出力)
数値的分解 発生す 最小 対角要素 含 へ インタ こ プ ン 数値的分解が imsl_c_lin_sol_posdef_coordinate 呼び出し 計算さ 場合
有効 す
IMSL_LARGEST_DIAGONAL_ELEMENT, float *large_element 。出力)
数値的分解 発生す 最大 対角要素 含 へ インタ こ プ ン 数値的分解が imsl_c_lin_sol_posdef_coordinate 呼び出し 計算さ 場合
有効 す
IMSL_NUM_NONZEROS_IN_FACTOR, int *num_nonzeros 。出力) 分解因子 非 ロ 合計数 含 へ インタ
IMSL_CSC_FORMAT, int *col_ptr, int *row_ind, float *values 。入力)
縮した疎 列形式。CSC) 係数行列 受け入 ます こ 格納方式 説明 イン ロ ダク ン 参照し く さい
説明
関数 imsl_c_lin_sol_posdef_coordinate 疎 ー 正定値係数行列A 持 線形方程式 解 ます こ 関数 フ 使用法 最初に係数行列 置換 ン ック分解が計算さ そ 数値的分解が計算さ ます 線形方程式 解 こ 数値的分解因子 使用し 求 こ が
ます
ン ック分解 ップ 最小次数 並び え 決定 コ キー因子L 疎 ータ構造体 設 定 ます こ ップ 疎 係数行列 パターン す わ 非 ロ要素 位置 け 必要
し 要素そ も 扱いませ 従 Imsl_c_sparse_elem 構造体 val欄 無視さ ます あ
アプ ー ンが 同 疎 パターン 持 要素 値が異 複数 疎 対称正定値係数行列 扱 う場合 IMSL_RETURN_SYMBOLIC_FACTOR IMSL_SUPPLY_SYMBOLIC_FACTOR 使用す こ によ ン ック分解 一度 け計算す よい す
ン ック分解因子によ コ キー因子L 疎 ータ構造体が与え 数値的分解 次式 が成 立 ように Lに ン 値 生成します：
PAP T = LL T ここ P 最小次数 並び えによ 決定さ 置換行列 す
数値的分解 2 方法 い 実行す こ が ます フ 疎 縮格納方式に基 い 標準 分解 行 います こ 方法 George Liu (1981年) に記述さ います プ ン
チフロンタ 法 使用す こ が ます チフロンタ 法 よ 多く 必要 します が 場合によ に よ 高速に ます チフロンタ 法 Liu (1987年) ーチンに基 い い
ます こ 方法 詳細 説明 Liu 1990年 よび Duff Reid 1983年 1984年 Ashcraft 1987
年 Ashcraftそ 他 1987年 Liu 1986年 1989年 参照し く さい
アプ ー ンが 同 係数行列 右辺が異 いく 線形方程式 解く必要があ 場合 プ ン IMSL_RETURN_NUMERIC_FACTOR IMSL_SUPPLY_NUMERIC_FACTOR 使 コ キー
分解因子 前処理 す こ が ます そ IMSL_SOLVE_ONLY プ ン 使い 全 方程
式 効率良く解くこ が ます
数値分解が与え そ 解 x 次によ 得 ます L y1
= Pb
LT y2= y
1x = PT y2
置換情報 P 数値因子構造体に保持さ ます
例題１
簡単 例 し 次 ー 正定値行列 考えます：
2 1 0
1 4 1 2
0 1 2 10
i
A i i
i
Ax = (
2 + 2i, 5 +15i, 36 + 28i)
T ように xT= (1 + i, 2 + 2i, 3 + 3i)
します A 下 角非 ロ 数 nz = 5 す
#include <imsl.h>
int main() {
Imsl_c_sparse_elem a[] = {0, 0, {2.0, 0.0}, 1, 1, {4.0, 0.0}, 2, 2, {10.0, 0.0}, 1, 0, {-1.0, -1.0}, 2, 1, {1.0, -2.0}};
f_complex b[] = {{-2.0, 2.0}, {5.0, 15.0}, {36.0, 28.0}};
int n = 3;
int nz = 5;
f_complex *x;
x = imsl_c_lin_sol_posdef_coordinate (n, nz, a, b, 0);
imsl_c_write_matrix ("Solution, x, of Ax = b", n, 1, x, 0);
imsl_free (x);
}
出力
Solution, x, of Ax = b 1 ( 1, 1) 2 ( 2, 2) 3 ( 3, 3)
例題
A = E(2500, 50) します 方程式 Ax = b1 解い 数値分解 返します そ こ 数値分解
使用し 方程式 Ax = b2 解 ます 実行時間 比がプ ン さ ます
#include <imsl.h>
#include <stdio.h>
int main() {
Imsl_c_sparse_elem *a;
Imsl_c_numeric_factor numeric_factor;
f_complex b_1[2500], b_2[2500], *x_1, *x_2;
int n, ic, nz, i, index;
double time_1, time_2;
float *rand_vec;
ic = 50;
n = ic*ic;
index = 0;
/* Generate two right hand sides */
rand_vec = imsl_f_random_uniform (4*n*sizeof(*rand_vec), 0);
for (i=0; i<n; i++) {
b_1[i].re = rand_vec[index++];
b_1[i].im = rand_vec[index++];
b_2[i].re = rand_vec[index++];
b_2[i].im = rand_vec[index++];
}
/* Build coefficient matrix a */
a = imsl_c_generate_test_coordinate (n, ic, &nz, IMSL_SYMMETRIC_STORAGE,
0);
/* Now solve Ax_1 = b_1 and return the numeric factorization */
time_1 = imsl_ctime ();
x_1 = imsl_c_lin_sol_posdef_coordinate (n, nz, a, b_1, IMSL_RETURN_NUMERIC_FACTOR, &numeric_factor,
0);
time_1 = imsl_ctime () - time_1;
/* Now solve Ax_2 = b_2 given the numeric factorization */
time_2 = imsl_ctime ();
x_2 = imsl_c_lin_sol_posdef_coordinate (n, nz, a, b_2, IMSL_SUPPLY_NUMERIC_FACTOR, &numeric_factor,
IMSL_SOLVE_ONLY, 0);
time_2 = imsl_ctime () - time_2;
printf("time_2/time_1 = %lf\n", time_2/time_1);
}
出力
time_2/time_1 = 0.096386