LAPACK++
LAPACK C++ API
Other auxiliary routines

Functions

void lapack::lacgv (int64_t n, std::complex< double > *x, int64_t incx)
 Conjugates a complex vector of length n. More...
 
void lapack::lacgv (int64_t n, std::complex< float > *X, int64_t incx)
 
double lapack::lapy2 (double x, double y)
 Returns \(\sqrt{ x^2 + y^2 },\) taking care not to cause unnecessary overflow. More...
 
float lapack::lapy2 (float x, float y)
 
double lapack::lapy3 (double x, double y, double z)
 Returns \(\sqrt{ x^2 + y^2 + z^2 },\) taking care not to cause unnecessary overflow. More...
 
float lapack::lapy3 (float x, float y, float z)
 
int64_t lapack::lascl (lapack::MatrixType matrixtype, int64_t kl, int64_t ku, double cfrom, double cto, int64_t m, int64_t n, std::complex< double > *A, int64_t lda)
 Multiplies the m-by-n complex matrix A by the real scalar cto / cfrom. More...
 
int64_t lapack::lascl (lapack::MatrixType type, int64_t kl, int64_t ku, double cfrom, double cto, int64_t m, int64_t n, double *A, int64_t lda)
 
int64_t lapack::lascl (lapack::MatrixType type, int64_t kl, int64_t ku, float cfrom, float cto, int64_t m, int64_t n, float *A, int64_t lda)
 
int64_t lapack::lascl (lapack::MatrixType type, int64_t kl, int64_t ku, float cfrom, float cto, int64_t m, int64_t n, std::complex< float > *A, int64_t lda)
 
void lapack::lassq (int64_t n, double const *X, int64_t incx, double *scale, double *sumsq)
 
void lapack::lassq (int64_t n, float const *X, int64_t incx, float *scale, float *sumsq)
 
void lapack::lassq (int64_t n, std::complex< double > const *x, int64_t incx, double *scale, double *sumsq)
 Compute sum-of-squares without unnecessary overflow. More...
 
void lapack::lassq (int64_t n, std::complex< float > const *X, int64_t incx, float *scale, float *sumsq)
 

Detailed Description

Function Documentation

◆ lacgv()

void lapack::lacgv ( int64_t  n,
std::complex< double > *  x,
int64_t  incx 
)

Conjugates a complex vector of length n.

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>. Real precisions are dummy inline functions that do nothing, to facilitate templating.

Parameters
[in]nThe length of the vector x. n >= 0.
[in,out]xThe vector x of length n, stored in an array of length 1+(n-1)*abs(incx). On entry, the vector of length n to be conjugated. On exit, x is overwritten with conj(x).
[in]incxThe spacing between successive elements of x.

◆ lapy2()

double lapack::lapy2 ( double  x,
double  y 
)

Returns \(\sqrt{ x^2 + y^2 },\) taking care not to cause unnecessary overflow.

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>.

Parameters
[in]x
[in]yx and y specify the values x and y.

◆ lapy3()

double lapack::lapy3 ( double  x,
double  y,
double  z 
)

Returns \(\sqrt{ x^2 + y^2 + z^2 },\) taking care not to cause unnecessary overflow.

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>.

Parameters
[in]x
[in]y
[in]zx, y and z specify the values x, y and z.

◆ lascl()

int64_t lapack::lascl ( lapack::MatrixType  matrixtype,
int64_t  kl,
int64_t  ku,
double  cfrom,
double  cto,
int64_t  m,
int64_t  n,
std::complex< double > *  A,
int64_t  lda 
)

Multiplies the m-by-n complex matrix A by the real scalar cto / cfrom.

This is done without over/underflow as long as the final result cto * A(i,j) / cfrom does not over/underflow. type specifies that A may be full, upper triangular, lower triangular, upper Hessenberg, or banded.

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>.

Parameters
[in]matrixtypetype indices the storage type of the input matrix.
  • lapack::MatrixType::General: A is a full matrix.
  • lapack::MatrixType::Lower: A is a lower triangular matrix.
  • lapack::MatrixType::Upper: A is an upper triangular matrix.
  • lapack::MatrixType::Hessenberg: A is an upper Hessenberg matrix.
  • lapack::MatrixType::LowerBand: A is a symmetric band matrix with lower bandwidth kl and upper bandwidth ku and with the only the lower half stored.
  • lapack::MatrixType::UpperBand: A is a symmetric band matrix with lower bandwidth kl and upper bandwidth ku and with the only the upper half stored.
  • lapack::MatrixType::Band: A is a band matrix with lower bandwidth kl and upper bandwidth ku. See lapack::gbtrf for storage details.
[in]klThe lower bandwidth of A. Referenced only if type = LowerBand, UpperBand, or Band.
[in]kuThe upper bandwidth of A. Referenced only if type = LowerBand, UpperBand, or Band.
[in]cfrom
[in]ctoThe matrix A is multiplied by cto/cfrom. A(i,j) is computed without over/underflow if the final result cto*A(i,j)/cfrom can be represented without over/underflow. cfrom must be nonzero.
[in]mThe number of rows of the matrix A. m >= 0.
[in]nThe number of columns of the matrix A. n >= 0.
[in,out]AThe m-by-n matrix A, stored in an lda-by-n array. The matrix to be multiplied by cto/cfrom. See matrixtype for the storage type.
[in]ldaThe leading dimension of the array A.
  • If matrixtype = General, Lower, Upper, or Hessenberg, lda >= max(1,m);
  • if matrixtype = LowerBand, lda >= kl+1;
  • if matrixtype = UpperBand, lda >= ku+1;
  • if matrixtype = Band, lda >= 2*kl+ku+1.
Returns
= 0: successful exit

◆ lassq()

void lapack::lassq ( int64_t  n,
std::complex< double > const *  x,
int64_t  incx,
double *  scale,
double *  sumsq 
)

Compute sum-of-squares without unnecessary overflow.

Returns the values scl and ssq such that

\[ scl^2 ssq = x_1^2 + \dots + x_n^2 + scale^2 sumsq, \]

where \(x_i = | x( 1 + ( i - 1 )*incx ) |, 1 \le i \le n.\) The value of sumsq is assumed to be at least unity and the value of ssq will then satisfy

1.0 <= ssq <= ( sumsq + 2*n ).

scale is assumed to be non-negative and scl returns the value

\[ scl = \max( scale, |real( x_i )|, |imag( x_i )| ). \]

scale and sumsq must be supplied in scale and sumsq respectively. scale and sumsq are overwritten by scl and ssq respectively.

The routine makes only one pass through the vector x.

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>.

Parameters
[in]nThe number of elements to be used from the vector x.
[in]xThe vector x of length 1+(n-1)*incx.
[in]incxThe increment between successive values of the vector x. incx > 0.
[in,out]scaleOn entry, the value scale in the equation above. On exit, scale is overwritten with the value scl.
[in,out]sumsqOn entry, the value sumsq in the equation above. On exit, sumsq is overwritten with the value ssq.