Skip Navigation Links
Numerical Libraries
Linear Algebra
Differential Equations
Optimization
Samples
Skip Navigation Links
Linear Algebra
CSLapack
CSBlas
   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 routine (version 3.1) --
  21:      /// Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
  22:      /// November 2006
  23:      /// Purpose
  24:      /// =======
  25:      /// 
  26:      /// DGBTRS solves a system of linear equations
  27:      /// A * X = B  or  A' * X = B
  28:      /// with a general band matrix A using the LU factorization computed
  29:      /// by DGBTRF.
  30:      /// 
  31:      ///</summary>
  32:      public class DGBTRS
  33:      {
  34:      
  35:   
  36:          #region Dependencies
  37:          
  38:          LSAME _lsame; DGEMV _dgemv; DGER _dger; DSWAP _dswap; DTBSV _dtbsv; XERBLA _xerbla; 
  39:   
  40:          #endregion
  41:   
  42:   
  43:          #region Fields
  44:          
  45:          const double ONE = 1.0E+0; bool LNOTI = false; bool NOTRAN = false; int I = 0; int J = 0; int KD = 0; int L = 0; 
  46:          int LM = 0;
  47:   
  48:          #endregion
  49:   
  50:          public DGBTRS(LSAME lsame, DGEMV dgemv, DGER dger, DSWAP dswap, DTBSV dtbsv, XERBLA xerbla)
  51:          {
  52:      
  53:   
  54:              #region Set Dependencies
  55:              
  56:              this._lsame = lsame; this._dgemv = dgemv; this._dger = dger; this._dswap = dswap; this._dtbsv = dtbsv; 
  57:              this._xerbla = xerbla;
  58:   
  59:              #endregion
  60:   
  61:          }
  62:      
  63:          public DGBTRS()
  64:          {
  65:      
  66:   
  67:              #region Dependencies (Initialization)
  68:              
  69:              LSAME lsame = new LSAME();
  70:              XERBLA xerbla = new XERBLA();
  71:              DSWAP dswap = new DSWAP();
  72:              DGEMV dgemv = new DGEMV(lsame, xerbla);
  73:              DGER dger = new DGER(xerbla);
  74:              DTBSV dtbsv = new DTBSV(lsame, xerbla);
  75:   
  76:              #endregion
  77:   
  78:   
  79:              #region Set Dependencies
  80:              
  81:              this._lsame = lsame; this._dgemv = dgemv; this._dger = dger; this._dswap = dswap; this._dtbsv = dtbsv; 
  82:              this._xerbla = xerbla;
  83:   
  84:              #endregion
  85:   
  86:          }
  87:          /// <summary>
  88:          /// Purpose
  89:          /// =======
  90:          /// 
  91:          /// DGBTRS solves a system of linear equations
  92:          /// A * X = B  or  A' * X = B
  93:          /// with a general band matrix A using the LU factorization computed
  94:          /// by DGBTRF.
  95:          /// 
  96:          ///</summary>
  97:          /// <param name="TRANS">
  98:          /// (input) CHARACTER*1
  99:          /// Specifies the form of the system of equations.
 100:          /// = 'N':  A * X = B  (No transpose)
 101:          /// = 'T':  A'* X = B  (Transpose)
 102:          /// = 'C':  A'* X = B  (Conjugate transpose = Transpose)
 103:          ///</param>
 104:          /// <param name="N">
 105:          /// (input) INTEGER
 106:          /// The order of the matrix A.  N .GE. 0.
 107:          ///</param>
 108:          /// <param name="KL">
 109:          /// (input) INTEGER
 110:          /// The number of subdiagonals within the band of A.  KL .GE. 0.
 111:          ///</param>
 112:          /// <param name="KU">
 113:          /// (input) INTEGER
 114:          /// The number of superdiagonals within the band of A.  KU .GE. 0.
 115:          ///</param>
 116:          /// <param name="NRHS">
 117:          /// (input) INTEGER
 118:          /// The number of right hand sides, i.e., the number of columns
 119:          /// of the matrix B.  NRHS .GE. 0.
 120:          ///</param>
 121:          /// <param name="AB">
 122:          /// (input) DOUBLE PRECISION array, dimension (LDAB,N)
 123:          /// Details of the LU factorization of the band matrix A, as
 124:          /// computed by DGBTRF.  U is stored as an upper triangular band
 125:          /// matrix with KL+KU superdiagonals in rows 1 to KL+KU+1, and
 126:          /// the multipliers used during the factorization are stored in
 127:          /// rows KL+KU+2 to 2*KL+KU+1.
 128:          ///</param>
 129:          /// <param name="LDAB">
 130:          /// (input) INTEGER
 131:          /// The leading dimension of the array AB.  LDAB .GE. 2*KL+KU+1.
 132:          ///</param>
 133:          /// <param name="IPIV">
 134:          /// (input) INTEGER array, dimension (N)
 135:          /// The pivot indices; for 1 .LE. i .LE. N, row i of the matrix was
 136:          /// interchanged with row IPIV(i).
 137:          ///</param>
 138:          /// <param name="B">
 139:          /// (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
 140:          /// On entry, the right hand side matrix B.
 141:          /// On exit, the solution matrix X.
 142:          ///</param>
 143:          /// <param name="LDB">
 144:          /// (input) INTEGER
 145:          /// The leading dimension of the array B.  LDB .GE. max(1,N).
 146:          ///</param>
 147:          /// <param name="INFO">
 148:          /// (output) INTEGER
 149:          /// = 0:  successful exit
 150:          /// .LT. 0: if INFO = -i, the i-th argument had an illegal value
 151:          ///</param>
 152:          public void Run(string TRANS, int N, int KL, int KU, int NRHS, double[] AB, int offset_ab
 153:                           , int LDAB, int[] IPIV, int offset_ipiv, ref double[] B, int offset_b, int LDB, ref int INFO)
 154:          {
 155:   
 156:              #region Array Index Correction
 157:              
 158:               int o_ab = -1 - LDAB + offset_ab;  int o_ipiv = -1 + offset_ipiv;  int o_b = -1 - LDB + offset_b; 
 159:   
 160:              #endregion
 161:   
 162:   
 163:              #region Strings
 164:              
 165:              TRANS = TRANS.Substring(0, 1);  
 166:   
 167:              #endregion
 168:   
 169:   
 170:              #region Prolog
 171:              
 172:              // *
 173:              // *  -- LAPACK routine (version 3.1) --
 174:              // *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
 175:              // *     November 2006
 176:              // *
 177:              // *     .. Scalar Arguments ..
 178:              // *     ..
 179:              // *     .. Array Arguments ..
 180:              // *     ..
 181:              // *
 182:              // *  Purpose
 183:              // *  =======
 184:              // *
 185:              // *  DGBTRS solves a system of linear equations
 186:              // *     A * X = B  or  A' * X = B
 187:              // *  with a general band matrix A using the LU factorization computed
 188:              // *  by DGBTRF.
 189:              // *
 190:              // *  Arguments
 191:              // *  =========
 192:              // *
 193:              // *  TRANS   (input) CHARACTER*1
 194:              // *          Specifies the form of the system of equations.
 195:              // *          = 'N':  A * X = B  (No transpose)
 196:              // *          = 'T':  A'* X = B  (Transpose)
 197:              // *          = 'C':  A'* X = B  (Conjugate transpose = Transpose)
 198:              // *
 199:              // *  N       (input) INTEGER
 200:              // *          The order of the matrix A.  N >= 0.
 201:              // *
 202:              // *  KL      (input) INTEGER
 203:              // *          The number of subdiagonals within the band of A.  KL >= 0.
 204:              // *
 205:              // *  KU      (input) INTEGER
 206:              // *          The number of superdiagonals within the band of A.  KU >= 0.
 207:              // *
 208:              // *  NRHS    (input) INTEGER
 209:              // *          The number of right hand sides, i.e., the number of columns
 210:              // *          of the matrix B.  NRHS >= 0.
 211:              // *
 212:              // *  AB      (input) DOUBLE PRECISION array, dimension (LDAB,N)
 213:              // *          Details of the LU factorization of the band matrix A, as
 214:              // *          computed by DGBTRF.  U is stored as an upper triangular band
 215:              // *          matrix with KL+KU superdiagonals in rows 1 to KL+KU+1, and
 216:              // *          the multipliers used during the factorization are stored in
 217:              // *          rows KL+KU+2 to 2*KL+KU+1.
 218:              // *
 219:              // *  LDAB    (input) INTEGER
 220:              // *          The leading dimension of the array AB.  LDAB >= 2*KL+KU+1.
 221:              // *
 222:              // *  IPIV    (input) INTEGER array, dimension (N)
 223:              // *          The pivot indices; for 1 <= i <= N, row i of the matrix was
 224:              // *          interchanged with row IPIV(i).
 225:              // *
 226:              // *  B       (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
 227:              // *          On entry, the right hand side matrix B.
 228:              // *          On exit, the solution matrix X.
 229:              // *
 230:              // *  LDB     (input) INTEGER
 231:              // *          The leading dimension of the array B.  LDB >= max(1,N).
 232:              // *
 233:              // *  INFO    (output) INTEGER
 234:              // *          = 0:  successful exit
 235:              // *          < 0: if INFO = -i, the i-th argument had an illegal value
 236:              // *
 237:              // *  =====================================================================
 238:              // *
 239:              // *     .. Parameters ..
 240:              // *     ..
 241:              // *     .. Local Scalars ..
 242:              // *     ..
 243:              // *     .. External Functions ..
 244:              // *     ..
 245:              // *     .. External Subroutines ..
 246:              // *     ..
 247:              // *     .. Intrinsic Functions ..
 248:              //      INTRINSIC          MAX, MIN;
 249:              // *     ..
 250:              // *     .. Executable Statements ..
 251:              // *
 252:              // *     Test the input parameters.
 253:              // *
 254:   
 255:              #endregion
 256:   
 257:   
 258:              #region Body
 259:              
 260:              INFO = 0;
 261:              NOTRAN = this._lsame.Run(TRANS, "N");
 262:              if (!NOTRAN && !this._lsame.Run(TRANS, "T") && !this._lsame.Run(TRANS, "C"))
 263:              {
 264:                  INFO =  - 1;
 265:              }
 266:              else
 267:              {
 268:                  if (N < 0)
 269:                  {
 270:                      INFO =  - 2;
 271:                  }
 272:                  else
 273:                  {
 274:                      if (KL < 0)
 275:                      {
 276:                          INFO =  - 3;
 277:                      }
 278:                      else
 279:                      {
 280:                          if (KU < 0)
 281:                          {
 282:                              INFO =  - 4;
 283:                          }
 284:                          else
 285:                          {
 286:                              if (NRHS < 0)
 287:                              {
 288:                                  INFO =  - 5;
 289:                              }
 290:                              else
 291:                              {
 292:                                  if (LDAB < (2 * KL + KU + 1))
 293:                                  {
 294:                                      INFO =  - 7;
 295:                                  }
 296:                                  else
 297:                                  {
 298:                                      if (LDB < Math.Max(1, N))
 299:                                      {
 300:                                          INFO =  - 10;
 301:                                      }
 302:                                  }
 303:                              }
 304:                          }
 305:                      }
 306:                  }
 307:              }
 308:              if (INFO != 0)
 309:              {
 310:                  this._xerbla.Run("DGBTRS",  - INFO);
 311:                  return;
 312:              }
 313:              // *
 314:              // *     Quick return if possible
 315:              // *
 316:              if (N == 0 || NRHS == 0) return;
 317:              // *
 318:              KD = KU + KL + 1;
 319:              LNOTI = KL > 0;
 320:              // *
 321:              if (NOTRAN)
 322:              {
 323:                  // *
 324:                  // *        Solve  A*X = B.
 325:                  // *
 326:                  // *        Solve L*X = B, overwriting B with X.
 327:                  // *
 328:                  // *        L is represented as a product of permutations and unit lower
 329:                  // *        triangular matrices L = P(1) * L(1) * ... * P(n-1) * L(n-1),
 330:                  // *        where each transformation L(i) is a rank-one modification of
 331:                  // *        the identity matrix.
 332:                  // *
 333:                  if (LNOTI)
 334:                  {
 335:                      for (J = 1; J <= N - 1; J++)
 336:                      {
 337:                          LM = Math.Min(KL, N - J);
 338:                          L = IPIV[J + o_ipiv];
 339:                          if (L != J) this._dswap.Run(NRHS, ref B, L+1 * LDB + o_b, LDB, ref B, J+1 * LDB + o_b, LDB);
 340:                          this._dger.Run(LM, NRHS,  - ONE, AB, KD + 1+J * LDAB + o_ab, 1, B, J+1 * LDB + o_b
 341:                                         , LDB, ref B, J + 1+1 * LDB + o_b, LDB);
 342:                      }
 343:                  }
 344:                  // *
 345:                  for (I = 1; I <= NRHS; I++)
 346:                  {
 347:                      // *
 348:                      // *           Solve U*X = B, overwriting B with X.
 349:                      // *
 350:                      this._dtbsv.Run("Upper", "No transpose", "Non-unit", N, KL + KU, AB, offset_ab
 351:                                      , LDAB, ref B, 1+I * LDB + o_b, 1);
 352:                  }
 353:                  // *
 354:              }
 355:              else
 356:              {
 357:                  // *
 358:                  // *        Solve A'*X = B.
 359:                  // *
 360:                  for (I = 1; I <= NRHS; I++)
 361:                  {
 362:                      // *
 363:                      // *           Solve U'*X = B, overwriting B with X.
 364:                      // *
 365:                      this._dtbsv.Run("Upper", "Transpose", "Non-unit", N, KL + KU, AB, offset_ab
 366:                                      , LDAB, ref B, 1+I * LDB + o_b, 1);
 367:                  }
 368:                  // *
 369:                  // *        Solve L'*X = B, overwriting B with X.
 370:                  // *
 371:                  if (LNOTI)
 372:                  {
 373:                      for (J = N - 1; J >= 1; J +=  - 1)
 374:                      {
 375:                          LM = Math.Min(KL, N - J);
 376:                          this._dgemv.Run("Transpose", LM, NRHS,  - ONE, B, J + 1+1 * LDB + o_b, LDB
 377:                                          , AB, KD + 1+J * LDAB + o_ab, 1, ONE, ref B, J+1 * LDB + o_b, LDB);
 378:                          L = IPIV[J + o_ipiv];
 379:                          if (L != J) this._dswap.Run(NRHS, ref B, L+1 * LDB + o_b, LDB, ref B, J+1 * LDB + o_b, LDB);
 380:                      }
 381:                  }
 382:              }
 383:              return;
 384:              // *
 385:              // *     End of DGBTRS
 386:              // *
 387:   
 388:              #endregion
 389:   
 390:          }
 391:      }
 392:  }