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: /// DLASD6 computes the SVD of an updated upper bidiagonal matrix B
27: /// obtained by merging two smaller ones by appending a row. This
28: /// routine is used only for the problem which requires all singular
29: /// values and optionally singular vector matrices in factored form.
30: /// B is an N-by-M matrix with N = NL + NR + 1 and M = N + SQRE.
31: /// A related subroutine, DLASD1, handles the case in which all singular
32: /// values and singular vectors of the bidiagonal matrix are desired.
33: ///
34: /// DLASD6 computes the SVD as follows:
35: ///
36: /// ( D1(in) 0 0 0 )
37: /// B = U(in) * ( Z1' a Z2' b ) * VT(in)
38: /// ( 0 0 D2(in) 0 )
39: ///
40: /// = U(out) * ( D(out) 0) * VT(out)
41: ///
42: /// where Z' = (Z1' a Z2' b) = u' VT', and u is a vector of dimension M
43: /// with ALPHA and BETA in the NL+1 and NL+2 th entries and zeros
44: /// elsewhere; and the entry b is empty if SQRE = 0.
45: ///
46: /// The singular values of B can be computed using D1, D2, the first
47: /// components of all the right singular vectors of the lower block, and
48: /// the last components of all the right singular vectors of the upper
49: /// block. These components are stored and updated in VF and VL,
50: /// respectively, in DLASD6. Hence U and VT are not explicitly
51: /// referenced.
52: ///
53: /// The singular values are stored in D. The algorithm consists of two
54: /// stages:
55: ///
56: /// The first stage consists of deflating the size of the problem
57: /// when there are multiple singular values or if there is a zero
58: /// in the Z vector. For each such occurence the dimension of the
59: /// secular equation problem is reduced by one. This stage is
60: /// performed by the routine DLASD7.
61: ///
62: /// The second stage consists of calculating the updated
63: /// singular values. This is done by finding the roots of the
64: /// secular equation via the routine DLASD4 (as called by DLASD8).
65: /// This routine also updates VF and VL and computes the distances
66: /// between the updated singular values and the old singular
67: /// values.
68: ///
69: /// DLASD6 is called from DLASDA.
70: ///
71: ///</summary>
72: public class DLASD6
73: {
74:
75:
76: #region Dependencies
77:
78: DCOPY _dcopy; DLAMRG _dlamrg; DLASCL _dlascl; DLASD7 _dlasd7; DLASD8 _dlasd8; XERBLA _xerbla;
79:
80: #endregion
81:
82:
83: #region Fields
84:
85: const double ONE = 1.0E+0; const double ZERO = 0.0E+0; int I = 0; int IDX = 0; int IDXC = 0; int IDXP = 0; int ISIGMA = 0;
86: int IVFW = 0;int IVLW = 0; int IW = 0; int M = 0; int N = 0; int N1 = 0; int N2 = 0; double ORGNRM = 0;
87:
88: #endregion
89:
90: public DLASD6(DCOPY dcopy, DLAMRG dlamrg, DLASCL dlascl, DLASD7 dlasd7, DLASD8 dlasd8, XERBLA xerbla)
91: {
92:
93:
94: #region Set Dependencies
95:
96: this._dcopy = dcopy; this._dlamrg = dlamrg; this._dlascl = dlascl; this._dlasd7 = dlasd7; this._dlasd8 = dlasd8;
97: this._xerbla = xerbla;
98:
99: #endregion
100:
101: }
102:
103: public DLASD6()
104: {
105:
106:
107: #region Dependencies (Initialization)
108:
109: DCOPY dcopy = new DCOPY();
110: DLAMRG dlamrg = new DLAMRG();
111: LSAME lsame = new LSAME();
112: DLAMC3 dlamc3 = new DLAMC3();
113: XERBLA xerbla = new XERBLA();
114: DROT drot = new DROT();
115: DLAPY2 dlapy2 = new DLAPY2();
116: DLASD5 dlasd5 = new DLASD5();
117: DDOT ddot = new DDOT();
118: DNRM2 dnrm2 = new DNRM2();
119: DLAMC1 dlamc1 = new DLAMC1(dlamc3);
120: DLAMC4 dlamc4 = new DLAMC4(dlamc3);
121: DLAMC5 dlamc5 = new DLAMC5(dlamc3);
122: DLAMC2 dlamc2 = new DLAMC2(dlamc3, dlamc1, dlamc4, dlamc5);
123: DLAMCH dlamch = new DLAMCH(lsame, dlamc2);
124: DLASCL dlascl = new DLASCL(lsame, dlamch, xerbla);
125: DLASD7 dlasd7 = new DLASD7(dcopy, dlamrg, drot, xerbla, dlamch, dlapy2);
126: DLAED6 dlaed6 = new DLAED6(dlamch);
127: DLASD4 dlasd4 = new DLASD4(dlaed6, dlasd5, dlamch);
128: DLASET dlaset = new DLASET(lsame);
129: DLASD8 dlasd8 = new DLASD8(dcopy, dlascl, dlasd4, dlaset, xerbla, ddot, dlamc3, dnrm2);
130:
131: #endregion
132:
133:
134: #region Set Dependencies
135:
136: this._dcopy = dcopy; this._dlamrg = dlamrg; this._dlascl = dlascl; this._dlasd7 = dlasd7; this._dlasd8 = dlasd8;
137: this._xerbla = xerbla;
138:
139: #endregion
140:
141: }
142: /// <summary>
143: /// Purpose
144: /// =======
145: ///
146: /// DLASD6 computes the SVD of an updated upper bidiagonal matrix B
147: /// obtained by merging two smaller ones by appending a row. This
148: /// routine is used only for the problem which requires all singular
149: /// values and optionally singular vector matrices in factored form.
150: /// B is an N-by-M matrix with N = NL + NR + 1 and M = N + SQRE.
151: /// A related subroutine, DLASD1, handles the case in which all singular
152: /// values and singular vectors of the bidiagonal matrix are desired.
153: ///
154: /// DLASD6 computes the SVD as follows:
155: ///
156: /// ( D1(in) 0 0 0 )
157: /// B = U(in) * ( Z1' a Z2' b ) * VT(in)
158: /// ( 0 0 D2(in) 0 )
159: ///
160: /// = U(out) * ( D(out) 0) * VT(out)
161: ///
162: /// where Z' = (Z1' a Z2' b) = u' VT', and u is a vector of dimension M
163: /// with ALPHA and BETA in the NL+1 and NL+2 th entries and zeros
164: /// elsewhere; and the entry b is empty if SQRE = 0.
165: ///
166: /// The singular values of B can be computed using D1, D2, the first
167: /// components of all the right singular vectors of the lower block, and
168: /// the last components of all the right singular vectors of the upper
169: /// block. These components are stored and updated in VF and VL,
170: /// respectively, in DLASD6. Hence U and VT are not explicitly
171: /// referenced.
172: ///
173: /// The singular values are stored in D. The algorithm consists of two
174: /// stages:
175: ///
176: /// The first stage consists of deflating the size of the problem
177: /// when there are multiple singular values or if there is a zero
178: /// in the Z vector. For each such occurence the dimension of the
179: /// secular equation problem is reduced by one. This stage is
180: /// performed by the routine DLASD7.
181: ///
182: /// The second stage consists of calculating the updated
183: /// singular values. This is done by finding the roots of the
184: /// secular equation via the routine DLASD4 (as called by DLASD8).
185: /// This routine also updates VF and VL and computes the distances
186: /// between the updated singular values and the old singular
187: /// values.
188: ///
189: /// DLASD6 is called from DLASDA.
190: ///
191: ///</summary>
192: /// <param name="ICOMPQ">
193: /// (input) INTEGER
194: /// Specifies whether singular vectors are to be computed in
195: /// factored form:
196: /// = 0: Compute singular values only.
197: /// = 1: Compute singular vectors in factored form as well.
198: ///</param>
199: /// <param name="NL">
200: /// (input) INTEGER
201: /// The row dimension of the upper block. NL .GE. 1.
202: ///</param>
203: /// <param name="NR">
204: /// (input) INTEGER
205: /// The row dimension of the lower block. NR .GE. 1.
206: ///</param>
207: /// <param name="SQRE">
208: /// (input) INTEGER
209: /// = 0: the lower block is an NR-by-NR square matrix.
210: /// = 1: the lower block is an NR-by-(NR+1) rectangular matrix.
211: ///
212: /// The bidiagonal matrix has row dimension N = NL + NR + 1,
213: /// and column dimension M = N + SQRE.
214: ///</param>
215: /// <param name="D">
216: /// (input/output) DOUBLE PRECISION array, dimension ( NL+NR+1 ).
217: /// On entry D(1:NL,1:NL) contains the singular values of the
218: /// upper block, and D(NL+2:N) contains the singular values
219: /// of the lower block. On exit D(1:N) contains the singular
220: /// values of the modified matrix.
221: ///</param>
222: /// <param name="VF">
223: /// (input/output) DOUBLE PRECISION array, dimension ( M )
224: /// On entry, VF(1:NL+1) contains the first components of all
225: /// right singular vectors of the upper block; and VF(NL+2:M)
226: /// contains the first components of all right singular vectors
227: /// of the lower block. On exit, VF contains the first components
228: /// of all right singular vectors of the bidiagonal matrix.
229: ///</param>
230: /// <param name="VL">
231: /// (input/output) DOUBLE PRECISION array, dimension ( M )
232: /// On entry, VL(1:NL+1) contains the last components of all
233: /// right singular vectors of the upper block; and VL(NL+2:M)
234: /// contains the last components of all right singular vectors of
235: /// the lower block. On exit, VL contains the last components of
236: /// all right singular vectors of the bidiagonal matrix.
237: ///</param>
238: /// <param name="ALPHA">
239: /// (input/output) DOUBLE PRECISION
240: /// Contains the diagonal element associated with the added row.
241: ///</param>
242: /// <param name="BETA">
243: /// (input/output) DOUBLE PRECISION
244: /// Contains the off-diagonal element associated with the added
245: /// row.
246: ///</param>
247: /// <param name="IDXQ">
248: /// (output) INTEGER array, dimension ( N )
249: /// This contains the permutation which will reintegrate the
250: /// subproblem just solved back into sorted order, i.e.
251: /// D( IDXQ( I = 1, N ) ) will be in ascending order.
252: ///</param>
253: /// <param name="PERM">
254: /// (output) INTEGER array, dimension ( N )
255: /// The permutations (from deflation and sorting) to be applied
256: /// to each block. Not referenced if ICOMPQ = 0.
257: ///</param>
258: /// <param name="GIVPTR">
259: /// (output) INTEGER
260: /// The number of Givens rotations which took place in this
261: /// subproblem. Not referenced if ICOMPQ = 0.
262: ///</param>
263: /// <param name="GIVCOL">
264: /// (output) INTEGER array, dimension ( LDGCOL, 2 )
265: /// Each pair of numbers indicates a pair of columns to take place
266: /// in a Givens rotation. Not referenced if ICOMPQ = 0.
267: ///</param>
268: /// <param name="LDGCOL">
269: /// (input) INTEGER
270: /// leading dimension of GIVCOL, must be at least N.
271: ///</param>
272: /// <param name="GIVNUM">
273: /// (output) DOUBLE PRECISION array, dimension ( LDGNUM, 2 )
274: /// Each number indicates the C or S value to be used in the
275: /// corresponding Givens rotation. Not referenced if ICOMPQ = 0.
276: ///</param>
277: /// <param name="LDGNUM">
278: /// (input) INTEGER
279: /// The leading dimension of GIVNUM and POLES, must be at least N.
280: ///</param>
281: /// <param name="POLES">
282: /// (output) DOUBLE PRECISION array, dimension ( LDGNUM, 2 )
283: /// On exit, POLES(1,*) is an array containing the new singular
284: /// values obtained from solving the secular equation, and
285: /// POLES(2,*) is an array containing the poles in the secular
286: /// equation. Not referenced if ICOMPQ = 0.
287: ///</param>
288: /// <param name="DIFL">
289: /// (output) DOUBLE PRECISION array, dimension ( N )
290: /// On exit, DIFL(I) is the distance between I-th updated
291: /// (undeflated) singular value and the I-th (undeflated) old
292: /// singular value.
293: ///</param>
294: /// <param name="DIFR">
295: /// (output) DOUBLE PRECISION array,
296: /// dimension ( LDGNUM, 2 ) if ICOMPQ = 1 and
297: /// dimension ( N ) if ICOMPQ = 0.
298: /// On exit, DIFR(I, 1) is the distance between I-th updated
299: /// (undeflated) singular value and the I+1-th (undeflated) old
300: /// singular value.
301: ///
302: /// If ICOMPQ = 1, DIFR(1:K,2) is an array containing the
303: /// normalizing factors for the right singular vector matrix.
304: ///
305: /// See DLASD8 for details on DIFL and DIFR.
306: ///</param>
307: /// <param name="Z">
308: /// (output) DOUBLE PRECISION array, dimension ( M )
309: /// The first elements of this array contain the components
310: /// of the deflation-adjusted updating row vector.
311: ///</param>
312: /// <param name="K">
313: /// (output) INTEGER
314: /// Contains the dimension of the non-deflated matrix,
315: /// This is the order of the related secular equation. 1 .LE. K .LE.N.
316: ///</param>
317: /// <param name="C">
318: /// (output) DOUBLE PRECISION
319: /// C contains garbage if SQRE =0 and the C-value of a Givens
320: /// rotation related to the right null space if SQRE = 1.
321: ///</param>
322: /// <param name="S">
323: /// (output) DOUBLE PRECISION
324: /// S contains garbage if SQRE =0 and the S-value of a Givens
325: /// rotation related to the right null space if SQRE = 1.
326: ///</param>
327: /// <param name="WORK">
328: /// (workspace) DOUBLE PRECISION array, dimension ( 4 * M )
329: ///</param>
330: /// <param name="IWORK">
331: /// (workspace) INTEGER array, dimension ( 3 * N )
332: ///</param>
333: /// <param name="INFO">
334: /// (output) INTEGER
335: /// = 0: successful exit.
336: /// .LT. 0: if INFO = -i, the i-th argument had an illegal value.
337: /// .GT. 0: if INFO = 1, an singular value did not converge
338: ///</param>
339: public void Run(int ICOMPQ, int NL, int NR, int SQRE, ref double[] D, int offset_d, ref double[] VF, int offset_vf
340: , ref double[] VL, int offset_vl, ref double ALPHA, ref double BETA, ref int[] IDXQ, int offset_idxq, ref int[] PERM, int offset_perm, ref int GIVPTR
341: , ref int[] GIVCOL, int offset_givcol, int LDGCOL, ref double[] GIVNUM, int offset_givnum, int LDGNUM, ref double[] POLES, int offset_poles, ref double[] DIFL, int offset_difl
342: , ref double[] DIFR, int offset_difr, ref double[] Z, int offset_z, ref int K, ref double C, ref double S, ref double[] WORK, int offset_work
343: , ref int[] IWORK, int offset_iwork, ref int INFO)
344: {
345:
346: #region Array Index Correction
347:
348: int o_d = -1 + offset_d; int o_vf = -1 + offset_vf; int o_vl = -1 + offset_vl; int o_idxq = -1 + offset_idxq;
349: int o_perm = -1 + offset_perm; int o_givcol = -1 - LDGCOL + offset_givcol;
350: int o_givnum = -1 - LDGNUM + offset_givnum; int o_poles = -1 - LDGNUM + offset_poles; int o_difl = -1 + offset_difl;
351: int o_difr = -1 + offset_difr; int o_z = -1 + offset_z; int o_work = -1 + offset_work;
352: int o_iwork = -1 + offset_iwork;
353:
354: #endregion
355:
356:
357: #region Prolog
358:
359: // *
360: // * -- LAPACK auxiliary routine (version 3.1) --
361: // * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
362: // * November 2006
363: // *
364: // * .. Scalar Arguments ..
365: // * ..
366: // * .. Array Arguments ..
367: // * ..
368: // *
369: // * Purpose
370: // * =======
371: // *
372: // * DLASD6 computes the SVD of an updated upper bidiagonal matrix B
373: // * obtained by merging two smaller ones by appending a row. This
374: // * routine is used only for the problem which requires all singular
375: // * values and optionally singular vector matrices in factored form.
376: // * B is an N-by-M matrix with N = NL + NR + 1 and M = N + SQRE.
377: // * A related subroutine, DLASD1, handles the case in which all singular
378: // * values and singular vectors of the bidiagonal matrix are desired.
379: // *
380: // * DLASD6 computes the SVD as follows:
381: // *
382: // * ( D1(in) 0 0 0 )
383: // * B = U(in) * ( Z1' a Z2' b ) * VT(in)
384: // * ( 0 0 D2(in) 0 )
385: // *
386: // * = U(out) * ( D(out) 0) * VT(out)
387: // *
388: // * where Z' = (Z1' a Z2' b) = u' VT', and u is a vector of dimension M
389: // * with ALPHA and BETA in the NL+1 and NL+2 th entries and zeros
390: // * elsewhere; and the entry b is empty if SQRE = 0.
391: // *
392: // * The singular values of B can be computed using D1, D2, the first
393: // * components of all the right singular vectors of the lower block, and
394: // * the last components of all the right singular vectors of the upper
395: // * block. These components are stored and updated in VF and VL,
396: // * respectively, in DLASD6. Hence U and VT are not explicitly
397: // * referenced.
398: // *
399: // * The singular values are stored in D. The algorithm consists of two
400: // * stages:
401: // *
402: // * The first stage consists of deflating the size of the problem
403: // * when there are multiple singular values or if there is a zero
404: // * in the Z vector. For each such occurence the dimension of the
405: // * secular equation problem is reduced by one. This stage is
406: // * performed by the routine DLASD7.
407: // *
408: // * The second stage consists of calculating the updated
409: // * singular values. This is done by finding the roots of the
410: // * secular equation via the routine DLASD4 (as called by DLASD8).
411: // * This routine also updates VF and VL and computes the distances
412: // * between the updated singular values and the old singular
413: // * values.
414: // *
415: // * DLASD6 is called from DLASDA.
416: // *
417: // * Arguments
418: // * =========
419: // *
420: // * ICOMPQ (input) INTEGER
421: // * Specifies whether singular vectors are to be computed in
422: // * factored form:
423: // * = 0: Compute singular values only.
424: // * = 1: Compute singular vectors in factored form as well.
425: // *
426: // * NL (input) INTEGER
427: // * The row dimension of the upper block. NL >= 1.
428: // *
429: // * NR (input) INTEGER
430: // * The row dimension of the lower block. NR >= 1.
431: // *
432: // * SQRE (input) INTEGER
433: // * = 0: the lower block is an NR-by-NR square matrix.
434: // * = 1: the lower block is an NR-by-(NR+1) rectangular matrix.
435: // *
436: // * The bidiagonal matrix has row dimension N = NL + NR + 1,
437: // * and column dimension M = N + SQRE.
438: // *
439: // * D (input/output) DOUBLE PRECISION array, dimension ( NL+NR+1 ).
440: // * On entry D(1:NL,1:NL) contains the singular values of the
441: // * upper block, and D(NL+2:N) contains the singular values
442: // * of the lower block. On exit D(1:N) contains the singular
443: // * values of the modified matrix.
444: // *
445: // * VF (input/output) DOUBLE PRECISION array, dimension ( M )
446: // * On entry, VF(1:NL+1) contains the first components of all
447: // * right singular vectors of the upper block; and VF(NL+2:M)
448: // * contains the first components of all right singular vectors
449: // * of the lower block. On exit, VF contains the first components
450: // * of all right singular vectors of the bidiagonal matrix.
451: // *
452: // * VL (input/output) DOUBLE PRECISION array, dimension ( M )
453: // * On entry, VL(1:NL+1) contains the last components of all
454: // * right singular vectors of the upper block; and VL(NL+2:M)
455: // * contains the last components of all right singular vectors of
456: // * the lower block. On exit, VL contains the last components of
457: // * all right singular vectors of the bidiagonal matrix.
458: // *
459: // * ALPHA (input/output) DOUBLE PRECISION
460: // * Contains the diagonal element associated with the added row.
461: // *
462: // * BETA (input/output) DOUBLE PRECISION
463: // * Contains the off-diagonal element associated with the added
464: // * row.
465: // *
466: // * IDXQ (output) INTEGER array, dimension ( N )
467: // * This contains the permutation which will reintegrate the
468: // * subproblem just solved back into sorted order, i.e.
469: // * D( IDXQ( I = 1, N ) ) will be in ascending order.
470: // *
471: // * PERM (output) INTEGER array, dimension ( N )
472: // * The permutations (from deflation and sorting) to be applied
473: // * to each block. Not referenced if ICOMPQ = 0.
474: // *
475: // * GIVPTR (output) INTEGER
476: // * The number of Givens rotations which took place in this
477: // * subproblem. Not referenced if ICOMPQ = 0.
478: // *
479: // * GIVCOL (output) INTEGER array, dimension ( LDGCOL, 2 )
480: // * Each pair of numbers indicates a pair of columns to take place
481: // * in a Givens rotation. Not referenced if ICOMPQ = 0.
482: // *
483: // * LDGCOL (input) INTEGER
484: // * leading dimension of GIVCOL, must be at least N.
485: // *
486: // * GIVNUM (output) DOUBLE PRECISION array, dimension ( LDGNUM, 2 )
487: // * Each number indicates the C or S value to be used in the
488: // * corresponding Givens rotation. Not referenced if ICOMPQ = 0.
489: // *
490: // * LDGNUM (input) INTEGER
491: // * The leading dimension of GIVNUM and POLES, must be at least N.
492: // *
493: // * POLES (output) DOUBLE PRECISION array, dimension ( LDGNUM, 2 )
494: // * On exit, POLES(1,*) is an array containing the new singular
495: // * values obtained from solving the secular equation, and
496: // * POLES(2,*) is an array containing the poles in the secular
497: // * equation. Not referenced if ICOMPQ = 0.
498: // *
499: // * DIFL (output) DOUBLE PRECISION array, dimension ( N )
500: // * On exit, DIFL(I) is the distance between I-th updated
501: // * (undeflated) singular value and the I-th (undeflated) old
502: // * singular value.
503: // *
504: // * DIFR (output) DOUBLE PRECISION array,
505: // * dimension ( LDGNUM, 2 ) if ICOMPQ = 1 and
506: // * dimension ( N ) if ICOMPQ = 0.
507: // * On exit, DIFR(I, 1) is the distance between I-th updated
508: // * (undeflated) singular value and the I+1-th (undeflated) old
509: // * singular value.
510: // *
511: // * If ICOMPQ = 1, DIFR(1:K,2) is an array containing the
512: // * normalizing factors for the right singular vector matrix.
513: // *
514: // * See DLASD8 for details on DIFL and DIFR.
515: // *
516: // * Z (output) DOUBLE PRECISION array, dimension ( M )
517: // * The first elements of this array contain the components
518: // * of the deflation-adjusted updating row vector.
519: // *
520: // * K (output) INTEGER
521: // * Contains the dimension of the non-deflated matrix,
522: // * This is the order of the related secular equation. 1 <= K <=N.
523: // *
524: // * C (output) DOUBLE PRECISION
525: // * C contains garbage if SQRE =0 and the C-value of a Givens
526: // * rotation related to the right null space if SQRE = 1.
527: // *
528: // * S (output) DOUBLE PRECISION
529: // * S contains garbage if SQRE =0 and the S-value of a Givens
530: // * rotation related to the right null space if SQRE = 1.
531: // *
532: // * WORK (workspace) DOUBLE PRECISION array, dimension ( 4 * M )
533: // *
534: // * IWORK (workspace) INTEGER array, dimension ( 3 * N )
535: // *
536: // * INFO (output) INTEGER
537: // * = 0: successful exit.
538: // * < 0: if INFO = -i, the i-th argument had an illegal value.
539: // * > 0: if INFO = 1, an singular value did not converge
540: // *
541: // * Further Details
542: // * ===============
543: // *
544: // * Based on contributions by
545: // * Ming Gu and Huan Ren, Computer Science Division, University of
546: // * California at Berkeley, USA
547: // *
548: // * =====================================================================
549: // *
550: // * .. Parameters ..
551: // * ..
552: // * .. Local Scalars ..
553: // * ..
554: // * .. External Subroutines ..
555: // * ..
556: // * .. Intrinsic Functions ..
557: // INTRINSIC ABS, MAX;
558: // * ..
559: // * .. Executable Statements ..
560: // *
561: // * Test the input parameters.
562: // *
563:
564: #endregion
565:
566:
567: #region Body
568:
569: INFO = 0;
570: N = NL + NR + 1;
571: M = N + SQRE;
572: // *
573: if ((ICOMPQ < 0) || (ICOMPQ > 1))
574: {
575: INFO = - 1;
576: }
577: else
578: {
579: if (NL < 1)
580: {
581: INFO = - 2;
582: }
583: else
584: {
585: if (NR < 1)
586: {
587: INFO = - 3;
588: }
589: else
590: {
591: if ((SQRE < 0) || (SQRE > 1))
592: {
593: INFO = - 4;
594: }
595: else
596: {
597: if (LDGCOL < N)
598: {
599: INFO = - 14;
600: }
601: else
602: {
603: if (LDGNUM < N)
604: {
605: INFO = - 16;
606: }
607: }
608: }
609: }
610: }
611: }
612: if (INFO != 0)
613: {
614: this._xerbla.Run("DLASD6", - INFO);
615: return;
616: }
617: // *
618: // * The following values are for bookkeeping purposes only. They are
619: // * integer pointers which indicate the portion of the workspace
620: // * used by a particular array in DLASD7 and DLASD8.
621: // *
622: ISIGMA = 1;
623: IW = ISIGMA + N;
624: IVFW = IW + M;
625: IVLW = IVFW + M;
626: // *
627: IDX = 1;
628: IDXC = IDX + N;
629: IDXP = IDXC + N;
630: // *
631: // * Scale.
632: // *
633: ORGNRM = Math.Max(Math.Abs(ALPHA), Math.Abs(BETA));
634: D[NL + 1 + o_d] = ZERO;
635: for (I = 1; I <= N; I++)
636: {
637: if (Math.Abs(D[I + o_d]) > ORGNRM)
638: {
639: ORGNRM = Math.Abs(D[I + o_d]);
640: }
641: }
642: this._dlascl.Run("G", 0, 0, ORGNRM, ONE, N
643: , 1, ref D, offset_d, N, ref INFO);
644: ALPHA = ALPHA / ORGNRM;
645: BETA = BETA / ORGNRM;
646: // *
647: // * Sort and Deflate singular values.
648: // *
649: this._dlasd7.Run(ICOMPQ, NL, NR, SQRE, ref K, ref D, offset_d
650: , ref Z, offset_z, ref WORK, IW + o_work, ref VF, offset_vf, ref WORK, IVFW + o_work, ref VL, offset_vl, ref WORK, IVLW + o_work
651: , ALPHA, BETA, ref WORK, ISIGMA + o_work, ref IWORK, IDX + o_iwork, ref IWORK, IDXP + o_iwork, ref IDXQ, offset_idxq
652: , ref PERM, offset_perm, ref GIVPTR, ref GIVCOL, offset_givcol, LDGCOL, ref GIVNUM, offset_givnum, LDGNUM
653: , ref C, ref S, ref INFO);
654: // *
655: // * Solve Secular Equation, compute DIFL, DIFR, and update VF, VL.
656: // *
657: this._dlasd8.Run(ICOMPQ, K, ref D, offset_d, ref Z, offset_z, ref VF, offset_vf, ref VL, offset_vl
658: , ref DIFL, offset_difl, ref DIFR, offset_difr, LDGNUM, ref WORK, ISIGMA + o_work, ref WORK, IW + o_work, ref INFO);
659: // *
660: // * Save the poles if ICOMPQ = 1.
661: // *
662: if (ICOMPQ == 1)
663: {
664: this._dcopy.Run(K, D, offset_d, 1, ref POLES, 1+1 * LDGNUM + o_poles, 1);
665: this._dcopy.Run(K, WORK, ISIGMA + o_work, 1, ref POLES, 1+2 * LDGNUM + o_poles, 1);
666: }
667: // *
668: // * Unscale.
669: // *
670: this._dlascl.Run("G", 0, 0, ONE, ORGNRM, N
671: , 1, ref D, offset_d, N, ref INFO);
672: // *
673: // * Prepare the IDXQ sorting permutation.
674: // *
675: N1 = K;
676: N2 = N - K;
677: this._dlamrg.Run(N1, N2, D, offset_d, 1, - 1, ref IDXQ, offset_idxq);
678: // *
679: return;
680: // *
681: // * End of DLASD6
682: // *
683:
684: #endregion
685:
686: }
687: }
688: }