Geek Answers Handbook

Calculate the number of inverse conditions using lapak (i.e. Rcond (x))

I want to do what rcond does in MATLAB / Octave using LAPACK from C. The MATLAB manual states that dgecon is used, and that uses a 1-based norm.

I wrote a simple test program for an extremely simple case; [1,1; 1.0] For this input, matlab and octave give 0.25 using rcond and 1 / cond (x, 1), but in the case of using LAPACK this sample program prints 0.0. For other cases, such as identification, it prints the correct value.

Since MATLAB really uses this procedure with success, what am I doing wrong? I am trying to decipher what Octave is doing, with little success, since its wrapped in

#include <stdio.h>

extern void dgecon_(const char *norm, const int *n, const double *a, 
     const int *lda, const double *anorm, double *rcond, double *work, 
     int *iwork, int *info, int len_norm);

int main()
{
    int i, info, n, lda;
    double anorm, rcond;

    double w[8] = { 0,0,0,0,0,0,0,0 };

    int iw[2] = { 0,0 };

    double x[4] = { 1, 1, 1, 0 };
    anorm = 2.0; /* maximum column sum, computed manually */
    n = 2;
    lda = 2;

    dgecon_("1", &n, x, &lda, &anorm, &rcond, w, iw, &info, 1);

    if (info != 0) fprintf(stderr, "failure with error %d\n", info);
    printf("%.5e\n", rcond);
    return 0;
}

Compiled with cc testdgecon.c -o testdgecon -llapack; ./ testdgecon

+5
1

.

LU- dgecon. , , . , .

- , LAPACK.

#include "stdio.h"

extern int dgecon_(const char *norm, const int *n, double *a, const int *lda, const double *anorm, double *rcond, double *work, int *iwork, int *info, int len);
extern int dgetrf_(const int *m, const int *n, double *a, const int *lda, int *lpiv, int *info);
extern double dlange_(const char *norm, const int *m, const int *n, const double *a, const int *lda, double *work, const int norm_len);

int main()
{
    int i, info, n, lda;
    double anorm, rcond;

    int iw[2];
    double w[8];
    double x[4] = {7,3,-9,2 };
    n = 2;
    lda = 2;

    /* Computes the norm of x */
    anorm = dlange_("1", &n, &n, x, &lda, w, 1);

    /* Modifies x in place with a LU decomposition */
    dgetrf_(&n, &n, x, &lda, iw, &info);
    if (info != 0) fprintf(stderr, "failure with error %d\n", info);

    /* Computes the reciprocal norm */
    dgecon_("1", &n, x, &lda, &anorm, &rcond, w, iw, &info, 1);
    if (info != 0) fprintf(stderr, "failure with error %d\n", info);

    printf("%.5e\n", rcond);
    return 0;
}
+9

More articles:


All Articles