Dsteqr.java
package org.mklab.sdpj.gpack.lapackwrap;
import org.mklab.nfc.matrix.ComplexNumericalMatrix;
import org.mklab.nfc.matrix.RealNumericalMatrix;
import org.mklab.nfc.scalar.ComplexNumericalScalar;
import org.mklab.nfc.scalar.RealNumericalScalar;
import org.mklab.sdpj.gpack.blaswrap.BLAS;
import org.mklab.sdpj.gpack.f2clibs.LibF77;
/**
* 陰的QL法またはQR法を使って対称3重対角行列のすべての固有値と固有ベクトルを求める
*
* @author takafumi
* @param <RS> type of real scalar
* @param <RM> type of real matrix
* @param <CS> type of complex scalar
* @param <CM> type of complex Matrix
*/
public class Dsteqr<RS extends RealNumericalScalar<RS, RM, CS, CM>, RM extends RealNumericalMatrix<RS, RM, CS, CM>, CS extends ComplexNumericalScalar<RS, RM, CS, CM>, CM extends ComplexNumericalMatrix<RS, RM, CS, CM>> {
/** */
private RS[] d;
/** */
private RS[] e;
/** */
private RS[] z;
/** */
private RS[] work;
/* -- LAPACK routine (version 3.0) --
Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
Courant Institute, Argonne National Lab, and Rice University
September 30, 1994
Purpose
=======
DSTEQR computes all eigenvalues and, optionally, eigenvectors of a
symmetric tridiagonal matrix using the implicit QL or QR method.
The eigenvectors of a full or band symmetric matrix can also be found
if DSYTRD or DSPTRD or DSBTRD has been used to reduce this matrix to
tridiagonal form.
Arguments
=========
COMPZ (input) CHARACTER*1
= 'N': Compute eigenvalues only.
= 'V': Compute eigenvalues and eigenvectors of the original
symmetric matrix. On entry, Z must contain the
orthogonal matrix used to reduce the original matrix
to tridiagonal form.
= 'I': Compute eigenvalues and eigenvectors of the
tridiagonal matrix. Z is initialized to the identity
matrix.
N (input) INTEGER
The order of the matrix. N >= 0.
D (input/output) DOUBLE PRECISION array, dimension (N)
On entry, the diagonal elements of the tridiagonal matrix.
On exit, if INFO = 0, the eigenvalues in ascending order.
RS (input/output) DOUBLE PRECISION array, dimension (N-1)
On entry, the (n-1) subdiagonal elements of the tridiagonal
matrix.
On exit, RS has been destroyed.
Z (input/output) DOUBLE PRECISION array, dimension (LDZ, N)
On entry, if COMPZ = 'V', then Z contains the orthogonal
matrix used in the reduction to tridiagonal form.
On exit, if INFO = 0, then if COMPZ = 'V', Z contains the
orthonormal eigenvectors of the original symmetric matrix,
and if COMPZ = 'I', Z contains the orthonormal eigenvectors
of the symmetric tridiagonal matrix.
If COMPZ = 'N', then Z is not referenced.
LDZ (input) INTEGER
The leading dimension of the array Z. LDZ >= 1, and if
eigenvectors are desired, then LDZ >= max(1,N).
WORK (workspace) DOUBLE PRECISION array, dimension (max(1,2*N-2))
If COMPZ = 'N', then WORK is not referenced.
INFO (output) INTEGER
= 0: successful exit
< 0: if INFO = -i, the i-th argument had an illegal value
> 0: the algorithm has failed to find all the eigenvalues in
a total of 30*N iterations; if INFO = i, then i
elements of RS have not converged to zero; on exit, D
and RS contain the elements of a symmetric tridiagonal
matrix which is orthogonally similar to the original
matrix.
=====================================================================
*/
/**
* @param compz compz
* @param n n
* @param di di
* @param ei ei
* @param zi zi
* @param ldz ldz
* @param worki worki
* @return result
*/
int execute(String compz, int n, RS[] di, RS[] ei, RS[] zi, int ldz, RS[] worki) {
final RS unit = di[0].createUnit();
this.d = di;
this.e = ei;
this.z = zi;
this.work = worki;
//E rt1 = unit.create(0);
//E rt2 = unit.create(0);
//E s = unit.create(0);
//E c = unit.create(0);
int pointer_d = -1;
int pointer_e = -1;
int z_dim1 = ldz;
int z_offset = 1 + z_dim1 * 1;
int pointer_z = -z_offset;
int pointer_work = -1;
int info = 0;
int icompz;
if (BLAS.lsame(compz, "N")) { //$NON-NLS-1$
icompz = 0;
} else if (BLAS.lsame(compz, "V")) { //$NON-NLS-1$
icompz = 1;
} else if (BLAS.lsame(compz, "I")) { //$NON-NLS-1$
icompz = 2;
} else {
icompz = -1;
}
if (icompz < 0) {
info = -1;
} else if (n < 0) {
info = -2;
} else if (ldz < 1 || icompz > 0 && ldz < Math.max(1, n)) {
info = -6;
}
if (info != 0) {
BLAS.xerbla("DSTEQR", -info); //$NON-NLS-1$
return info;
}
/* Quick return if possible */
if (n == 0) {
return info;
}
if (n == 1) {
if (icompz == 2) {
this.z[z_dim1 + 1 + pointer_z] = unit.create(1);
}
return info;
}
/* Determine the unit roundoff and over/underflow thresholds. */
RS eps = Clapack.dlamch("E",unit); //$NON-NLS-1$
RS eps2 = eps.multiply(eps);
RS safmin = Clapack.dlamch("S",unit); //$NON-NLS-1$
RS safmax = unit.create(1).divide(safmin);
RS ssfmax = (safmax.sqrt()).divide(3);
RS ssfmin = safmin.sqrt().divide(eps2);
/*
* Compute the eigenvalues and eigenvectors of the tridiagonal matrix.
*/
if (icompz == 2) {
RS c_b9 = unit.create(0);
RS c_b10 = unit.create(1);
dlaset("Full", n, n, c_b9, c_b10, pointer_z + z_offset, ldz); //$NON-NLS-1$
}
int nmaxit = n * 30;
int jtot = 0;
/*
* Determine where the matrix splits and choose QL or QR iteration for each
* block, according to whether top or bottom diagonal element is smaller.
*/
int l;
int lend;
int iscale;
int lendsv;
int lsv;
int l1 = 1;
int nm1 = n - 1;
RS anorm;
do {
do {
do {
if (l1 > n) {
info = gotoL160(pointer_d, pointer_z, n, icompz, z_dim1);
return info;
}
if (l1 > 1) {
this.e[pointer_e + l1 - 1] = unit.create(0);
}
int m = 0;
boolean isGotoL30 = false;
if (l1 <= nm1) {
for (m = l1; m <= nm1; ++m) {
RS d1 = this.e[pointer_e + m];
RS tst = d1.abs();
if (tst.isZero()) {
isGotoL30 = true;
break;
}
d1 = this.d[pointer_d + m];
RS d2 = this.d[pointer_d + m + 1];
if (tst.isLessThanOrEquals((d1.abs().sqrt()).multiply((d2.abs().sqrt()).multiply(eps)))) {
this.e[pointer_e + m] = unit.create(0);
isGotoL30 = true;
break;
}
/* L20: */
}
}
if (!isGotoL30) {
m = n;
}
// L30:
l = l1;
lsv = l;
lend = m;
lendsv = lend;
l1 = m + 1;
} while (lend == l);
/* Scale submatrix in rows and columns L to LEND */
anorm = dlanst("I", lend - l + 1, pointer_d + l, pointer_e + l); //$NON-NLS-1$
iscale = 0;
} while (anorm.isZero());
if (anorm.isGreaterThan(ssfmax)) {
iscale = 1;
info = dlascl_d("G", 0, 0, anorm, ssfmax, lend - l + 1, 1, pointer_d + l, n); //$NON-NLS-1$
info = dlascl_e("G", 0, 0, anorm, ssfmax, lend - l, 1, pointer_e + l, n); //$NON-NLS-1$
} else if (anorm.isLessThan(ssfmin)) {
iscale = 2;
info = dlascl_d("G", 0, 0, anorm, ssfmin, lend - l + 1, 1, pointer_d + l, n); //$NON-NLS-1$
info = dlascl_e("G", 0, 0, anorm, ssfmin, lend - l, 1, pointer_e + l, n); //$NON-NLS-1$
}
/* Choose between QL and QR iteration */
RS d1 = this.d[pointer_d + lend];
RS d2 = this.d[pointer_d + l];
if ((d1.abs()).isLessThan(d2.abs())) {
lend = lsv;
l = lendsv;
}
if (lend > l) {
/*
* QL Iteration
* Look for small subdiagonal element.
*/
int m = 0;
boolean isGotoL40 = false;
do {
// L40:
isGotoL40 = false;
boolean isGotoL60 = false;
if (l != lend) {
int lendm1 = lend - 1;
for (m = l; m <= lendm1; ++m) {
d1 = this.e[pointer_e + m];
d2 = d1.abs();
RS tst = d2.multiply(d2);
d1 = this.d[pointer_d + m];
d2 = this.d[pointer_d + m + 1];
if (tst.isLessThanOrEquals((eps2.multiply(d1.abs()).multiply(d2.abs())).add(safmin))) {
isGotoL60 = true;
break;
}
/* L50: */
}
}
if (!isGotoL60) {
m = lend;
}
// L60:
if (m < lend) {
this.e[pointer_e + m] = unit.create(0);
}
RS p = this.d[pointer_d + l];
if (m == l) {
// Eigenvalue found.
// L80:
this.d[pointer_d + l] = p;
++l;
if (l <= lend) {
isGotoL40 = true;
continue;
}
info = gotoL140(pointer_d, pointer_e, n, iscale, lendsv, lsv, anorm, ssfmin, ssfmax);
break;
}
/*
* If remaining matrix is 2-by-2, use DLAE2 or SLAEV2 to compute its eigensystem.
*/
if (m == l + 1) {
RS rt1;
RS rt2;
if (icompz > 0) {
RS[] ans = Clapack.dlaev2(this.d[pointer_d + l], this.e[pointer_e + l], this.d[pointer_d + l + 1]);
rt1 = ans[0];
rt2 = ans[1];
RS c = ans[2];
RS s = ans[3];
this.work[pointer_work + l] = c;
this.work[pointer_work + n - 1 + l] = s;
dlasr("R", "V", "B", n, 2, pointer_work + l, pointer_work + n - 1 + l, (l) * z_dim1 + 1 + pointer_z, ldz); //$NON-NLS-1$ //$NON-NLS-2$ //$NON-NLS-3$
} else {
RS[] temp = Clapack.dlae2(this.d[pointer_d + l], this.e[pointer_e + l], this.d[pointer_d + l + 1]);
rt1 = temp[0];
rt2 = temp[1];
}
this.d[pointer_d + l] = rt1;
this.d[pointer_d + l + 1] = rt2;
this.e[pointer_e + l] = unit.create(0);
l += 2;
if (l <= lend) {
isGotoL40 = true;
continue;
}
info = gotoL140(pointer_d, pointer_e, n, iscale, lendsv, lsv, anorm, ssfmin, ssfmax);
break;
}
if (jtot == nmaxit) {
info = gotoL140(pointer_d, pointer_e, n, iscale, lendsv, lsv, anorm, ssfmin, ssfmax);
break;
}
++jtot;
/* Form shift. */
RS g = (this.d[pointer_d + l + 1].subtract(p)).divide(this.e[pointer_e + l].multiply(2));
RS c_b10 = unit.create(1);
RS r = Clapack.dlapy2(g, c_b10);
g = this.d[pointer_d + m].subtract(p).add(this.e[pointer_e + l].divide(g.add(LibF77.d_sign(r, g))));
RS s = unit.create(1);
RS c = unit.create(1);
p = unit.create(0);
for (int i = m - 1; i >= l; --i) {
RS f = s.multiply(this.e[pointer_e + i]);
RS b = c.multiply(this.e[pointer_e + i]);
RS[] ans = Clapack.dlartg(g, f);
c = ans[0];
s = ans[1];
r = ans[2];
if (i != m - 1) {
this.e[pointer_e + i + 1] = r;
}
g = this.d[pointer_d + i + 1].subtract(p);
r = ((this.d[pointer_d + i].subtract(g)).multiply(s)).add(c.multiply(2).multiply(b));
p = s.multiply(r);
this.d[pointer_d + i + 1] = g.add(p);
g = (c.multiply(r)).subtract(b);
/* If eigenvectors are desired, then save rotations. */
if (icompz > 0) {
this.work[pointer_work + i] = c;
this.work[pointer_work + n - 1 + i] = s.unaryMinus();
}
/* L70: */
}
/* If eigenvectors are desired, then apply saved rotations. */
if (icompz > 0) {
int mm = m - l + 1;
dlasr("R", "V", "B", n, mm, pointer_work + l, pointer_work + n - 1 + l, (l) * z_dim1 + 1 + pointer_z, ldz); //$NON-NLS-1$ //$NON-NLS-2$ //$NON-NLS-3$
}
this.d[pointer_d + l] = this.d[pointer_d + l].subtract(p);
this.e[pointer_e + l] = g;
isGotoL40 = true;
continue;
/* Eigenvalue found. */
// L80:
// d__[pointer_d__+l] = p;
//
// ++l;
// if (l <= lend) {
// goto L40;
// isGotoL40 = true;
// continue;
// }
// goto L140;
// gotoL140(pointer_d__,pointer_e,n,info);
// break;
} while (isGotoL40);
} else {
/*
* QR Iteration
* Look for small superdiagonal element.
*/
int m = 0;
boolean isGotoL90 = false;
do {
isGotoL90 = false;
// L90:
boolean isGotoL110 = false;
if (l != lend) {
int lendp1 = lend + 1;
for (m = l; m >= lendp1; --m) {
d1 = this.e[pointer_e + m - 1];
d2 = d1.abs();
RS tst = d2.multiply(d2);
d1 = this.d[pointer_d + m];
d2 = this.d[pointer_d + m - 1];
if (tst.isLessThanOrEquals((eps2.multiply(d1.abs()).multiply(d2.abs())).add(safmin))) {
isGotoL110 = true;
break;
}
/* L100: */
}
}
if (!isGotoL110) {
m = lend;
}
//L110:
if (m > lend) {
this.e[pointer_e + m - 1] = unit.create(0);
}
RS p = this.d[pointer_d + l];
if (m == l) {
// L130:
this.d[pointer_d + l] = p;
--l;
if (l >= lend) {
isGotoL90 = true;
continue;
}
info = gotoL140(pointer_d, pointer_e, n, iscale, lendsv, lsv, anorm, ssfmin, ssfmax);
break;//TODO 間違ってるかも
}
/*
* If remaining matrix is 2-by-2, use DLAE2 or SLAEV2 to compute its eigensystem.
*/
if (m == l - 1) {
RS rt1;
RS rt2;
if (icompz > 0) {
RS[] ans = Clapack.dlaev2(this.d[pointer_d + l - 1], this.e[pointer_e + l - 1], this.d[pointer_d + l]);
rt1 = ans[0];
rt2 = ans[1];
RS c = ans[2];
RS s = ans[3];
this.work[pointer_work + m] = c;
this.work[pointer_work + n - 1 + m] = s;
dlasr("R", "V", "F", n, 2, pointer_work + m, pointer_work + n - 1 + m, (l - 1) * z_dim1 + 1 + pointer_z, ldz); //$NON-NLS-1$ //$NON-NLS-2$ //$NON-NLS-3$
} else {
RS[] temp = Clapack.dlae2(this.d[pointer_d + l - 1], this.e[pointer_e + l - 1], this.d[pointer_d + l]);
rt1 = temp[0];
rt2 = temp[1];
}
this.d[pointer_d + l - 1] = rt1;
this.d[pointer_d + l] = rt2;
this.e[pointer_e + l - 1] = unit.create(0);
l += -2;
if (l >= lend) {
isGotoL90 = true;
continue;
}
info = gotoL140(pointer_d, pointer_e, n, iscale, lendsv, lsv, anorm, ssfmin, ssfmax);
break;
}
if (jtot == nmaxit) {
info = gotoL140(pointer_d, pointer_e, n, iscale, lendsv, lsv, anorm, ssfmin, ssfmax);
break;
}
++jtot;
/* Form shift. */
RS g = (this.d[pointer_d + l - 1].subtract(p)).divide((this.e[pointer_e + l - 1].multiply(2)));
RS c_b10 = unit.create(1);
RS r = Clapack.dlapy2(g, c_b10);
g = this.d[pointer_d + m].subtract(p).add(this.e[pointer_e + l - 1].divide(g.add(LibF77.d_sign(r, g))));
RS s = unit.create(1);
RS c = unit.create(1);
p = unit.create(0);
/* Inner loop */
int lm1 = l - 1;
for (int i = m; i <= lm1; ++i) {
RS f = s.multiply(this.e[pointer_e + i]);
RS b = c.multiply(this.e[pointer_e + i]);
RS[] ans = Clapack.dlartg(g, f);
c = ans[0];
s = ans[1];
r = ans[2];
if (i != m) {
this.e[pointer_e + i - 1] = r;
}
g = this.d[pointer_d + i].subtract(p);
r = ((this.d[pointer_d + i + 1].subtract(g)).multiply(s)).add(c.multiply(2).multiply(b));
p = s.multiply(r);
this.d[pointer_d + i] = g.add(p);
g = (c.multiply(r)).subtract(b);
/* If eigenvectors are desired, then save rotations. */
if (icompz > 0) {
this.work[pointer_work + i] = c;
this.work[pointer_work + n - 1 + i] = s;
}
/* L120: */
}
/* If eigenvectors are desired, then apply saved rotations. */
if (icompz > 0) {
int mm = l - m + 1;
dlasr("R", "V", "F", n, mm, pointer_work + m, pointer_work + n - 1 + m, (m) * z_dim1 + 1 + pointer_z, ldz); //$NON-NLS-1$ //$NON-NLS-2$ //$NON-NLS-3$
}
this.d[pointer_d + l] = this.d[pointer_d + l].subtract(p);
this.e[pointer_e + lm1] = g;
isGotoL90 = true;
continue;
/* Eigenvalue found. */
// L130:
// d__[pointer_d__+l] = p;
//
// --l;
// if (l >= lend) {
// goto L90;
// isGotoL90 = true;
// }
// goto L140;
// gotoL140(pointer_d__,pointer_e,n,info);
// break;
} while (isGotoL90);
}
/* Undo scaling if necessary */
/*
* Check for no convergence to an eigenvalue after a total of NMAXIT
* iterations.
*/
} while (jtot < nmaxit);
for (int i = 1; i <= n - 1; ++i) {
if (!this.e[pointer_e + i].isZero()) {
info = info + 1;
}
/* L150: */
}
return info;
// /* Order eigenvalues and eigenvectors. */
//
// L160:
// if (icompz == 0) {
//
// /* Use Quick Sort */
//
// dlasrt("I", n, pointer_d__+1, info);
//
// } else {
//
// /* Use Selection Sort to minimize swaps of eigenvectors */
//
// i__1 = n;
// for (ii = 2; ii <= i__1; ++ii) {
// i__ = ii - 1;
// k = i__;
// p = d__[pointer_d__+i__];
// i__2 = n;
// for (j = ii; j <= i__2; ++j) {
// if (d__[pointer_d__+j].isLessThan( p)) {
// k = j;
// p = d__[pointer_d__+j];
// }
// /* L170: */
// }
// if (k != i__) {
// d__[pointer_d__+k] = d__[pointer_d__+i__];
// d__[pointer_d__+i__] = p;
// dswap(n, (i__)*z_dim1 + 1+pointer_z__, c__1, (k)*z_dim1 + 1+pointer_z__,
// c__1);
// }
// /* L180: */
// }
// }
//
// L190:
// return 0;
//
/* End of DSTEQR */
}
/**
* @param pointer_d pointer_d
* @param pointer_e pointer_e
* @param n n
* @param iscale iscale
* @param lendsv lendsv
* @param lsv lsv
* @param anorm anorm
* @param ssfmin ssfmin
* @param ssfmax ssfmax
* @return result
*/
private int gotoL140(int pointer_d, int pointer_e, int n, int iscale, int lendsv, int lsv, RS anorm, RS ssfmin, RS ssfmax) {
/* Undo scaling if necessary */
//L140:
int info = 0;
if (iscale == 1) {
info = dlascl_d("G", 0, 0, ssfmax, anorm, lendsv - lsv + 1, 1, pointer_d + lsv, n); //$NON-NLS-1$
info = dlascl_e("G", 0, 0, ssfmax, anorm, lendsv - lsv, 1, pointer_e + lsv, n); //$NON-NLS-1$
} else if (iscale == 2) {
info = dlascl_d("G", 0, 0, ssfmin, anorm, lendsv - lsv + 1, 1, pointer_d + lsv, n); //$NON-NLS-1$
info = dlascl_e("G", 0, 0, ssfmin, anorm, lendsv - lsv, 1, pointer_e + lsv, n); //$NON-NLS-1$
}
return info;
}
/**
* @param pointer_d pointer_d
* @param pointer_z pointer_z
* @param n n
* @param icompz icompz
* @param z_dim1 z_dim1
* @return result
*/
private int gotoL160(int pointer_d, int pointer_z, int n, int icompz, int z_dim1) {
/* Order eigenvalues and eigenvectors. */
int info = 0;
if (icompz == 0) {
/* Use Quick Sort */
info = dlasrt("I", n, pointer_d + 1); //$NON-NLS-1$
} else {
/* Use Selection Sort to minimize swaps of eigenvectors */
for (int ii = 2; ii <= n; ++ii) {
int i = ii - 1;
int k = i;
RS p = this.d[pointer_d + i];
for (int j = ii; j <= n; ++j) {
if (this.d[pointer_d + j].isLessThan(p)) {
k = j;
p = this.d[pointer_d + j];
}
}
if (k != i) {
this.d[pointer_d + k] = this.d[pointer_d + i];
this.d[pointer_d + i] = p;
dswap(n, (i) * z_dim1 + 1 + pointer_z, 1, (k) * z_dim1 + 1 + pointer_z, 1);
}
}
}
return info;
}
/**
* @param side side
* @param pivot pivot
* @param direct direct
* @param m m
* @param n n
* @param cIndex work cIndex work
* @param sIndex work sndex work
* @param aIndex z__ (output)
* @param lda lda
*/
private void dlasr(String side, String pivot, String direct, int m, int n, int cIndex, int sIndex, int aIndex, int lda) {
final RS unit = this.d[0].createUnit();
RS[] cTemp = unit.createArray(this.work.length - cIndex);
RS[] sTemp = unit.createArray(this.work.length - sIndex);
RS[] aTemp = unit.createArray(this.z.length - aIndex);
System.arraycopy(this.work, cIndex, cTemp, 0, cTemp.length);
System.arraycopy(this.work, sIndex, sTemp, 0, sTemp.length);
System.arraycopy(this.z, aIndex, aTemp, 0, aTemp.length);
Clapack.dlasr(side, pivot, direct, m, n, cTemp, sTemp, aTemp, lda);
System.arraycopy(aTemp, 0, this.z, aIndex, aTemp.length);
}
/**
* @param n n
* @param dxIndex dxIndex
* @param incx incx
* @param dyIndex dyIndex
* @param incy incy
*/
private void dswap(int n, int dxIndex, int incx, int dyIndex, int incy) {
final RS unit = this.d[0].createUnit();
RS[] dx = unit.createArray(this.z.length - dxIndex);
RS[] dy = unit.createArray(this.z.length - dyIndex);
System.arraycopy(this.z, dxIndex, dx, 0, dx.length);
System.arraycopy(this.z, dyIndex, dy, 0, dy.length);
BLAS.dswap(n, dx, incx, dy, incy);
System.arraycopy(dx, 0, this.z, dxIndex, dx.length);
System.arraycopy(dy, 0, this.z, dyIndex, dy.length);
}
/**
* @param id id
* @param n n
* @param dIndex dIndex
* @return result
*/
private int dlasrt(String id, int n, int dIndex) {
final RS unit = this.d[0].createUnit();
RS[] dTemp = unit.createArray(this.d.length - dIndex);
System.arraycopy(this.d, dIndex, dTemp, 0, dTemp.length);
int info = Clapack.dlasrt(id, n, dTemp);
System.arraycopy(dTemp, 0, this.d, dIndex, dTemp.length);
return info;
}
/**
* @param norm input
* @param n input
* @param dIndex d input
* @param eIndex e input
* @return ans output
*/
private RS dlanst(String norm, int n, int dIndex, int eIndex) {
final RS unit = this.d[0].createUnit();
RS[] dTemp = unit.createArray(this.d.length - dIndex);
RS[] eTemp = unit.createArray(this.e.length - eIndex);
System.arraycopy(this.d, dIndex, dTemp, 0, dTemp.length);
System.arraycopy(this.e, eIndex, eTemp, 0, eTemp.length);
return Clapack.dlanst(norm, n, dTemp, eTemp);
}
/**
* @param uplo uplo
* @param m m
* @param n2 n2
* @param alpha alpha
* @param beta beta
* @param aIndex aIndex
* @param lda lda
*/
private void dlaset(String uplo, int m, int n2, RS alpha, RS beta, int aIndex, int lda) {
final RS unit = alpha.createUnit();
RS[] aTemp = unit.createArray(this.z.length - aIndex);
System.arraycopy(this.z, aIndex, aTemp, 0, aTemp.length);
Clapack.dlaset(uplo, m, n2, alpha, beta, aTemp, lda);
System.arraycopy(aTemp, 0, this.z, aIndex, aTemp.length);
}
/**
* @param type type
* @param kl kl
* @param ku ku
* @param cfrom cfrom
* @param cto cto
* @param m m
* @param n2 n2
* @param aIndex d
* @param lda lda
* @return result
*/
private int dlascl_d(String type, int kl, int ku, RS cfrom, RS cto, int m, int n2, int aIndex, int lda) {
final RS unit = cfrom.createUnit();
RS[] aTemp = unit.createArray(this.d.length - aIndex);
System.arraycopy(this.d, aIndex, aTemp, 0, aTemp.length);
int info = Clapack.dlascl(type, kl, ku, cfrom, cto, m, n2, aTemp, lda);
System.arraycopy(aTemp, 0, this.d, aIndex, aTemp.length);
return info;
}
/**
* @param type type
* @param kl kl
* @param ku ku
* @param cfrom cfrom
* @param cto cto
* @param m m
* @param n2 n2
* @param aIndex e
* @param lda lda
* @return result
*/
private int dlascl_e(String type, int kl, int ku, RS cfrom, RS cto, int m, int n2, int aIndex, int lda) {
final RS unit = cfrom.createUnit();
RS[] aTemp = unit.createArray(this.e.length - aIndex);
System.arraycopy(this.e, aIndex, aTemp, 0, aTemp.length);
int info = Clapack.dlascl(type, kl, ku, cfrom, cto, m, n2, aTemp, lda);
System.arraycopy(aTemp, 0, this.e, aIndex, aTemp.length);
return info;
}
}