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:      /// Purpose
21:      /// =======
22:      ///
23:      /// DTRMV  performs one of the matrix-vector operations
24:      ///
25:      /// x := A*x,   or   x := A'*x,
26:      ///
27:      /// where x is an n element vector and  A is an n by n unit, or non-unit,
28:      /// upper or lower triangular matrix.
29:      ///
30:      /// Parameters
31:      /// ==========
32:      ///
33:      /// UPLO   - CHARACTER*1.
34:      /// On entry, UPLO specifies whether the matrix is an upper or
35:      /// lower triangular matrix as follows:
36:      ///
37:      /// UPLO = 'U' or 'u'   A is an upper triangular matrix.
38:      ///
39:      /// UPLO = 'L' or 'l'   A is a lower triangular matrix.
40:      ///
41:      /// Unchanged on exit.
42:      ///
43:      /// TRANS  - CHARACTER*1.
44:      /// On entry, TRANS specifies the operation to be performed as
45:      /// follows:
46:      ///
47:      /// TRANS = 'N' or 'n'   x := A*x.
48:      ///
49:      /// TRANS = 'T' or 't'   x := A'*x.
50:      ///
51:      /// TRANS = 'C' or 'c'   x := A'*x.
52:      ///
53:      /// Unchanged on exit.
54:      ///
55:      /// DIAG   - CHARACTER*1.
56:      /// On entry, DIAG specifies whether or not A is unit
57:      /// triangular as follows:
58:      ///
59:      /// DIAG = 'U' or 'u'   A is assumed to be unit triangular.
60:      ///
61:      /// DIAG = 'N' or 'n'   A is not assumed to be unit
62:      /// triangular.
63:      ///
64:      /// Unchanged on exit.
65:      ///
66:      /// N      - INTEGER.
67:      /// On entry, N specifies the order of the matrix A.
68:      /// N must be at least zero.
69:      /// Unchanged on exit.
70:      ///
71:      /// A      - DOUBLE PRECISION array of DIMENSION ( LDA, n ).
72:      /// Before entry with  UPLO = 'U' or 'u', the leading n by n
73:      /// upper triangular part of the array A must contain the upper
74:      /// triangular matrix and the strictly lower triangular part of
75:      /// A is not referenced.
76:      /// Before entry with UPLO = 'L' or 'l', the leading n by n
77:      /// lower triangular part of the array A must contain the lower
78:      /// triangular matrix and the strictly upper triangular part of
79:      /// A is not referenced.
80:      /// Note that when  DIAG = 'U' or 'u', the diagonal elements of
81:      /// A are not referenced either, but are assumed to be unity.
82:      /// Unchanged on exit.
83:      ///
84:      /// LDA    - INTEGER.
85:      /// On entry, LDA specifies the first dimension of A as declared
86:      /// in the calling (sub) program. LDA must be at least
87:      /// max( 1, n ).
88:      /// Unchanged on exit.
89:      ///
90:      /// X      - DOUBLE PRECISION array of dimension at least
91:      /// ( 1 + ( n - 1 )*abs( INCX ) ).
92:      /// Before entry, the incremented array X must contain the n
93:      /// element vector x. On exit, X is overwritten with the
94:      /// tranformed vector x.
95:      ///
96:      /// INCX   - INTEGER.
97:      /// On entry, INCX specifies the increment for the elements of
98:      /// X. INCX must not be zero.
99:      /// Unchanged on exit.
100:      ///</summary>
101:      public class DTRMV
102:      {
103:
104:
105:          #region Dependencies
106:
107:          LSAME _lsame; XERBLA _xerbla;
108:
109:          #endregion
110:
111:
112:          #region Fields
113:
114:          const double ZERO = 0.0E+0; double TEMP = 0; int I = 0; int INFO = 0; int IX = 0; int J = 0; int JX = 0; int KX = 0;
115:          bool NOUNIT = false;
116:
117:          #endregion
118:
119:          public DTRMV(LSAME lsame, XERBLA xerbla)
120:          {
121:
122:
123:              #region Set Dependencies
124:
125:              this._lsame = lsame; this._xerbla = xerbla;
126:
127:              #endregion
128:
129:          }
130:
131:          public DTRMV()
132:          {
133:
134:
135:              #region Dependencies (Initialization)
136:
137:              LSAME lsame = new LSAME();
138:              XERBLA xerbla = new XERBLA();
139:
140:              #endregion
141:
142:
143:              #region Set Dependencies
144:
145:              this._lsame = lsame; this._xerbla = xerbla;
146:
147:              #endregion
148:
149:          }
150:          /// <summary>
151:          /// Purpose
152:          /// =======
153:          ///
154:          /// DTRMV  performs one of the matrix-vector operations
155:          ///
156:          /// x := A*x,   or   x := A'*x,
157:          ///
158:          /// where x is an n element vector and  A is an n by n unit, or non-unit,
159:          /// upper or lower triangular matrix.
160:          ///
161:          /// Parameters
162:          /// ==========
163:          ///
164:          /// UPLO   - CHARACTER*1.
165:          /// On entry, UPLO specifies whether the matrix is an upper or
166:          /// lower triangular matrix as follows:
167:          ///
168:          /// UPLO = 'U' or 'u'   A is an upper triangular matrix.
169:          ///
170:          /// UPLO = 'L' or 'l'   A is a lower triangular matrix.
171:          ///
172:          /// Unchanged on exit.
173:          ///
174:          /// TRANS  - CHARACTER*1.
175:          /// On entry, TRANS specifies the operation to be performed as
176:          /// follows:
177:          ///
178:          /// TRANS = 'N' or 'n'   x := A*x.
179:          ///
180:          /// TRANS = 'T' or 't'   x := A'*x.
181:          ///
182:          /// TRANS = 'C' or 'c'   x := A'*x.
183:          ///
184:          /// Unchanged on exit.
185:          ///
186:          /// DIAG   - CHARACTER*1.
187:          /// On entry, DIAG specifies whether or not A is unit
188:          /// triangular as follows:
189:          ///
190:          /// DIAG = 'U' or 'u'   A is assumed to be unit triangular.
191:          ///
192:          /// DIAG = 'N' or 'n'   A is not assumed to be unit
193:          /// triangular.
194:          ///
195:          /// Unchanged on exit.
196:          ///
197:          /// N      - INTEGER.
198:          /// On entry, N specifies the order of the matrix A.
199:          /// N must be at least zero.
200:          /// Unchanged on exit.
201:          ///
202:          /// A      - DOUBLE PRECISION array of DIMENSION ( LDA, n ).
203:          /// Before entry with  UPLO = 'U' or 'u', the leading n by n
204:          /// upper triangular part of the array A must contain the upper
205:          /// triangular matrix and the strictly lower triangular part of
206:          /// A is not referenced.
207:          /// Before entry with UPLO = 'L' or 'l', the leading n by n
208:          /// lower triangular part of the array A must contain the lower
209:          /// triangular matrix and the strictly upper triangular part of
210:          /// A is not referenced.
211:          /// Note that when  DIAG = 'U' or 'u', the diagonal elements of
212:          /// A are not referenced either, but are assumed to be unity.
213:          /// Unchanged on exit.
214:          ///
215:          /// LDA    - INTEGER.
216:          /// On entry, LDA specifies the first dimension of A as declared
217:          /// in the calling (sub) program. LDA must be at least
218:          /// max( 1, n ).
219:          /// Unchanged on exit.
220:          ///
221:          /// X      - DOUBLE PRECISION array of dimension at least
222:          /// ( 1 + ( n - 1 )*abs( INCX ) ).
223:          /// Before entry, the incremented array X must contain the n
224:          /// element vector x. On exit, X is overwritten with the
225:          /// tranformed vector x.
226:          ///
227:          /// INCX   - INTEGER.
228:          /// On entry, INCX specifies the increment for the elements of
229:          /// X. INCX must not be zero.
230:          /// Unchanged on exit.
231:          ///</summary>
232:          /// <param name="UPLO">
233:          /// - CHARACTER*1.
234:          /// On entry, UPLO specifies whether the matrix is an upper or
235:          /// lower triangular matrix as follows:
236:          ///
237:          /// UPLO = 'U' or 'u'   A is an upper triangular matrix.
238:          ///
239:          /// UPLO = 'L' or 'l'   A is a lower triangular matrix.
240:          ///
241:          /// Unchanged on exit.
242:          ///</param>
243:          /// <param name="TRANS">
244:          /// - CHARACTER*1.
245:          /// On entry, TRANS specifies the operation to be performed as
246:          /// follows:
247:          ///
248:          /// TRANS = 'N' or 'n'   x := A*x.
249:          ///
250:          /// TRANS = 'T' or 't'   x := A'*x.
251:          ///
252:          /// TRANS = 'C' or 'c'   x := A'*x.
253:          ///
254:          /// Unchanged on exit.
255:          ///</param>
256:          /// <param name="DIAG">
257:          /// - CHARACTER*1.
258:          /// On entry, DIAG specifies whether or not A is unit
259:          /// triangular as follows:
260:          ///
261:          /// DIAG = 'U' or 'u'   A is assumed to be unit triangular.
262:          ///
263:          /// DIAG = 'N' or 'n'   A is not assumed to be unit
264:          /// triangular.
265:          ///
266:          /// Unchanged on exit.
267:          ///</param>
268:          /// <param name="N">
269:          /// - INTEGER.
270:          /// On entry, N specifies the order of the matrix A.
271:          /// N must be at least zero.
272:          /// Unchanged on exit.
273:          ///</param>
274:          /// <param name="A">
275:          /// - DOUBLE PRECISION array of DIMENSION ( LDA, n ).
276:          /// Before entry with  UPLO = 'U' or 'u', the leading n by n
277:          /// upper triangular part of the array A must contain the upper
278:          /// triangular matrix and the strictly lower triangular part of
279:          /// A is not referenced.
280:          /// Before entry with UPLO = 'L' or 'l', the leading n by n
281:          /// lower triangular part of the array A must contain the lower
282:          /// triangular matrix and the strictly upper triangular part of
283:          /// A is not referenced.
284:          /// Note that when  DIAG = 'U' or 'u', the diagonal elements of
285:          /// A are not referenced either, but are assumed to be unity.
286:          /// Unchanged on exit.
287:          ///</param>
288:          /// <param name="LDA">
289:          /// - INTEGER.
290:          /// On entry, LDA specifies the first dimension of A as declared
291:          /// in the calling (sub) program. LDA must be at least
292:          /// max( 1, n ).
293:          /// Unchanged on exit.
294:          ///</param>
295:          /// <param name="X">
296:          /// := A*x,   or   x := A'*x,
297:          ///</param>
298:          /// <param name="INCX">
299:          /// - INTEGER.
300:          /// On entry, INCX specifies the increment for the elements of
301:          /// X. INCX must not be zero.
302:          /// Unchanged on exit.
303:          ///
304:          ///</param>
305:          public void Run(string UPLO, string TRANS, string DIAG, int N, double[] A, int offset_a, int LDA
306:                           , ref double[] X, int offset_x, int INCX)
307:          {
308:
309:              #region Array Index Correction
310:
311:               int o_a = -1 - LDA + offset_a;  int o_x = -1 + offset_x;
312:
313:              #endregion
314:
315:
316:              #region Strings
317:
318:              UPLO = UPLO.Substring(0, 1);  TRANS = TRANS.Substring(0, 1);  DIAG = DIAG.Substring(0, 1);
319:
320:              #endregion
321:
322:
323:              #region Prolog
324:
325:              // *     .. Scalar Arguments ..
326:              // *     .. Array Arguments ..
327:              // *     ..
328:              // *
329:              // *  Purpose
330:              // *  =======
331:              // *
332:              // *  DTRMV  performs one of the matrix-vector operations
333:              // *
334:              // *     x := A*x,   or   x := A'*x,
335:              // *
336:              // *  where x is an n element vector and  A is an n by n unit, or non-unit,
337:              // *  upper or lower triangular matrix.
338:              // *
339:              // *  Parameters
340:              // *  ==========
341:              // *
342:              // *  UPLO   - CHARACTER*1.
343:              // *           On entry, UPLO specifies whether the matrix is an upper or
344:              // *           lower triangular matrix as follows:
345:              // *
346:              // *              UPLO = 'U' or 'u'   A is an upper triangular matrix.
347:              // *
348:              // *              UPLO = 'L' or 'l'   A is a lower triangular matrix.
349:              // *
350:              // *           Unchanged on exit.
351:              // *
352:              // *  TRANS  - CHARACTER*1.
353:              // *           On entry, TRANS specifies the operation to be performed as
354:              // *           follows:
355:              // *
356:              // *              TRANS = 'N' or 'n'   x := A*x.
357:              // *
358:              // *              TRANS = 'T' or 't'   x := A'*x.
359:              // *
360:              // *              TRANS = 'C' or 'c'   x := A'*x.
361:              // *
362:              // *           Unchanged on exit.
363:              // *
364:              // *  DIAG   - CHARACTER*1.
365:              // *           On entry, DIAG specifies whether or not A is unit
366:              // *           triangular as follows:
367:              // *
368:              // *              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
369:              // *
370:              // *              DIAG = 'N' or 'n'   A is not assumed to be unit
371:              // *                                  triangular.
372:              // *
373:              // *           Unchanged on exit.
374:              // *
375:              // *  N      - INTEGER.
376:              // *           On entry, N specifies the order of the matrix A.
377:              // *           N must be at least zero.
378:              // *           Unchanged on exit.
379:              // *
380:              // *  A      - DOUBLE PRECISION array of DIMENSION ( LDA, n ).
381:              // *           Before entry with  UPLO = 'U' or 'u', the leading n by n
382:              // *           upper triangular part of the array A must contain the upper
383:              // *           triangular matrix and the strictly lower triangular part of
384:              // *           A is not referenced.
385:              // *           Before entry with UPLO = 'L' or 'l', the leading n by n
386:              // *           lower triangular part of the array A must contain the lower
387:              // *           triangular matrix and the strictly upper triangular part of
388:              // *           A is not referenced.
389:              // *           Note that when  DIAG = 'U' or 'u', the diagonal elements of
390:              // *           A are not referenced either, but are assumed to be unity.
391:              // *           Unchanged on exit.
392:              // *
393:              // *  LDA    - INTEGER.
394:              // *           On entry, LDA specifies the first dimension of A as declared
395:              // *           in the calling (sub) program. LDA must be at least
396:              // *           max( 1, n ).
397:              // *           Unchanged on exit.
398:              // *
399:              // *  X      - DOUBLE PRECISION array of dimension at least
400:              // *           ( 1 + ( n - 1 )*abs( INCX ) ).
401:              // *           Before entry, the incremented array X must contain the n
402:              // *           element vector x. On exit, X is overwritten with the
403:              // *           tranformed vector x.
404:              // *
405:              // *  INCX   - INTEGER.
406:              // *           On entry, INCX specifies the increment for the elements of
407:              // *           X. INCX must not be zero.
408:              // *           Unchanged on exit.
409:              // *
410:              // *
411:              // *  Level 2 Blas routine.
412:              // *
413:              // *  -- Written on 22-October-1986.
414:              // *     Jack Dongarra, Argonne National Lab.
415:              // *     Jeremy Du Croz, Nag Central Office.
416:              // *     Sven Hammarling, Nag Central Office.
417:              // *     Richard Hanson, Sandia National Labs.
418:              // *
419:              // *
420:              // *     .. Parameters ..
421:              // *     .. Local Scalars ..
422:              // *     .. External Functions ..
423:              // *     .. External Subroutines ..
424:              // *     .. Intrinsic Functions ..
425:              //      INTRINSIC          MAX;
426:              // *     ..
427:              // *     .. Executable Statements ..
428:              // *
429:              // *     Test the input parameters.
430:              // *
431:
432:              #endregion
433:
434:
435:              #region Body
436:
437:              INFO = 0;
438:              if (!this._lsame.Run(UPLO, "U") && !this._lsame.Run(UPLO, "L"))
439:              {
440:                  INFO = 1;
441:              }
442:              else
443:              {
444:                  if (!this._lsame.Run(TRANS, "N") && !this._lsame.Run(TRANS, "T") && !this._lsame.Run(TRANS, "C"))
445:                  {
446:                      INFO = 2;
447:                  }
448:                  else
449:                  {
450:                      if (!this._lsame.Run(DIAG, "U") && !this._lsame.Run(DIAG, "N"))
451:                      {
452:                          INFO = 3;
453:                      }
454:                      else
455:                      {
456:                          if (N < 0)
457:                          {
458:                              INFO = 4;
459:                          }
460:                          else
461:                          {
462:                              if (LDA < Math.Max(1, N))
463:                              {
464:                                  INFO = 6;
465:                              }
466:                              else
467:                              {
468:                                  if (INCX == 0)
469:                                  {
470:                                      INFO = 8;
471:                                  }
472:                              }
473:                          }
474:                      }
475:                  }
476:              }
477:              if (INFO != 0)
478:              {
479:                  this._xerbla.Run("DTRMV ", INFO);
480:                  return;
481:              }
482:              // *
483:              // *     Quick return if possible.
484:              // *
485:              if (N == 0) return;
486:              // *
487:              NOUNIT = this._lsame.Run(DIAG, "N");
488:              // *
489:              // *     Set up the start point in X if the increment is not unity. This
490:              // *     will be  ( N - 1 )*INCX  too small for descending loops.
491:              // *
492:              if (INCX <= 0)
493:              {
494:                  KX = 1 - (N - 1) * INCX;
495:              }
496:              else
497:              {
498:                  if (INCX != 1)
499:                  {
500:                      KX = 1;
501:                  }
502:              }
503:              // *
504:              // *     Start the operations. In this version the elements of A are
505:              // *     accessed sequentially with one pass through A.
506:              // *
507:              if (this._lsame.Run(TRANS, "N"))
508:              {
509:                  // *
510:                  // *        Form  x := A*x.
511:                  // *
512:                  if (this._lsame.Run(UPLO, "U"))
513:                  {
514:                      if (INCX == 1)
515:                      {
516:                          for (J = 1; J <= N; J++)
517:                          {
518:                              if (X[J + o_x] != ZERO)
519:                              {
520:                                  TEMP = X[J + o_x];
521:                                  for (I = 1; I <= J - 1; I++)
522:                                  {
523:                                      X[I + o_x] = X[I + o_x] + TEMP * A[I+J * LDA + o_a];
524:                                  }
525:                                  if (NOUNIT) X[J + o_x] = X[J + o_x] * A[J+J * LDA + o_a];
526:                              }
527:                          }
528:                      }
529:                      else
530:                      {
531:                          JX = KX;
532:                          for (J = 1; J <= N; J++)
533:                          {
534:                              if (X[JX + o_x] != ZERO)
535:                              {
536:                                  TEMP = X[JX + o_x];
537:                                  IX = KX;
538:                                  for (I = 1; I <= J - 1; I++)
539:                                  {
540:                                      X[IX + o_x] = X[IX + o_x] + TEMP * A[I+J * LDA + o_a];
541:                                      IX = IX + INCX;
542:                                  }
543:                                  if (NOUNIT) X[JX + o_x] = X[JX + o_x] * A[J+J * LDA + o_a];
544:                              }
545:                              JX = JX + INCX;
546:                          }
547:                      }
548:                  }
549:                  else
550:                  {
551:                      if (INCX == 1)
552:                      {
553:                          for (J = N; J >= 1; J +=  - 1)
554:                          {
555:                              if (X[J + o_x] != ZERO)
556:                              {
557:                                  TEMP = X[J + o_x];
558:                                  for (I = N; I >= J + 1; I +=  - 1)
559:                                  {
560:                                      X[I + o_x] = X[I + o_x] + TEMP * A[I+J * LDA + o_a];
561:                                  }
562:                                  if (NOUNIT) X[J + o_x] = X[J + o_x] * A[J+J * LDA + o_a];
563:                              }
564:                          }
565:                      }
566:                      else
567:                      {
568:                          KX = KX + (N - 1) * INCX;
569:                          JX = KX;
570:                          for (J = N; J >= 1; J +=  - 1)
571:                          {
572:                              if (X[JX + o_x] != ZERO)
573:                              {
574:                                  TEMP = X[JX + o_x];
575:                                  IX = KX;
576:                                  for (I = N; I >= J + 1; I +=  - 1)
577:                                  {
578:                                      X[IX + o_x] = X[IX + o_x] + TEMP * A[I+J * LDA + o_a];
579:                                      IX = IX - INCX;
580:                                  }
581:                                  if (NOUNIT) X[JX + o_x] = X[JX + o_x] * A[J+J * LDA + o_a];
582:                              }
583:                              JX = JX - INCX;
584:                          }
585:                      }
586:                  }
587:              }
588:              else
589:              {
590:                  // *
591:                  // *        Form  x := A'*x.
592:                  // *
593:                  if (this._lsame.Run(UPLO, "U"))
594:                  {
595:                      if (INCX == 1)
596:                      {
597:                          for (J = N; J >= 1; J +=  - 1)
598:                          {
599:                              TEMP = X[J + o_x];
600:                              if (NOUNIT) TEMP = TEMP * A[J+J * LDA + o_a];
601:                              for (I = J - 1; I >= 1; I +=  - 1)
602:                              {
603:                                  TEMP = TEMP + A[I+J * LDA + o_a] * X[I + o_x];
604:                              }
605:                              X[J + o_x] = TEMP;
606:                          }
607:                      }
608:                      else
609:                      {
610:                          JX = KX + (N - 1) * INCX;
611:                          for (J = N; J >= 1; J +=  - 1)
612:                          {
613:                              TEMP = X[JX + o_x];
614:                              IX = JX;
615:                              if (NOUNIT) TEMP = TEMP * A[J+J * LDA + o_a];
616:                              for (I = J - 1; I >= 1; I +=  - 1)
617:                              {
618:                                  IX = IX - INCX;
619:                                  TEMP = TEMP + A[I+J * LDA + o_a] * X[IX + o_x];
620:                              }
621:                              X[JX + o_x] = TEMP;
622:                              JX = JX - INCX;
623:                          }
624:                      }
625:                  }
626:                  else
627:                  {
628:                      if (INCX == 1)
629:                      {
630:                          for (J = 1; J <= N; J++)
631:                          {
632:                              TEMP = X[J + o_x];
633:                              if (NOUNIT) TEMP = TEMP * A[J+J * LDA + o_a];
634:                              for (I = J + 1; I <= N; I++)
635:                              {
636:                                  TEMP = TEMP + A[I+J * LDA + o_a] * X[I + o_x];
637:                              }
638:                              X[J + o_x] = TEMP;
639:                          }
640:                      }
641:                      else
642:                      {
643:                          JX = KX;
644:                          for (J = 1; J <= N; J++)
645:                          {
646:                              TEMP = X[JX + o_x];
647:                              IX = JX;
648:                              if (NOUNIT) TEMP = TEMP * A[J+J * LDA + o_a];
649:                              for (I = J + 1; I <= N; I++)
650:                              {
651:                                  IX = IX + INCX;
652:                                  TEMP = TEMP + A[I+J * LDA + o_a] * X[IX + o_x];
653:                              }
654:                              X[JX + o_x] = TEMP;
655:                              JX = JX + INCX;
656:                          }
657:                      }
658:                  }
659:              }
660:              // *
661:              return;
662:              // *
663:              // *     End of DTRMV .
664:              // *
665:
666:              #endregion
667:
668:          }
669:      }
670:  }