Keyword: スパース行列, LU分解
概要
本サンプルは実スパース行列のLU分解を行うC言語によるサンプルプログラムです。 本サンプルは以下に示される実スパース行列をLU分解し、結果を出力します。
※本サンプルはnAG Cライブラリに含まれる関数 nag_superlu_lu_factorize() のExampleコードです。本サンプル及び関数の詳細情報は nag_superlu_lu_factorize のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで
入力データ
（本関数の詳細はnag_superlu_lu_factorize のマニュアルページを参照）1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
| このデータをダウンロード |
nag_superlu_lu_factorize (f11mec) Example Program Data 5 n 1 3 5 7 9 12 icolzp(i) i=0..n 2. 1 4. 3 1. 1 -2. 5 1. 2 1. 3 -1. 2 1. 4 1. 3 2. 4 3. 5 a(i), irowix(i) i=0..nnz-1
- １行目はタイトル行で読み飛ばされます。
- ２行目に行列Aの次数(n)を指定しています。
- ３～８行目に行列Aの各列の最初の非ゼロ要素のインデックス（行列内の位置）(icolzp)を指定しています。
- ９～１９行目に行列Aの非ゼロ要素とその行インデックス(irowix)を指定しています。
出力結果
（本関数の詳細はnag_superlu_lu_factorize のマニュアルページを参照）1 2 3 4 5 6 7 8 9 10 11
| この出力例をダウンロード |
nag_superlu_lu_factorize (f11mec) Example Program Results
Number of nonzeros in factors (excluding unit diagonal)
14
Factor elements in LVAL
-2.00 -0.50 4.00 0.50 1.00 -0.50 1.00 -0.50 2.00 2.00
Factor elements in LVAL
1.00 -1.00 1.00 3.00
- ５行目にLU分解後の非ゼロ要素（単位対角要素を除く）の数が出力されています。
- ９行目にLU分解された行列Lの非ゼロ値と行列Uの非ゼロ値の一部が出力されています。
- １３行目にLU分解された行列Uの非ゼロ値の一部が出力されています。
ソースコード
（本関数の詳細はnag_superlu_lu_factorize のマニュアルページを参照）
※本サンプルソースコードはnAG数値計算ライブラリ（Windows, Linux, MAC等に対応）の関数を呼び出します。
サンプルのコンパイル及び実行方法
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138
| このソースコードをダウンロード |
/* nag_superlu_lu_factorize (f11mec) Example Program.
*
* CLL6I261D/CLL6I261DL Version.
*
* Copyright 2017 Numerical Algorithms Group.
*
* Mark 26.1, 2017.
*/
#include <stdio.h>
#include <nag.h>
#include <nagx04.h>
#include <nag_stdlib.h>
#include <nagf11.h>
/* Table of constant values */
int main(void)
{
double flop, thresh;
Integer exit_status = 0, i, n, nnz, nnzl, nnzu, nzlmx, nzlumx, nzumx;
double *a = 0, *lval = 0, *uval = 0;
Integer *icolzp = 0, *il = 0, *iprm = 0, *irowix = 0;
Integer *iu = 0;
/* Nag types */
NagError fail;
Nag_ColumnPermutationType ispec;
Nag_OrderType order = Nag_ColMajor;
Nag_MatrixType matrix = Nag_GeneralMatrix;
Nag_DiagType diag = Nag_NonUnitDiag;
INIT_FAIL(fail);
/* nag_superlu_lu_factorize (f11mec).
* LU factorization of real sparse matrix
*/
printf("nag_superlu_lu_factorize (f11mec) Example Program Results\n\n");
/* Skip heading in data file */
scanf("%*[^\n] ");
/* Read order of matrix */
scanf("%ld%*[^\n] ", &n);
/* Read the matrix A */
/* Allocate memory */
if (!(icolzp = nAG_ALLOC(n + 1, Integer)))
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
for (i = 1; i <= n + 1; ++i)
scanf("%ld%*[^\n] ", &icolzp[i - 1]);
nnz = icolzp[n] - 1;
/* Allocate memory */
if (!(irowix = nAG_ALLOC(nnz, Integer)) ||
!(a = nAG_ALLOC(nnz, double)) ||
!(il = nAG_ALLOC(7 * n + 8 * nnz + 4, Integer)) ||
!(iu = nAG_ALLOC(2 * n + 8 * nnz + 1, Integer)) ||
!(uval = nAG_ALLOC(8 * nnz, double)) ||
!(lval = nAG_ALLOC(8 * nnz, double)) ||
!(iprm = nAG_ALLOC(7 * n, Integer)))
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
for (i = 1; i <= nnz; ++i)
scanf("%lf%ld%*[^\n] ", &a[i - 1], &irowix[i - 1]);
/* Calculate COLAMD permutation */
ispec = Nag_Sparse_Colamd;
/* nag_superlu_column_permutation (f11mdc).
* Real sparse nonsymmetric linear systems, setup for
* nag_superlu_lu_factorize (f11mec)
*/
nag_superlu_column_permutation(ispec, n, icolzp, irowix, iprm, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_superlu_column_permutation (f11mdc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
/* Factorise */
thresh = 1.;
nzlmx = 8 * nnz;
nzlumx = 8 * nnz;
nzumx = 8 * nnz;
/* nag_superlu_lu_factorize (f11mec), see above. */
nag_superlu_lu_factorize(n, irowix, a, iprm, thresh, nzlmx, &nzlumx, nzumx,
il, lval, iu, uval, &nnzl, &nnzu, &flop, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_superlu_lu_factorize (f11mec).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
/* Output results */
printf("\n");
printf("Number of nonzeros in factors (excluding unit diagonal)\n");
printf("%8ld\n\n", nnzl + nnzu - n);
/* nag_gen_real_mat_print_comp (x04cbc).
* Print real general matrix (comprehensive)
*/
fflush(stdout);
nag_gen_real_mat_print_comp(order, matrix, diag, 1, 10, lval, 1, "%6.2f",
"Factor elements in LVAL", Nag_NoLabels, NULL,
Nag_NoLabels, NULL, 80, 1, 0, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_gen_real_mat_print_comp (x04cbc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
printf("\n");
/* nag_gen_real_mat_print_comp (x04cbc), see above. */
fflush(stdout);
nag_gen_real_mat_print_comp(order, matrix, diag, 1, 4, uval, 1, "%6.2f",
"Factor elements in LVAL", Nag_NoLabels, NULL,
Nag_NoLabels, NULL, 80, 1, 0, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_gen_real_mat_print_comp (x04cbc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
END:
nAG_FREE(a);
nAG_FREE(lval);
nAG_FREE(uval);
nAG_FREE(icolzp);
nAG_FREE(il);
nAG_FREE(iprm);
nAG_FREE(irowix);
nAG_FREE(iu);
return exit_status;
}