LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Files Functions Typedefs Macros
example_DGESV_colmajor.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 (col-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 /* Includes */
35 #include <stdlib.h>
36 #include <stdio.h>
37 #include <string.h>
38 #include "lapacke.h"
39 #include "lapacke_example_aux.h"
40 
41 /* Main program */
42 int main(int argc, char **argv) {
43 
44  /* Locals */
45  lapack_int n, nrhs, lda, ldb, info;
46  int i, j;
47  double normr, normb;
48  /* Local arrays */
49  double *A, *b, *Acopy, *bcopy;
50  lapack_int *ipiv;
51 
52  /* Default Value */
53  n = 5; nrhs = 1;
54 
55  /* Arguments */
56  for( i = 1; i < argc; i++ ) {
57  if( strcmp( argv[i], "-n" ) == 0 ) {
58  n = atoi(argv[i+1]);
59  i++;
60  }
61  if( strcmp( argv[i], "-nrhs" ) == 0 ) {
62  nrhs = atoi(argv[i+1]);
63  i++;
64  }
65  }
66 
67  /* Initialization */
68  lda=n, ldb=n;
69  A = (double *)malloc(n*n*sizeof(double)) ;
70  if (A==NULL){ printf("error of memory allocation\n"); exit(0); }
71  b = (double *)malloc(n*nrhs*sizeof(double)) ;
72  if (b==NULL){ printf("error of memory allocation\n"); exit(0); }
73  ipiv = (lapack_int *)malloc(n*sizeof(lapack_int)) ;
74  if (ipiv==NULL){ printf("error of memory allocation\n"); exit(0); }
75 
76  for( i = 0; i < n; i++ ) {
77  for( j = 0; j < n; j++ ) A[i+j*lda] = ((double) rand()) / ((double) RAND_MAX) - 0.5;
78  }
79 
80  for(i=0;i<n*nrhs;i++)
81  b[i] = ((double) rand()) / ((double) RAND_MAX) - 0.5;
82 
83  /* Print Entry Matrix */
84  print_matrix_colmajor( "Entry Matrix A", n, n, A, lda );
85  /* Print Right Rand Side */
86  print_matrix_colmajor( "Right Rand Side b", n, nrhs, b, ldb );
87  printf( "\n" );
88 
89  /* Executable statements */
90  printf( "LAPACKE_dgesv (row-major, high-level) Example Program Results\n" );
91  /* Solve the equations A*X = B */
92  info = LAPACKE_dgesv( LAPACK_COL_MAJOR, n, nrhs, A, lda, ipiv,
93  b, ldb );
94 
95  /* Check for the exact singularity */
96  if( info > 0 ) {
97  printf( "The diagonal element of the triangular factor of A,\n" );
98  printf( "U(%i,%i) is zero, so that A is singular;\n", info, info );
99  printf( "the solution could not be computed.\n" );
100  exit( 1 );
101  }
102  if (info <0) exit( 1 );
103  /* Print solution */
104  print_matrix_colmajor( "Solution", n, nrhs, b, ldb );
105  /* Print details of LU factorization */
106  print_matrix_colmajor( "Details of LU factorization", n, n, A, lda );
107  /* Print pivot indices */
108  print_vector( "Pivot indices", n, ipiv );
109  exit( 0 );
110 } /* End of LAPACKE_dgesv Example */
111 
#define LAPACK_COL_MAJOR
Definition: lapacke.h:120
int main(int argc, char **argv)
void print_matrix_colmajor(char *desc, lapack_int m, lapack_int n, double *mat, lapack_int ldm)
void print_vector(char *desc, lapack_int n, lapack_int *vec)
#define lapack_int
Definition: lapacke.h:47
lapack_int LAPACKE_dgesv(int matrix_order, lapack_int n, lapack_int nrhs, double *a, lapack_int lda, lapack_int *ipiv, double *b, lapack_int ldb)
Definition: lapacke_dgesv.c:36