`   1:  #region Translated by Jose Antonio De Santiago-Castillo.`
`   2:   `
`   3:  //Translated by Jose Antonio De Santiago-Castillo. `
`   4:  //E-mail:JAntonioDeSantiago@gmail.com`
`   5:  //Web: www.DotNumerics.com`
`   6:  //`
`   7:  //Fortran to C# Translation.`
`   8:  //Translated by:`
`   9:  //F2CSharp Version 0.71 (November 10, 2009)`
`  10:  //Code Optimizations: None`
`  11:  //`
`  12:  #endregion`
`  13:   `
`  14:  using System;`
`  15:  using DotNumerics.FortranLibrary;`
`  16:   `
`  17:  namespace DotNumerics.CSLapack`
`  18:  {`
`  19:      /// <summary>`
`  20:      /// -- LAPACK driver routine (version 3.1) --`
`  21:      /// Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..`
`  22:      /// November 2006`
`  23:      /// Purpose`
`  24:      /// =======`
`  25:      /// `
`  26:      /// DSYEVD computes all eigenvalues and, optionally, eigenvectors of a`
`  27:      /// real symmetric matrix A. If eigenvectors are desired, it uses a`
`  28:      /// divide and conquer algorithm.`
`  29:      /// `
`  30:      /// The divide and conquer algorithm makes very mild assumptions about`
`  31:      /// floating point arithmetic. It will work on machines with a guard`
`  32:      /// digit in add/subtract, or on those binary machines without guard`
`  33:      /// digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or`
`  34:      /// Cray-2. It could conceivably fail on hexadecimal or decimal machines`
`  35:      /// without guard digits, but we know of none.`
`  36:      /// `
`  37:      /// Because of large use of BLAS of level 3, DSYEVD needs N**2 more`
`  38:      /// workspace than DSYEVX.`
`  39:      /// `
`  40:      ///</summary>`
`  41:      public class DSYEVD`
`  42:      {`
`  43:      `
`  44:   `
`  45:          #region Dependencies`
`  46:          `
`  47:          LSAME _lsame; DLAMCH _dlamch; DLANSY _dlansy; ILAENV _ilaenv; DLACPY _dlacpy; DLASCL _dlascl; DORMTR _dormtr; `
`  48:          DSCAL _dscal;DSTEDC _dstedc; DSTERF _dsterf; DSYTRD _dsytrd; XERBLA _xerbla; `
`  49:   `
`  50:          #endregion`
`  51:   `
`  52:   `
`  53:          #region Fields`
`  54:          `
`  55:          const double ZERO = 0.0E+0; const double ONE = 1.0E+0; bool LOWER = false; bool LQUERY = false; bool WANTZ = false; `
`  56:          int IINFO = 0;int INDE = 0; int INDTAU = 0; int INDWK2 = 0; int INDWRK = 0; int ISCALE = 0; int LIOPT = 0; int LIWMIN = 0; `
`  57:          int LLWORK = 0;int LLWRK2 = 0; int LOPT = 0; int LWMIN = 0; double ANRM = 0; double BIGNUM = 0; double EPS = 0; `
`  58:          double RMAX = 0;double RMIN = 0; double SAFMIN = 0; double SIGMA = 0; double SMLNUM = 0; `
`  59:   `
`  60:          #endregion`
`  61:   `
`  62:          public DSYEVD(LSAME lsame, DLAMCH dlamch, DLANSY dlansy, ILAENV ilaenv, DLACPY dlacpy, DLASCL dlascl, DORMTR dormtr, DSCAL dscal, DSTEDC dstedc, DSTERF dsterf`
`  63:                        , DSYTRD dsytrd, XERBLA xerbla)`
`  64:          {`
`  65:      `
`  66:   `
`  67:              #region Set Dependencies`
`  68:              `
`  69:              this._lsame = lsame; this._dlamch = dlamch; this._dlansy = dlansy; this._ilaenv = ilaenv; this._dlacpy = dlacpy; `
`  70:              this._dlascl = dlascl;this._dormtr = dormtr; this._dscal = dscal; this._dstedc = dstedc; this._dsterf = dsterf; `
`  71:              this._dsytrd = dsytrd;this._xerbla = xerbla; `
`  72:   `
`  73:              #endregion`
`  74:   `
`  75:          }`
`  76:      `
`  77:          public DSYEVD()`
`  78:          {`
`  79:      `
`  80:   `
`  81:              #region Dependencies (Initialization)`
`  82:              `
`  83:              LSAME lsame = new LSAME();`
`  84:              DLAMC3 dlamc3 = new DLAMC3();`
`  85:              DLASSQ dlassq = new DLASSQ();`
`  86:              IEEECK ieeeck = new IEEECK();`
`  87:              IPARMQ iparmq = new IPARMQ();`
`  88:              XERBLA xerbla = new XERBLA();`
`  89:              DCOPY dcopy = new DCOPY();`
`  90:              DSCAL dscal = new DSCAL();`
`  91:              IDAMAX idamax = new IDAMAX();`
`  92:              DLAPY2 dlapy2 = new DLAPY2();`
`  93:              DLAMRG dlamrg = new DLAMRG();`
`  94:              DROT drot = new DROT();`
`  95:              DNRM2 dnrm2 = new DNRM2();`
`  96:              DLAED5 dlaed5 = new DLAED5();`
`  97:              DLAE2 dlae2 = new DLAE2();`
`  98:              DLAEV2 dlaev2 = new DLAEV2();`
`  99:              DSWAP dswap = new DSWAP();`
` 100:              DAXPY daxpy = new DAXPY();`
` 101:              DDOT ddot = new DDOT();`
` 102:              DLAMC1 dlamc1 = new DLAMC1(dlamc3);`
` 103:              DLAMC4 dlamc4 = new DLAMC4(dlamc3);`
` 104:              DLAMC5 dlamc5 = new DLAMC5(dlamc3);`
` 105:              DLAMC2 dlamc2 = new DLAMC2(dlamc3, dlamc1, dlamc4, dlamc5);`
` 106:              DLAMCH dlamch = new DLAMCH(lsame, dlamc2);`
` 107:              DLANSY dlansy = new DLANSY(dlassq, lsame);`
` 108:              ILAENV ilaenv = new ILAENV(ieeeck, iparmq);`
` 109:              DLACPY dlacpy = new DLACPY(lsame);`
` 110:              DLASCL dlascl = new DLASCL(lsame, dlamch, xerbla);`
` 111:              DGEMM dgemm = new DGEMM(lsame, xerbla);`
` 112:              DTRMM dtrmm = new DTRMM(lsame, xerbla);`
` 113:              DLARFB dlarfb = new DLARFB(lsame, dcopy, dgemm, dtrmm);`
` 114:              DGEMV dgemv = new DGEMV(lsame, xerbla);`
` 115:              DTRMV dtrmv = new DTRMV(lsame, xerbla);`
` 116:              DLARFT dlarft = new DLARFT(dgemv, dtrmv, lsame);`
` 117:              DGER dger = new DGER(xerbla);`
` 118:              DLARF dlarf = new DLARF(dgemv, dger, lsame);`
` 119:              DORM2L dorm2l = new DORM2L(lsame, dlarf, xerbla);`
` 120:              DORMQL dormql = new DORMQL(lsame, ilaenv, dlarfb, dlarft, dorm2l, xerbla);`
` 121:              DORM2R dorm2r = new DORM2R(lsame, dlarf, xerbla);`
` 122:              DORMQR dormqr = new DORMQR(lsame, ilaenv, dlarfb, dlarft, dorm2r, xerbla);`
` 123:              DORMTR dormtr = new DORMTR(lsame, ilaenv, dormql, dormqr, xerbla);`
` 124:              DLANST dlanst = new DLANST(lsame, dlassq);`
` 125:              DLAED2 dlaed2 = new DLAED2(idamax, dlamch, dlapy2, dcopy, dlacpy, dlamrg, drot, dscal, xerbla);`
` 126:              DLAED6 dlaed6 = new DLAED6(dlamch);`
` 127:              DLAED4 dlaed4 = new DLAED4(dlamch, dlaed5, dlaed6);`
` 128:              DLASET dlaset = new DLASET(lsame);`
` 129:              DLAED3 dlaed3 = new DLAED3(dlamc3, dnrm2, dcopy, dgemm, dlacpy, dlaed4, dlaset, xerbla);`
` 130:              DLAED1 dlaed1 = new DLAED1(dcopy, dlaed2, dlaed3, dlamrg, xerbla);`
` 131:              DLAED8 dlaed8 = new DLAED8(idamax, dlamch, dlapy2, dcopy, dlacpy, dlamrg, drot, dscal, xerbla);`
` 132:              DLAED9 dlaed9 = new DLAED9(dlamc3, dnrm2, dcopy, dlaed4, xerbla);`
` 133:              DLAEDA dlaeda = new DLAEDA(dcopy, dgemv, drot, xerbla);`
` 134:              DLAED7 dlaed7 = new DLAED7(dgemm, dlaed8, dlaed9, dlaeda, dlamrg, xerbla);`
` 135:              DLARTG dlartg = new DLARTG(dlamch);`
` 136:              DLASR dlasr = new DLASR(lsame, xerbla);`
` 137:              DLASRT dlasrt = new DLASRT(lsame, xerbla);`
` 138:              DSTEQR dsteqr = new DSTEQR(lsame, dlamch, dlanst, dlapy2, dlae2, dlaev2, dlartg, dlascl, dlaset, dlasr`
` 139:                                         , dlasrt, dswap, xerbla);`
` 140:              DLAED0 dlaed0 = new DLAED0(dcopy, dgemm, dlacpy, dlaed1, dlaed7, dsteqr, xerbla, ilaenv);`
` 141:              DSTERF dsterf = new DSTERF(dlamch, dlanst, dlapy2, dlae2, dlascl, dlasrt, xerbla);`
` 142:              DSTEDC dstedc = new DSTEDC(lsame, ilaenv, dlamch, dlanst, dgemm, dlacpy, dlaed0, dlascl, dlaset, dlasrt`
` 143:                                         , dsteqr, dsterf, dswap, xerbla);`
` 144:              DLARFG dlarfg = new DLARFG(dlamch, dlapy2, dnrm2, dscal);`
` 145:              DSYMV dsymv = new DSYMV(lsame, xerbla);`
` 146:              DLATRD dlatrd = new DLATRD(daxpy, dgemv, dlarfg, dscal, dsymv, lsame, ddot);`
` 147:              DSYR2K dsyr2k = new DSYR2K(lsame, xerbla);`
` 148:              DSYR2 dsyr2 = new DSYR2(lsame, xerbla);`
` 149:              DSYTD2 dsytd2 = new DSYTD2(daxpy, dlarfg, dsymv, dsyr2, xerbla, lsame, ddot);`
` 150:              DSYTRD dsytrd = new DSYTRD(dlatrd, dsyr2k, dsytd2, xerbla, lsame, ilaenv);`
` 151:   `
` 152:              #endregion`
` 153:   `
` 154:   `
` 155:              #region Set Dependencies`
` 156:              `
` 157:              this._lsame = lsame; this._dlamch = dlamch; this._dlansy = dlansy; this._ilaenv = ilaenv; this._dlacpy = dlacpy; `
` 158:              this._dlascl = dlascl;this._dormtr = dormtr; this._dscal = dscal; this._dstedc = dstedc; this._dsterf = dsterf; `
` 159:              this._dsytrd = dsytrd;this._xerbla = xerbla; `
` 160:   `
` 161:              #endregion`
` 162:   `
` 163:          }`
` 164:          /// <summary>`
` 165:          /// Purpose`
` 166:          /// =======`
` 167:          /// `
` 168:          /// DSYEVD computes all eigenvalues and, optionally, eigenvectors of a`
` 169:          /// real symmetric matrix A. If eigenvectors are desired, it uses a`
` 170:          /// divide and conquer algorithm.`
` 171:          /// `
` 172:          /// The divide and conquer algorithm makes very mild assumptions about`
` 173:          /// floating point arithmetic. It will work on machines with a guard`
` 174:          /// digit in add/subtract, or on those binary machines without guard`
` 175:          /// digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or`
` 176:          /// Cray-2. It could conceivably fail on hexadecimal or decimal machines`
` 177:          /// without guard digits, but we know of none.`
` 178:          /// `
` 179:          /// Because of large use of BLAS of level 3, DSYEVD needs N**2 more`
` 180:          /// workspace than DSYEVX.`
` 181:          /// `
` 182:          ///</summary>`
` 183:          /// <param name="JOBZ">`
` 184:          /// (input) CHARACTER*1`
` 185:          /// = 'N':  Compute eigenvalues only;`
` 186:          /// = 'V':  Compute eigenvalues and eigenvectors.`
` 187:          ///</param>`
` 188:          /// <param name="UPLO">`
` 189:          /// (input) CHARACTER*1`
` 190:          /// = 'U':  Upper triangle of A is stored;`
` 191:          /// = 'L':  Lower triangle of A is stored.`
` 192:          ///</param>`
` 193:          /// <param name="N">`
` 194:          /// (input) INTEGER`
` 195:          /// The order of the matrix A.  N .GE. 0.`
` 196:          ///</param>`
` 197:          /// <param name="A">`
` 198:          /// (input/output) DOUBLE PRECISION array, dimension (LDA, N)`
` 199:          /// On entry, the symmetric matrix A.  If UPLO = 'U', the`
` 200:          /// leading N-by-N upper triangular part of A contains the`
` 201:          /// upper triangular part of the matrix A.  If UPLO = 'L',`
` 202:          /// the leading N-by-N lower triangular part of A contains`
` 203:          /// the lower triangular part of the matrix A.`
` 204:          /// On exit, if JOBZ = 'V', then if INFO = 0, A contains the`
` 205:          /// orthonormal eigenvectors of the matrix A.`
` 206:          /// If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')`
` 207:          /// or the upper triangle (if UPLO='U') of A, including the`
` 208:          /// diagonal, is destroyed.`
` 209:          ///</param>`
` 210:          /// <param name="LDA">`
` 211:          /// (input) INTEGER`
` 212:          /// The leading dimension of the array A.  LDA .GE. max(1,N).`
` 213:          ///</param>`
` 214:          /// <param name="W">`
` 215:          /// (output) DOUBLE PRECISION array, dimension (N)`
` 216:          /// If INFO = 0, the eigenvalues in ascending order.`
` 217:          ///</param>`
` 218:          /// <param name="WORK">`
` 219:          /// (workspace/output) DOUBLE PRECISION array,`
` 220:          /// dimension (LWORK)`
` 221:          /// On exit, if INFO = 0, WORK(1) returns the optimal LWORK.`
` 222:          ///</param>`
` 223:          /// <param name="LWORK">`
` 224:          /// (input) INTEGER`
` 225:          /// The dimension of the array WORK.`
` 226:          /// If N .LE. 1,               LWORK must be at least 1.`
` 227:          /// If JOBZ = 'N' and N .GT. 1, LWORK must be at least 2*N+1.`
` 228:          /// If JOBZ = 'V' and N .GT. 1, LWORK must be at least`
` 229:          /// 1 + 6*N + 2*N**2.`
` 230:          /// `
` 231:          /// If LWORK = -1, then a workspace query is assumed; the routine`
` 232:          /// only calculates the optimal sizes of the WORK and IWORK`
` 233:          /// arrays, returns these values as the first entries of the WORK`
` 234:          /// and IWORK arrays, and no error message related to LWORK or`
` 235:          /// LIWORK is issued by XERBLA.`
` 236:          ///</param>`
` 237:          /// <param name="IWORK">`
` 238:          /// (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))`
` 239:          /// On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.`
` 240:          ///</param>`
` 241:          /// <param name="LIWORK">`
` 242:          /// (input) INTEGER`
` 243:          /// The dimension of the array IWORK.`
` 244:          /// If N .LE. 1,                LIWORK must be at least 1.`
` 245:          /// If JOBZ  = 'N' and N .GT. 1, LIWORK must be at least 1.`
` 246:          /// If JOBZ  = 'V' and N .GT. 1, LIWORK must be at least 3 + 5*N.`
` 247:          /// `
` 248:          /// If LIWORK = -1, then a workspace query is assumed; the`
` 249:          /// routine only calculates the optimal sizes of the WORK and`
` 250:          /// IWORK arrays, returns these values as the first entries of`
` 251:          /// the WORK and IWORK arrays, and no error message related to`
` 252:          /// LWORK or LIWORK is issued by XERBLA.`
` 253:          ///</param>`
` 254:          /// <param name="INFO">`
` 255:          /// (output) INTEGER`
` 256:          /// = 0:  successful exit`
` 257:          /// .LT. 0:  if INFO = -i, the i-th argument had an illegal value`
` 258:          /// .GT. 0:  if INFO = i and JOBZ = 'N', then the algorithm failed`
` 259:          /// to converge; i off-diagonal elements of an intermediate`
` 260:          /// tridiagonal form did not converge to zero;`
` 261:          /// if INFO = i and JOBZ = 'V', then the algorithm failed`
` 262:          /// to compute an eigenvalue while working on the submatrix`
` 263:          /// lying in rows and columns INFO/(N+1) through`
` 264:          /// mod(INFO,N+1).`
` 265:          ///</param>`
` 266:          public void Run(string JOBZ, string UPLO, int N, ref double[] A, int offset_a, int LDA, ref double[] W, int offset_w`
` 267:                           , ref double[] WORK, int offset_work, int LWORK, ref int[] IWORK, int offset_iwork, int LIWORK, ref int INFO)`
` 268:          {`
` 269:   `
` 270:              #region Array Index Correction`
` 271:              `
` 272:               int o_a = -1 - LDA + offset_a;  int o_w = -1 + offset_w;  int o_work = -1 + offset_work; `
` 273:               int o_iwork = -1 + offset_iwork;`
` 274:   `
` 275:              #endregion`
` 276:   `
` 277:   `
` 278:              #region Strings`
` 279:              `
` 280:              JOBZ = JOBZ.Substring(0, 1);  UPLO = UPLO.Substring(0, 1);  `
` 281:   `
` 282:              #endregion`
` 283:   `
` 284:   `
` 285:              #region Prolog`
` 286:              `
` 287:              // *`
` 288:              // *  -- LAPACK driver routine (version 3.1) --`
` 289:              // *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..`
` 290:              // *     November 2006`
` 291:              // *`
` 292:              // *     .. Scalar Arguments ..`
` 293:              // *     ..`
` 294:              // *     .. Array Arguments ..`
` 295:              // *     ..`
` 296:              // *`
` 297:              // *  Purpose`
` 298:              // *  =======`
` 299:              // *`
` 300:              // *  DSYEVD computes all eigenvalues and, optionally, eigenvectors of a`
` 301:              // *  real symmetric matrix A. If eigenvectors are desired, it uses a`
` 302:              // *  divide and conquer algorithm.`
` 303:              // *`
` 304:              // *  The divide and conquer algorithm makes very mild assumptions about`
` 305:              // *  floating point arithmetic. It will work on machines with a guard`
` 306:              // *  digit in add/subtract, or on those binary machines without guard`
` 307:              // *  digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or`
` 308:              // *  Cray-2. It could conceivably fail on hexadecimal or decimal machines`
` 309:              // *  without guard digits, but we know of none.`
` 310:              // *`
` 311:              // *  Because of large use of BLAS of level 3, DSYEVD needs N**2 more`
` 312:              // *  workspace than DSYEVX.`
` 313:              // *`
` 314:              // *  Arguments`
` 315:              // *  =========`
` 316:              // *`
` 317:              // *  JOBZ    (input) CHARACTER*1`
` 318:              // *          = 'N':  Compute eigenvalues only;`
` 319:              // *          = 'V':  Compute eigenvalues and eigenvectors.`
` 320:              // *`
` 321:              // *  UPLO    (input) CHARACTER*1`
` 322:              // *          = 'U':  Upper triangle of A is stored;`
` 323:              // *          = 'L':  Lower triangle of A is stored.`
` 324:              // *`
` 325:              // *  N       (input) INTEGER`
` 326:              // *          The order of the matrix A.  N >= 0.`
` 327:              // *`
` 328:              // *  A       (input/output) DOUBLE PRECISION array, dimension (LDA, N)`
` 329:              // *          On entry, the symmetric matrix A.  If UPLO = 'U', the`
` 330:              // *          leading N-by-N upper triangular part of A contains the`
` 331:              // *          upper triangular part of the matrix A.  If UPLO = 'L',`
` 332:              // *          the leading N-by-N lower triangular part of A contains`
` 333:              // *          the lower triangular part of the matrix A.`
` 334:              // *          On exit, if JOBZ = 'V', then if INFO = 0, A contains the`
` 335:              // *          orthonormal eigenvectors of the matrix A.`
` 336:              // *          If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')`
` 337:              // *          or the upper triangle (if UPLO='U') of A, including the`
` 338:              // *          diagonal, is destroyed.`
` 339:              // *`
` 340:              // *  LDA     (input) INTEGER`
` 341:              // *          The leading dimension of the array A.  LDA >= max(1,N).`
` 342:              // *`
` 343:              // *  W       (output) DOUBLE PRECISION array, dimension (N)`
` 344:              // *          If INFO = 0, the eigenvalues in ascending order.`
` 345:              // *`
` 346:              // *  WORK    (workspace/output) DOUBLE PRECISION array,`
` 347:              // *                                         dimension (LWORK)`
` 348:              // *          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.`
` 349:              // *`
` 350:              // *  LWORK   (input) INTEGER`
` 351:              // *          The dimension of the array WORK.`
` 352:              // *          If N <= 1,               LWORK must be at least 1.`
` 353:              // *          If JOBZ = 'N' and N > 1, LWORK must be at least 2*N+1.`
` 354:              // *          If JOBZ = 'V' and N > 1, LWORK must be at least`
` 355:              // *                                                1 + 6*N + 2*N**2.`
` 356:              // *`
` 357:              // *          If LWORK = -1, then a workspace query is assumed; the routine`
` 358:              // *          only calculates the optimal sizes of the WORK and IWORK`
` 359:              // *          arrays, returns these values as the first entries of the WORK`
` 360:              // *          and IWORK arrays, and no error message related to LWORK or`
` 361:              // *          LIWORK is issued by XERBLA.`
` 362:              // *`
` 363:              // *  IWORK   (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))`
` 364:              // *          On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.`
` 365:              // *`
` 366:              // *  LIWORK  (input) INTEGER`
` 367:              // *          The dimension of the array IWORK.`
` 368:              // *          If N <= 1,                LIWORK must be at least 1.`
` 369:              // *          If JOBZ  = 'N' and N > 1, LIWORK must be at least 1.`
` 370:              // *          If JOBZ  = 'V' and N > 1, LIWORK must be at least 3 + 5*N.`
` 371:              // *`
` 372:              // *          If LIWORK = -1, then a workspace query is assumed; the`
` 373:              // *          routine only calculates the optimal sizes of the WORK and`
` 374:              // *          IWORK arrays, returns these values as the first entries of`
` 375:              // *          the WORK and IWORK arrays, and no error message related to`
` 376:              // *          LWORK or LIWORK is issued by XERBLA.`
` 377:              // *`
` 378:              // *  INFO    (output) INTEGER`
` 379:              // *          = 0:  successful exit`
` 380:              // *          < 0:  if INFO = -i, the i-th argument had an illegal value`
` 381:              // *          > 0:  if INFO = i and JOBZ = 'N', then the algorithm failed`
` 382:              // *                to converge; i off-diagonal elements of an intermediate`
` 383:              // *                tridiagonal form did not converge to zero;`
` 384:              // *                if INFO = i and JOBZ = 'V', then the algorithm failed`
` 385:              // *                to compute an eigenvalue while working on the submatrix`
` 386:              // *                lying in rows and columns INFO/(N+1) through`
` 387:              // *                mod(INFO,N+1).`
` 388:              // *`
` 389:              // *  Further Details`
` 390:              // *  ===============`
` 391:              // *`
` 392:              // *  Based on contributions by`
` 393:              // *     Jeff Rutter, Computer Science Division, University of California`
` 394:              // *     at Berkeley, USA`
` 395:              // *  Modified by Francoise Tisseur, University of Tennessee.`
` 396:              // *`
` 397:              // *  Modified description of INFO. Sven, 16 Feb 05.`
` 398:              // *  =====================================================================`
` 399:              // *`
` 400:              // *     .. Parameters ..`
` 401:              // *     ..`
` 402:              // *     .. Local Scalars ..`
` 403:              // *`
` 404:              // *     ..`
` 405:              // *     .. External Functions ..`
` 406:              // *     ..`
` 407:              // *     .. External Subroutines ..`
` 408:              // *     ..`
` 409:              // *     .. Intrinsic Functions ..`
` 410:              //      INTRINSIC          MAX, SQRT;`
` 411:              // *     ..`
` 412:              // *     .. Executable Statements ..`
` 413:              // *`
` 414:              // *     Test the input parameters.`
` 415:              // *`
` 416:   `
` 417:              #endregion`
` 418:   `
` 419:   `
` 420:              #region Body`
` 421:              `
` 422:              WANTZ = this._lsame.Run(JOBZ, "V");`
` 423:              LOWER = this._lsame.Run(UPLO, "L");`
` 424:              LQUERY = (LWORK ==  - 1 || LIWORK ==  - 1);`
` 425:              // *`
` 426:              INFO = 0;`
` 427:              if (!(WANTZ || this._lsame.Run(JOBZ, "N")))`
` 428:              {`
` 429:                  INFO =  - 1;`
` 430:              }`
` 431:              else`
` 432:              {`
` 433:                  if (!(LOWER || this._lsame.Run(UPLO, "U")))`
` 434:                  {`
` 435:                      INFO =  - 2;`
` 436:                  }`
` 437:                  else`
` 438:                  {`
` 439:                      if (N < 0)`
` 440:                      {`
` 441:                          INFO =  - 3;`
` 442:                      }`
` 443:                      else`
` 444:                      {`
` 445:                          if (LDA < Math.Max(1, N))`
` 446:                          {`
` 447:                              INFO =  - 5;`
` 448:                          }`
` 449:                      }`
` 450:                  }`
` 451:              }`
` 452:              // *`
` 453:              if (INFO == 0)`
` 454:              {`
` 455:                  if (N <= 1)`
` 456:                  {`
` 457:                      LIWMIN = 1;`
` 458:                      LWMIN = 1;`
` 459:                      LOPT = LWMIN;`
` 460:                      LIOPT = LIWMIN;`
` 461:                  }`
` 462:                  else`
` 463:                  {`
` 464:                      if (WANTZ)`
` 465:                      {`
` 466:                          LIWMIN = 3 + 5 * N;`
` 467:                          LWMIN = 1 + 6 * N + 2 * (int)Math.Pow(N, 2);`
` 468:                      }`
` 469:                      else`
` 470:                      {`
` 471:                          LIWMIN = 1;`
` 472:                          LWMIN = 2 * N + 1;`
` 473:                      }`
` 474:                      LOPT = Math.Max(LWMIN, 2 * N + this._ilaenv.Run(1, "DSYTRD", UPLO, N,  - 1,  - 1,  - 1));`
` 475:                      LIOPT = LIWMIN;`
` 476:                  }`
` 477:                  WORK[1 + o_work] = LOPT;`
` 478:                  IWORK[1 + o_iwork] = LIOPT;`
` 479:                  // *`
` 480:                  if (LWORK < LWMIN && !LQUERY)`
` 481:                  {`
` 482:                      INFO =  - 8;`
` 483:                  }`
` 484:                  else`
` 485:                  {`
` 486:                      if (LIWORK < LIWMIN && !LQUERY)`
` 487:                      {`
` 488:                          INFO =  - 10;`
` 489:                      }`
` 490:                  }`
` 491:              }`
` 492:              // *`
` 493:              if (INFO != 0)`
` 494:              {`
` 495:                  this._xerbla.Run("DSYEVD",  - INFO);`
` 496:                  return;`
` 497:              }`
` 498:              else`
` 499:              {`
` 500:                  if (LQUERY)`
` 501:                  {`
` 502:                      return;`
` 503:                  }`
` 504:              }`
` 505:              // *`
` 506:              // *     Quick return if possible`
` 507:              // *`
` 508:              if (N == 0) return;`
` 509:              // *`
` 510:              if (N == 1)`
` 511:              {`
` 512:                  W[1 + o_w] = A[1+1 * LDA + o_a];`
` 513:                  if (WANTZ) A[1+1 * LDA + o_a] = ONE;`
` 514:                  return;`
` 515:              }`
` 516:              // *`
` 517:              // *     Get machine constants.`
` 518:              // *`
` 519:              SAFMIN = this._dlamch.Run("Safe minimum");`
` 520:              EPS = this._dlamch.Run("Precision");`
` 521:              SMLNUM = SAFMIN / EPS;`
` 522:              BIGNUM = ONE / SMLNUM;`
` 523:              RMIN = Math.Sqrt(SMLNUM);`
` 524:              RMAX = Math.Sqrt(BIGNUM);`
` 525:              // *`
` 526:              // *     Scale matrix to allowable range, if necessary.`
` 527:              // *`
` 528:              ANRM = this._dlansy.Run("M", UPLO, N, A, offset_a, LDA, ref WORK, offset_work);`
` 529:              ISCALE = 0;`
` 530:              if (ANRM > ZERO && ANRM < RMIN)`
` 531:              {`
` 532:                  ISCALE = 1;`
` 533:                  SIGMA = RMIN / ANRM;`
` 534:              }`
` 535:              else`
` 536:              {`
` 537:                  if (ANRM > RMAX)`
` 538:                  {`
` 539:                      ISCALE = 1;`
` 540:                      SIGMA = RMAX / ANRM;`
` 541:                  }`
` 542:              }`
` 543:              if (ISCALE == 1)`
` 544:              {`
` 545:                  this._dlascl.Run(UPLO, 0, 0, ONE, SIGMA, N`
` 546:                                   , N, ref A, offset_a, LDA, ref INFO);`
` 547:              }`
` 548:              // *`
` 549:              // *     Call DSYTRD to reduce symmetric matrix to tridiagonal form.`
` 550:              // *`
` 551:              INDE = 1;`
` 552:              INDTAU = INDE + N;`
` 553:              INDWRK = INDTAU + N;`
` 554:              LLWORK = LWORK - INDWRK + 1;`
` 555:              INDWK2 = INDWRK + N * N;`
` 556:              LLWRK2 = LWORK - INDWK2 + 1;`
` 557:              // *`
` 558:              this._dsytrd.Run(UPLO, N, ref A, offset_a, LDA, ref W, offset_w, ref WORK, INDE + o_work`
` 559:                               , ref WORK, INDTAU + o_work, ref WORK, INDWRK + o_work, LLWORK, ref IINFO);`
` 560:              LOPT = 2 * N + (int)WORK[INDWRK + o_work];`
` 561:              // *`
` 562:              // *     For eigenvalues only, call DSTERF.  For eigenvectors, first call`
` 563:              // *     DSTEDC to generate the eigenvector matrix, WORK(INDWRK), of the`
` 564:              // *     tridiagonal matrix, then call DORMTR to multiply it by the`
` 565:              // *     Householder transformations stored in A.`
` 566:              // *`
` 567:              if (!WANTZ)`
` 568:              {`
` 569:                  this._dsterf.Run(N, ref W, offset_w, ref WORK, INDE + o_work, ref INFO);`
` 570:              }`
` 571:              else`
` 572:              {`
` 573:                  this._dstedc.Run("I", N, ref W, offset_w, ref WORK, INDE + o_work, ref WORK, INDWRK + o_work, N`
` 574:                                   , ref WORK, INDWK2 + o_work, LLWRK2, ref IWORK, offset_iwork, LIWORK, ref INFO);`
` 575:                  this._dormtr.Run("L", UPLO, "N", N, N, ref A, offset_a`
` 576:                                   , LDA, WORK, INDTAU + o_work, ref WORK, INDWRK + o_work, N, ref WORK, INDWK2 + o_work, LLWRK2`
` 577:                                   , ref IINFO);`
` 578:                  this._dlacpy.Run("A", N, N, WORK, INDWRK + o_work, N, ref A, offset_a`
` 579:                                   , LDA);`
` 580:                  LOPT = (int)Math.Max(LOPT, 1 + 6 * N + 2 * Math.Pow(N, 2));`
` 581:              }`
` 582:              // *`
` 583:              // *     If matrix was scaled, then rescale eigenvalues appropriately.`
` 584:              // *`
` 585:              if (ISCALE == 1) this._dscal.Run(N, ONE / SIGMA, ref W, offset_w, 1);`
` 586:              // *`
` 587:              WORK[1 + o_work] = LOPT;`
` 588:              IWORK[1 + o_iwork] = LIOPT;`
` 589:              // *`
` 590:              return;`
` 591:              // *`
` 592:              // *     End of DSYEVD`
` 593:              // *`
` 594:   `
` 595:              #endregion`
` 596:   `
` 597:          }`
` 598:      }`
` 599:  }`