Riemann Zeta
function can be defined in the complex plane. Why beside its complications, in C++17
, GSL
and Boost C++
the argument can be only real(double)
?
Riemann zeta function, function useful in number theory for investigating properties of prime numbers. Written as ζ(x), it was originally defined as the infinite series ζ(x) = 1 + 2−x + 3−x + 4−x + ⋯. When x = 1, this series is called the harmonic series, which increases without bound—i.e., its sum is infinite.
The Riemann hypothesis, a formula related to the distribution of prime numbers, has remained unsolved for more than a century.
So tell me, how exactly are values for the Riemann Zeta Function computed? Actually, this formula: ξ(s)=π−s/2Γ(s2)ζ(s). The ξ(s) satisfies ξ(s)=ξ(1−s) and is defined everywhere except s=1.
The expression states that the sum of the zeta function is equal to the product of the reciprocal of one minus the reciprocal of primes to the power s. This astonishing connection laid the foundation for modern prime number theory, which from this point on used the zeta function ζ(s) as a way of studying primes.
Here is an Implementation in C++
const long double LOWER_THRESHOLD = 1.0e-6;
const long double UPPER_BOUND = 1.0e+4;
const int MAXNUM = 100;
std::complex<long double> zeta(const std::complex<long double>& s)
{
std::complex<long double> a_arr[MAXNUM + 1];
std::complex<long double> half(0.5, 0.0);
std::complex<long double> one(1.0, 0.0);
std::complex<long double> two(2.0, 0.0);
std::complex<long double> rev(-1.0, 0.0);
std::complex<long double> sum(0.0, 0.0);
std::complex<long double> prev(1.0e+20, 0.0);
a_arr[0] = half / (one - std::pow(two, (one - s))); //initialize with a_0 = 0.5 / (1 - 2^(1-s))
sum += a_arr[0];
for (int n = 1; n <= MAXNUM; n++)
{
std::complex<long double> nCplx(n, 0.0); //complex index
for (int k = 0; k < n; k++)
{
std::complex<long double> kCplx(k, 0.0); //complex index
a_arr[k] *= half * (nCplx / (nCplx - kCplx));
sum += a_arr[k];
}
a_arr[n] = (rev * a_arr[n - 1] * std::pow((nCplx / (nCplx + one)), s) / nCplx);
sum += a_arr[n];
if (std::abs(prev - sum) < LOWER_THRESHOLD)//If the difference is less than or equal to the threshold value, it is considered to be convergent and the calculation is terminated.
break;
if (std::abs(sum) > UPPER_BOUND)//doesn't work for large values, so it gets terminated when it exceeds UPPER_BOUND
break;
prev = sum;
}
return sum;
}
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With