LAPACK Archives

[Lapack] LAPACKE_dlange and LAPACKE_dgecon

Dear Scott,
DGECON estimates the reciprocal of the condition number of a general real 
matrix A, in either the 1-norm or the infinity-norm, using the LU factorization 
computed by DGETRF.
This means that you need to call DGETRF before calling DGECON.

        aNorm = LAPACKE_dlange(LAPACK_ROW_MAJOR, ONE_NORM, NROWS, NCOLS, a, 
LEADING_DIMENSION_A);
                   info = LAPACKE_dgetrf(LAPACK_ROW_MAJOR, NROWS, NCOLS, a, 
LEADING_DIMENSION_A, ipiv);
        info = LAPACKE_dgecon(LAPACK_ROW_MAJOR, ONE_NORM, n, a, 
LEADING_DIMENSION_A, aNorm, &rcond); // aNorm should be 35.019999999999996

will give you the correct result.
Regards,
Julie

On Mar 9, 2012, at 8:49 AM, Lee, Scott wrote:

Julie,
 
Thanks for finding my silly mistake.  I have another question though.  Now 
when I run the example I get:
 
LAPACKE_dlange / One-norm of A = 35.020000
LAPACKE_dgecon / RCOND of A    = 0.000221
 
I compared these results with the one-norm and RCOND of matrix A using Matlab 
(below).  The one-norm from LAPACKE agrees with Matlab, but the reciprocal 
condition number is different by two orders of magnitude.   Am I still doing 
something wrong?  For reference, I?ve included the updated LAPACKE code at 
the end.
 
Thanks.
-Scott
 
Matlab results:
-------------------
A
A =
  Columns 1 through 4
   6.800000000000000  -6.050000000000000  -0.450000000000000   
8.320000000000000
  -2.110000000000000  -3.300000000000000   2.580000000000000   
2.710000000000000
   5.660000000000000   5.360000000000000  -2.700000000000000   
4.350000000000000
   5.970000000000000  -4.440000000000000   0.270000000000000  
-7.170000000000000
   8.230000000000000   1.080000000000000   9.039999999999999   
2.140000000000000
  Column 5
  -9.670000000000000
  -5.140000000000000
  -7.260000000000000
   6.080000000000000
  -6.870000000000000
norm(A,1)
ans =
  35.019999999999996
rcond(A)
ans =
   0.035419346601492
 
LAPACKE example:
------------------------
#include <stdlib.h>
#include <stdio.h>
#include "lapacke.h"
 
/* Auxiliary routines prototypes */
extern void print_matrix( char* desc, lapack_int m, lapack_int n, double* a, 
lapack_int lda );
extern void print_int_vector( char* desc, lapack_int n, lapack_int* a );
 
/* Parameters */
#define N 5
#define NRHS 3
#define LDA N
#define LDB NRHS
 
/* Main program */
int main() {
        /* Locals */
        lapack_int n = N, nrhs = NRHS, lda = LDA, ldb = LDB, info;
        /* Local arrays */
        lapack_int ipiv[N];
        double a[LDA*N] = {
            6.80, -6.05, -0.45,  8.32, -9.67,
           -2.11, -3.30,  2.58,  2.71, -5.14,
            5.66, 5.36, -2.70,  4.35, -7.26,
            5.97, -4.44,  0.27, -7.17, 6.08,
            8.23, 1.08,  9.04,  2.14, -6.87
        };
        double b[LDB*N] = {
            4.02, -1.56, 9.81,
            6.19,  4.00, -4.09,
           -8.22, -8.67, -4.57,
           -7.57,  1.75, -8.61,
           -3.03,  2.86, 8.99
        };
 
       double aNorm;
       double rcond;
       char ONE_NORM = '1';
       lapack_int NROWS = n;
       lapack_int NCOLS = n;
       lapack_int LEADING_DIMENSION_A = n;
 
              /* Print Entry Matrix */
        print_matrix( "Entry Matrix A", n, n, a, lda );
        /* Print Right Rand Side */
        print_matrix( "Right Rand Side", n, nrhs, b, ldb );
        printf( "\n" );
        /* Executable statements */
        printf( "LAPACKE_dgecon Example Program Results\n" );
       aNorm = LAPACKE_dlange(LAPACK_ROW_MAJOR, ONE_NORM, NROWS, NCOLS, a, 
LEADING_DIMENSION_A);
        info = LAPACKE_dgecon(LAPACK_ROW_MAJOR, ONE_NORM, n, a, 
LEADING_DIMENSION_A, aNorm, &rcond); // aNorm should be 35.019999999999996
        /* Check for the exact singularity */
              if (info == 0)
              {
                     printf("LAPACKE_dgecon completed SUCCESSFULLY...\n");
              }
              else if ( info < 0 )
              {
            printf( "Element %d of A had an illegal value\n", -info );
            exit( 1 );
        }
              else
              {
            printf( "Unrecognized value of INFO = %d\n", info );
            exit( 1 );
              }
 
        /* Print solution */
       printf("LAPACKE_dlange / One-norm of A = %lf\n", aNorm);
        printf("LAPACKE_dgecon / RCOND of A    = %lf\n", rcond);
        exit( 0 );
} /* End of LAPACKE_dgesv Example */
 
/* Auxiliary routine: printing a matrix */
void print_matrix( char* desc, lapack_int m, lapack_int n, double* a, 
lapack_int lda ) {
        lapack_int i, j;
        printf( "\n %s\n", desc );
        for( i = 0; i < m; i++ ) {
                for( j = 0; j < n; j++ ) printf( " %6.2f", a[i*lda+j] );
                printf( "\n" );
        }
}
 
/* Auxiliary routine: printing a vector of integers */
void print_int_vector( char* desc, lapack_int n, lapack_int* a ) {
        lapack_int j;
        printf( "\n %s\n", desc );
        for( j = 0; j < n; j++ ) printf( " %6i", a[j] );
        printf( "\n" );
}
 
 
From: julie langou [mailto:julie@Domain.Removed] 
Sent: Thursday, March 08, 2012 6:14 PM
To: Lee, Scott
Cc: lapack@Domain.Removed
Subject: Re: [Lapack] LAPACKE_dlange and LAPACKE_dgecon
 
Dear Scott,
here is your program working..
First. it is better to use lapack_int than just int
I just had to correct one line 
        aNorm = LAPACKE_dlange( LAPACK_ROW_MAJOR, ONE_NORM, NROWS, NCOLS, a, 
n); //, WORK);
Also you can just pass numbers directly instead of declaring variable.
The following works too:
        aNorm = LAPACKE_dlange( LAPACK_ROW_MAJOR, '1', n, n, a, n); //, WORK);
        info = LAPACKE_dgecon( LAPACK_ROW_MAJOR, '1', n, a, n, aNorm, rcond); 
// aNorm should be 35.019999999999996
Let me know if you have any question
Julie
 
 
#include <stdlib.h>
#include <stdio.h>
#include "lapacke.h"
 
/* Auxiliary routines prototypes */
extern void print_matrix( char* desc, lapack_int m, lapack_int n, double* a, 
lapack_int lda );
extern void print_int_vector( char* desc, lapack_int n, lapack_int* a );
 
/* Parameters */
#define N 5
#define NRHS 3
#define LDA N
#define LDB NRHS
 
/* Main program */
int main() {
        /* Locals */
        lapack_int n = N, nrhs = NRHS, lda = LDA, ldb = LDB, info;
        /* Local arrays */
        lapack_int ipiv[N];
        double a[LDA*N] = {
            6.80, -6.05, -0.45,  8.32, -9.67,
           -2.11, -3.30,  2.58,  2.71, -5.14,
            5.66, 5.36, -2.70,  4.35, -7.26,
            5.97, -4.44,  0.27, -7.17, 6.08,
            8.23, 1.08,  9.04,  2.14, -6.87
        };
        double b[LDB*N] = {
            4.02, -1.56, 9.81,
            6.19,  4.00, -4.09,
           -8.22, -8.67, -4.57,
           -7.57,  1.75, -8.61,
           -3.03,  2.86, 8.99
        };
 
       double aNorm;
       double rcond[1] = { 0.0 };
       char ONE_NORM = '1';
       lapack_int MATRIX_ORDER = n;
       lapack_int NROWS = n;
       lapack_int NCOLS = n;
       lapack_int LEADING_DIMENSION_A = n;
       double WORK[N];
 
              /* Print Entry Matrix */
        print_matrix( "Entry Matrix A", n, n, a, lda );
        /* Print Right Rand Side */
        print_matrix( "Right Rand Side", n, nrhs, b, ldb );
        printf( "\n" );
        /* Executable statements */
        printf( "LAPACKE_dgecon Example Program Results\n" );
              // DLANGE( NORM, M, N, A, LDA, WORK )
              // LEFT OFF...When running this example, I get "Wrong parameter 
1 in LAPACKE_dlange
        aNorm = LAPACKE_dlange( LAPACK_ROW_MAJOR, ONE_NORM, NROWS, NCOLS, a, 
n); //, WORK);
        info = LAPACKE_dgecon( LAPACK_ROW_MAJOR, ONE_NORM, MATRIX_ORDER, a, 
n, aNorm, rcond); // aNorm should be 35.019999999999996
        /* Check for the exact singularity */
              if (info == 0)
              {
                     printf("LAPACKE_dgecon completed SUCCESSFULLY...\n");
              }
              else if ( info < 0 )
              {
            printf( "Element %d of A had an illegal value\n", -info );
            exit( 1 );
        }
              else
              {
            printf( "Unrecognized value of INFO = %d\n", info );
            exit( 1 );
              }
 
        /* Print solution */
       printf("LAPACKE_dlange / One-norm of A = %lf\n", aNorm);
        printf("LAPACKE_dgecon / RCOND of A    = %e", rcond[1] );
       exit( 0 );
} /* End of LAPACKE_dgesv Example */
 
/* Auxiliary routine: printing a matrix */
void print_matrix( char* desc, lapack_int m, lapack_int n, double* a, 
lapack_int lda ) {
        lapack_int i, j;
        printf( "\n %s\n", desc );
        for( i = 0; i < m; i++ ) {
                for( j = 0; j < n; j++ ) printf( " %6.2f", a[i*lda+j] );
                printf( "\n" );
        }
}
 
/* Auxiliary routine: printing a vector of integers */
void print_int_vector( char* desc, lapack_int n, lapack_int* a ) {
        lapack_int j;
        printf( "\n %s\n", desc );
        for( j = 0; j < n; j++ ) printf( " %6i", a[j] );
        printf( "\n" );
}
 
 
On Mar 8, 2012, at 2:09 PM, Lee, Scott wrote:


Hi,
 
I modified one of the LAPACKE_examples to compute the matrix condition number 
using LAPACKE_dlange and LAPACKE_dgecon.  See changes to the example 
highlighted in yellow, at end of the e-mail.
 
When I run the modified example, I get the following output:
 
Entry Matrix A
   6.80  -6.05  -0.45   8.32  -9.67
  -2.11  -3.30   2.58   2.71  -5.14
   5.66   5.36  -2.70   4.35  -7.26
   5.97  -4.44   0.27  -7.17   6.08
   8.23   1.08   9.04   2.14  -6.87
 
Right Rand Side
   4.02  -1.56   9.81
   6.19   4.00  -4.09
  -8.22  -8.67  -4.57
  -7.57   1.75  -8.61
  -3.03   2.86   8.99
 
LAPACKE_dgecon Example Program Results
Wrong parameter 1 in LAPACKE_dlange
** On entry to DGECON parameter number  5 had an illegal value
 
I?ve tried several things but I can?t get past this error.  Any ideas?
 
Thanks.
-Scott Lee
 
 
Modified example:
---------------------------
#include <stdlib.h>
#include <stdio.h>
#include "lapacke.h"
 
/* Auxiliary routines prototypes */
extern void print_matrix( char* desc, lapack_int m, lapack_int n, double* a, 
lapack_int lda );
extern void print_int_vector( char* desc, lapack_int n, lapack_int* a );
 
/* Parameters */
#define N 5
#define NRHS 3
#define LDA N
#define LDB NRHS
 
/* Main program */
int main() {
        /* Locals */
        lapack_int n = N, nrhs = NRHS, lda = LDA, ldb = LDB, info;
        /* Local arrays */
        lapack_int ipiv[N];
        double a[LDA*N] = {
            6.80, -6.05, -0.45,  8.32, -9.67,
           -2.11, -3.30,  2.58,  2.71, -5.14,
            5.66, 5.36, -2.70,  4.35, -7.26,
            5.97, -4.44,  0.27, -7.17, 6.08,
            8.23, 1.08,  9.04,  2.14, -6.87
        };
        double b[LDB*N] = {
            4.02, -1.56, 9.81,
            6.19,  4.00, -4.09,
           -8.22, -8.67, -4.57,
           -7.57,  1.75, -8.61,
           -3.03,  2.86, 8.99
        };
 
       double aNorm;
      double rcond[1] = { 0.0 };
       char ONE_NORM = '1';
       int MATRIX_ORDER = n;
       int NROWS = n;
       int NCOLS = n;
       int LEADING_DIMENSION_A = n;
       double WORK[N];
 
              /* Print Entry Matrix */
        print_matrix( "Entry Matrix A", n, n, a, lda );
        /* Print Right Rand Side */
        print_matrix( "Right Rand Side", n, nrhs, b, ldb );
        printf( "\n" );
        /* Executable statements */
        printf( "LAPACKE_dgecon Example Program Results\n" );
              // DLANGE( NORM, M, N, A, LDA, WORK )
              // LEFT OFF...When running this example, I get "Wrong parameter 
1 in LAPACKE_dlange
       aNorm = LAPACKE_dlange(MATRIX_ORDER, ONE_NORM, NROWS, NCOLS, a, 
LEADING_DIMENSION_A); //, WORK);
        info = LAPACKE_dgecon( LAPACK_ROW_MAJOR, ONE_NORM, MATRIX_ORDER, a, 
LEADING_DIMENSION_A, aNorm, rcond); // aNorm should be 35.019999999999996
        /* Check for the exact singularity */
              if (info == 0)
              {
                     printf("LAPACKE_dgecon completed SUCCESSFULLY...\n");
              }
              else if ( info < 0 )
              {
            printf( "Element %d of A had an illegal value\n", -info );
            exit( 1 );
        }
              else
              {
            printf( "Unrecognized value of INFO = %d\n", info );
            exit( 1 );
              }
 
        /* Print solution */
       printf("LAPACKE_dlange / One-norm of A = %lf\n", aNorm);
        printf("LAPACKE_dgecon / RCOND of A    = %e", rcond[1] );
       exit( 0 );
} /* End of LAPACKE_dgesv Example */
 
/* Auxiliary routine: printing a matrix */
void print_matrix( char* desc, lapack_int m, lapack_int n, double* a, 
lapack_int lda ) {
        lapack_int i, j;
        printf( "\n %s\n", desc );
        for( i = 0; i < m; i++ ) {
                for( j = 0; j < n; j++ ) printf( " %6.2f", a[i*lda+j] );
                printf( "\n" );
        }
}
 
/* Auxiliary routine: printing a vector of integers */
void print_int_vector( char* desc, lapack_int n, lapack_int* a ) {
        lapack_int j;
        printf( "\n %s\n", desc );
        for( j = 0; j < n; j++ ) printf( " %6i", a[j] );
        printf( "\n" );
}
 
 
_______________________________________________
Lapack mailing list
Lapack@Domain.Removed
http://lists.eecs.utk.edu/mailman/listinfo/lapack
 

-------------- next part --------------
An HTML attachment was scrubbed...
URL: 
http://lists.eecs.utk.edu/mailman/private/lapack/attachments/20120309/639cbe43/attachment-0001.html
 

<Prev in Thread] Current Thread [Next in Thread>


For additional information you may use the LAPACK/ScaLAPACK Forum.
Or one of the mailing lists, or