1
编辑：预期的输出结束时发布。Java-比较itr（）和itr.next（）以获得双精度值
我写了一个程序，该程序通过Riemann-Siegel formula计算Riemann Zeta Function的零。我很有兴趣调整程序中的一种方法来观察所谓的Lehmer's Phenomenon。
在程序中，我期待调整findRoots()方法的结束。
完整的程序引用
/**************************************************************************
**
** Riemann-Siegel Formula for roots of Zeta(s) on critical line.
**
**************************************************************************
** Axion004
** 07/31/2015
**
** This program finds the roots of Zeta(s) using the well known Riemann-
** Siegel formula. The Riemann–Siegel theta function is approximated
** using Stirling's approximation. It also uses an interpolation method to
** locate zeroes. The coefficients for R(t) are handled by the Taylor
** Series approximation originally listed by Haselgrove in 1960. It is
** necessary to use these coefficients in order to increase computational
** speed.
**************************************************************************/
import java.util.Iterator;
import java.util.LinkedHashSet;
//These two imports are from the Apache Commons Math library
import org.apache.commons.math3.analysis.UnivariateFunction;
import org.apache.commons.math3.analysis.solvers.BracketingNthOrderBrentSolver;
public class RiemannSiegelTwo{
public static void main(String[] args){
SiegelMain();
}
// Main method
public static void SiegelMain() {
System.out.println("Zeroes inside the critical line for " +
"Zeta(1/2 + it). The t values are referenced below.");
System.out.println();
findRoots();
}
/**
* The sign of a calculated double value.
* @param x - the double value.
* @return the sign in -1, 1, or 0 format.
*/
private static int sign(double x) {
if (x < 0.0)
return -1;
else if (x > 0.0)
return 1;
else
return 0;
}
/**
* Finds the roots of a specified function through the
* BracketingNthOrderBrentSolver from the Apache Commons Math library.
* See http://commons.apache.org/proper/commons-math/apidocs/org/
* apache/commons/math3/analysis/solvers/BracketingNthOrderBrentSolver
* .html
* The zeroes inside the interval of 0 < t < 10000 are printed from
* a TreeSet.
*/
public static void findRoots() {
BracketingNthOrderBrentSolver f = new
BracketingNthOrderBrentSolver();
UnivariateFunction func = (double x) -> RiemennZ(x, 4);
LinkedHashSet<Double> set = new LinkedHashSet<>();
double i = 1.0;
while (i < 1000) {
i+= 0.1;
if(sign(func.value(i)) != sign(func.value(i+0.1))) {
set.add(f.solve(1000, func, i, i+0.1));
}
}
Iterator<Double> itr = set.iterator();
while(itr.hasNext()) {
System.out.println(itr.next());
}
}
/**
* Riemann-Siegel theta function using the approximation by the
* Stirling series.
* @param t - the value of t inside the Z(t) function.
* @return Stirling's approximation for theta(t).
*/
public static double theta (double t) {
return (t/2.0 * Math.log(t/(2.0*Math.PI)) - t/2.0 - Math.PI/8.0
+ 1.0/(48.0*Math.pow(t, 1)) + 7.0/(5760*Math.pow(t, 3)));
}
/**
* Computes Math.Floor of the absolute value term passed in as t.
* @param t - the value of t inside the Z(t) function.
* @return Math.floor of the absolute value of t.
*/
public static double fAbs(double t) {
return Math.floor(Math.abs(t));
}
/**
* Riemann-Siegel Z(t) function implemented per the Riemenn Siegel
* formula. See http://mathworld.wolfram.com/Riemann-SiegelFormula.html
* for details
* @param t - the value of t inside the Z(t) function.
* @param r - referenced for calculating the remainder terms by the
* Taylor series approximations.
* @return the approximate value of Z(t) through the Riemann-Siegel
* formula
*/
public static double RiemennZ(double t, int r) {
double twopi = Math.PI * 2.0;
double val = Math.sqrt(t/twopi);
double n = fAbs(val);
double sum = 0.0;
for (int i = 1; i <= n; i++) {
sum += (Math.cos(theta(t) - t * Math.log(i)))/Math.sqrt(i);
}
sum = 2.0 * sum;
double remainder;
double frac = val - n;
int k = 0;
double R = 0.0;
// Necessary to individually calculate each remainder term by using
// Taylor Series co-efficients. These coefficients are defined below.
while (k <= r) {
R = R + C(k, 2.0*frac-1.0) * Math.pow(t/twopi,
((double) k) * -0.5);
k++;
}
remainder = Math.pow(-1, (int)n-1) * Math.pow(t/twopi, -0.25) * R;
return sum + remainder;
}
/**
* C terms for the Riemann-Siegel formula. See
* https://web.viu.ca/pughg/thesis.d/masters.thesis.pdf for details.
* Calculates the Taylor Series coefficients for C0, C1, C2, C3,
* and C4.
* @param n - the number of coefficient terms to use.
* @param z - referenced per the Taylor series calculations.
* @return the Taylor series approximation of the remainder terms.
*/
public static double C (int n, double z) {
if (n==0)
return(.38268343236508977173 * Math.pow(z, 0.0)
+.43724046807752044936 * Math.pow(z, 2.0)
+.13237657548034352332 * Math.pow(z, 4.0)
-.01360502604767418865 * Math.pow(z, 6.0)
-.01356762197010358089 * Math.pow(z, 8.0)
-.00162372532314446528 * Math.pow(z,10.0)
+.00029705353733379691 * Math.pow(z,12.0)
+.00007943300879521470 * Math.pow(z,14.0)
+.00000046556124614505 * Math.pow(z,16.0)
-.00000143272516309551 * Math.pow(z,18.0)
-.00000010354847112313 * Math.pow(z,20.0)
+.000000* Math.pow(z,22.0)
+.00000000178810838580 * Math.pow(z,24.0)
-.00000000003391414390 * Math.pow(z,26.0)
-.00000000001632663390 * Math.pow(z,28.0)
-.00000000000037851093 * Math.pow(z,30.0)
+.00000000000009327423 * Math.pow(z,32.0)
+.00000000000000522184 * Math.pow(z,34.0)
-.00000000000000033507 * Math.pow(z,36.0)
-.00000000000000003412 * Math.pow(z,38.0)
+.00000000000000000058 * Math.pow(z,40.0)
+.00000000000000000015 * Math.pow(z,42.0));
else if (n==1)
return(-.02682510262837534703 * Math.pow(z, 1.0)
+.01378477342635185305 * Math.pow(z, 3.0)
+.03849125048223508223 * Math.pow(z, 5.0)
+.00987106629906207647 * Math.pow(z, 7.0)
-.00331075976085840433 * Math.pow(z, 9.0)
-.00146478085779541508 * Math.pow(z,11.0)
-.00001320794062487696 * Math.pow(z,13.0)
+.00005922748701847141 * Math.pow(z,15.0)
+.00000598024258537345 * Math.pow(z,17.0)
-.00000096413224561698 * Math.pow(z,19.0)
-.00000018334733722714 * Math.pow(z,21.0)
+.00000000446708756272 * Math.pow(z,23.0)
+.00000000270963508218 * Math.pow(z,25.0)
+.00000000007785288654 * Math.pow(z,27.0)
-.00000000002343762601 * Math.pow(z,29.0)
-.00000000000158301728 * Math.pow(z,31.0)
+.00000000000012119942 * Math.pow(z,33.0)
+.00000000000001458378 * Math.pow(z,35.0)
-.00000000000000028786 * Math.pow(z,37.0)
-.00000000000000008663 * Math.pow(z,39.0)
-.00000000000000000084 * Math.pow(z,41.0)
+.00000000000000000036 * Math.pow(z,43.0)
+.00000000000000000001 * Math.pow(z,45.0));
else if (n==2)
return(+.00518854283029316849 * Math.pow(z, 0.0)
+.00030946583880634746 * Math.pow(z, 2.0)
-.01133594107822937338 * Math.pow(z, 4.0)
+.00223304574195814477 * Math.pow(z, 6.0)
+.00519663740886233021 * Math.pow(z, 8.0)
+.00034399144076208337 * Math.pow(z,10.0)
-.00059106484274705828 * Math.pow(z,12.0)
-.00010229972547935857 * Math.pow(z,14.0)
+.00002088839221699276 * Math.pow(z,16.0)
+.00000592766549309654 * Math.pow(z,18.0)
-.00000016423838362436 * Math.pow(z,20.0)
-.00000015161199700941 * Math.pow(z,22.0)
-.00000000590780369821 * Math.pow(z,24.0)
+.00000000209115148595 * Math.pow(z,26.0)
+.00000000017815649583 * Math.pow(z,28.0)
-.00000000001616407246 * Math.pow(z,30.0)
-.00000000000238069625 * Math.pow(z,32.0)
+.00000000000005398265 * Math.pow(z,34.0)
+.00000000000001975014 * Math.pow(z,36.0)
+.00000000000000023333 * Math.pow(z,38.0)
-.00000000000000011188 * Math.pow(z,40.0)
-.00000000000000000416 * Math.pow(z,42.0)
+.00000000000000000044 * Math.pow(z,44.0)
+.00000000000000000003 * Math.pow(z,46.0));
else if (n==3)
return(-.00133971609071945690 * Math.pow(z, 1.0)
+.00374421513637939370 * Math.pow(z, 3.0)
-.00133031789193214681 * Math.pow(z, 5.0)
-.00226546607654717871 * Math.pow(z, 7.0)
+.00095484999985067304 * Math.pow(z, 9.0)
+.00060100384589636039 * Math.pow(z,11.0)
-.00010128858286776622 * Math.pow(z,13.0)
-.00006865733449299826 * Math.pow(z,15.0)
+.00000059853667915386 * Math.pow(z,17.0)
+.00000333165985123995 * Math.pow(z,19.0)
+.00000021919289102435 * Math.pow(z,21.0)
-.00000007890884245681 * Math.pow(z,23.0)
-.00000000941468508130 * Math.pow(z,25.0)
+.00000000095701162109 * Math.pow(z,27.0)
+.00000000018763137453 * Math.pow(z,29.0)
-.00000000000443783768 * Math.pow(z,31.0)
-.00000000000224267385 * Math.pow(z,33.0)
-.00000000000003627687 * Math.pow(z,35.0)
+.00000000000001763981 * Math.pow(z,37.0)
+.00000000000000079608 * Math.pow(z,39.0)
-.00000000000000009420 * Math.pow(z,41.0)
-.00000000000000000713 * Math.pow(z,43.0)
+.00000000000000000033 * Math.pow(z,45.0)
+.00000000000000000004 * Math.pow(z,47.0));
else
return(+.00046483389361763382 * Math.pow(z, 0.0)
-.00100566073653404708 * Math.pow(z, 2.0)
+.00024044856573725793 * Math.pow(z, 4.0)
+.00102830861497023219 * Math.pow(z, 6.0)
-.00076578610717556442 * Math.pow(z, 8.0)
-.00020365286803084818 * Math.pow(z,10.0)
+.00023212290491068728 * Math.pow(z,12.0)
+.00003260214424386520 * Math.pow(z,14.0)
-.00002557906251794953 * Math.pow(z,16.0)
-.00000410746443891574 * Math.pow(z,18.0)
+.00000117811136403713 * Math.pow(z,20.0)
+.00000024456561422485 * Math.pow(z,22.0)
-.00000002391582476734 * Math.pow(z,24.0)
-.00000000750521420704 * Math.pow(z,26.0)
+.00000000013312279416 * Math.pow(z,28.0)
+.00000000013440626754 * Math.pow(z,30.0)
+.00000000000351377004 * Math.pow(z,32.0)
-.00000000000151915445 * Math.pow(z,34.0)
-.00000000000008915418 * Math.pow(z,36.0)
+.00000000000001119589 * Math.pow(z,38.0)
+.00000000000000105160 * Math.pow(z,40.0)
-.00000000000000005179 * Math.pow(z,42.0)
-.00000000000000000807 * Math.pow(z,44.0)
+.00000000000000000011 * Math.pow(z,46.0)
+.00000000000000000004 * Math.pow(z,48.0));
}
}
其上印有零包括
14.134728277620736
21.022037047686375
25.010858848857314
30.42487545102874
32.93506194020564
37.5861781599404
40.91871891395697
43.32707339344611
48.0051507135618
49.773832659813834
52.97032138572673
56.446247773482625
59.34704389126199
60.83177861205533
65.11254400260545
67.07981057423619
69.54640168845049
72.067157682421
75.70469070082693
77.14484006152179
79.33737502795745
82.91038084149456
84.73549299832806
87.42527458836074
88.809111229007
92.49189925490495
94.65134407077679
95.87063420239032
98.83119422921574
101.31785099451231
103.72553805583749
105.44662303529545
107.1686111943489
111.02953552739187
111.87465919250342
114.32022090805314
116.22668032564269
118.79078286343636
121.37012500478542
122.9468292915669
124.25681855486741
127.51668388015422
129.5787041982112
131.08768853318705
133.49773719981772
134.7565097562393
138.1160420512064
139.73620895846443
141.1237073980395
143.11184581101554
146.00098248274014
147.4227653471707
150.05352041306668
150.92525762000432
153.02469380843965
156.11290929776234
157.5975918115767
158.8499881762274
161.18896413483722
163.0307096895234
165.5370691851963
167.18443998174794
169.09451541024382
169.91197648339372
173.41153651776858
174.7541915257848
176.441434295885
178.37740777747973
179.91648401935672
182.2070784847937
184.87446784765135
185.59878367831283
187.22892258339354
189.41615865594937
192.02665636108685
193.07972660281158
195.26539668016733
196.87648183995336
198.0153096770279
201.26475194277046
202.49359451564763
204.18967180130994
205.39469720336314
207.90625888685815
209.57650971791733
211.69086259398296
213.34791936192354
214.5470447814087
216.16953850912452
219.06759634795878
220.71491884220737
221.43070555208672
224.0070002561315
224.98332466649444
...等等。观察莱默现象，我要找调整以下
Iterator<Double> itr = set.iterator();
while(itr.hasNext()) {
System.out.println(itr.next());
}
类似的东西
// Not workable code
Iterator<Double> itr = set.iterator();
while(itr.hasNext()) {
double EPSILON = .001;
if(itr.nextVal - itr.currentVal < EPSILON)
System.out.println(itr.next());
}
因此，如果这两个值接近零一起只打印。我可以用Iterator来做这件事吗？我需要使用ListIterator吗？
我尝试了一个增强的for循环，并遇到了同样的问题。是否有更好或更有效的方法来比较从LinkedHashSet集打印的值？
（注：我最初创建一个TreeSet，而不是一个LinkedHashSet我后来改变了这一个LinkedHashSet当我意识到LinkedHashSet比一个TreeSet快得多）
预期输出（改变findRoots内i < 50000（）方法）
public static void findRoots() {
BracketingNthOrderBrentSolver f = new
BracketingNthOrderBrentSolver();
UnivariateFunction func = (double x) -> RiemennZ(x, 4);
LinkedHashSet<Double> set = new LinkedHashSet<>();
double i = 1.0;
while (i < 50000) {
i+= 0.1;
if(sign(func.value(i)) != sign(func.value(i+0.1))) {
set.add(f.solve(1000, func, i, i+0.1));
}
}
double EPSILON = .05;
Iterator<Double> itr = set.iterator();
Double prevVal = null;
while(itr.hasNext()) {
Double currentVal = itr.next();
if (prevVal != null) {
if(currentVal - prevVal < EPSILON) {
System.out.println(prevVal);
System.out.println(currentVal);
}
}
prevVal = currentVal;
}
}
输出。两个零都在EPSILON内（相互间为0.05）。这个重要原因处理Z（t）函数背后的数学问题。 EPSILON可以调整以找到更接近的零点。
Zeroes inside the critical line for Zeta(1/2 + it). The t values are referenced below.
5229.198557200015
5229.2418112597425
7005.062866174953
7005.100564672575
17143.786536183896
17143.82184350522
33179.36529436468
33179.40157549113
42525.79593689609
42525.835168267266
为什么你需要设置？可以重复吗？ – ilj
该设置用于删除重复项。这将很难通过重复数据进行分类。 – Axion004
你可以发布预期的输出。 –