LAPACK
3.5.0
LAPACK: Linear Algebra PACKage
Main Page
Data Types List
Files
File List
File Members
All
Classes
Files
Functions
Variables
Typedefs
Macros
example_DGESV_rowmajor.c
Go to the documentation of this file.
1
/*
2
LAPACKE_dgesv Example
3
=====================
4
5
The program computes the solution to the system of linear
6
equations with a square matrix A and multiple
7
right-hand sides B, where A is the coefficient matrix
8
and b is the right-hand side matrix:
9
10
Description
11
===========
12
13
The routine solves for X the system of linear equations A*X = B,
14
where A is an n-by-n matrix, the columns of matrix B are individual
15
right-hand sides, and the columns of X are the corresponding
16
solutions.
17
18
The LU decomposition with partial pivoting and row interchanges is
19
used to factor A as A = P*L*U, where P is a permutation matrix, L
20
is unit lower triangular, and U is upper triangular. The factored
21
form of A is then used to solve the system of equations A*X = B.
22
23
LAPACKE Interface
24
=================
25
26
LAPACKE_dgesv (row-major, high-level) Example Program Results
27
28
-- LAPACKE Example routine (version 3.5.0) --
29
-- LAPACK is a software package provided by Univ. of Tennessee, --
30
-- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
31
February 2012
32
33
*/
34
#include <stdlib.h>
35
#include <stdio.h>
36
#include <string.h>
37
#include <
lapacke.h
>
38
#include "
lapacke_example_aux.h
"
39
40
/* Main program */
41
int
main
(
int
argc,
char
**argv) {
42
43
/* Locals */
44
lapack_int
n, nrhs, lda, ldb, info;
45
int
i,
j
;
46
double
normr, normb;
47
/* Local arrays */
48
double
*
A
, *
b
, *Acopy, *bcopy;
49
lapack_int
*ipiv;
50
51
/* Default Value */
52
n = 5; nrhs = 1;
53
54
/* Arguments */
55
for
( i = 1; i < argc; i++ ) {
56
if
( strcmp( argv[i],
"-n"
) == 0 ) {
57
n = atoi(argv[i+1]);
58
i++;
59
}
60
if
( strcmp( argv[i],
"-nrhs"
) == 0 ) {
61
nrhs = atoi(argv[i+1]);
62
i++;
63
}
64
}
65
66
/* Initialization */
67
lda=n, ldb=nrhs;
68
A = (
double
*)malloc(n*n*
sizeof
(
double
)) ;
69
if
(A==NULL){ printf(
"error of memory allocation\n"
); exit(0); }
70
b = (
double
*)malloc(n*nrhs*
sizeof
(
double
)) ;
71
if
(b==NULL){ printf(
"error of memory allocation\n"
); exit(0); }
72
ipiv = (
lapack_int
*)malloc(n*
sizeof
(
lapack_int
)) ;
73
if
(ipiv==NULL){ printf(
"error of memory allocation\n"
); exit(0); }
74
75
for
( i = 0; i < n; i++ ) {
76
for
( j = 0; j < n; j++ ) A[i*lda+j] = ((
double
) rand()) / ((double) RAND_MAX) - 0.5;
77
}
78
79
for
(i=0;i<n*nrhs;i++)
80
b[i] = ((
double
) rand()) / ((
double
) RAND_MAX) - 0.5;
81
82
/* Print Entry Matrix */
83
print_matrix_rowmajor
(
"Entry Matrix A"
, n, n, A, lda );
84
/* Print Right Rand Side */
85
print_matrix_rowmajor
(
"Right Rand Side b"
, n, nrhs, b, ldb );
86
printf(
"\n"
);
87
/* Executable statements */
88
printf(
"LAPACKE_dgesv (row-major, high-level) Example Program Results\n"
);
89
/* Solve the equations A*X = B */
90
info =
LAPACKE_dgesv
(
LAPACK_ROW_MAJOR
, n, nrhs, A, lda, ipiv,
91
b, ldb );
92
/* Check for the exact singularity */
93
if
( info > 0 ) {
94
printf(
"The diagonal element of the triangular factor of A,\n"
);
95
printf(
"U(%i,%i) is zero, so that A is singular;\n"
, info, info );
96
printf(
"the solution could not be computed.\n"
);
97
exit( 1 );
98
}
99
if
(info <0) exit( 1 );
100
/* Print solution */
101
print_matrix_rowmajor
(
"Solution"
, n, nrhs, b, ldb );
102
/* Print details of LU factorization */
103
print_matrix_rowmajor
(
"Details of LU factorization"
, n, n, A, lda );
104
/* Print pivot indices */
105
print_vector
(
"Pivot indices"
, n, ipiv );
106
exit( 0 );
107
}
/* End of LAPACKE_dgesv Example */
108
lapacke
example
example_DGESV_rowmajor.c
Generated on Mon Dec 30 2013 16:09:14 for LAPACK by
1.8.1.2