SPSOLVE 通过一系列步骤来计算稀疏方程组的解：初始化、排序以简化填充、符号因式分解、数值分解和三角求解。用户代码可以调用个体例程或使用一次调用的接口来执行这些步骤。
下表列出了 SPSOLVE 中用户可访问的例程及其用途。
|
可以按如下所示的顺序调用各个 SPSOLVE 例程，对结构相同但数值不同的矩阵进行求解：
call dgssin() ! initialization, input coefficient matrix structure
call dgssor() ! fill-reducing ordering, symbolic factorization
! (or call dgssuo() to specify a user ordering,
! and perform symbolic factorization)
do m = 1, number_of_structurally_identical_matrices
call dgssfa() ! input coefficient matrix values, numeric ! factorization
do r = 1, number_of_right_hand_sides
call dgsssl() ! triangular solve
enddo
enddo
一次调用接口不像标准接口那样灵活，但它涵盖了对单一矩阵进行因式分解以及对一些右侧数字求解的最常见情况。可以通过另外调用 dgsssl() 来对其他右侧进行求解，如下例所示。
call dgssfs() ! initialization, input coefficient matrix structure
! fill-reducing ordering, symbolic factorization
! input coefficient matrix values, numeric factorization
! triangular solve
do r = 1, number_of_right_hand_sides
call dgsssl() ! triangular solve
enddo
要使用 SPSOLVE，必须按如下所示的顺序调用其例程：
一次调用接口：用于对单一矩阵求解
DGSSFS() - 初始化、排序、分解因子、求解
DGSSSL() - 其他求解（可选）：根据需要重复 DGSSSL()
DGSSDA() - 取消分配工作存储
标准接口：用于对结构相同的多个矩阵进行求解
DGSSIN() - 初始化
DGSSOR() 或 DGSSUO() - 排序并以符号方式分解因子
DGSSFA() - 分解因子
DGSSSL() - 求解：根据需要重复 DGSSFA() 或 DGSSSL()
DGSSDA() - 取消分配工作存储
以下示例显示了使用一次调用接口对对称方程组求解，并使用标准接口对对称方程组求解。
在示例 1中，使用了一次调用接口对对称方程组求解，在示例 2中，调用了个体例程来对对称方程组求解。
示例 5显示了如何从 C 程序调用 Fortran SPSOLVE 接口。有关如何从 C 程序调用 Fortran 例程的更多信息，请参见《Oracle Solaris Studio 12.4: Fortran Programming Guide》。
示例 1 对对称方程组求解－一次调用接口my_system% cat example_1call.f
program example_1call
c
c This program is an example driver that calls the sparse solver.
c It factors and solves a symmetric system, by calling the
c one-call interface.
c
implicit none
integer neqns, ier, msglvl, outunt, ldrhs, nrhs
character mtxtyp*2, pivot*1, ordmthd*3
double precision handle(150)
integer colstr(6), rowind(9)
double precision values(9), rhs(5), xexpct(5)
integer i
c
c Sparse matrix structure and value arrays. From George and Liu,
c page 3.
c Ax = b, (solve for x) where:
c
c 4.0 1.0 2.0 0.5 2.0 2.0 7.0
c 1.0 0.5 0.0 0.0 0.0 2.0 3.0
c A = 2.0 0.0 3.0 0.0 0.0 x = 1.0 b = 7.0
c 0.5 0.0 0.0 0.625 0.0 -8.0 -4.0
c 2.0 0.0 0.0 0.0 16.0 -0.5 -4.0
c
data colstr / 1, 6, 7, 8, 9, 10 /
data rowind / 1, 2, 3, 4, 5, 2, 3, 4, 5 /
data values / 4.0d0, 1.0d0, 2.0d0, 0.5d0, 2.0d0, 0.5d0, 3.0d0,
& 0.625d0, 16.0d0 /
data rhs / 7.0d0, 3.0d0, 7.0d0, -4.0d0, -4.0d0 /
data xexpct / 2.0d0, 2.0d0, 1.0d0, -8.0d0, -0.5d0 /
c
c set calling parameters
c
mtxtyp= 'ss'
pivot = 'n'
neqns = 5
nrhs = 1
ldrhs = 5
outunt = 6
msglvl = 0
ordmthd = 'mmd'
c
c call single call interface
c
call dgssfs ( mtxtyp, pivot, neqns , colstr, rowind,
& values, nrhs , rhs, ldrhs , ordmthd,
& outunt, msglvl, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c deallocate sparse solver storage
c
call dgssda ( handle, ier )
if ( ier .ne. 0 ) goto 110
c
c print values of sol
c
write(6,200) 'i', 'rhs(i)', 'expected rhs(i)', 'error'
do i = 1, neqns
write(6,300) i, rhs(i), xexpct(i), (rhs(i)-xexpct(i))
enddo
stop
110 continue
c
c call to sparse solver returns an error
c
write ( 6 , 400 )
& ' example: FAILED sparse solver error number = ', ier
stop
200 format(a5,3a20)
300 format(i5,3d20.12) ! i/sol/xexpct values
400 format(a60,i20) ! fail message, sparse solver error number
end
my_system% f95 -dalign example_1call.f -library=sunperf
my_sytem% a.out
i rhs(i) expected rhs(i) error
1 0.200000000000D+01 0.200000000000D+01 -0.528466159722D-13
2 0.200000000000D+01 0.200000000000D+01 0.105249142734D-12
3 0.100000000000D+01 0.100000000000D+01 0.350830475782D-13
4 -0.800000000000D+01 -0.800000000000D+01 0.426325641456D-13
5 -0.500000000000D+00 -0.500000000000D+00 0.660582699652D-14
示例 2 对对称方程组求解－标准接口
my_system% cat example_ss.f
program example_ss
c
c This program is an example driver that calls the sparse solver.
c It factors and solves a symmetric system.
implicit none
integer neqns, ier, msglvl, outunt, ldrhs, nrhs
character mtxtyp*2, pivot*1, ordmthd*3
double precision handle(150)
integer colstr(6), rowind(9)
double precision values(9), rhs(5), xexpct(5)
integer i
c
c Sparse matrix structure and value arrays. From George and Liu,
c page 3.
c Ax = b, (solve for x) where:
c
c 4.0 1.0 2.0 0.5 2.0 2.0 7.0
c 1.0 0.5 0.0 0.0 0.0 2.0 3.0
c A = 2.0 0.0 3.0 0.0 0.0 x = 1.0 b = 7.0
c 0.5 0.0 0.0 0.625 0.0 -8.0 -4.0
c 2.0 0.0 0.0 0.0 16.0 -0.5 -4.0
c
data colstr / 1, 6, 7, 8, 9, 10 /
data rowind / 1, 2, 3, 4, 5, 2, 3, 4, 5 /
data values / 4.0d0, 1.0d0, 2.0d0, 0.5d0, 2.0d0, 0.5d0,
& 3.0d0, 0.625d0, 16.0d0 /
data rhs / 7.0d0, 3.0d0, 7.0d0, -4.0d0, -4.0d0 /
data xexpct / 2.0d0, 2.0d0, 1.0d0, -8.0d0, -0.5d0 /
c
c initialize solver
c
mtxtyp= 'ss'
pivot = 'n'
neqns = 5
outunt = 6
msglvl = 0
c
c call regular interface
c
call dgssin ( mtxtyp, pivot, neqns , colstr, rowind,
& outunt, msglvl, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c ordering and symbolic factorization
c
ordmthd = 'mmd'
call dgssor ( ordmthd, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c numeric factorization
c
call dgssfa ( neqns, colstr, rowind, values, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c solution
c
nrhs = 1
ldrhs = 5
call dgsssl ( nrhs, rhs, ldrhs, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c deallocate sparse solver storage
c
call dgssda ( handle, ier )
if ( ier .ne. 0 ) goto 110
c
c print values of sol
c
write(6,200) 'i', 'rhs(i)', 'expected rhs(i)', 'error'
do i = 1, neqns
write(6,300) i, rhs(i), xexpct(i), (rhs(i)-xexpct(i))
enddo
stop
110 continue
c
c call to sparse solver returns an error
c
write ( 6 , 400 )
& ' example: FAILED sparse solver error number = ', ier
stop
200 format(a5,3a20)
300 format(i5,3d20.12) ! i/sol/xexpct values
400 format(a60,i20) ! fail message, sparse solver error number
end
my_system% f95 -dalign example_ss.f -library=sunperf
my_sytem% a.out
i rhs(i) expected rhs(i) error
1 0.200000000000D+01 0.200000000000D+01 -0.528466159722D-13
2 0.200000000000D+01 0.200000000000D+01 0.105249142734D-12
3 0.100000000000D+01 0.100000000000D+01 0.350830475782D-13
4 -0.800000000000D+01 -0.800000000000D+01 0.426325641456D-13
5 -0.500000000000D+00 -0.500000000000D+00 0.660582699652D-14
示例 3 对具有非对称值的结构对称方程组求解－标准接口
my_system% cat example_su.f
program example_su
c
c This program is an example driver that calls the sparse solver.
c It factors and solves a structurally symmetric system
c (w/unsymmetric values).
c
implicit none
integer neqns, ier, msglvl, outunt, ldrhs, nrhs
character mtxtyp*2, pivot*1, ordmthd*3
double precision handle(150)
integer colstr(5), rowind(8)
double precision values(8), rhs(4), xexpct(4)
integer i
c
c Sparse matrix structure and value arrays. Coefficient matrix
c has a symmetric structure and unsymmetric values.
c Ax = b, (solve for x) where:
c
c 1.0 3.0 0.0 0.0 1.0 7.0
c 2.0 4.0 0.0 7.0 2.0 38.0
c A = 0.0 0.0 6.0 0.0 x = 3.0 b = 18.0
c 0.0 5.0 0.0 8.0 4.0 42.0
c
data colstr / 1, 3, 6, 7, 9 /
data rowind / 1, 2, 1, 2, 4, 3, 2, 4 /
data values / 1.0d0, 2.0d0, 3.0d0, 4.0d0, 5.0d0, 6.0d0, 7.0d0,
& 8.0d0 /
data rhs / 7.0d0, 38.0d0, 18.0d0, 42.0d0 /
data xexpct / 1.0d0, 2.0d0, 3.0d0, 4.0d0 /
c
c initialize solver
c
mtxtyp= 'su'
pivot = 'n'
neqns = 4
outunt = 6
msglvl = 0
c
c call regular interface
c
call dgssin ( mtxtyp, pivot, neqns , colstr, rowind,
& outunt, msglvl, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c ordering and symbolic factorization
c
ordmthd = 'mmd'
call dgssor ( ordmthd, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c numeric factorization
c
call dgssfa ( neqns, colstr, rowind, values, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c solution
c
nrhs = 1
ldrhs = 4
call dgsssl ( nrhs, rhs, ldrhs, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c deallocate sparse solver storage
c
call dgssda ( handle, ier )
if ( ier .ne. 0 ) goto 110
c
c print values of sol
c
write(6,200) 'i', 'rhs(i)', 'expected rhs(i)', 'error'
do i = 1, neqns
write(6,300) i, rhs(i), xexpct(i), (rhs(i)-xexpct(i))
enddo
stop
110 continue
c
c call to sparse solver returns an error
c
write ( 6 , 400 )
& ' example: FAILED sparse solver error number = ', ier
stop
200 format(a5,3a20)
300 format(i5,3d20.12) ! i/sol/xexpct values
400 format(a60,i20) ! fail message, sparse solver error number
end
my_system% f95 -dalign example_su.f -library=sunperf
my_system% a.out
i rhs(i) expected rhs(i) error
1 0.100000000000D+01 0.100000000000D+01 0.000000000000D+00
2 0.200000000000D+01 0.200000000000D+01 0.000000000000D+00
3 0.300000000000D+01 0.300000000000D+01 0.000000000000D+00
4 0.400000000000D+01 0.400000000000D+01 0.000000000000D+00
示例 4 对非对称方程组求解－标准接口
my_system% cat example_uu.f
program example_uu
c
c This program is an example driver that calls the sparse solver.
c It factors and solves an unsymmetric system.
c
implicit none
integer neqns, ier, msglvl, outunt, ldrhs, nrhs
character mtxtyp*2, pivot*1, ordmthd*3
double precision handle(150)
integer colstr(6), rowind(10)
double precision values(10), rhs(5), xexpct(5)
integer i
c
c Sparse matrix structure and value arrays. Unsummetric matrix A.
c Ax = b, (solve for x) where:
c
c 1.0 0.0 0.0 0.0 0.0 1.0 1.0
c 2.0 6.0 0.0 0.0 9.0 2.0 59.0
c A = 3.0 0.0 7.0 0.0 0.0 x = 3.0 b = 24.0
c 4.0 0.0 0.0 8.0 0.0 4.0 36.0
c 5.0 0.0 0.0 0.0 10.0 5.0 55.0
c
data colstr / 1, 6, 7, 8, 9, 11 /
data rowind / 1, 2, 3, 4, 5, 2, 3, 4, 2, 5 /
data values / 1.0d0, 2.0d0, 3.0d0, 4.0d0, 5.0d0, 6.0d0, 7.0d0,
& 8.0d0, 9.0d0, 10.0d0 /
data rhs / 1.0d0, 59.0d0, 24.0d0, 36.0d0, 55.0d0 /
data xexpct / 1.0d0, 2.0d0, 3.0d0, 4.0d0, 5.0d0 /
c
c initialize solver
c
mtxtyp= 'uu'
pivot = 'n'
neqns = 5
outunt = 6
msglvl = 3
call dgssin ( mtxtyp, pivot, neqns , colstr, rowind,
& outunt, msglvl, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c ordering and symbolic factorization
c
ordmthd = 'mmd'
call dgssor ( ordmthd, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c numeric factorization
c
call dgssfa ( neqns, colstr, rowind, values, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c solution
c
nrhs = 1
ldrhs = 5
call dgsssl ( nrhs, rhs, ldrhs, handle, ier )
if ( ier .ne. 0 ) goto 110
c
c deallocate sparse solver storage
c
call dgssda ( handle, ier )
if ( ier .ne. 0 ) goto 110
c
c print values of sol
c
write(6,200) 'i', 'rhs(i)', 'expected rhs(i)', 'error'
do i = 1, neqns
write(6,300) i, rhs(i), xexpct(i), (rhs(i)-xexpct(i))
enddo
stop
110 continue
c
c call to sparse solver returns an error
c
write ( 6 , 400 )
& ' example: FAILED sparse solver error number = ', ier
stop
200 format(a5,3a20)
300 format(i5,3d20.12) ! i/sol/xexpct values
400 format(a60,i20) ! fail message, sparse solver error number
end
my_system% f95 -dalign example_uu.f -library=sunperf
my_system% a.out
i rhs(i) expected rhs(i) error
1 0.100000000000D+01 0.100000000000D+01 0.000000000000D+00
2 0.200000000000D+01 0.200000000000D+01 0.000000000000D+00
3 0.300000000000D+01 0.300000000000D+01 0.000000000000D+00
4 0.400000000000D+01 0.400000000000D+01 0.000000000000D+00
5 0.500000000000D+01 0.500000000000D+01 0.000000000000D+00
示例 5 从 C 调用 SPSOLVE 例程
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <sys/time.h>
#include <sunperf.h>
int main() {
/*
Sparse matrix structure and value arrays. Coefficient matrix
is a general unsymmetric sparse matrix.
Ax = b, (solve for x) where:
1.0 0.0 7.0 9.0 0.0 1.0 17.0
2.0 4.0 0.0 0.0 0.0 1.0 6.0
A = 0.0 5.0 8.0 0.0 0.0 x = 1.0 b = 13.0
0.0 0.0 0.0 10.0 11.0 1.0 21.0
3.0 6.0 0.0 0.0 12.0 1.0 21.0
*/
/* Array indices must be one-based for calling SPSOLVE routines */
int colstr[] = {1, 4, 7, 9, 11, 13};
int rowind[] = {1, 2, 5, 2, 3, 5, 1, 3, 1, 4, 4, 5};
double values[] = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0,
7.0, 8.0, 9.0, 10.0, 11.0, 12.0};
double rhs[] = {17.0, 6.0, 13.0, 21.0, 21.0};
double xexpct[] = {1.0, 1.0, 1.0, 1.0, 1.0};
int n = 5, nnz = 12, nrhs = 1, msglvl = 0, outunt = 6, ierr,
i,j,k, int_ierr;
double t[4], handle[150];
char type[] = "uu", piv = 'n';
/* Last two parameters in argument list indicate lengths of
* character arguments type and piv
*/
dgssin_(type, &piv, &n, colstr, rowind, &outunt, &msglvl,
handle, &ierr,2,1);
if (ierr != 0) {
int_ierr = ierr;
printf("dgssin err = %d\n", int_ierr);
return -1;
}
char ordmth[] = "mmd";
dgssor_(ordmth, handle, &ierr, 3);
if (ierr != 0) {
int_ierr = ierr;
printf("dgssor err = %d\n", int_ierr);
return -1;
}
dgssfa_(&n, colstr, rowind, values, handle, &ierr);
if (ierr != 0) {
int_ierr = ierr;
printf("dgssfa err = %d\n", int_ierr);
return -1;
}
dgsssl_(&nrhs, rhs, &n, handle, &ierr);
if (ierr != 0) {
int_ierr = ierr;
printf("dgsssl err = %d\n", int_ierr);
return -1;
}
printf("i computed solution expected solution\n");
for (i=0; i<n; i++)
printf("%d %lf %lf\n", i,rhs[i], 1.0);
}
my_system% cc -m32 -xmemalign=8s dr.c -library=sunperf
my_system% ./a.out
i computed solution expected solution
0 1.000000 1.000000
1 1.000000 1.000000
2 1.000000 1.000000
3 1.000000 1.000000
4 1.000000 1.000000