EN
注册
我要评分
文档获取效率
文档正确性
内容完整性
文档易理解
在线提单
论坛求助
鲲鹏小智

?gelqf

计算矩阵的LQ分解,即A = L * Q。

接口定义

C Interface:

void sgelqf_(const int *M, const int *N, float *A, const int *LDA, float *TAU, float *WORK, const int *LWORK, int *INFO);

void dgelqf_(const int *M, const int *N, double *A, const int *LDA, double *TAU, double *WORK, const int *LWORK, int *INFO);

void cgelqf_(const int *M, const int *N, float _Complex *A, const int *LDA, float _Complex *TAU, float _Complex *WORK, const int *LWORK, int *INFO);

void zgelqf_(const int *M, const int *N, double _Complex *A, const int *LDA, double _Complex *TAU, double _Complex *WORK, const int *LWORK, int *INFO);

Fortran Interface:

DGELQF(m, n, a, lda, tau, work, lwork, info);

SGELQF(m, n, a, lda, tau, work, lwork, info);

CGELQF(m, n, a, lda, tau, work, lwork, info);

ZGELQF(m, n, a, lda, tau, work, lwork, info);

参数

参数名

类型

描述

输入/输出

m

整数型

矩阵A的行数。

输入

n

整数型

矩阵A的列数。

输入

a

  • 在sgelqf中为单精度浮点型数组。
  • 在dgelqf中为双精度浮点型数组。
  • 在cgelqf中为单精度复数型数组。
  • 在zgelqf中为双精度复数型数组。
  • 调用前保存待分解的矩阵A。
  • 调用后对角线以及对角线以上的部分保存大小为min(m,n)*n的矩阵R(当m≥n时L为下三角矩阵),对角线以下元素和tau共同表示正交矩阵Q(参见说明)。

输入/输出

lda

整数型

A的leading dimension大小,要求lda≥max(1, m)。

输出

tau

  • 在sgelqf中为单精度浮点型数组。
  • 在dgelqf中为双精度浮点型数组。
  • 在cgelqf中为单精度复数型数组。
  • 在zgelqf中为双精度复数型数组。

初等反射的系数,长度为min(m,n)(参见说明)。

输出

work

  • 在sgelqf中为单精度浮点型数组。
  • 在dgelqf中为双精度浮点型数组。
  • 在cgelqf中为单精度复数型数组。
  • 在zgelqf中为双精度复数型数组。

临时存储空间,使用lwork=-1调用后work[0]为最优的lwork值。

输出

lwork

整数型

work数组的长度。

lwork=-1时查询最优work大小,结果保存在work[0]中,否则要求lwork≥n。

输入

info

整数型

执行结果:

  • 等于0:成功。
  • 小于0:第-info个参数值不合法。

输出

分解结果矩阵Q通过一系列初等反射乘积表示:Q=H(1)*H(2)*...*H(k), k=min(m,n)。H(i)=I-tau*v*v'。tau为标量,v为向量且前i-1个元素为0,第i个元素为1,剩余元素保存在a的第i列中(a的下三角部分)。

依赖

#include "klapack.h"

示例

C Interface:

    int m = 6; 
    int n = 4; 
    int lda = 6; 
    int info = 0; 
    double tau[4]; 
    double *work = NULL; 
    double qwork; 
    int lwork = -1; 
    /* 
     * A (6x4, stored in column-major): 
     *  2.229  1.273  0.087  0.035 
     *  8.667  4.317  4.091  3.609 
     *  0.205  7.810  2.553  6.507 
     *  2.758  2.911  8.791  5.051 
     *  8.103  1.396  1.317  4.738 
     *  8.859  3.161  0.808  5.972 
     */ 
    double a[] = {2.229, 8.667, 0.205, 2.758, 8.103, 8.859, 
                    1.273, 4.317, 7.810, 2.911, 1.396, 3.161, 
                    0.087, 4.091, 2.553, 8.791, 1.317, 0.808, 
                    0.035, 3.609, 6.507, 5.051, 4.738, 5.972}; 
    /* Query optimal work size */ 
    dgelqf_(&m, &n, a, &lda, tau, &qwork, &lwork, &info); 
    if (info != 0) { 
        return ERROR; 
    } 
    lwork = (int)qwork; 
    work = (double *)malloc(sizeof(double) * lwork); 
    /* Calculate LQ */ 
    dgelqf_(&m, &n, a, &lda, tau, work, &lwork, &info); 
    free(work); 
    /* 
     * Output: 
     * tau 
     *   1.867784        1.115696        1.413422        0.000000 
     * A output (stored in column-major) 
     *  -2.568611       -9.848324       -4.223656       -4.202616
     *  -7.832678       -9.363028       0.265340        5.150248
     *  5.402586        9.567440        4.194706        4.480431
     *  0.018134        -0.653528       -7.929047       -1.156667
     *  1.133611        -0.463361       0.007295        -0.604570
     *  0.644209        -2.887734       3.399866        4.103192
     */

Fortran Interface:

        PARAMETER (m = 6) 
        PARAMETER (n = 4) 
        PARAMETER (lda = 6) 
        INTEGER :: info = 0 
        REAL(8) :: tau(4) 
        REAL(8) :: qwork(1) 
        INTEGER :: lwork = -1 
        REAL(8), ALLOCATABLE :: work(:) 
* 
*       A (6x4, stored in column-major): 
*         2.229  1.273  0.087  0.035 
*         8.667  4.317  4.091  3.609 
*         0.205  7.810  2.553  6.507 
*         2.758  2.911  8.791  5.051 
*         8.103  1.396  1.317  4.738 
*         8.859  3.161  0.808  5.972 
* 
        REAL(8) :: a(m, n) 
        DATA a / 2.229, 8.667, 0.205, 2.758, 8.103, 8.859, 
     $           1.273, 4.317, 7.810, 2.911, 1.396, 3.161, 
     $           0.087, 4.091, 2.553, 8.791, 1.317, 0.808, 
     $           0.035, 3.609, 6.507, 5.051, 4.738, 5.972 / 
 
        EXTERNAL DGELQF 
*       Query optimal work size 
        CALL DGELQF(m, n, a, lda, tau, qwork, lwork, info) 
        IF (info.NE.0) THEN 
            CALL EXIT(1) 
        END IF 
        lwork = INT(qwork(1)) 
        ALLOCATE(work(lwork)) 
*       Calculate LQ 
        CALL DGELQF(m, n, a, lda, tau, work, lwork, info) 
        DEALLOCATE(work) 
 
*       Output: 
*       tau 
*         1.867784        1.115696        1.413422        0.000000 
*       A output (stored in column-major) 
*        -2.568611       -9.848324       -4.223656       -4.202616
*        -7.832678       -9.363028       0.265340        5.150248
*        5.402586        9.567440        4.194706        4.480431
*        0.018134        -0.653528       -7.929047       -1.156667
*        1.133611        -0.463361       0.007295        -0.604570
*        0.644209        -2.887734       3.399866        4.103192
搜索结果
找到“0”个结果

当前产品无相关内容

未找到相关内容,请尝试其他搜索词