Changeset 2563 for trunk/sources/ALGLIB/bdsvd.cs
- Timestamp:
- 12/17/09 17:05:22 (14 years ago)
- File:
-
- 1 edited
-
trunk/sources/ALGLIB/bdsvd.cs (modified) (34 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/sources/ALGLIB/bdsvd.cs
r2430 r2563 154 154 155 155 156 /*************************************************************************157 Obsolete 1-based subroutine. See RMatrixBDSVD for 0-based replacement.158 159 History:160 * 31 March, 2007.161 changed MAXITR from 6 to 12.162 163 -- LAPACK routine (version 3.0) --164 Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,165 Courant Institute, Argonne National Lab, and Rice University166 October 31, 1999.167 *************************************************************************/168 156 public static bool bidiagonalsvddecomposition(ref double[] d, 169 157 double[] e, … … 277 265 if( n==1 ) 278 266 { 279 if( d[1]<0)267 if( (double)(d[1])<(double)(0) ) 280 268 { 281 269 d[1] = -d[1]; … … 384 372 } 385 373 sminl = 0; 386 if( tol>=0)374 if( (double)(tol)>=(double)(0) ) 387 375 { 388 376 … … 391 379 // 392 380 sminoa = Math.Abs(d[1]); 393 if( sminoa!=0)381 if( (double)(sminoa)!=(double)(0) ) 394 382 { 395 383 mu = sminoa; … … 398 386 mu = Math.Abs(d[i])*(mu/(mu+Math.Abs(e[i-1]))); 399 387 sminoa = Math.Min(sminoa, mu); 400 if( sminoa==0)388 if( (double)(sminoa)==(double)(0) ) 401 389 { 402 390 break; … … 453 441 // Find diagonal block of matrix to work on 454 442 // 455 if( tol<0 & Math.Abs(d[m])<=thresh)443 if( (double)(tol)<(double)(0) & (double)(Math.Abs(d[m]))<=(double)(thresh) ) 456 444 { 457 445 d[m] = 0; … … 465 453 abss = Math.Abs(d[ll]); 466 454 abse = Math.Abs(e[ll]); 467 if( tol<0 & abss<=thresh)455 if( (double)(tol)<(double)(0) & (double)(abss)<=(double)(thresh) ) 468 456 { 469 457 d[ll] = 0; 470 458 } 471 if( abse<=thresh)459 if( (double)(abse)<=(double)(thresh) ) 472 460 { 473 461 matrixsplitflag = true; … … 606 594 // 607 595 bchangedir = false; 608 if( idir==1 & Math.Abs(d[ll])<1.0E-3*Math.Abs(d[m]) )596 if( idir==1 & (double)(Math.Abs(d[ll]))<(double)(1.0E-3*Math.Abs(d[m])) ) 609 597 { 610 598 bchangedir = true; 611 599 } 612 if( idir==2 & Math.Abs(d[m])<1.0E-3*Math.Abs(d[ll]) )600 if( idir==2 & (double)(Math.Abs(d[m]))<(double)(1.0E-3*Math.Abs(d[ll])) ) 613 601 { 614 602 bchangedir = true; … … 616 604 if( ll!=oldll | m!=oldm | bchangedir ) 617 605 { 618 if( Math.Abs(d[ll])>=Math.Abs(d[m]) )606 if( (double)(Math.Abs(d[ll]))>=(double)(Math.Abs(d[m])) ) 619 607 { 620 608 … … 644 632 // First apply standard test to bottom of matrix 645 633 // 646 if( Math.Abs(e[m-1])<=Math.Abs(tol)*Math.Abs(d[m]) | tol<0 & Math.Abs(e[m-1])<=thresh)634 if( (double)(Math.Abs(e[m-1]))<=(double)(Math.Abs(tol)*Math.Abs(d[m])) | (double)(tol)<(double)(0) & (double)(Math.Abs(e[m-1]))<=(double)(thresh) ) 647 635 { 648 636 e[m-1] = 0; 649 637 continue; 650 638 } 651 if( tol>=0)639 if( (double)(tol)>=(double)(0) ) 652 640 { 653 641 … … 661 649 for(lll=ll; lll<=m-1; lll++) 662 650 { 663 if( Math.Abs(e[lll])<=tol*mu)651 if( (double)(Math.Abs(e[lll]))<=(double)(tol*mu) ) 664 652 { 665 653 e[lll] = 0; … … 684 672 // First apply standard test to top of matrix 685 673 // 686 if( Math.Abs(e[ll])<=Math.Abs(tol)*Math.Abs(d[ll]) | tol<0 & Math.Abs(e[ll])<=thresh)674 if( (double)(Math.Abs(e[ll]))<=(double)(Math.Abs(tol)*Math.Abs(d[ll])) | (double)(tol)<(double)(0) & (double)(Math.Abs(e[ll]))<=(double)(thresh) ) 687 675 { 688 676 e[ll] = 0; 689 677 continue; 690 678 } 691 if( tol>=0)679 if( (double)(tol)>=(double)(0) ) 692 680 { 693 681 … … 701 689 for(lll=m-1; lll>=ll; lll--) 702 690 { 703 if( Math.Abs(e[lll])<=tol*mu)691 if( (double)(Math.Abs(e[lll]))<=(double)(tol*mu) ) 704 692 { 705 693 e[lll] = 0; … … 724 712 // accuracy, and if so set the shift to zero. 725 713 // 726 if( tol>=0 & n*tol*(sminl/smax)<=Math.Max(eps, 0.01*tol) )714 if( (double)(tol)>=(double)(0) & (double)(n*tol*(sminl/smax))<=(double)(Math.Max(eps, 0.01*tol)) ) 727 715 { 728 716 … … 752 740 // Test if shift negligible, and if so set to zero 753 741 // 754 if( sll>0)755 { 756 if( AP.Math.Sqr(shift/sll)<eps)742 if( (double)(sll)>(double)(0) ) 743 { 744 if( (double)(AP.Math.Sqr(shift/sll))<(double)(eps) ) 757 745 { 758 746 shift = 0; … … 769 757 // If SHIFT = 0, do simplified QR iteration 770 758 // 771 if( shift==0)759 if( (double)(shift)==(double)(0) ) 772 760 { 773 761 if( idir==1 ) … … 817 805 // Test convergence 818 806 // 819 if( Math.Abs(e[m-1])<=thresh)807 if( (double)(Math.Abs(e[m-1]))<=(double)(thresh) ) 820 808 { 821 809 e[m-1] = 0; … … 868 856 // Test convergence 869 857 // 870 if( Math.Abs(e[ll])<=thresh)858 if( (double)(Math.Abs(e[ll]))<=(double)(thresh) ) 871 859 { 872 860 e[ll] = 0; … … 935 923 // Test convergence 936 924 // 937 if( Math.Abs(e[m-1])<=thresh)925 if( (double)(Math.Abs(e[m-1]))<=(double)(thresh) ) 938 926 { 939 927 e[m-1] = 0; … … 979 967 // Test convergence 980 968 // 981 if( Math.Abs(e[ll])<=thresh)969 if( (double)(Math.Abs(e[ll]))<=(double)(thresh) ) 982 970 { 983 971 e[ll] = 0; … … 1013 1001 for(i=1; i<=n; i++) 1014 1002 { 1015 if( d[i]<0)1003 if( (double)(d[i])<(double)(0) ) 1016 1004 { 1017 1005 d[i] = -d[i]; … … 1044 1032 for(j=2; j<=n+1-i; j++) 1045 1033 { 1046 if( d[j]<=smin)1034 if( (double)(d[j])<=(double)(smin) ) 1047 1035 { 1048 1036 isub = j; … … 1117 1105 double result = 0; 1118 1106 1119 if( b>=0)1107 if( (double)(b)>=(double)(0) ) 1120 1108 { 1121 1109 result = Math.Abs(a); … … 1150 1138 fhmn = Math.Min(fa, ha); 1151 1139 fhmx = Math.Max(fa, ha); 1152 if( fhmn==0)1140 if( (double)(fhmn)==(double)(0) ) 1153 1141 { 1154 1142 ssmin = 0; 1155 if( fhmx==0)1143 if( (double)(fhmx)==(double)(0) ) 1156 1144 { 1157 1145 ssmax = ga; … … 1164 1152 else 1165 1153 { 1166 if( ga<fhmx)1154 if( (double)(ga)<(double)(fhmx) ) 1167 1155 { 1168 1156 aas = 1+fhmn/fhmx; … … 1176 1164 { 1177 1165 au = fhmx/ga; 1178 if( au==0)1166 if( (double)(au)==(double)(0) ) 1179 1167 { 1180 1168 … … 1249 1237 // 1250 1238 pmax = 1; 1251 swp = ha>fa;1239 swp = (double)(ha)>(double)(fa); 1252 1240 if( swp ) 1253 1241 { … … 1266 1254 gt = g; 1267 1255 ga = Math.Abs(gt); 1268 if( ga==0)1256 if( (double)(ga)==(double)(0) ) 1269 1257 { 1270 1258 … … 1282 1270 { 1283 1271 gasmal = true; 1284 if( ga>fa)1272 if( (double)(ga)>(double)(fa) ) 1285 1273 { 1286 1274 pmax = 2; 1287 if( fa/ga<AP.Math.MachineEpsilon)1275 if( (double)(fa/ga)<(double)(AP.Math.MachineEpsilon) ) 1288 1276 { 1289 1277 … … 1293 1281 gasmal = false; 1294 1282 ssmax = ga; 1295 if( ha>1)1283 if( (double)(ha)>(double)(1) ) 1296 1284 { 1297 1285 v = ga/ha; … … 1316 1304 // 1317 1305 d = fa-ha; 1318 if( d==fa)1306 if( (double)(d)==(double)(fa) ) 1319 1307 { 1320 1308 l = 1; … … 1329 1317 tt = t*t; 1330 1318 s = Math.Sqrt(tt+mm); 1331 if( l==0)1319 if( (double)(l)==(double)(0) ) 1332 1320 { 1333 1321 r = Math.Abs(m); … … 1340 1328 ssmin = ha/a; 1341 1329 ssmax = fa*a; 1342 if( mm==0)1330 if( (double)(mm)==(double)(0) ) 1343 1331 { 1344 1332 … … 1346 1334 // Note that M is very tiny 1347 1335 // 1348 if( l==0)1336 if( (double)(l)==(double)(0) ) 1349 1337 { 1350 1338 t = extsignbdsqr(2, ft)*extsignbdsqr(1, gt);
Note: See TracChangeset
for help on using the changeset viewer.
