superlu_smp

プン引数

IMSL_EQUILIBRATE, int equilibrate (入力)

入力行列 A 因子分解前に平衡させに指定す equilibrate 説明

0 平衡化しい

1 平衡化す

フ : equilibrate = 0

IMSL_COLUMN_ORDERING_METHOD, Imsl_col_ordering method (入力)

分解前に疎保ようにす列並びえ方法次一 method に設定し並びえ法選択す

method 説明

IMSL_NATURAL 自然順序すわ入力行列列順

IMSL_MMD_ATA A^TA ^{構造体最小次} ^順 IMSL_MMD_AT_PLUS_A A^T + A ^構造体 ^最小次 ^順

IMSL_COLAMD 列近似最小次順

IMSL_PERMC プン引数 IMSL_COLPERM_VECTOR ー

が入力した置換ベク permc 与え順ベク permc 数 0,1,…,n-1 置換す

フ : method = IMSL_COLAMD

IMSL_COLPERM_VECTOR, int permc[] (入力)

後並びえ(postordering) 前置換行列 P 定義す長さ n 配列こ引数

IMSL_COLUMN_ORDERING_METHOD が method = IMSL_PERMC 用い必須あそ以外無視さ

IMSL_TRANSPOSE, int transpose (入力)

転置問題 A^Tx = b ^{解く場合に指定す} ^こ ^プ ^ン ^因子分解 ^指定す ^プ

ン一緒に使うこが

transpose 説明

0 Ax = b 解く

1 A^Tx = b ^解く.

フ : transpose = 0

IMSL_ITERATIVE_REFINEMENT, int refine (入力) 反復改良すう指定す

refine 説明

0 反復改良しい

1 反復改良す

フ : refine = 1

IMSL_FACTOR_SOLVE, int factsol (入力)

LU 分解け線形方程式解くけ両方行う指定す factsol 説明

0 A LU 分解方程式 Ax = b ^解も求

1 LU 分解け行う LU 分解プン引数

IMSL_RETURN_SPARSE_LU_FACTOR 引渡さ入力引数 b 無視さ

factsol 説明

2 与えたLU分解 Ax = b ^解く LU分解プン引

数IMSL_SUPPLY_SPARSE_LU_FACTOR 渡さ

い入力引数 A 反復改良条件数計算あいッ増大因子逆数計算以外無視さ

フ : factsol = 0

IMSL_DIAG_PIVOT_THRESH, double diag_pivot_thresh (入力) 対角項がッし扱わう閾値指定す 0.0



diag_pivot_thresh



1.0.

フ : diag_pivot_thresh = 1.0

IMSL_SNODE_PREDICTION, int snode_prediction (入力)

L ーパーーにあ非ロ数予測方法示す

snode_prediction 説明

0 静的方法使う

1 動的方法使う

フ : snode_prediction = 0

IMSL_PERFORMANCE_TUNING, int sp_ienv[] (入力)

長さ 8 配列ーが因子分解ア性能チーニングす際正パータす

i sp_ienv[i] 説明

0 パネイ

フ : sp_ienv[0] = 10

1 ーパーー融合制御す緩和パータ

フ : sp_ienv[1] = 5

2 ーパーー最大許容イ

フ : sp_ienv[2] = 100

3 2D ロッキングに使わ最小行次元

フ : sp_ienv[3] = 200

4 2D ロッキングに使わ最小列次元

フ : sp_ienv[4] = 40

5 L ーパーー値格納す配列 nzval 大さ負場合充

填増加因子すわそ大さ絶対値元行列A 非ロ数積が格納領域アローすに用います正場合保存領域がアローさ非ロ数表す

配列 sp_ienv こ要素 L ーパーーにあ非ロ数予測

が動的方法場合すわ snode_prediction = 1 けに用い

フ : sp_ienv[5] = -20

6 列 U に格納す配列 rowind nzval 大さ負場合充填増加因子すわそ大さ絶対値元行列A 非ロ数積が格納領域アローすに用います正場合格納領域がアローさ非ロ数表す

フ : sp_ienv[6] = -20

7 L ーパーー添え字格納す配列 rowind 大さ負場合

充填増加因子すわそ大さ絶対値元行列A 非ロ数積が格納領域アローすに用います正場合格納領域がアローさ非ロ数表す

フ : sp_ienv[7] = -10

IMSL_CSC_FORMAT, int HB_col_ptr[], int HB_row_ind[], float HB_values[] (入力)

縮さた疎列ー (CSC) 係数行列受け入こー説明インロダ

クン参照

IMSL_SUPPLY_SPARSE_LU_FACTOR, Imsl_f_super_lu_factor lu_factor_supplied (入力) プン IMSL_RETURN_SPARSE_LU_FACTOR 計算さた入力行列 LU 分解含

Imsl_f_super_lu_factor 型構造体こ構造体にい説明参照こ構造体

内に割当た領域解放すに関数 imsl_f_superlu_factor_free 使う IMSL_RETURN_SPARSE_LU_FACTOR, Imsl_f_super_lu_smp_factor *lu_factor_returned (出力)

入力行列 LU 分解含 Imsl_f_super_lu_smp_factor型構造体アこ構造体

にい説明参照こ構造体内に割当た領域解放すに関数 imsl_f_superlu_smp_factor_free 使う

IMSL_CONDITION, float *condition (出力) 平衡化後行列条件数逆数推定

IMSL_PIVOT_GROWTH_FACTOR, float *recip_pivot_growth (出力) ッ増大因子逆数

 

min (j P D AD Pr r c c j) _/Uj_

recip_pivot_growth が 1 よ非常に小さい場合 LU 分解安定性良くい IMSL_FORWARD_ERROR_BOUND, float *ferr (出力)

解ベク x に対す推定さ前方誤差範囲こプン引数 IMSL_ITERATIVE_REFINEMENT 1 に設定すこが必要す IMSL_BACKWARD_ERROR, float *berr (出力)

解ベク x 成分毎相対後方誤差こプン引数 IMSL_ITERATIVE_REFINEMENT 1 に設定すこが必要す IMSL_RETURN_USER, float x[] (出力)

線形方程式解 x 含ー割当長さ n 配列

説明

imsl_f_superlu_smp 線型方程式解法逐次版 imsl_f_superlu 同す

関数 imsl_f_superlu_smp 行列 A LU分解ーパーーダ格納方式用います逐次版比べ L LU 因子連続す列ーパーー連続いこがあます従

各列あいーパーー始インタく各列あいーパーー終わインタも必要ます因子分解構造体 Imsl_f_super_lu_smp_factor そ構造体 Imsl_f_hbp_format Imsl_f_scp_format に含まます

こ構造体以下に記述します

表 1.1 構造体 Imsl_f_hbp_format

パータータ型説明

nnz int 行列非ロ数

nzval float * 列パックさた非ロ値配列

rowind int * 非ロ行番号配列

colbeg int * 大さ ncol+1 配列 colbeg[j] nzval[]

rowind[] に列 j が始ま位置格納す要素

colbeg[ncol] 配列 nzval[] rowind[] 最初フー位置指す

colend int * 大さ ncol 配列 colend[j] nzval[]

rowind[] に列 j 最後要素一先位置格納す

表 1.2 構造体 Imsl_f_scp_format

パータータ型説明

nnz int ーパーーダ行列非ロ数

nsuper int ーパーー数 1 引いたも

nzval float * 列パックさた非ロ値配列

nzval_colbeg int * 大さ ncol+1 配列 nzval_colbeg[j] -nzval[]

列 j 始指すン nzval_colbeg[ncol]

nzval[] 最初フー位置指す

nzval_colend int * 大さ ncol 配列 nzval_colend[j] nzval[]

列 j 最後要素一先位置指す

rowind int * 矩形ーパーー縮さた行インック配列

rowind_colbeg int * 大さ ncol+1 配列 rowind_colbeg[j] rowind[]

列 j 始指す要素 rowind_colbeg[ncol]

rowind[] 最初フー位置指す

rowind_colend int * 大さ ncol 配列 rowind_colend[j] rowind[]

列 j 最後要素一先位置指す

col_to_sup int * 大さ ncol+1 配列 col_to_sup[j] 列 j が所属すーパーー番号 col_to_sup[] 最初 ncol 個ンけが定義さ

sup_to_colbeg int * 大さ ncol+1 配列 sup_to_colbeg[s] s 番目ーパーー最初列指すこ配列最初 nsuper+1 位置けが用い

sup_to_colend int * 大さ ncol 配列 sup_to_colend[s] s 番目ーパーー最後列一先指すこ配列最初 nsuper+1 位置けが用い

表 1.3 構造体 Imsl_f_super_lu_smp_factor

パータータ型説明

nrow int 行列 A 行数

ncol int 行列 A 列数

equilibration_

method int A 平衡化す方法：

0 – 平衡化しい

1 – 行平衡化

2 – 列平衡化

3 – 行列平衡化

rowscale float * 大さ nrow 配列 A に対す行ー因子含

columnscale float * 大さ ncol 配列 A に対す列ー因子含

rowperm int * 大さ nrow 行置換配列行置換行列 Pr 表す

colperm int * 大さ ncol 列置換配列列置換行列 Pc 表す

U Imsl_f_hbp

_format * ーパーーダロック外 A U 因子部分

Harwell-Boeing フーッ格納さ

L Imsl_f_scp

_format * A L 因子 ロック下角行列しーパーーダ

フーッ格納さ

構造体 Imsl_d_super_lu_smp_factor そ構造体も同様にそ定義 float double に

Imsl_f_hbp_format Imsl_d_hbp_format に Imsl_f_scp_format Imsl_d_scp_format に置換えたもす

逐次版比べ LU 分解こが並列処理さます逐次版よう動的拡大並列コー

に実装す困難ーが L 配列 rowind U 配列 rowind nzval ( 表構造

体 Imsl_f_scp_format Imsl_f_hbp_format 参照) 大さ推定値性能チーニング配列 sp_ienv

要素 6 7 提供せませ

L ーパーー各々が連続したに格納さようにすたに L ーパーー大さにい静的あい動的予測すこがます静的バージン関数 imsl_f_superlu_smp

がフ使う方法すが PA = LU 行置換 P に対し L 非ロ構造 Householder パー

QR 分解 A = QR Householder行列 H 非ロ構造に含まいこ使いますさに L 基

本ーパーー各々常に H 基本ーパーーに含まいこ示すこがます従

L ーパーー L 因子配列 nzval H ーパーー大さに基い分解前に推

定アローすこがますーパーー分割 H ーパーー大さ計算量行列 A 非ロ数にほ線型す

実際に静的予測方法多く問題十分あませしし H 非ロ数が L 非ロ数よ大い場合プン引数 IMSL_SNODE_PREDICTION 1 にッし動的予測試こがますこ方法も H ーパーー分割使いますが L ーパーーダグフ動的に探査し要求さーよ厳密限得こがます動的方法使うにー性能チーニング配列 sp_ienv 要素 5 L 因子配列 nzval 大さ決

ませ

並列ア詳細 Demmel et al. (1999c) 参照

Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:

(1) Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.

(2) Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.

(3) Neither the name of Lawrence Berkeley National Laboratory, U.S. Dept. of Energy nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS “AS IS” AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE

DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

例題 1

疎 6×6 行列

LU 分解します

y = (1,2,3,4,5,6)^T す b1 := Ay = (10,7,45,33,-34,31)^T b2 := A^Ty = (-9,8,39,13,1,21)^T A LU 分解疎線型方程式 Ax = b1 よび A^Tx = b2 解くに用います

#include <imsl.h>

int main(){

Imsl_f_sparse_elem a[] = { 0, 0, 10.0, 1, 1, 10.0,

1, 2, -3.0, 1, 3, -1.0, 2, 2, 15.0, 3, 0, -2.0, 3, 3, 10.0, 3, 4, -1.0, 4, 0, -1.0, 4, 3, -5.0, 4, 4, 1.0, 4, 5, -3.0, 5, 0, -1.0, 5, 1, -2.0, 5, 5, 6.0};

float b1[] = {10.0, 7.0, 45.0, 33.0, -34.0, 31.0};

float b2[] = { -9.0, 8.0, 39.0, 13.0, 1.0, 21.0 };

int n = 6, nz = 15;

float *x = NULL;

x = imsl_f_superlu_smp (n, nz, a, b1, 0);

imsl_f_write_matrix ("solution to A*x = b1", 1, n, x, 0);

imsl_free (x);

x = imsl_f_superlu_smp (n, nz, a, b2, IMSL_TRANSPOSE, 1, 0);

imsl_f_write_matrix ("solution to A^T*x = b2", 1, n, x, 0);

imsl_free (x);

}

出力

solution to A*x = b1

1 2 3 4 5 6 1 2 3 4 5 6 solution to A^T*x = b2

1 2 3 4 5 6 1 2 3 4 5 6

例題 2

ここ行列 A = E(1000,10) 用い LU 分解が異右辺持線型方程式解くに用い

こ示します最大絶対誤差プンします計算後 LU 分解たにアローさた領域関数 imsl_f_superlu_smp_factor_free 解放します

#include <imsl.h>

#include <stdlib.h>

#include <stdio.h>

int main(){

Imsl_f_sparse_elem *a = NULL;

Imsl_f_super_lu_smp_factor lu_factor;

float *b = NULL, *x = NULL, *mod_five = NULL, *mod_ten = NULL;

float error_factor_solve, error_solve;

int n = 1000, c = 10;

int i, nz, index;

/* Get the coefficient matrix */

a = imsl_f_generate_test_coordinate (n, c, &nz, 0);

/* Set two different predetermined solutions */

mod_five = (float*) malloc (n*sizeof(*mod_five));

mod_ten = (float*) malloc (n*sizeof(*mod_ten));

for (i=0; i<n; i++) {

mod_five[i] = (float) (i % 5);

mod_ten[i] = (float) (i % 10);

}

/* Choose b so that x will approximate mod_five */

b = (float *) imsl_f_mat_mul_rect_coordinate ("A*x", IMSL_A_MATRIX, n, n, nz, a,

IMSL_X_VECTOR, n, mod_five, 0);

/* Solve Ax = b */

x = imsl_f_superlu_smp (n, nz, a, b,

IMSL_RETURN_SPARSE_LU_FACTOR, &lu_factor, 0);

/* Compute max absolute error */

error_factor_solve = imsl_f_vector_norm (n, x, IMSL_SECOND_VECTOR, mod_five,

IMSL_INF_NORM, &index, 0);

free (mod_five);

imsl_free (b);

imsl_free (x);

/* Get new right hand side -- b = A * mod_ten */

b = (float *) imsl_f_mat_mul_rect_coordinate ("A*x", IMSL_A_MATRIX, n, n, nz, a,

IMSL_X_VECTOR, n, mod_ten, 0);

/* Use the previously computed factorization to solve Ax = b */

x = imsl_f_superlu_smp (n, nz, a, b,

IMSL_SUPPLY_SPARSE_LU_FACTOR, &lu_factor, IMSL_FACTOR_SOLVE, 2,

0);

error_solve = imsl_f_vector_norm (n, x, IMSL_SECOND_VECTOR, mod_ten,

IMSL_INF_NORM, &index, 0);

free (mod_ten);

imsl_free (b);

imsl_free (x);

imsl_free (a);

/* Free sparse LU structure */

imsl_f_superlu_smp_factor_free (&lu_factor);

/* Print errors */

printf ("absolute error (factor/solve) = %e\n", error_factor_solve);

printf ("absolute error (solve) = %e\n", error_solve);

}

出力

absolute error (factor/solve) = 1.096725e-005 absolute error (solve) = 5.435944e-005

警告ー

IMSL_ILL_CONDITIONED 入力行列非常に悪条件あそ L1 条件数逆

数推定値 “rcond” = # あ解正確いもしい

重大ー

IMSL_SINGULAR_MATRIX 入力行列特異あ

superlu_smp ^複素数

左方コ法一般複素疎行列 LU 分解 OpenMP並列行い複素数疎線形方程式 Ax = b 解ます

梗概

#include <imsl.h>

f_complex *imsl_c_superlu_smp (int n, int nz, Imsl_c_sparse_elem a[], f_complex b[],…,0) void imsl_c_superlu_smp_factor_free (Imsl_c_super_lu_smp_factor *factor)

d_complex 型関数 imsl_z_superlu_smp imsl_z_superlu_smp_factor_free

必須引数

int n (入力)

入力行列次数 int nz (入力)

行列非ロ数 Imsl_c_sparse_elem a[] (入力)

長さ nz 配列行列非ロン位置値含 Imsl_c_sparse_elem構造体

説明本ニアインロダクン参照し下さい f_complex b[] (入力)

右辺含長さ n 配列

戻値

疎線形方程式 Ax = b ^解 x へインタこ領域解放すに imsl_free 使う

解が得い場合 NULL が返さます

プン引数梗概

#include <imsl.h>

IMSL_RETURN_USER

, f_complex

x[]

,

)

プン引数

IMSL_EQUILIBRATE, int equilibrate (入力)

入力行列 A 因子分解前に平衡させに指定す equilibrate 説明

0 平衡化しい

1 平衡化す

フ : equilibrate = 0

IMSL_COLUMN_ORDERING_METHOD, Imsl_col_ordering method (入力)

分解前に疎保ようにす列並びえ方法次一 method に設定し並びえ法選択す

method 説明

IMSL_NATURAL 自然順序すわ入力行列列順

IMSL_MMD_ATA A^TA ^{構造体最小次} ^順 IMSL_MMD_AT_PLUS_A A^T + A ^構造体 ^最小次 ^順

IMSL_COLAMD 列近似最小次順

IMSL_PERMC プン引数 IMSL_COLPERM_VECTOR ーが入力した置換ベク permc 与え

順ベク permc 数 0,1,…,n-1 置換す

フ : method = IMSL_COLAMD

IMSL_COLPERM_VECTOR, int permc[] (入力)

後並びえ(postordering) 前置換行列 P 定義す長さ n 配列こ引数

IMSL_COLUMN_ORDERING_METHOD が method = IMSL_PERMC 用い必須あそ以外無視さ

IMSL_TRANSPOSE, int transpose (入力)

転置問題 A^Tx = b ^あ ^い A^Hx = b ^{解く場合に指定す} ^こ ^プ ^ン ^因子分解

指定すプン一緒に使うこが transpose 説明

0 Ax = b 解く

1 A^Tx = b ^解く. 2 A^Hx = b ^解く.

フ : transpose = 0

IMSL_ITERATIVE_REFINEMENT, int refine (入力) 反復改良すう指定す

refine 説明

0 反復改良しい

1 反復改良す

フ : refine = 1

ドキュメント内 imslc2016 1 math IMSL Cライブラリ v2016.1 ユーザーズガイド : Math (ページ 107-126)

プ ン引数





 

説明

例題 1

出力

例題 2

出力

警告 ー

重大 ー

superlu_smp 複素数

梗概

必須 引数

戻 値

プ ン引数 梗概

(int

, int

, Imsl_c_sparse_elem

, f_complex

,

, int

,

, Imsl_col_ordering

,

, int

,

, int

,

, int

,

, int

,

, float

,

, int

,

, int

,

, int

, int

, f_complex

,

,

,

,

, float

,

, float

, float

,

, float

,

, f_complex

,

)

プ ン引数

プン引数

警告ー

重大ー

superlu_smp ^複素数

必須引数

戻値

プン引数梗概

プン引数