Changeset 2563
- Timestamp:
- 12/17/09 17:05:22 (15 years ago)
- Location:
- trunk/sources/ALGLIB
- Files:
-
- 105 added
- 2 deleted
- 17 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/sources/ALGLIB/ALGLIB.csproj
r2445 r2563 80 80 </ItemGroup> 81 81 <ItemGroup> 82 <Compile Include="airyf.cs" /> 82 83 <Compile Include="ap.cs" /> 84 <Compile Include="autogk.cs" /> 85 <Compile Include="bdss.cs" /> 83 86 <Compile Include="bdsvd.cs" /> 87 <Compile Include="bessel.cs" /> 88 <Compile Include="betaf.cs" /> 84 89 <Compile Include="bidiagonal.cs" /> 90 <Compile Include="binomialdistr.cs" /> 85 91 <Compile Include="blas.cs" /> 92 <Compile Include="cblas.cs" /> 93 <Compile Include="cdet.cs" /> 94 <Compile Include="chebyshev.cs" /> 95 <Compile Include="chisquaredistr.cs" /> 96 <Compile Include="cholesky.cs" /> 97 <Compile Include="cinverse.cs" /> 98 <Compile Include="clu.cs" /> 99 <Compile Include="conv.cs" /> 100 <Compile Include="corr.cs" /> 86 101 <Compile Include="correlation.cs" /> 87 102 <Compile Include="AlglibPlugin.cs"> 88 103 <CopyToOutputDirectory>Always</CopyToOutputDirectory> 89 104 </Compile> 105 <Compile Include="correlationtests.cs" /> 106 <Compile Include="cqr.cs" /> 107 <Compile Include="crcond.cs" /> 108 <Compile Include="creflections.cs" /> 109 <Compile Include="csolve.cs" /> 110 <Compile Include="ctrinverse.cs" /> 111 <Compile Include="ctrlinsolve.cs" /> 112 <Compile Include="dawson.cs" /> 113 <Compile Include="densesolver.cs" /> 90 114 <Compile Include="descriptivestatistics.cs" /> 91 <Compile Include="gammaf.cs" /> 115 <Compile Include="det.cs" /> 116 <Compile Include="dforest.cs" /> 117 <Compile Include="elliptic.cs" /> 118 <Compile Include="estnorm.cs" /> 119 <Compile Include="expintegrals.cs" /> 120 <Compile Include="fdistr.cs" /> 121 <Compile Include="fft.cs" /> 122 <Compile Include="fht.cs" /> 123 <Compile Include="fresnel.cs" /> 124 <Compile Include="ftbase.cs" /> 125 <Compile Include="gammafunc.cs" /> 126 <Compile Include="gkq.cs" /> 127 <Compile Include="gq.cs" /> 128 <Compile Include="hbisinv.cs" /> 129 <Compile Include="hblas.cs" /> 130 <Compile Include="hcholesky.cs" /> 131 <Compile Include="hermite.cs" /> 132 <Compile Include="hessenberg.cs" /> 133 <Compile Include="hevd.cs" /> 134 <Compile Include="hqrnd.cs" /> 135 <Compile Include="hsschur.cs" /> 136 <Compile Include="htridiagonal.cs" /> 137 <Compile Include="ibetaf.cs" /> 92 138 <Compile Include="igammaf.cs" /> 139 <Compile Include="inv.cs" /> 140 <Compile Include="inverseupdate.cs" /> 141 <Compile Include="jacobianelliptic.cs" /> 142 <Compile Include="jarquebera.cs" /> 143 <Compile Include="kmeans.cs" /> 144 <Compile Include="laguerre.cs" /> 145 <Compile Include="lbfgs.cs" /> 146 <Compile Include="lda.cs" /> 147 <Compile Include="ldlt.cs" /> 93 148 <Compile Include="leastsquares.cs" /> 149 <Compile Include="legendre.cs" /> 94 150 <Compile Include="linreg.cs" /> 151 <Compile Include="logit.cs" /> 95 152 <Compile Include="lq.cs" /> 153 <Compile Include="lsfit.cs" /> 154 <Compile Include="lu.cs" /> 155 <Compile Include="mannwhitneyu.cs" /> 156 <Compile Include="matgen.cs" /> 157 <Compile Include="minlm.cs" /> 158 <Compile Include="mlpbase.cs" /> 159 <Compile Include="mlpe.cs" /> 160 <Compile Include="mlptrain.cs" /> 161 <Compile Include="nearunityunit.cs" /> 96 162 <Compile Include="normaldistr.cs" /> 163 <Compile Include="nsevd.cs" /> 164 <Compile Include="odesolver.cs" /> 165 <Compile Include="pca.cs" /> 166 <Compile Include="poissondistr.cs" /> 167 <Compile Include="polint.cs" /> 97 168 <Compile Include="Properties\AssemblyInfo.cs" /> 169 <Compile Include="psif.cs" /> 98 170 <Compile Include="qr.cs" /> 171 <Compile Include="ratint.cs" /> 172 <Compile Include="ratinterpolation.cs" /> 173 <Compile Include="rcond.cs" /> 99 174 <Compile Include="reflections.cs" /> 100 175 <Compile Include="rotations.cs" /> 176 <Compile Include="sbisinv.cs" /> 177 <Compile Include="sblas.cs" /> 178 <Compile Include="schur.cs" /> 179 <Compile Include="sdet.cs" /> 180 <Compile Include="sevd.cs" /> 181 <Compile Include="sinverse.cs" /> 182 <Compile Include="spddet.cs" /> 183 <Compile Include="spdgevd.cs" /> 184 <Compile Include="spdinverse.cs" /> 185 <Compile Include="spdrcond.cs" /> 186 <Compile Include="spdsolve.cs" /> 187 <Compile Include="spline1d.cs" /> 188 <Compile Include="spline2d.cs" /> 101 189 <Compile Include="spline3.cs" /> 190 <Compile Include="srcond.cs" /> 191 <Compile Include="ssolve.cs" /> 192 <Compile Include="stest.cs" /> 193 <Compile Include="studenttdistr.cs" /> 194 <Compile Include="studentttests.cs" /> 102 195 <Compile Include="svd.cs" /> 196 <Compile Include="taskgen.cs" /> 197 <Compile Include="tdbisinv.cs" /> 198 <Compile Include="tdevd.cs" /> 199 <Compile Include="tridiagonal.cs" /> 200 <Compile Include="trigintegrals.cs" /> 201 <Compile Include="trinverse.cs" /> 202 <Compile Include="trlinsolve.cs" /> 203 <Compile Include="tsort.cs" /> 204 <Compile Include="variancetests.cs" /> 205 <Compile Include="wsr.cs" /> 206 <Compile Include="xblas.cs" /> 103 207 </ItemGroup> 104 208 <ItemGroup> -
trunk/sources/ALGLIB/ap.cs
r2430 r2563 45 45 public static bool operator==(Complex lhs, Complex rhs) 46 46 { 47 return ( lhs.x==rhs.x) & (lhs.y==rhs.y);47 return ((double)lhs.x==(double)rhs.x) & ((double)lhs.y==(double)rhs.y); 48 48 } 49 49 public static bool operator!=(Complex lhs, Complex rhs) 50 50 { 51 return ( lhs.x!=rhs.x) | (lhs.y!=rhs.y);51 return ((double)lhs.x!=(double)rhs.x) | ((double)lhs.y!=(double)rhs.y); 52 52 } 53 53 public static Complex operator+(Complex lhs) … … 92 92 return result; 93 93 } 94 } 94 public override int GetHashCode() 95 { 96 return x.GetHashCode() ^ y.GetHashCode(); 97 } 98 public override bool Equals(object obj) 99 { 100 if( obj is byte) 101 return Equals(new Complex((byte)obj)); 102 if( obj is sbyte) 103 return Equals(new Complex((sbyte)obj)); 104 if( obj is short) 105 return Equals(new Complex((short)obj)); 106 if( obj is ushort) 107 return Equals(new Complex((ushort)obj)); 108 if( obj is int) 109 return Equals(new Complex((int)obj)); 110 if( obj is uint) 111 return Equals(new Complex((uint)obj)); 112 if( obj is long) 113 return Equals(new Complex((long)obj)); 114 if( obj is ulong) 115 return Equals(new Complex((ulong)obj)); 116 if( obj is float) 117 return Equals(new Complex((float)obj)); 118 if( obj is double) 119 return Equals(new Complex((double)obj)); 120 if( obj is decimal) 121 return Equals(new Complex((double)(decimal)obj)); 122 return base.Equals(obj); 123 } 124 } 95 125 96 126 /******************************************************************** -
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); -
trunk/sources/ALGLIB/bidiagonal.cs
r2430 r2563 775 775 776 776 777 /*************************************************************************778 Obsolete 1-based subroutine.779 See RMatrixBD for 0-based replacement.780 *************************************************************************/781 777 public static void tobidiagonal(ref double[,] a, 782 778 int m, … … 938 934 939 935 940 /*************************************************************************941 Obsolete 1-based subroutine.942 See RMatrixBDUnpackQ for 0-based replacement.943 *************************************************************************/944 936 public static void unpackqfrombidiagonal(ref double[,] qp, 945 937 int m, … … 1020 1012 1021 1013 1022 /*************************************************************************1023 Obsolete 1-based subroutine.1024 See RMatrixBDMultiplyByQ for 0-based replacement.1025 *************************************************************************/1026 1014 public static void multiplybyqfrombidiagonal(ref double[,] qp, 1027 1015 int m, … … 1169 1157 1170 1158 1171 /*************************************************************************1172 Obsolete 1-based subroutine.1173 See RMatrixBDUnpackPT for 0-based replacement.1174 *************************************************************************/1175 1159 public static void unpackptfrombidiagonal(ref double[,] qp, 1176 1160 int m, … … 1251 1235 1252 1236 1253 /*************************************************************************1254 Obsolete 1-based subroutine.1255 See RMatrixBDMultiplyByP for 0-based replacement.1256 *************************************************************************/1257 1237 public static void multiplybypfrombidiagonal(ref double[,] qp, 1258 1238 int m, … … 1402 1382 1403 1383 1404 /*************************************************************************1405 Obsolete 1-based subroutine.1406 See RMatrixBDUnpackDiagonals for 0-based replacement.1407 *************************************************************************/1408 1384 public static void unpackdiagonalsfrombidiagonal(ref double[,] b, 1409 1385 int m, -
trunk/sources/ALGLIB/blas.cs
r2430 r2563 51 51 for(ix=i1; ix<=i2; ix++) 52 52 { 53 if( x[ix]!=0)53 if( (double)(x[ix])!=(double)(0) ) 54 54 { 55 55 absxi = Math.Abs(x[ix]); 56 if( scl<absxi)56 if( (double)(scl)<(double)(absxi) ) 57 57 { 58 58 ssq = 1+ssq*AP.Math.Sqr(scl/absxi); … … 82 82 for(i=i1+1; i<=i2; i++) 83 83 { 84 if( Math.Abs(x[i])>Math.Abs(x[result]) )84 if( (double)(Math.Abs(x[i]))>(double)(Math.Abs(x[result])) ) 85 85 { 86 86 result = i; … … 104 104 for(i=i1+1; i<=i2; i++) 105 105 { 106 if( Math.Abs(x[i,j])>Math.Abs(x[result,j]) )106 if( (double)(Math.Abs(x[i,j]))>(double)(Math.Abs(x[result,j])) ) 107 107 { 108 108 result = i; … … 126 126 for(j=j1+1; j<=j2; j++) 127 127 { 128 if( Math.Abs(x[i,j])>Math.Abs(x[i,result]) )128 if( (double)(Math.Abs(x[i,j]))>(double)(Math.Abs(x[i,result])) ) 129 129 { 130 130 result = j; … … 316 316 // beta*y 317 317 // 318 if( beta==0)318 if( (double)(beta)==(double)(0) ) 319 319 { 320 320 for(i=iy1; i<=iy2; i++) … … 361 361 // beta*y 362 362 // 363 if( beta==0)363 if( (double)(beta)==(double)(0) ) 364 364 { 365 365 for(i=iy1; i<=iy2; i++) … … 405 405 w = Math.Max(xabs, yabs); 406 406 z = Math.Min(xabs, yabs); 407 if( z==0)407 if( (double)(z)==(double)(0) ) 408 408 { 409 409 result = w; … … 497 497 // Prepare C 498 498 // 499 if( beta==0)499 if( (double)(beta)==(double)(0) ) 500 500 { 501 501 for(i=ci1; i<=ci2; i++) -
trunk/sources/ALGLIB/correlation.cs
r2430 r2563 88 88 s = s+t1*t2; 89 89 } 90 if( xv==0 | yv==0)90 if( (double)(xv)==(double)(0) | (double)(yv)==(double)(0) ) 91 91 { 92 92 result = 0; … … 174 174 { 175 175 k = t/2; 176 if( r[k-1]>=r[t-1])176 if( (double)(r[k-1])>=(double)(r[t-1]) ) 177 177 { 178 178 t = 1; … … 213 213 if( k<i ) 214 214 { 215 if( r[k]>r[k-1])215 if( (double)(r[k])>(double)(r[k-1]) ) 216 216 { 217 217 k = k+1; 218 218 } 219 219 } 220 if( r[t-1]>=r[k-1])220 if( (double)(r[t-1])>=(double)(r[k-1]) ) 221 221 { 222 222 t = 0; … … 248 248 while( j<=n-1 ) 249 249 { 250 if( r[j]!=r[i])250 if( (double)(r[j])!=(double)(r[i]) ) 251 251 { 252 252 break; -
trunk/sources/ALGLIB/descriptivestatistics.cs
r2445 r2563 90 90 v2 = AP.Math.Sqr(v2)/n; 91 91 variance = (v1-v2)/(n-1); 92 if( variance<0)92 if( (double)(variance)<(double)(0) ) 93 93 { 94 94 variance = 0; … … 100 100 // Skewness and kurtosis 101 101 // 102 if( stddev!=0)102 if( (double)(stddev)!=(double)(0) ) 103 103 { 104 104 for(i=0; i<=n-1; i++) … … 226 226 // 1 or 2 elements in partition 227 227 // 228 if( ir==l+1 & x[ir]<x[l])228 if( ir==l+1 & (double)(x[ir])<(double)(x[l]) ) 229 229 { 230 230 tval = x[l]; … … 240 240 x[midp] = x[l+1]; 241 241 x[l+1] = tval; 242 if( x[l]>x[ir])242 if( (double)(x[l])>(double)(x[ir]) ) 243 243 { 244 244 tval = x[l]; … … 246 246 x[ir] = tval; 247 247 } 248 if( x[l+1]>x[ir])248 if( (double)(x[l+1])>(double)(x[ir]) ) 249 249 { 250 250 tval = x[l+1]; … … 252 252 x[ir] = tval; 253 253 } 254 if( x[l]>x[l+1])254 if( (double)(x[l])>(double)(x[l+1]) ) 255 255 { 256 256 tval = x[l]; … … 267 267 i = i+1; 268 268 } 269 while( x[i]<a);269 while( (double)(x[i])<(double)(a) ); 270 270 do 271 271 { 272 272 j = j-1; 273 273 } 274 while( x[j]>a);274 while( (double)(x[j])>(double)(a) ); 275 275 if( j<i ) 276 276 { … … 305 305 for(i=k+1; i<=n-1; i++) 306 306 { 307 if( x[i]<a)307 if( (double)(x[i])<(double)(a) ) 308 308 { 309 309 a = x[i]; … … 339 339 340 340 System.Diagnostics.Debug.Assert(n>1, "CalculatePercentile: N<=1!"); 341 System.Diagnostics.Debug.Assert( p>=0 & p<=1, "CalculatePercentile: incorrect P!");341 System.Diagnostics.Debug.Assert((double)(p)>=(double)(0) & (double)(p)<=(double)(1), "CalculatePercentile: incorrect P!"); 342 342 internalstatheapsort(ref x, n); 343 if( p==0)343 if( (double)(p)==(double)(0) ) 344 344 { 345 345 v = x[0]; 346 346 return; 347 347 } 348 if( p==1)348 if( (double)(p)==(double)(1) ) 349 349 { 350 350 v = x[n-1]; … … 378 378 { 379 379 k = t/2; 380 if( arr[k-1]>=arr[t-1])380 if( (double)(arr[k-1])>=(double)(arr[t-1]) ) 381 381 { 382 382 t = 1; … … 411 411 if( k<i ) 412 412 { 413 if( arr[k]>arr[k-1])413 if( (double)(arr[k])>(double)(arr[k-1]) ) 414 414 { 415 415 k = k+1; 416 416 } 417 417 } 418 if( arr[t-1]>=arr[k-1])418 if( (double)(arr[t-1])>=(double)(arr[k-1]) ) 419 419 { 420 420 t = 0; -
trunk/sources/ALGLIB/igammaf.cs
r2445 r2563 73 73 74 74 igammaepsilon = 0.000000000000001; 75 if( x<=0 | a<=0)75 if( (double)(x)<=(double)(0) | (double)(a)<=(double)(0) ) 76 76 { 77 77 result = 0; 78 78 return result; 79 79 } 80 if( x>1 & x>a)80 if( (double)(x)>(double)(1) & (double)(x)>(double)(a) ) 81 81 { 82 82 result = 1-incompletegammac(a, x); 83 83 return result; 84 84 } 85 ax = a*Math.Log(x)-x-gammaf .lngamma(a, ref tmp);86 if( ax<-709.78271289338399)85 ax = a*Math.Log(x)-x-gammafunc.lngamma(a, ref tmp); 86 if( (double)(ax)<(double)(-709.78271289338399) ) 87 87 { 88 88 result = 0; … … 99 99 ans = ans+c; 100 100 } 101 while( c/ans>igammaepsilon);101 while( (double)(c/ans)>(double)(igammaepsilon) ); 102 102 result = ans*ax/a; 103 103 return result; … … 164 164 igammabignumber = 4503599627370496.0; 165 165 igammabignumberinv = 2.22044604925031308085*0.0000000000000001; 166 if( x<=0 | a<=0)166 if( (double)(x)<=(double)(0) | (double)(a)<=(double)(0) ) 167 167 { 168 168 result = 1; 169 169 return result; 170 170 } 171 if( x<1 | x<a)171 if( (double)(x)<(double)(1) | (double)(x)<(double)(a) ) 172 172 { 173 173 result = 1-incompletegamma(a, x); 174 174 return result; 175 175 } 176 ax = a*Math.Log(x)-x-gammaf .lngamma(a, ref tmp);177 if( ax<-709.78271289338399)176 ax = a*Math.Log(x)-x-gammafunc.lngamma(a, ref tmp); 177 if( (double)(ax)<(double)(-709.78271289338399) ) 178 178 { 179 179 result = 0; … … 197 197 pk = pkm1*z-pkm2*yc; 198 198 qk = qkm1*z-qkm2*yc; 199 if( qk!=0)199 if( (double)(qk)!=(double)(0) ) 200 200 { 201 201 r = pk/qk; … … 211 211 qkm2 = qkm1; 212 212 qkm1 = qk; 213 if( Math.Abs(pk)>igammabignumber)213 if( (double)(Math.Abs(pk))>(double)(igammabignumber) ) 214 214 { 215 215 pkm2 = pkm2*igammabignumberinv; … … 219 219 } 220 220 } 221 while( t>igammaepsilon);221 while( (double)(t)>(double)(igammaepsilon) ); 222 222 result = ans*ax; 223 223 return result; … … 290 290 y = 1-d-normaldistr.invnormaldistribution(y0)*Math.Sqrt(d); 291 291 x = a*y*y*y; 292 lgm = gammaf .lngamma(a, ref tmp);292 lgm = gammafunc.lngamma(a, ref tmp); 293 293 i = 0; 294 294 while( i<10 ) 295 295 { 296 if( x>x0 | x<x1)296 if( (double)(x)>(double)(x0) | (double)(x)<(double)(x1) ) 297 297 { 298 298 d = 0.0625; … … 300 300 } 301 301 y = incompletegammac(a, x); 302 if( y<yl | y>yh)302 if( (double)(y)<(double)(yl) | (double)(y)>(double)(yh) ) 303 303 { 304 304 d = 0.0625; 305 305 break; 306 306 } 307 if( y<y0)307 if( (double)(y)<(double)(y0) ) 308 308 { 309 309 x0 = x; … … 316 316 } 317 317 d = (a-1)*Math.Log(x)-x-lgm; 318 if( d<-709.78271289338399)318 if( (double)(d)<(double)(-709.78271289338399) ) 319 319 { 320 320 d = 0.0625; … … 323 323 d = -Math.Exp(d); 324 324 d = (y-y0)/d; 325 if( Math.Abs(d/x)<igammaepsilon)325 if( (double)(Math.Abs(d/x))<(double)(igammaepsilon) ) 326 326 { 327 327 result = x; … … 331 331 i = i+1; 332 332 } 333 if( x0==iinvgammabignumber)334 { 335 if( x<=0)333 if( (double)(x0)==(double)(iinvgammabignumber) ) 334 { 335 if( (double)(x)<=(double)(0) ) 336 336 { 337 337 x = 1; 338 338 } 339 while( x0==iinvgammabignumber)339 while( (double)(x0)==(double)(iinvgammabignumber) ) 340 340 { 341 341 x = (1+d)*x; 342 342 y = incompletegammac(a, x); 343 if( y<y0)343 if( (double)(y)<(double)(y0) ) 344 344 { 345 345 x0 = x; … … 358 358 y = incompletegammac(a, x); 359 359 lgm = (x0-x1)/(x1+x0); 360 if( Math.Abs(lgm)<dithresh)360 if( (double)(Math.Abs(lgm))<(double)(dithresh) ) 361 361 { 362 362 break; 363 363 } 364 364 lgm = (y-y0)/y0; 365 if( Math.Abs(lgm)<dithresh)366 { 367 break; 368 } 369 if( x<=0.0)370 { 371 break; 372 } 373 if( y>=y0)365 if( (double)(Math.Abs(lgm))<(double)(dithresh) ) 366 { 367 break; 368 } 369 if( (double)(x)<=(double)(0.0) ) 370 { 371 break; 372 } 373 if( (double)(y)>=(double)(y0) ) 374 374 { 375 375 x1 = x; -
trunk/sources/ALGLIB/leastsquares.cs
r2430 r2563 190 190 // B2 := (d^(-1) * B2')' 191 191 // 192 if( d[1]!=0)192 if( (double)(d[1])!=(double)(0) ) 193 193 { 194 194 for(i=1; i<=ni; i++) 195 195 { 196 if( d[i]>AP.Math.MachineEpsilon*10*Math.Sqrt(ni)*d[1])196 if( (double)(d[i])>(double)(AP.Math.MachineEpsilon*10*Math.Sqrt(ni)*d[1]) ) 197 197 { 198 198 b2[1,i] = b2[1,i]/d[i]; … … 286 286 d1 = AP.Math.Sqr(c)*pp+AP.Math.Sqr(s)*qq-2*s*c*pq; 287 287 d2 = AP.Math.Sqr(s)*pp+AP.Math.Sqr(c)*qq+2*s*c*pq; 288 if( Math.Abs(d1)>Math.Abs(d2) )288 if( (double)(Math.Abs(d1))>(double)(Math.Abs(d2)) ) 289 289 { 290 290 m = Math.Abs(d1); … … 296 296 t1 = c*b1-s*b2; 297 297 t2 = s*b1+c*b2; 298 if( Math.Abs(d1)>m*AP.Math.MachineEpsilon*1000)298 if( (double)(Math.Abs(d1))>(double)(m*AP.Math.MachineEpsilon*1000) ) 299 299 { 300 300 t1 = t1/d1; … … 304 304 t1 = 0; 305 305 } 306 if( Math.Abs(d2)>m*AP.Math.MachineEpsilon*1000)306 if( (double)(Math.Abs(d2))>(double)(m*AP.Math.MachineEpsilon*1000) ) 307 307 { 308 308 t2 = t2/d2; … … 493 493 for(i=1; i<=ni; i++) 494 494 { 495 if( d[i]>AP.Math.MachineEpsilon*10*Math.Sqrt(ni)*d[1])495 if( (double)(d[i])>(double)(AP.Math.MachineEpsilon*10*Math.Sqrt(ni)*d[1]) ) 496 496 { 497 497 b2[1,i] = b2[1,i]/d[i]; … … 585 585 for(i=1; i<=n-1; i++) 586 586 { 587 if( x[i]>maxx)587 if( (double)(x[i])>(double)(maxx) ) 588 588 { 589 589 maxx = x[i]; 590 590 } 591 if( x[i]<minx)591 if( (double)(x[i])<(double)(minx) ) 592 592 { 593 593 minx = x[i]; 594 594 } 595 595 } 596 if( minx==maxx)596 if( (double)(minx)==(double)(maxx) ) 597 597 { 598 598 minx = minx-0.5; … … 869 869 for(i=1; i<=ni; i++) 870 870 { 871 if( d[i]>AP.Math.MachineEpsilon*10*Math.Sqrt(ni)*d[1])871 if( (double)(d[i])>(double)(AP.Math.MachineEpsilon*10*Math.Sqrt(ni)*d[1]) ) 872 872 { 873 873 b2[1,i] = b2[1,i]/d[i]; … … 1133 1133 return result; 1134 1134 } 1135 if( ws[1]==0 | ws[nc]<=AP.Math.MachineEpsilon*10*Math.Sqrt(nc)*ws[1])1135 if( (double)(ws[1])==(double)(0) | (double)(ws[nc])<=(double)(AP.Math.MachineEpsilon*10*Math.Sqrt(nc)*ws[1]) ) 1136 1136 { 1137 1137 result = false; … … 1215 1215 for(i=1; i<=reducedsize; i++) 1216 1216 { 1217 if( ws[i]!=0 & ws[i]>AP.Math.MachineEpsilon*10*Math.Sqrt(nc)*ws[1])1217 if( (double)(ws[i])!=(double)(0) & (double)(ws[i])>(double)(AP.Math.MachineEpsilon*10*Math.Sqrt(nc)*ws[1]) ) 1218 1218 { 1219 1219 tmp[i] = tmp[i]/ws[i]; -
trunk/sources/ALGLIB/linreg.cs
r2445 r2563 232 232 means[j] = mean; 233 233 sigmas[j] = Math.Sqrt(variance); 234 if( sigmas[j]==0)234 if( (double)(sigmas[j])==(double)(0) ) 235 235 { 236 236 sigmas[j] = 1; … … 357 357 } 358 358 descriptivestatistics.calculatemoments(ref x, npoints, ref mean, ref variance, ref skewness, ref kurtosis); 359 if( Math.Abs(mean)>Math.Sqrt(variance) )359 if( (double)(Math.Abs(mean))>(double)(Math.Sqrt(variance)) ) 360 360 { 361 361 … … 372 372 // variation is large, it is better to bring variance to 1 373 373 // 374 if( variance==0)374 if( (double)(variance)==(double)(0) ) 375 375 { 376 376 variance = 1; … … 698 698 for(i=0; i<=npoints-1; i++) 699 699 { 700 if( xy[i,nvars]!=0)700 if( (double)(xy[i,nvars])!=(double)(0) ) 701 701 { 702 702 i1_ = (offs)-(0); … … 806 806 807 807 808 /*************************************************************************809 Obsolete subroutine, use LRBuildS810 811 -- ALGLIB --812 Copyright 26.04.2008 by Bochkanov Sergey813 814 References:815 1. Numerical Recipes in C, "15.2 Fitting Data to a Straight Line"816 *************************************************************************/817 808 public static void lrlines(ref double[,] xy, 818 809 ref double[] s, … … 845 836 for(i=0; i<=n-1; i++) 846 837 { 847 if( s[i]<=0)838 if( (double)(s[i])<=(double)(0) ) 848 839 { 849 840 info = -2; … … 875 866 e1 = 0.5*(ss+sxx+t); 876 867 e2 = 0.5*(ss+sxx-t); 877 if( Math.Min(e1, e2)<=1000*AP.Math.MachineEpsilon*Math.Max(e1, e2) )868 if( (double)(Math.Min(e1, e2))<=(double)(1000*AP.Math.MachineEpsilon*Math.Max(e1, e2)) ) 878 869 { 879 870 info = -3; … … 923 914 924 915 925 /*************************************************************************926 Obsolete subroutine, use LRBuild927 928 -- ALGLIB --929 Copyright 02.08.2008 by Bochkanov Sergey930 *************************************************************************/931 916 public static void lrline(ref double[,] xy, 932 917 int n, … … 1005 990 for(i=0; i<=npoints-1; i++) 1006 991 { 1007 if( s[i]<=0)992 if( (double)(s[i])<=(double)(0) ) 1008 993 { 1009 994 info = -2; … … 1062 1047 return; 1063 1048 } 1064 if( sv[0]<=0)1049 if( (double)(sv[0])<=(double)(0) ) 1065 1050 { 1066 1051 … … 1090 1075 return; 1091 1076 } 1092 if( sv[nvars-1]<=epstol*AP.Math.MachineEpsilon*sv[0])1077 if( (double)(sv[nvars-1])<=(double)(epstol*AP.Math.MachineEpsilon*sv[0]) ) 1093 1078 { 1094 1079 … … 1103 1088 for(k=nvars; k>=1; k--) 1104 1089 { 1105 if( sv[k-1]>epstol*AP.Math.MachineEpsilon*sv[0])1090 if( (double)(sv[k-1])>(double)(epstol*AP.Math.MachineEpsilon*sv[0]) ) 1106 1091 { 1107 1092 … … 1173 1158 for(i=0; i<=nvars-1; i++) 1174 1159 { 1175 if( sv[i]>epstol*AP.Math.MachineEpsilon*sv[0])1160 if( (double)(sv[i])>(double)(epstol*AP.Math.MachineEpsilon*sv[0]) ) 1176 1161 { 1177 1162 svi[i] = 1/sv[i]; … … 1292 1277 ar.rmserror = ar.rmserror+AP.Math.Sqr(r-xy[i,nvars]); 1293 1278 ar.avgerror = ar.avgerror+Math.Abs(r-xy[i,nvars]); 1294 if( xy[i,nvars]!=0)1279 if( (double)(xy[i,nvars])!=(double)(0) ) 1295 1280 { 1296 1281 ar.avgrelerror = ar.avgrelerror+Math.Abs((r-xy[i,nvars])/xy[i,nvars]); … … 1306 1291 p += u[i,i_]*u[i,i_]; 1307 1292 } 1308 if( p>1-epstol*AP.Math.MachineEpsilon)1293 if( (double)(p)>(double)(1-epstol*AP.Math.MachineEpsilon) ) 1309 1294 { 1310 1295 ar.cvdefects[ar.ncvdefects] = i; … … 1315 1300 ar.cvrmserror = ar.cvrmserror+AP.Math.Sqr(r-xy[i,nvars]); 1316 1301 ar.cvavgerror = ar.cvavgerror+Math.Abs(r-xy[i,nvars]); 1317 if( xy[i,nvars]!=0)1302 if( (double)(xy[i,nvars])!=(double)(0) ) 1318 1303 { 1319 1304 ar.cvavgrelerror = ar.cvavgrelerror+Math.Abs((r-xy[i,nvars])/xy[i,nvars]); -
trunk/sources/ALGLIB/lq.cs
r2430 r2563 247 247 248 248 249 /*************************************************************************250 Obsolete 1-based subroutine251 See RMatrixLQ for 0-based replacement.252 *************************************************************************/253 249 public static void lqdecomposition(ref double[,] a, 254 250 int m, … … 309 305 310 306 311 /*************************************************************************312 Obsolete 1-based subroutine313 See RMatrixLQUnpackQ for 0-based replacement.314 *************************************************************************/315 307 public static void unpackqfromlq(ref double[,] a, 316 308 int m, … … 380 372 381 373 382 /*************************************************************************383 Obsolete 1-based subroutine384 *************************************************************************/385 374 public static void lqdecompositionunpacked(double[,] a, 386 375 int m, -
trunk/sources/ALGLIB/normaldistr.cs
r2445 r2563 68 68 s = Math.Sign(x); 69 69 x = Math.Abs(x); 70 if( x<0.5)70 if( (double)(x)<(double)(0.5) ) 71 71 { 72 72 xsq = x*x; … … 88 88 return result; 89 89 } 90 if( x>=10)90 if( (double)(x)>=(double)(10) ) 91 91 { 92 92 result = s; … … 131 131 double q = 0; 132 132 133 if( x<0)133 if( (double)(x)<(double)(0) ) 134 134 { 135 135 result = 2-erfc(-x); 136 136 return result; 137 137 } 138 if( x<0.5)138 if( (double)(x)<(double)(0.5) ) 139 139 { 140 140 result = 1.0-erf(x); 141 141 return result; 142 142 } 143 if( x>=10)143 if( (double)(x)>=(double)(10) ) 144 144 { 145 145 result = 0; … … 269 269 expm2 = 0.13533528323661269189; 270 270 s2pi = 2.50662827463100050242; 271 if( y0<=0)271 if( (double)(y0)<=(double)(0) ) 272 272 { 273 273 result = -AP.Math.MaxRealNumber; 274 274 return result; 275 275 } 276 if( y0>=1)276 if( (double)(y0)>=(double)(1) ) 277 277 { 278 278 result = AP.Math.MaxRealNumber; … … 281 281 code = 1; 282 282 y = y0; 283 if( y>1.0-expm2)283 if( (double)(y)>(double)(1.0-expm2) ) 284 284 { 285 285 y = 1.0-y; 286 286 code = 0; 287 287 } 288 if( y>expm2)288 if( (double)(y)>(double)(expm2) ) 289 289 { 290 290 y = y-0.5; … … 312 312 x0 = x-Math.Log(x)/x; 313 313 z = 1.0/x; 314 if( x<8.0)314 if( (double)(x)<(double)(8.0) ) 315 315 { 316 316 p1 = 4.05544892305962419923; -
trunk/sources/ALGLIB/qr.cs
r2430 r2563 263 263 264 264 265 /*************************************************************************266 Obsolete 1-based subroutine. See RMatrixQR for 0-based replacement.267 *************************************************************************/268 265 public static void qrdecomposition(ref double[,] a, 269 266 int m, … … 322 319 323 320 324 /*************************************************************************325 Obsolete 1-based subroutine. See RMatrixQRUnpackQ for 0-based replacement.326 *************************************************************************/327 321 public static void unpackqfromqr(ref double[,] a, 328 322 int m, … … 392 386 393 387 394 /*************************************************************************395 Obsolete 1-based subroutine. See RMatrixQR for 0-based replacement.396 *************************************************************************/397 388 public static void qrdecompositionunpacked(double[,] a, 398 389 int m, -
trunk/sources/ALGLIB/reflections.cs
r2430 r2563 98 98 } 99 99 s = 1; 100 if( mx!=0)101 { 102 if( mx<=AP.Math.MinRealNumber/AP.Math.MachineEpsilon)100 if( (double)(mx)!=(double)(0) ) 101 { 102 if( (double)(mx)<=(double)(AP.Math.MinRealNumber/AP.Math.MachineEpsilon) ) 103 103 { 104 104 s = AP.Math.MinRealNumber/AP.Math.MachineEpsilon; … … 112 112 else 113 113 { 114 if( mx>=AP.Math.MaxRealNumber*AP.Math.MachineEpsilon)114 if( (double)(mx)>=(double)(AP.Math.MaxRealNumber*AP.Math.MachineEpsilon) ) 115 115 { 116 116 s = AP.Math.MaxRealNumber*AP.Math.MachineEpsilon; … … 130 130 alpha = x[1]; 131 131 xnorm = 0; 132 if( mx!=0)132 if( (double)(mx)!=(double)(0) ) 133 133 { 134 134 for(j=2; j<=n; j++) … … 138 138 xnorm = Math.Sqrt(xnorm)*mx; 139 139 } 140 if( xnorm==0)140 if( (double)(xnorm)==(double)(0) ) 141 141 { 142 142 … … 154 154 mx = Math.Max(Math.Abs(alpha), Math.Abs(xnorm)); 155 155 beta = -(mx*Math.Sqrt(AP.Math.Sqr(alpha/mx)+AP.Math.Sqr(xnorm/mx))); 156 if( alpha<0)156 if( (double)(alpha)<(double)(0) ) 157 157 { 158 158 beta = -beta; … … 215 215 int i_ = 0; 216 216 217 if( tau==0| n1>n2 | m1>m2 )217 if( (double)(tau)==(double)(0) | n1>n2 | m1>m2 ) 218 218 { 219 219 return; … … 294 294 int i1_ = 0; 295 295 296 if( tau==0| n1>n2 | m1>m2 )296 if( (double)(tau)==(double)(0) | n1>n2 | m1>m2 ) 297 297 { 298 298 return; -
trunk/sources/ALGLIB/rotations.cs
r2430 r2563 93 93 ctemp = c[j-m1+1]; 94 94 stemp = s[j-m1+1]; 95 if( ctemp!=1 | stemp!=0)95 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 96 96 { 97 97 jp1 = j+1; … … 129 129 ctemp = c[j-m1+1]; 130 130 stemp = s[j-m1+1]; 131 if( ctemp!=1 | stemp!=0)131 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 132 132 { 133 133 temp = a[j+1,n1]; … … 150 150 ctemp = c[j-m1+1]; 151 151 stemp = s[j-m1+1]; 152 if( ctemp!=1 | stemp!=0)152 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 153 153 { 154 154 jp1 = j+1; … … 186 186 ctemp = c[j-m1+1]; 187 187 stemp = s[j-m1+1]; 188 if( ctemp!=1 | stemp!=0)188 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 189 189 { 190 190 temp = a[j+1,n1]; … … 256 256 ctemp = c[j-n1+1]; 257 257 stemp = s[j-n1+1]; 258 if( ctemp!=1 | stemp!=0)258 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 259 259 { 260 260 jp1 = j+1; … … 292 292 ctemp = c[j-n1+1]; 293 293 stemp = s[j-n1+1]; 294 if( ctemp!=1 | stemp!=0)294 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 295 295 { 296 296 temp = a[m1,j+1]; … … 313 313 ctemp = c[j-n1+1]; 314 314 stemp = s[j-n1+1]; 315 if( ctemp!=1 | stemp!=0)315 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 316 316 { 317 317 jp1 = j+1; … … 349 349 ctemp = c[j-n1+1]; 350 350 stemp = s[j-n1+1]; 351 if( ctemp!=1 | stemp!=0)351 if( (double)(ctemp)!=(double)(1) | (double)(stemp)!=(double)(0) ) 352 352 { 353 353 temp = a[m1,j+1]; … … 378 378 double g1 = 0; 379 379 380 if( g==0)380 if( (double)(g)==(double)(0) ) 381 381 { 382 382 cs = 1; … … 386 386 else 387 387 { 388 if( f==0)388 if( (double)(f)==(double)(0) ) 389 389 { 390 390 cs = 0; … … 399 399 cs = f1/r; 400 400 sn = g1/r; 401 if( Math.Abs(f)>Math.Abs(g) & cs<0)401 if( (double)(Math.Abs(f))>(double)(Math.Abs(g)) & (double)(cs)<(double)(0) ) 402 402 { 403 403 cs = -cs; … … 408 408 } 409 409 } 410 411 412 private static void testrotations()413 {414 double[,] al1 = new double[0,0];415 double[,] al2 = new double[0,0];416 double[,] ar1 = new double[0,0];417 double[,] ar2 = new double[0,0];418 double[] cl = new double[0];419 double[] sl = new double[0];420 double[] cr = new double[0];421 double[] sr = new double[0];422 double[] w = new double[0];423 int m = 0;424 int n = 0;425 int maxmn = 0;426 double t = 0;427 int pass = 0;428 int passcount = 0;429 int i = 0;430 int j = 0;431 double err = 0;432 double maxerr = 0;433 bool isforward = new bool();434 435 passcount = 1000;436 maxerr = 0;437 for(pass=1; pass<=passcount; pass++)438 {439 440 //441 // settings442 //443 m = 2+AP.Math.RandomInteger(50);444 n = 2+AP.Math.RandomInteger(50);445 isforward = AP.Math.RandomReal()>0.5;446 maxmn = Math.Max(m, n);447 al1 = new double[m+1, n+1];448 al2 = new double[m+1, n+1];449 ar1 = new double[m+1, n+1];450 ar2 = new double[m+1, n+1];451 cl = new double[m-1+1];452 sl = new double[m-1+1];453 cr = new double[n-1+1];454 sr = new double[n-1+1];455 w = new double[maxmn+1];456 457 //458 // matrices and rotaions459 //460 for(i=1; i<=m; i++)461 {462 for(j=1; j<=n; j++)463 {464 al1[i,j] = 2*AP.Math.RandomReal()-1;465 al2[i,j] = al1[i,j];466 ar1[i,j] = al1[i,j];467 ar2[i,j] = al1[i,j];468 }469 }470 for(i=1; i<=m-1; i++)471 {472 t = 2*Math.PI*AP.Math.RandomReal();473 cl[i] = Math.Cos(t);474 sl[i] = Math.Sin(t);475 }476 for(j=1; j<=n-1; j++)477 {478 t = 2*Math.PI*AP.Math.RandomReal();479 cr[j] = Math.Cos(t);480 sr[j] = Math.Sin(t);481 }482 483 //484 // Test left485 //486 applyrotationsfromtheleft(isforward, 1, m, 1, n, ref cl, ref sl, ref al1, ref w);487 for(j=1; j<=n; j++)488 {489 applyrotationsfromtheleft(isforward, 1, m, j, j, ref cl, ref sl, ref al2, ref w);490 }491 err = 0;492 for(i=1; i<=m; i++)493 {494 for(j=1; j<=n; j++)495 {496 err = Math.Max(err, Math.Abs(al1[i,j]-al2[i,j]));497 }498 }499 maxerr = Math.Max(err, maxerr);500 501 //502 // Test right503 //504 applyrotationsfromtheright(isforward, 1, m, 1, n, ref cr, ref sr, ref ar1, ref w);505 for(i=1; i<=m; i++)506 {507 applyrotationsfromtheright(isforward, i, i, 1, n, ref cr, ref sr, ref ar2, ref w);508 }509 err = 0;510 for(i=1; i<=m; i++)511 {512 for(j=1; j<=n; j++)513 {514 err = Math.Max(err, Math.Abs(ar1[i,j]-ar2[i,j]));515 }516 }517 maxerr = Math.Max(err, maxerr);518 }519 System.Console.Write("TESTING ROTATIONS");520 System.Console.WriteLine();521 System.Console.Write("Pass count ");522 System.Console.Write("{0,0:d}",passcount);523 System.Console.WriteLine();524 System.Console.Write("Error is ");525 System.Console.Write("{0,5:E3}",maxerr);526 System.Console.WriteLine();527 }528 410 } 529 411 } -
trunk/sources/ALGLIB/spline3.cs
r2430 r2563 25 25 public class spline3 26 26 { 27 /*************************************************************************28 This subroutine builds linear spline coefficients table.29 30 Input parameters:31 X - spline nodes, array[0..N-1]32 Y - function values, array[0..N-1]33 N - points count, N>=234 35 Output parameters:36 C - coefficients table. Used by SplineInterpolation and other37 subroutines from this file.38 39 -- ALGLIB PROJECT --40 Copyright 24.06.2007 by Bochkanov Sergey41 *************************************************************************/42 27 public static void buildlinearspline(double[] x, 43 28 double[] y, … … 86 71 87 72 88 /*************************************************************************89 This subroutine builds cubic spline coefficients table.90 91 Input parameters:92 X - spline nodes, array[0..N-1]93 Y - function values, array[0..N-1]94 N - points count, N>=295 BoundLType - boundary condition type for the left boundary96 BoundL - left boundary condition (first or second derivative,97 depending on the BoundLType)98 BoundRType - boundary condition type for the right boundary99 BoundR - right boundary condition (first or second derivative,100 depending on the BoundRType)101 102 Output parameters:103 C - coefficients table. Used by SplineInterpolation and104 other subroutines from this file.105 106 The BoundLType/BoundRType parameters can have the following values:107 * 0, which corresponds to the parabolically terminated spline108 (BoundL/BoundR are ignored).109 * 1, which corresponds to the first derivative boundary condition110 * 2, which corresponds to the second derivative boundary condition111 112 -- ALGLIB PROJECT --113 Copyright 23.06.2007 by Bochkanov Sergey114 *************************************************************************/115 73 public static void buildcubicspline(double[] x, 116 74 double[] y, … … 240 198 241 199 242 /*************************************************************************243 This subroutine builds cubic Hermite spline coefficients table.244 245 Input parameters:246 X - spline nodes, array[0..N-1]247 Y - function values, array[0..N-1]248 D - derivatives, array[0..N-1]249 N - points count, N>=2250 251 Output parameters:252 C - coefficients table. Used by SplineInterpolation and253 other subroutines from this file.254 255 -- ALGLIB PROJECT --256 Copyright 23.06.2007 by Bochkanov Sergey257 *************************************************************************/258 200 public static void buildhermitespline(double[] x, 259 201 double[] y, … … 310 252 311 253 312 /*************************************************************************313 This subroutine builds Akima spline coefficients table.314 315 Input parameters:316 X - spline nodes, array[0..N-1]317 Y - function values, array[0..N-1]318 N - points count, N>=5319 320 Output parameters:321 C - coefficients table. Used by SplineInterpolation and322 other subroutines from this file.323 324 -- ALGLIB PROJECT --325 Copyright 24.06.2007 by Bochkanov Sergey326 *************************************************************************/327 254 public static void buildakimaspline(double[] x, 328 255 double[] y, … … 365 292 for(i=2; i<=n-3; i++) 366 293 { 367 if( Math.Abs(w[i-1])+Math.Abs(w[i+1])!=0)294 if( (double)(Math.Abs(w[i-1])+Math.Abs(w[i+1]))!=(double)(0) ) 368 295 { 369 296 d[i] = (w[i+1]*diff[i-1]+w[i-1]*diff[i])/(w[i+1]+w[i-1]); … … 386 313 387 314 388 /*************************************************************************389 This subroutine calculates the value of the spline at the given point X.390 391 Input parameters:392 C - coefficients table. Built by BuildLinearSpline,393 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.394 X - point395 396 Result:397 S(x)398 399 -- ALGLIB PROJECT --400 Copyright 23.06.2007 by Bochkanov Sergey401 *************************************************************************/402 315 public static double splineinterpolation(ref double[] c, 403 316 double x) … … 420 333 { 421 334 m = (l+r)/2; 422 if( c[m]>=x)335 if( (double)(c[m])>=(double)(x) ) 423 336 { 424 337 r = m; … … 440 353 441 354 442 /*************************************************************************443 This subroutine differentiates the spline.444 445 Input parameters:446 C - coefficients table. Built by BuildLinearSpline,447 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.448 X - point449 450 Result:451 S - S(x)452 DS - S'(x)453 D2S - S''(x)454 455 -- ALGLIB PROJECT --456 Copyright 24.06.2007 by Bochkanov Sergey457 *************************************************************************/458 355 public static void splinedifferentiation(ref double[] c, 459 356 double x, … … 478 375 { 479 376 m = (l+r)/2; 480 if( c[m]>=x)377 if( (double)(c[m])>=(double)(x) ) 481 378 { 482 379 r = m; … … 499 396 500 397 501 /*************************************************************************502 This subroutine makes the copy of the spline.503 504 Input parameters:505 C - coefficients table. Built by BuildLinearSpline,506 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.507 508 Result:509 CC - spline copy510 511 -- ALGLIB PROJECT --512 Copyright 29.06.2007 by Bochkanov Sergey513 *************************************************************************/514 398 public static void splinecopy(ref double[] c, 515 399 ref double[] cc) … … 527 411 528 412 529 /*************************************************************************530 This subroutine unpacks the spline into the coefficients table.531 532 Input parameters:533 C - coefficients table. Built by BuildLinearSpline,534 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.535 X - point536 537 Result:538 Tbl - coefficients table, unpacked format, array[0..N-2, 0..5].539 For I = 0...N-2:540 Tbl[I,0] = X[i]541 Tbl[I,1] = X[i+1]542 Tbl[I,2] = C0543 Tbl[I,3] = C1544 Tbl[I,4] = C2545 Tbl[I,5] = C3546 On [x[i], x[i+1]] spline is equals to:547 S(x) = C0 + C1*t + C2*t^2 + C3*t^3548 t = x-x[i]549 550 -- ALGLIB PROJECT --551 Copyright 29.06.2007 by Bochkanov Sergey552 *************************************************************************/553 413 public static void splineunpack(ref double[] c, 554 414 ref int n, … … 576 436 577 437 578 /*************************************************************************579 This subroutine performs linear transformation of the spline argument.580 581 Input parameters:582 C - coefficients table. Built by BuildLinearSpline,583 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.584 A, B- transformation coefficients: x = A*t + B585 Result:586 C - transformed spline587 588 -- ALGLIB PROJECT --589 Copyright 30.06.2007 by Bochkanov Sergey590 *************************************************************************/591 438 public static void splinelintransx(ref double[] c, 592 439 double a, … … 608 455 // Special case: A=0 609 456 // 610 if( a==0)457 if( (double)(a)==(double)(0) ) 611 458 { 612 459 v = splineinterpolation(ref c, b); … … 641 488 642 489 643 /*************************************************************************644 This subroutine performs linear transformation of the spline.645 646 Input parameters:647 C - coefficients table. Built by BuildLinearSpline,648 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.649 A, B- transformation coefficients: S2(x) = A*S(x) + B650 Result:651 C - transformed spline652 653 -- ALGLIB PROJECT --654 Copyright 30.06.2007 by Bochkanov Sergey655 *************************************************************************/656 490 public static void splinelintransy(ref double[] c, 657 491 double a, … … 683 517 684 518 685 /*************************************************************************686 This subroutine integrates the spline.687 688 Input parameters:689 C - coefficients table. Built by BuildLinearSpline,690 BuildHermiteSpline, BuildCubicSpline, BuildAkimaSpline.691 X - right bound of the integration interval [a, x]692 Result:693 integral(S(t)dt,a,x)694 695 -- ALGLIB PROJECT --696 Copyright 23.06.2007 by Bochkanov Sergey697 *************************************************************************/698 519 public static double splineintegration(ref double[] c, 699 520 double x) … … 718 539 { 719 540 m = (l+r)/2; 720 if( c[m]>=x)541 if( (double)(c[m])>=(double)(x) ) 721 542 { 722 543 r = m; … … 751 572 752 573 753 /*************************************************************************754 Obsolete subroutine, left for backward compatibility.755 *************************************************************************/756 574 public static void spline3buildtable(int n, 757 575 int diffn, … … 798 616 do 799 617 { 800 if( x[j]<=x[j+g])618 if( (double)(x[j])<=(double)(x[j+g]) ) 801 619 { 802 620 c = false; … … 892 710 893 711 894 /*************************************************************************895 Obsolete subroutine, left for backward compatibility.896 *************************************************************************/897 712 public static double spline3interpolate(int n, 898 713 ref double[,] c, … … 913 728 half = l/2; 914 729 middle = first+half; 915 if( c[4,middle]<x)730 if( (double)(c[4,middle])<(double)(x) ) 916 731 { 917 732 first = middle+1; … … 933 748 934 749 935 /*************************************************************************936 Internal subroutine. Heap sort.937 *************************************************************************/938 750 private static void heapsortpoints(ref double[] x, 939 751 ref double[] y, … … 956 768 for(i=1; i<=n-1; i++) 957 769 { 958 isascending = isascending & x[i]>x[i-1];959 isdescending = isdescending & x[i]<x[i-1];770 isascending = isascending & (double)(x[i])>(double)(x[i-1]); 771 isdescending = isdescending & (double)(x[i])<(double)(x[i-1]); 960 772 } 961 773 if( isascending ) … … 1000 812 { 1001 813 k = t/2; 1002 if( x[k-1]>=x[t-1])814 if( (double)(x[k-1])>=(double)(x[t-1]) ) 1003 815 { 1004 816 t = 1; … … 1039 851 if( k<i ) 1040 852 { 1041 if( x[k]>x[k-1])853 if( (double)(x[k])>(double)(x[k-1]) ) 1042 854 { 1043 855 k = k+1; 1044 856 } 1045 857 } 1046 if( x[t-1]>=x[k-1])858 if( (double)(x[t-1])>=(double)(x[k-1]) ) 1047 859 { 1048 860 t = 0; … … 1066 878 1067 879 1068 /*************************************************************************1069 Internal subroutine. Heap sort.1070 *************************************************************************/1071 880 private static void heapsortdpoints(ref double[] x, 1072 881 ref double[] y, … … 1090 899 for(i=1; i<=n-1; i++) 1091 900 { 1092 isascending = isascending & x[i]>x[i-1];1093 isdescending = isdescending & x[i]<x[i-1];901 isascending = isascending & (double)(x[i])>(double)(x[i-1]); 902 isdescending = isdescending & (double)(x[i])<(double)(x[i-1]); 1094 903 } 1095 904 if( isascending ) … … 1137 946 { 1138 947 k = t/2; 1139 if( x[k-1]>=x[t-1])948 if( (double)(x[k-1])>=(double)(x[t-1]) ) 1140 949 { 1141 950 t = 1; … … 1182 991 if( k<i ) 1183 992 { 1184 if( x[k]>x[k-1])993 if( (double)(x[k])>(double)(x[k-1]) ) 1185 994 { 1186 995 k = k+1; 1187 996 } 1188 997 } 1189 if( x[t-1]>=x[k-1])998 if( (double)(x[t-1])>=(double)(x[k-1]) ) 1190 999 { 1191 1000 t = 0; … … 1212 1021 1213 1022 1214 /*************************************************************************1215 Internal subroutine. Tridiagonal solver.1216 *************************************************************************/1217 1023 private static void solvetridiagonal(double[] a, 1218 1024 double[] b, … … 1247 1053 1248 1054 1249 /*************************************************************************1250 Internal subroutine. Three-point differentiation1251 *************************************************************************/1252 1055 private static double diffthreepoint(double t, 1253 1056 double x0, -
trunk/sources/ALGLIB/svd.cs
r2430 r2563 154 154 // Use bidiagonal reduction with QR-decomposition 155 155 // 156 if( m>1.6*n)156 if( (double)(m)>(double)(1.6*n) ) 157 157 { 158 158 if( uneeded==0 ) … … 224 224 // Use bidiagonal reduction with LQ-decomposition 225 225 // 226 if( n>1.6*m)226 if( (double)(n)>(double)(1.6*m) ) 227 227 { 228 228 if( vtneeded==0 ) … … 341 341 342 342 343 /*************************************************************************344 Obsolete 1-based subroutine.345 See RMatrixSVD for 0-based replacement.346 *************************************************************************/347 343 public static bool svddecomposition(double[,] a, 348 344 int m, … … 422 418 // Use bidiagonal reduction with QR-decomposition 423 419 // 424 if( m>1.6*n)420 if( (double)(m)>(double)(1.6*n) ) 425 421 { 426 422 if( uneeded==0 ) … … 492 488 // Use bidiagonal reduction with LQ-decomposition 493 489 // 494 if( n>1.6*m)490 if( (double)(n)>(double)(1.6*m) ) 495 491 { 496 492 if( vtneeded==0 )
Note: See TracChangeset
for help on using the changeset viewer.