?(or,un)mbr

计算C=Q*C或Q**T*C或C*Q或C*Q**T或C=P*C或P**T*C或C*P或C*P**T，其中正交矩阵Q和P**T由函数GEBRD计算得到，即A=Q*B*P**T。

接口定义

C Interface：

sormbr_(const char *vect, const char *side, const char *trans, const int *m, const int *n, const int *k, const float *a, const int *lda, const float *tau, float *c, const int *ldc, float *work, const int *lwork, int *info);

dormbr_(const char *vect, const char *side, const char *trans, const int *m, const int *n, const int *k, const double *a, const int *lda, const double *tau, double *c, const int *ldc, double *work, const int *lwork, int *info);

cunmbr_(const char *vect, const char *side, const char *trans, const int *m, const int *n, const int *k, const float _Complex *a, const int *lda, const float _Complex *tau, float _Complex *c, const int *ldc, float _Complex *work, const int *lwork, int *info);

zunmbr_(const char *vect, const char *side, const char *trans, const int *m, const int *n, const int *k, const double _Complex *a, const int *lda, const double _Complex *tau, double _Complex *c, const int *ldc, double _Complex *work, const int *lwork, int *info);

Fortran Interface：

SORMBR(VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO);

DORMBR(VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO);

CUNMBR(VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO);

ZUNMBR(VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO);

参数

参数

类型

说明

输入/输出

vect

字符型
  • 'Q'：则使用Q计算。
  • 'P'：则使用P计算。

输入

side

字符型
  • 'L'：则Q，P左乘C。
  • 'R'：则Q，P右乘C。

输入

trans

字符型
  • 'N'：Q，P不转置。
  • 'T'：Q，P转置。

输入

m

整数型

矩阵C的行数，m≥0。

输入

n

整数型

矩阵C的列数，n≥0。

输入

k

整数型

Q的列数或P的行数，k≥0。

输入

a
  • 在sormbr中为单精度浮点型数组。
  • 在dormbr中为双精度浮点型数组。
  • 在cunmbr中为单精度复数型数组。
  • 在zunmbr中为双精度复数型数组。
  • 若vect='Q'，大小为lda*min(nq,k)。
  • 若vect='P'，大小为lda*nq。

该数组表示了反射器H和G，可通过两者的乘积计算出Q和P。

输入

lda

整数型

矩阵A的主维。

  • 若vect='Q'，则lda≥max(1,nq)。
  • 若vect='P'，则lda≥(1, min(nq, k))。

输入

tau
  • 在sormbr中为单精度浮点型数组。
  • 在dormbr中为双精度浮点型数组。
  • 在cunmbr中为单精度复数型数组。
  • 在zunmbr中为双精度复数型数组。

大小为min(nq,k)。

该数组表示H和G的标量因子。

输入

c
  • 在sormbr中为单精度浮点型数组。
  • 在dormbr中为双精度浮点型数组。
  • 在cunmbr中为单精度复数型数组。
  • 在zunmbr中为双精度复数型数组。

大小为ldc*n。

  • 进入时：为m*n的矩阵。
  • 退出时：为Q*C或Q**H*C或C*Q**H或C*Q或P*C或P**H*C或C*P或C*P**H。

输入，输出

ldc

整数型

矩阵c的主维，ldc≥max(1, m)。

输入

work
  • 在sormbr中为单精度浮点型数组。
  • 在dormbr中为双精度浮点型数组。
  • 在cunmbr中为单精度复数型数组。
  • 在zunmbr中为双精度复数型数组。

工作数组，大小为max(1, lwork)。

info=0时，work(0)返回最优lwork大小。

输出

lwork

整数型

work的大小。

  • 若side='L'， lwork≥max(1,n)。
  • 若side='R'，lwork≥max(1,m)。
  • 若lwork=-1，则返回最优的大小。

输入

Info

整数型
  • 等于0：表示成功。
  • 小于0：info=-i，表示第i个参数非法。

输出

依赖

#include "klapack.h"

示例

C Interface：

const char vect = 'Q';
const char side = 'L';
const char trans = 'N';
const int n = 4;
const int m = 4;
const int k = 4;
const int lda = n;
const int ldc = m;

int info = 0;
double *work = NULL;
double qwork;
int lwork = -1;

double *d = (double*)malloc(n * sizeof(double));
double *e = (double*)malloc((n - 1) * sizeof(double));
double *tauq = (double*)malloc(n * sizeof(double));
double *taup = (double*)malloc(n * sizeof(double));

double a[] = {7.027, 8.710, 1.015, 6.929,
             8.710, 0.839, 2.469, 3.850,
             1.015, 2.469, 1.930, 6.761,
             6.929, 3.850, 6.761, 4.344};

double c[] = {7.027, 8.710, 1.015, 6.929,
             8.710, 0.839, 2.469, 3.850,
             1.015, 2.469, 1.930, 6.761,
             6.929, 3.850, 6.761, 4.344};

/* Query optimal work size */
dgebrd_(&m, &n, a, &lda, d, e, tauq, taup, &qwork, &lwork, &info);
if (info != 0) {
    printf("ERROR, info = %d\n", info);
    return info;
}
lwork = (int)qwork;
work = (double *)malloc(sizeof(double) * lwork);
/* Calculate Q */
dgebrd_(&m, &n, a, &lda, d, e, tauq, taup, work, &lwork, &info);
if (info != 0) {
    printf("ERROR, info = %d\n", info);
    return info;
}

dormbr_(&vect, &side, &trans, &m, &n, &k, a, &lda, tauq, c, &ldc, work, &lwork, &info);
if (info != 0) {
    printf("ERROR, info = %d\n", info);
    return info;
}
/* output */
* a
* -13.201670      13.064517       0.286648        0.441152
* 0.430577        -8.159976       2.308820        -0.069390
* 0.050176        -0.414492       5.045717        -1.018082
* 0.342534        -0.018475       0.799162        -3.517170

Fortran Interface:

CHARACTER :: vect = "Q"
CHARACTER :: side = "L"
CHARACTER :: tarns = "N"
PARAMETER (n = 4) 
PARAMETER (m = 4)
PARAMETER (lda = 4) 
PARAMETER (ldc = 4)  
PARAMETER (k = 4) 

INTEGER :: info = 0 
REAL(8) :: a(lda, n) 
REAL(8) :: c(ldc, n)
REAL(8) :: d(n) 
REAL(8) :: e(n-1) 
REAL(8) :: tauq(n) 
REAL(8) :: taup(n) 
REAL(8), ALLOCATABLE :: work(:)
REAL(8) :: qwork(1)
INTEGER :: lwork = -1
  
DATA a / 7.027, 8.710, 1.015, 6.929,
&        8.710, 0.839, 2.469, 3.850,
&        1.015, 2.469, 1.930, 6.761,
&        6.929, 3.850, 6.761, 4.344 / 

DATA c / 7.027, 8.710, 1.015, 6.929,
&        8.710, 0.839, 2.469, 3.850,
&        1.015, 2.469, 1.930, 6.761,
&        6.929, 3.850, 6.761, 4.344 /
EXTERNAL DGEBRD,DORMBR
CALL DGEBRD(m, n, a, lda, d, e, tauq, taup, qwork, lwork, info)
IF (info.NE.0) THEN 
    CALL EXIT(1) 
END IF
lwork = INT(qwork(1))  
ALLOCATE(work(lwork))
CALL DGEBRD(m, n, a, lda, d, e, tauq, taup, qwork, lwork, info)
IF (info.NE.0) THEN 
    CALL EXIT(1) 
END IF
CALL DORMBR(vect, side, trans, m, n, k, a, lda, tauq, c, ldc, work, lwork, info)
DEALLOCATE(work);
* 
* Output: 
* a
* -13.201670      13.064517       0.286648        0.441152
* 0.430577        -8.159976       2.308820        -0.069390
* 0.050176        -0.414492       5.045717        -1.018082
* 0.342534        -0.018475       0.799162        -3.517170
父主题： 其他函数