`   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 auxiliary routine (version 3.1) --`
`  21:      /// Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..`
`  22:      /// November 2006`
`  23:      /// Purpose`
`  24:      /// =======`
`  25:      /// `
`  26:      /// DLASY2 solves for the N1 by N2 matrix X, 1 .LE. N1,N2 .LE. 2, in`
`  27:      /// `
`  28:      /// op(TL)*X + ISGN*X*op(TR) = SCALE*B,`
`  29:      /// `
`  30:      /// where TL is N1 by N1, TR is N2 by N2, B is N1 by N2, and ISGN = 1 or`
`  31:      /// -1.  op(T) = T or T', where T' denotes the transpose of T.`
`  32:      /// `
`  33:      ///</summary>`
`  34:      public class DLASY2`
`  35:      {`
`  36:      `
`  37:   `
`  38:          #region Dependencies`
`  39:          `
`  40:          IDAMAX _idamax; DLAMCH _dlamch; DCOPY _dcopy; DSWAP _dswap; `
`  41:   `
`  42:          #endregion`
`  43:   `
`  44:   `
`  45:          #region Fields`
`  46:          `
`  47:          const double ZERO = 0.0E+0; const double ONE = 1.0E+0; const double TWO = 2.0E+0; const double HALF = 0.5E+0; `
`  48:          const double EIGHT = 8.0E+0;bool BSWAP = false; bool XSWAP = false; int I = 0; int IP = 0; int IPIV = 0; int IPSV = 0; `
`  49:          int J = 0;int JP = 0; int JPSV = 0; int K = 0; double BET = 0; double EPS = 0; double GAM = 0; double L21 = 0; `
`  50:          double SGN = 0;double SMIN = 0; double SMLNUM = 0; double TAU1 = 0; double TEMP = 0; double U11 = 0; double U12 = 0; `
`  51:          double U22 = 0;double XMAX = 0; bool[] BSWPIV = new bool[4]; int o_bswpiv = -1; `
`  52:          bool[] XSWPIV = new bool[4]; int o_xswpiv = -1;int[] JPIV = new int[4]; int o_jpiv = -1; `
`  53:          int[] LOCL21 = new int[4];  int o_locl21 = -1;int[] LOCU12 = new int[4];  int o_locu12 = -1; `
`  54:          int[] LOCU22 = new int[4];  int o_locu22 = -1;double[] BTMP = new double[4]; int offset_btmp = 0; int o_btmp = -1; `
`  55:          double[] T16 = new double[4 * 4]; int offset_t16 = 0; int o_t16 = -5;double[] TMP = new double[4]; int offset_tmp = 0; int o_tmp = -1; `
`  56:          double[] X2 = new double[2];  int o_x2 = -1;`
`  57:   `
`  58:          #endregion`
`  59:   `
`  60:          public DLASY2(IDAMAX idamax, DLAMCH dlamch, DCOPY dcopy, DSWAP dswap)`
`  61:          {`
`  62:      `
`  63:   `
`  64:              #region Set Dependencies`
`  65:              `
`  66:              this._idamax = idamax; this._dlamch = dlamch; this._dcopy = dcopy; this._dswap = dswap; `
`  67:   `
`  68:              #endregion`
`  69:   `
`  70:   `
`  71:              #region Data Inicializacion`
`  72:              `
`  73:              //LOCU12/3,4,1,2`
`  74:              LOCU12[1 + o_locu12] = 3;`
`  75:              LOCU12[2 + o_locu12] = 4;`
`  76:              LOCU12[3 + o_locu12] = 1;`
`  77:              LOCU12[4 + o_locu12] = 2;`
`  78:              //LOCL21/2,1,4,3`
`  79:              LOCL21[1 + o_locl21] = 2;`
`  80:              LOCL21[2 + o_locl21] = 1;`
`  81:              LOCL21[3 + o_locl21] = 4;`
`  82:              LOCL21[4 + o_locl21] = 3;`
`  83:              //LOCU22/4,3,2,1`
`  84:              LOCU22[1 + o_locu22] = 4;`
`  85:              LOCU22[2 + o_locu22] = 3;`
`  86:              LOCU22[3 + o_locu22] = 2;`
`  87:              LOCU22[4 + o_locu22] = 1;`
`  88:              //XSWPIV/.FALSE.,.FALSE.,.TRUE.,.TRUE.`
`  89:              XSWPIV[1 + o_xswpiv] = false;`
`  90:              XSWPIV[2 + o_xswpiv] = false;`
`  91:              XSWPIV[3 + o_xswpiv] = true;`
`  92:              XSWPIV[4 + o_xswpiv] = true;`
`  93:              //BSWPIV/.FALSE.,.TRUE.,.FALSE.,.TRUE.`
`  94:              BSWPIV[1 + o_bswpiv] = false;`
`  95:              BSWPIV[2 + o_bswpiv] = true;`
`  96:              BSWPIV[3 + o_bswpiv] = false;`
`  97:              BSWPIV[4 + o_bswpiv] = true;`
`  98:   `
`  99:              #endregion`
` 100:   `
` 101:          }`
` 102:      `
` 103:          public DLASY2()`
` 104:          {`
` 105:      `
` 106:   `
` 107:              #region Dependencies (Initialization)`
` 108:              `
` 109:              IDAMAX idamax = new IDAMAX();`
` 110:              LSAME lsame = new LSAME();`
` 111:              DLAMC3 dlamc3 = new DLAMC3();`
` 112:              DCOPY dcopy = new DCOPY();`
` 113:              DSWAP dswap = new DSWAP();`
` 114:              DLAMC1 dlamc1 = new DLAMC1(dlamc3);`
` 115:              DLAMC4 dlamc4 = new DLAMC4(dlamc3);`
` 116:              DLAMC5 dlamc5 = new DLAMC5(dlamc3);`
` 117:              DLAMC2 dlamc2 = new DLAMC2(dlamc3, dlamc1, dlamc4, dlamc5);`
` 118:              DLAMCH dlamch = new DLAMCH(lsame, dlamc2);`
` 119:   `
` 120:              #endregion`
` 121:   `
` 122:   `
` 123:              #region Set Dependencies`
` 124:              `
` 125:              this._idamax = idamax; this._dlamch = dlamch; this._dcopy = dcopy; this._dswap = dswap; `
` 126:   `
` 127:              #endregion`
` 128:   `
` 129:   `
` 130:              #region Data Inicializacion`
` 131:              `
` 132:              //LOCU12/3,4,1,2`
` 133:              LOCU12[1 + o_locu12] = 3;`
` 134:              LOCU12[2 + o_locu12] = 4;`
` 135:              LOCU12[3 + o_locu12] = 1;`
` 136:              LOCU12[4 + o_locu12] = 2;`
` 137:              //LOCL21/2,1,4,3`
` 138:              LOCL21[1 + o_locl21] = 2;`
` 139:              LOCL21[2 + o_locl21] = 1;`
` 140:              LOCL21[3 + o_locl21] = 4;`
` 141:              LOCL21[4 + o_locl21] = 3;`
` 142:              //LOCU22/4,3,2,1`
` 143:              LOCU22[1 + o_locu22] = 4;`
` 144:              LOCU22[2 + o_locu22] = 3;`
` 145:              LOCU22[3 + o_locu22] = 2;`
` 146:              LOCU22[4 + o_locu22] = 1;`
` 147:              //XSWPIV/.FALSE.,.FALSE.,.TRUE.,.TRUE.`
` 148:              XSWPIV[1 + o_xswpiv] = false;`
` 149:              XSWPIV[2 + o_xswpiv] = false;`
` 150:              XSWPIV[3 + o_xswpiv] = true;`
` 151:              XSWPIV[4 + o_xswpiv] = true;`
` 152:              //BSWPIV/.FALSE.,.TRUE.,.FALSE.,.TRUE.`
` 153:              BSWPIV[1 + o_bswpiv] = false;`
` 154:              BSWPIV[2 + o_bswpiv] = true;`
` 155:              BSWPIV[3 + o_bswpiv] = false;`
` 156:              BSWPIV[4 + o_bswpiv] = true;`
` 157:   `
` 158:              #endregion`
` 159:   `
` 160:          }`
` 161:          /// <summary>`
` 162:          /// Purpose`
` 163:          /// =======`
` 164:          /// `
` 165:          /// DLASY2 solves for the N1 by N2 matrix X, 1 .LE. N1,N2 .LE. 2, in`
` 166:          /// `
` 167:          /// op(TL)*X + ISGN*X*op(TR) = SCALE*B,`
` 168:          /// `
` 169:          /// where TL is N1 by N1, TR is N2 by N2, B is N1 by N2, and ISGN = 1 or`
` 170:          /// -1.  op(T) = T or T', where T' denotes the transpose of T.`
` 171:          /// `
` 172:          ///</summary>`
` 173:          /// <param name="LTRANL">`
` 174:          /// (input) LOGICAL`
` 175:          /// On entry, LTRANL specifies the op(TL):`
` 176:          /// = .FALSE., op(TL) = TL,`
` 177:          /// = .TRUE., op(TL) = TL'.`
` 178:          ///</param>`
` 179:          /// <param name="LTRANR">`
` 180:          /// (input) LOGICAL`
` 181:          /// On entry, LTRANR specifies the op(TR):`
` 182:          /// = .FALSE., op(TR) = TR,`
` 183:          /// = .TRUE., op(TR) = TR'.`
` 184:          ///</param>`
` 185:          /// <param name="ISGN">`
` 186:          /// (input) INTEGER`
` 187:          /// On entry, ISGN specifies the sign of the equation`
` 188:          /// as described before. ISGN may only be 1 or -1.`
` 189:          ///</param>`
` 190:          /// <param name="N1">`
` 191:          /// (input) INTEGER`
` 192:          /// On entry, N1 specifies the order of matrix TL.`
` 193:          /// N1 may only be 0, 1 or 2.`
` 194:          ///</param>`
` 195:          /// <param name="N2">`
` 196:          /// (input) INTEGER`
` 197:          /// On entry, N2 specifies the order of matrix TR.`
` 198:          /// N2 may only be 0, 1 or 2.`
` 199:          ///</param>`
` 200:          /// <param name="TL">`
` 201:          /// (input) DOUBLE PRECISION array, dimension (LDTL,2)`
` 202:          /// On entry, TL contains an N1 by N1 matrix.`
` 203:          ///</param>`
` 204:          /// <param name="LDTL">`
` 205:          /// (input) INTEGER`
` 206:          /// The leading dimension of the matrix TL. LDTL .GE. max(1,N1).`
` 207:          ///</param>`
` 208:          /// <param name="TR">`
` 209:          /// (input) DOUBLE PRECISION array, dimension (LDTR,2)`
` 210:          /// On entry, TR contains an N2 by N2 matrix.`
` 211:          ///</param>`
` 212:          /// <param name="LDTR">`
` 213:          /// (input) INTEGER`
` 214:          /// The leading dimension of the matrix TR. LDTR .GE. max(1,N2).`
` 215:          ///</param>`
` 216:          /// <param name="B">`
` 217:          /// (input) DOUBLE PRECISION array, dimension (LDB,2)`
` 218:          /// On entry, the N1 by N2 matrix B contains the right-hand`
` 219:          /// side of the equation.`
` 220:          ///</param>`
` 221:          /// <param name="LDB">`
` 222:          /// (input) INTEGER`
` 223:          /// The leading dimension of the matrix B. LDB .GE. max(1,N1).`
` 224:          ///</param>`
` 225:          /// <param name="SCALE">`
` 226:          /// (output) DOUBLE PRECISION`
` 227:          /// On exit, SCALE contains the scale factor. SCALE is chosen`
` 228:          /// less than or equal to 1 to prevent the solution overflowing.`
` 229:          ///</param>`
` 230:          /// <param name="X">`
` 231:          /// (output) DOUBLE PRECISION array, dimension (LDX,2)`
` 232:          /// On exit, X contains the N1 by N2 solution.`
` 233:          ///</param>`
` 234:          /// <param name="LDX">`
` 235:          /// (input) INTEGER`
` 236:          /// The leading dimension of the matrix X. LDX .GE. max(1,N1).`
` 237:          ///</param>`
` 238:          /// <param name="XNORM">`
` 239:          /// (output) DOUBLE PRECISION`
` 240:          /// On exit, XNORM is the infinity-norm of the solution.`
` 241:          ///</param>`
` 242:          /// <param name="INFO">`
` 243:          /// (output) INTEGER`
` 244:          /// On exit, INFO is set to`
` 245:          /// 0: successful exit.`
` 246:          /// 1: TL and TR have too close eigenvalues, so TL or`
` 247:          /// TR is perturbed to get a nonsingular equation.`
` 248:          /// NOTE: In the interests of speed, this routine does not`
` 249:          /// check the inputs for errors.`
` 250:          ///</param>`
` 251:          public void Run(bool LTRANL, bool LTRANR, int ISGN, int N1, int N2, double[] TL, int offset_tl`
` 252:                           , int LDTL, double[] TR, int offset_tr, int LDTR, double[] B, int offset_b, int LDB, ref double SCALE`
` 253:                           , ref double[] X, int offset_x, int LDX, ref double XNORM, ref int INFO)`
` 254:          {`
` 255:   `
` 256:              #region Array Index Correction`
` 257:              `
` 258:               int o_tl = -1 - LDTL + offset_tl;  int o_tr = -1 - LDTR + offset_tr;  int o_b = -1 - LDB + offset_b; `
` 259:               int o_x = -1 - LDX + offset_x;`
` 260:   `
` 261:              #endregion`
` 262:   `
` 263:   `
` 264:              #region Prolog`
` 265:              `
` 266:              // *`
` 267:              // *  -- LAPACK auxiliary routine (version 3.1) --`
` 268:              // *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..`
` 269:              // *     November 2006`
` 270:              // *`
` 271:              // *     .. Scalar Arguments ..`
` 272:              // *     ..`
` 273:              // *     .. Array Arguments ..`
` 274:              // *     ..`
` 275:              // *`
` 276:              // *  Purpose`
` 277:              // *  =======`
` 278:              // *`
` 279:              // *  DLASY2 solves for the N1 by N2 matrix X, 1 <= N1,N2 <= 2, in`
` 280:              // *`
` 281:              // *         op(TL)*X + ISGN*X*op(TR) = SCALE*B,`
` 282:              // *`
` 283:              // *  where TL is N1 by N1, TR is N2 by N2, B is N1 by N2, and ISGN = 1 or`
` 284:              // *  -1.  op(T) = T or T', where T' denotes the transpose of T.`
` 285:              // *`
` 286:              // *  Arguments`
` 287:              // *  =========`
` 288:              // *`
` 289:              // *  LTRANL  (input) LOGICAL`
` 290:              // *          On entry, LTRANL specifies the op(TL):`
` 291:              // *             = .FALSE., op(TL) = TL,`
` 292:              // *             = .TRUE., op(TL) = TL'.`
` 293:              // *`
` 294:              // *  LTRANR  (input) LOGICAL`
` 295:              // *          On entry, LTRANR specifies the op(TR):`
` 296:              // *            = .FALSE., op(TR) = TR,`
` 297:              // *            = .TRUE., op(TR) = TR'.`
` 298:              // *`
` 299:              // *  ISGN    (input) INTEGER`
` 300:              // *          On entry, ISGN specifies the sign of the equation`
` 301:              // *          as described before. ISGN may only be 1 or -1.`
` 302:              // *`
` 303:              // *  N1      (input) INTEGER`
` 304:              // *          On entry, N1 specifies the order of matrix TL.`
` 305:              // *          N1 may only be 0, 1 or 2.`
` 306:              // *`
` 307:              // *  N2      (input) INTEGER`
` 308:              // *          On entry, N2 specifies the order of matrix TR.`
` 309:              // *          N2 may only be 0, 1 or 2.`
` 310:              // *`
` 311:              // *  TL      (input) DOUBLE PRECISION array, dimension (LDTL,2)`
` 312:              // *          On entry, TL contains an N1 by N1 matrix.`
` 313:              // *`
` 314:              // *  LDTL    (input) INTEGER`
` 315:              // *          The leading dimension of the matrix TL. LDTL >= max(1,N1).`
` 316:              // *`
` 317:              // *  TR      (input) DOUBLE PRECISION array, dimension (LDTR,2)`
` 318:              // *          On entry, TR contains an N2 by N2 matrix.`
` 319:              // *`
` 320:              // *  LDTR    (input) INTEGER`
` 321:              // *          The leading dimension of the matrix TR. LDTR >= max(1,N2).`
` 322:              // *`
` 323:              // *  B       (input) DOUBLE PRECISION array, dimension (LDB,2)`
` 324:              // *          On entry, the N1 by N2 matrix B contains the right-hand`
` 325:              // *          side of the equation.`
` 326:              // *`
` 327:              // *  LDB     (input) INTEGER`
` 328:              // *          The leading dimension of the matrix B. LDB >= max(1,N1).`
` 329:              // *`
` 330:              // *  SCALE   (output) DOUBLE PRECISION`
` 331:              // *          On exit, SCALE contains the scale factor. SCALE is chosen`
` 332:              // *          less than or equal to 1 to prevent the solution overflowing.`
` 333:              // *`
` 334:              // *  X       (output) DOUBLE PRECISION array, dimension (LDX,2)`
` 335:              // *          On exit, X contains the N1 by N2 solution.`
` 336:              // *`
` 337:              // *  LDX     (input) INTEGER`
` 338:              // *          The leading dimension of the matrix X. LDX >= max(1,N1).`
` 339:              // *`
` 340:              // *  XNORM   (output) DOUBLE PRECISION`
` 341:              // *          On exit, XNORM is the infinity-norm of the solution.`
` 342:              // *`
` 343:              // *  INFO    (output) INTEGER`
` 344:              // *          On exit, INFO is set to`
` 345:              // *             0: successful exit.`
` 346:              // *             1: TL and TR have too close eigenvalues, so TL or`
` 347:              // *                TR is perturbed to get a nonsingular equation.`
` 348:              // *          NOTE: In the interests of speed, this routine does not`
` 349:              // *                check the inputs for errors.`
` 350:              // *`
` 351:              // * =====================================================================`
` 352:              // *`
` 353:              // *     .. Parameters ..`
` 354:              // *     ..`
` 355:              // *     .. Local Scalars ..`
` 356:              // *     ..`
` 357:              // *     .. Local Arrays ..`
` 358:              // *     ..`
` 359:              // *     .. External Functions ..`
` 360:              // *     ..`
` 361:              // *     .. External Subroutines ..`
` 362:              // *     ..`
` 363:              // *     .. Intrinsic Functions ..`
` 364:              //      INTRINSIC          ABS, MAX;`
` 365:              // *     ..`
` 366:              // *     .. Data statements ..`
` 367:              // *     ..`
` 368:              // *     .. Executable Statements ..`
` 369:              // *`
` 370:              // *     Do not check the input parameters for errors`
` 371:              // *`
` 372:   `
` 373:              #endregion`
` 374:   `
` 375:   `
` 376:              #region Body`
` 377:              `
` 378:              INFO = 0;`
` 379:              // *`
` 380:              // *     Quick return if possible`
` 381:              // *`
` 382:              if (N1 == 0 || N2 == 0) return;`
` 383:              // *`
` 384:              // *     Set constants to control overflow`
` 385:              // *`
` 386:              EPS = this._dlamch.Run("P");`
` 387:              SMLNUM = this._dlamch.Run("S") / EPS;`
` 388:              SGN = ISGN;`
` 389:              // *`
` 390:              K = N1 + N1 + N2 - 2;`
` 391:              switch (K)`
` 392:              {`
` 393:                  case 1: goto LABEL10;`
` 394:                  case 2: goto LABEL20;`
` 395:                  case 3: goto LABEL30;`
` 396:                  case 4: goto LABEL50;`
` 397:              }`
` 398:              // *`
` 399:              // *     1 by 1: TL11*X + SGN*X*TR11 = B11`
` 400:              // *`
` 401:          LABEL10:;`
` 402:              TAU1 = TL[1+1 * LDTL + o_tl] + SGN * TR[1+1 * LDTR + o_tr];`
` 403:              BET = Math.Abs(TAU1);`
` 404:              if (BET <= SMLNUM)`
` 405:              {`
` 406:                  TAU1 = SMLNUM;`
` 407:                  BET = SMLNUM;`
` 408:                  INFO = 1;`
` 409:              }`
` 410:              // *`
` 411:              SCALE = ONE;`
` 412:              GAM = Math.Abs(B[1+1 * LDB + o_b]);`
` 413:              if (SMLNUM * GAM > BET) SCALE = ONE / GAM;`
` 414:              // *`
` 415:              X[1+1 * LDX + o_x] = (B[1+1 * LDB + o_b] * SCALE) / TAU1;`
` 416:              XNORM = Math.Abs(X[1+1 * LDX + o_x]);`
` 417:              return;`
` 418:              // *`
` 419:              // *     1 by 2:`
` 420:              // *     TL11*[X11 X12] + ISGN*[X11 X12]*op[TR11 TR12]  = [B11 B12]`
` 421:              // *                                       [TR21 TR22]`
` 422:              // *`
` 423:          LABEL20:;`
` 424:              // *`
` 425:              SMIN = Math.Max(EPS * Math.Max(Math.Abs(TL[1+1 * LDTL + o_tl]), Math.Max(Math.Abs(TR[1+1 * LDTR + o_tr]), Math.Max(Math.Abs(TR[1+2 * LDTR + o_tr]), Math.Max(Math.Abs(TR[2+1 * LDTR + o_tr]), Math.Abs(TR[2+2 * LDTR + o_tr]))))), SMLNUM);`
` 426:              TMP[1 + o_tmp] = TL[1+1 * LDTL + o_tl] + SGN * TR[1+1 * LDTR + o_tr];`
` 427:              TMP[4 + o_tmp] = TL[1+1 * LDTL + o_tl] + SGN * TR[2+2 * LDTR + o_tr];`
` 428:              if (LTRANR)`
` 429:              {`
` 430:                  TMP[2 + o_tmp] = SGN * TR[2+1 * LDTR + o_tr];`
` 431:                  TMP[3 + o_tmp] = SGN * TR[1+2 * LDTR + o_tr];`
` 432:              }`
` 433:              else`
` 434:              {`
` 435:                  TMP[2 + o_tmp] = SGN * TR[1+2 * LDTR + o_tr];`
` 436:                  TMP[3 + o_tmp] = SGN * TR[2+1 * LDTR + o_tr];`
` 437:              }`
` 438:              BTMP[1 + o_btmp] = B[1+1 * LDB + o_b];`
` 439:              BTMP[2 + o_btmp] = B[1+2 * LDB + o_b];`
` 440:              goto LABEL40;`
` 441:              // *`
` 442:              // *     2 by 1:`
` 443:              // *          op[TL11 TL12]*[X11] + ISGN* [X11]*TR11  = [B11]`
` 444:              // *            [TL21 TL22] [X21]         [X21]         [B21]`
` 445:              // *`
` 446:          LABEL30:;`
` 447:              SMIN = Math.Max(EPS * Math.Max(Math.Abs(TR[1+1 * LDTR + o_tr]), Math.Max(Math.Abs(TL[1+1 * LDTL + o_tl]), Math.Max(Math.Abs(TL[1+2 * LDTL + o_tl]), Math.Max(Math.Abs(TL[2+1 * LDTL + o_tl]), Math.Abs(TL[2+2 * LDTL + o_tl]))))), SMLNUM);`
` 448:              TMP[1 + o_tmp] = TL[1+1 * LDTL + o_tl] + SGN * TR[1+1 * LDTR + o_tr];`
` 449:              TMP[4 + o_tmp] = TL[2+2 * LDTL + o_tl] + SGN * TR[1+1 * LDTR + o_tr];`
` 450:              if (LTRANL)`
` 451:              {`
` 452:                  TMP[2 + o_tmp] = TL[1+2 * LDTL + o_tl];`
` 453:                  TMP[3 + o_tmp] = TL[2+1 * LDTL + o_tl];`
` 454:              }`
` 455:              else`
` 456:              {`
` 457:                  TMP[2 + o_tmp] = TL[2+1 * LDTL + o_tl];`
` 458:                  TMP[3 + o_tmp] = TL[1+2 * LDTL + o_tl];`
` 459:              }`
` 460:              BTMP[1 + o_btmp] = B[1+1 * LDB + o_b];`
` 461:              BTMP[2 + o_btmp] = B[2+1 * LDB + o_b];`
` 462:          LABEL40:;`
` 463:              // *`
` 464:              // *     Solve 2 by 2 system using complete pivoting.`
` 465:              // *     Set pivots less than SMIN to SMIN.`
` 466:              // *`
` 467:              IPIV = this._idamax.Run(4, TMP, offset_tmp, 1);`
` 468:              U11 = TMP[IPIV + o_tmp];`
` 469:              if (Math.Abs(U11) <= SMIN)`
` 470:              {`
` 471:                  INFO = 1;`
` 472:                  U11 = SMIN;`
` 473:              }`
` 474:              U12 = TMP[LOCU12[IPIV + o_locu12] + o_tmp];`
` 475:              L21 = TMP[LOCL21[IPIV + o_locl21] + o_tmp] / U11;`
` 476:              U22 = TMP[LOCU22[IPIV + o_locu22] + o_tmp] - U12 * L21;`
` 477:              XSWAP = XSWPIV[IPIV + o_xswpiv];`
` 478:              BSWAP = BSWPIV[IPIV + o_bswpiv];`
` 479:              if (Math.Abs(U22) <= SMIN)`
` 480:              {`
` 481:                  INFO = 1;`
` 482:                  U22 = SMIN;`
` 483:              }`
` 484:              if (BSWAP)`
` 485:              {`
` 486:                  TEMP = BTMP[2 + o_btmp];`
` 487:                  BTMP[2 + o_btmp] = BTMP[1 + o_btmp] - L21 * TEMP;`
` 488:                  BTMP[1 + o_btmp] = TEMP;`
` 489:              }`
` 490:              else`
` 491:              {`
` 492:                  BTMP[2 + o_btmp] = BTMP[2 + o_btmp] - L21 * BTMP[1 + o_btmp];`
` 493:              }`
` 494:              SCALE = ONE;`
` 495:              if ((TWO * SMLNUM) * Math.Abs(BTMP[2 + o_btmp]) > Math.Abs(U22) || (TWO * SMLNUM) * Math.Abs(BTMP[1 + o_btmp]) > Math.Abs(U11))`
` 496:              {`
` 497:                  SCALE = HALF / Math.Max(Math.Abs(BTMP[1 + o_btmp]), Math.Abs(BTMP[2 + o_btmp]));`
` 498:                  BTMP[1 + o_btmp] = BTMP[1 + o_btmp] * SCALE;`
` 499:                  BTMP[2 + o_btmp] = BTMP[2 + o_btmp] * SCALE;`
` 500:              }`
` 501:              X2[2 + o_x2] = BTMP[2 + o_btmp] / U22;`
` 502:              X2[1 + o_x2] = BTMP[1 + o_btmp] / U11 - (U12 / U11) * X2[2 + o_x2];`
` 503:              if (XSWAP)`
` 504:              {`
` 505:                  TEMP = X2[2 + o_x2];`
` 506:                  X2[2 + o_x2] = X2[1 + o_x2];`
` 507:                  X2[1 + o_x2] = TEMP;`
` 508:              }`
` 509:              X[1+1 * LDX + o_x] = X2[1 + o_x2];`
` 510:              if (N1 == 1)`
` 511:              {`
` 512:                  X[1+2 * LDX + o_x] = X2[2 + o_x2];`
` 513:                  XNORM = Math.Abs(X[1+1 * LDX + o_x]) + Math.Abs(X[1+2 * LDX + o_x]);`
` 514:              }`
` 515:              else`
` 516:              {`
` 517:                  X[2+1 * LDX + o_x] = X2[2 + o_x2];`
` 518:                  XNORM = Math.Max(Math.Abs(X[1+1 * LDX + o_x]), Math.Abs(X[2+1 * LDX + o_x]));`
` 519:              }`
` 520:              return;`
` 521:              // *`
` 522:              // *     2 by 2:`
` 523:              // *     op[TL11 TL12]*[X11 X12] +ISGN* [X11 X12]*op[TR11 TR12] = [B11 B12]`
` 524:              // *       [TL21 TL22] [X21 X22]        [X21 X22]   [TR21 TR22]   [B21 B22]`
` 525:              // *`
` 526:              // *     Solve equivalent 4 by 4 system using complete pivoting.`
` 527:              // *     Set pivots less than SMIN to SMIN.`
` 528:              // *`
` 529:          LABEL50:;`
` 530:              SMIN = Math.Max(Math.Abs(TR[1+1 * LDTR + o_tr]), Math.Max(Math.Abs(TR[1+2 * LDTR + o_tr]), Math.Max(Math.Abs(TR[2+1 * LDTR + o_tr]), Math.Abs(TR[2+2 * LDTR + o_tr]))));`
` 531:              SMIN = Math.Max(SMIN, Math.Max(Math.Abs(TL[1+1 * LDTL + o_tl]), Math.Max(Math.Abs(TL[1+2 * LDTL + o_tl]), Math.Max(Math.Abs(TL[2+1 * LDTL + o_tl]), Math.Abs(TL[2+2 * LDTL + o_tl])))));`
` 532:              SMIN = Math.Max(EPS * SMIN, SMLNUM);`
` 533:              BTMP[1 + o_btmp] = ZERO;`
` 534:              this._dcopy.Run(16, BTMP, offset_btmp, 0, ref T16, offset_t16, 1);`
` 535:              T16[1+1 * 4 + o_t16] = TL[1+1 * LDTL + o_tl] + SGN * TR[1+1 * LDTR + o_tr];`
` 536:              T16[2+2 * 4 + o_t16] = TL[2+2 * LDTL + o_tl] + SGN * TR[1+1 * LDTR + o_tr];`
` 537:              T16[3+3 * 4 + o_t16] = TL[1+1 * LDTL + o_tl] + SGN * TR[2+2 * LDTR + o_tr];`
` 538:              T16[4+4 * 4 + o_t16] = TL[2+2 * LDTL + o_tl] + SGN * TR[2+2 * LDTR + o_tr];`
` 539:              if (LTRANL)`
` 540:              {`
` 541:                  T16[1+2 * 4 + o_t16] = TL[2+1 * LDTL + o_tl];`
` 542:                  T16[2+1 * 4 + o_t16] = TL[1+2 * LDTL + o_tl];`
` 543:                  T16[3+4 * 4 + o_t16] = TL[2+1 * LDTL + o_tl];`
` 544:                  T16[4+3 * 4 + o_t16] = TL[1+2 * LDTL + o_tl];`
` 545:              }`
` 546:              else`
` 547:              {`
` 548:                  T16[1+2 * 4 + o_t16] = TL[1+2 * LDTL + o_tl];`
` 549:                  T16[2+1 * 4 + o_t16] = TL[2+1 * LDTL + o_tl];`
` 550:                  T16[3+4 * 4 + o_t16] = TL[1+2 * LDTL + o_tl];`
` 551:                  T16[4+3 * 4 + o_t16] = TL[2+1 * LDTL + o_tl];`
` 552:              }`
` 553:              if (LTRANR)`
` 554:              {`
` 555:                  T16[1+3 * 4 + o_t16] = SGN * TR[1+2 * LDTR + o_tr];`
` 556:                  T16[2+4 * 4 + o_t16] = SGN * TR[1+2 * LDTR + o_tr];`
` 557:                  T16[3+1 * 4 + o_t16] = SGN * TR[2+1 * LDTR + o_tr];`
` 558:                  T16[4+2 * 4 + o_t16] = SGN * TR[2+1 * LDTR + o_tr];`
` 559:              }`
` 560:              else`
` 561:              {`
` 562:                  T16[1+3 * 4 + o_t16] = SGN * TR[2+1 * LDTR + o_tr];`
` 563:                  T16[2+4 * 4 + o_t16] = SGN * TR[2+1 * LDTR + o_tr];`
` 564:                  T16[3+1 * 4 + o_t16] = SGN * TR[1+2 * LDTR + o_tr];`
` 565:                  T16[4+2 * 4 + o_t16] = SGN * TR[1+2 * LDTR + o_tr];`
` 566:              }`
` 567:              BTMP[1 + o_btmp] = B[1+1 * LDB + o_b];`
` 568:              BTMP[2 + o_btmp] = B[2+1 * LDB + o_b];`
` 569:              BTMP[3 + o_btmp] = B[1+2 * LDB + o_b];`
` 570:              BTMP[4 + o_btmp] = B[2+2 * LDB + o_b];`
` 571:              // *`
` 572:              // *     Perform elimination`
` 573:              // *`
` 574:              for (I = 1; I <= 3; I++)`
` 575:              {`
` 576:                  XMAX = ZERO;`
` 577:                  for (IP = I; IP <= 4; IP++)`
` 578:                  {`
` 579:                      for (JP = I; JP <= 4; JP++)`
` 580:                      {`
` 581:                          if (Math.Abs(T16[IP+JP * 4 + o_t16]) >= XMAX)`
` 582:                          {`
` 583:                              XMAX = Math.Abs(T16[IP+JP * 4 + o_t16]);`
` 584:                              IPSV = IP;`
` 585:                              JPSV = JP;`
` 586:                          }`
` 587:                      }`
` 588:                  }`
` 589:                  if (IPSV != I)`
` 590:                  {`
` 591:                      this._dswap.Run(4, ref T16, IPSV+1 * 4 + o_t16, 4, ref T16, I+1 * 4 + o_t16, 4);`
` 592:                      TEMP = BTMP[I + o_btmp];`
` 593:                      BTMP[I + o_btmp] = BTMP[IPSV + o_btmp];`
` 594:                      BTMP[IPSV + o_btmp] = TEMP;`
` 595:                  }`
` 596:                  if (JPSV != I) this._dswap.Run(4, ref T16, 1+JPSV * 4 + o_t16, 1, ref T16, 1+I * 4 + o_t16, 1);`
` 597:                  JPIV[I + o_jpiv] = JPSV;`
` 598:                  if (Math.Abs(T16[I+I * 4 + o_t16]) < SMIN)`
` 599:                  {`
` 600:                      INFO = 1;`
` 601:                      T16[I+I * 4 + o_t16] = SMIN;`
` 602:                  }`
` 603:                  for (J = I + 1; J <= 4; J++)`
` 604:                  {`
` 605:                      T16[J+I * 4 + o_t16] = T16[J+I * 4 + o_t16] / T16[I+I * 4 + o_t16];`
` 606:                      BTMP[J + o_btmp] = BTMP[J + o_btmp] - T16[J+I * 4 + o_t16] * BTMP[I + o_btmp];`
` 607:                      for (K = I + 1; K <= 4; K++)`
` 608:                      {`
` 609:                          T16[J+K * 4 + o_t16] = T16[J+K * 4 + o_t16] - T16[J+I * 4 + o_t16] * T16[I+K * 4 + o_t16];`
` 610:                      }`
` 611:                  }`
` 612:              }`
` 613:              if (Math.Abs(T16[4+4 * 4 + o_t16]) < SMIN) T16[4+4 * 4 + o_t16] = SMIN;`
` 614:              SCALE = ONE;`
` 615:              if ((EIGHT * SMLNUM) * Math.Abs(BTMP[1 + o_btmp]) > Math.Abs(T16[1+1 * 4 + o_t16]) || (EIGHT * SMLNUM) * Math.Abs(BTMP[2 + o_btmp]) > Math.Abs(T16[2+2 * 4 + o_t16]) || (EIGHT * SMLNUM) * Math.Abs(BTMP[3 + o_btmp]) > Math.Abs(T16[3+3 * 4 + o_t16]) || (EIGHT * SMLNUM) * Math.Abs(BTMP[4 + o_btmp]) > Math.Abs(T16[4+4 * 4 + o_t16]))`
` 616:              {`
` 617:                  SCALE = (ONE / EIGHT) / Math.Max(Math.Abs(BTMP[1 + o_btmp]), Math.Max(Math.Abs(BTMP[2 + o_btmp]), Math.Max(Math.Abs(BTMP[3 + o_btmp]), Math.Abs(BTMP[4 + o_btmp]))));`
` 618:                  BTMP[1 + o_btmp] = BTMP[1 + o_btmp] * SCALE;`
` 619:                  BTMP[2 + o_btmp] = BTMP[2 + o_btmp] * SCALE;`
` 620:                  BTMP[3 + o_btmp] = BTMP[3 + o_btmp] * SCALE;`
` 621:                  BTMP[4 + o_btmp] = BTMP[4 + o_btmp] * SCALE;`
` 622:              }`
` 623:              for (I = 1; I <= 4; I++)`
` 624:              {`
` 625:                  K = 5 - I;`
` 626:                  TEMP = ONE / T16[K+K * 4 + o_t16];`
` 627:                  TMP[K + o_tmp] = BTMP[K + o_btmp] * TEMP;`
` 628:                  for (J = K + 1; J <= 4; J++)`
` 629:                  {`
` 630:                      TMP[K + o_tmp] = TMP[K + o_tmp] - (TEMP * T16[K+J * 4 + o_t16]) * TMP[J + o_tmp];`
` 631:                  }`
` 632:              }`
` 633:              for (I = 1; I <= 3; I++)`
` 634:              {`
` 635:                  if (JPIV[4 - I + o_jpiv] != 4 - I)`
` 636:                  {`
` 637:                      TEMP = TMP[4 - I + o_tmp];`
` 638:                      TMP[4 - I + o_tmp] = TMP[JPIV[4 - I + o_jpiv] + o_tmp];`
` 639:                      TMP[JPIV[4 - I + o_jpiv] + o_tmp] = TEMP;`
` 640:                  }`
` 641:              }`
` 642:              X[1+1 * LDX + o_x] = TMP[1 + o_tmp];`
` 643:              X[2+1 * LDX + o_x] = TMP[2 + o_tmp];`
` 644:              X[1+2 * LDX + o_x] = TMP[3 + o_tmp];`
` 645:              X[2+2 * LDX + o_x] = TMP[4 + o_tmp];`
` 646:              XNORM = Math.Max(Math.Abs(TMP[1 + o_tmp]) + Math.Abs(TMP[3 + o_tmp]), Math.Abs(TMP[2 + o_tmp]) + Math.Abs(TMP[4 + o_tmp]));`
` 647:              return;`
` 648:              // *`
` 649:              // *     End of DLASY2`
` 650:              // *`
` 651:   `
` 652:              #endregion`
` 653:   `
` 654:          }`
` 655:      }`
` 656:  }`