Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion mega/__init__.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
from .op.function import Gamma, Haversine, JordanTotient
from .op.arithmethic import SigmaZ, EulerPhi
from .op.arithmethic import SigmaZ, EulerPhi, Chebyshev
from .op.tensor import Tensor

# utils folders
Expand Down
96 changes: 95 additions & 1 deletion mega/op/arithmethic.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@
# cython: language_level=3
# mega/op/arithmetic.pyx

from libc.math cimport pow, sqrt
from libc.math cimport pow, sqrt, log

cdef class SigmaZ:
"""
Expand Down Expand Up @@ -267,4 +267,98 @@ cdef class EulerPhi:
def __repr__(self):
return f"EulerPhi({self.n})"

cdef class Chebyshev:
"""
compute chebyshev function ϑ(x), with formula:

ϑ(x) = Σ_{p ≤ x} log(p)

Attribute:
x (double): upper bound of summation

Methods:
compute(): compute ϑ(x) using trial division
__repr__(): representation of chebyshev function
__getitem__(key): retrieve cache value or compute it
__setitem__(key, value): manually cache a compute result

Example:
>>> cheb = Chebyshev(10.0)
>>> cheb.compute()
5.347107530637821
"""
cdef public double x
cdef dict _cache

def __cinit__(self, double x):
"""
initialize chebyshev instance with input validation

Parameter:
x (double): must be >= 2
"""
if x < 2:
raise ValueError("x must be >= 2")
self.x = x
self._cache = {}

def __dealloc__(self):
self._cache.clear()

cpdef double compute(self):
"""
compute the first chebyshev function

this method using trial divison to check primality
then accumulate log(p) for each prime <= floor(x)

Return:
(double): result of chebyshev functon
"""
cdef double result = 0.0
cdef int i = 2
cdef int limit = <int>self.x
cdef int j, is_prime

while i <= limit:
is_prime = 1
for j in range(2, <int>sqrt(i) + 1):
if i % j == 0:
is_prime = 0
break
if is_prime:
result += log(i)
i += 1
return result

def __getitem__(self, double key):
"""
retrieve chebyshev function from internal cache, or compute
and store it

Parameter:
key (double): positive number >= 2

Return:
(double): Chebyshev(key)
"""
if key in self._cache:
return self._cache[key]
temp = Chebyshev(key)
result = temp.compute()
self._cache[key] = result
return result

def __setitem__(self, double key, double value):
"""
manually set cached value for chebyshev(key)

Parameter:
key(double): positive number >= 2
value(doube): precomputed value of chebyshev(key)
"""
self._cache[key] = value


def __repr__(self):
return f"Chebyshev(x={self.x})"
65 changes: 65 additions & 0 deletions test/arithmetic_test.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
import mega
import pytest
import math

SIGMA_VALUES: list = [
(1, 0, 1),
Expand Down Expand Up @@ -45,6 +46,16 @@
(105, 48),
]

CHEBYSHEV_VALUE: list = [
(2.0, math.log(2)),
(3.0, math.log(2) + math.log(3)),
(5.0, math.log(2) + math.log(3) + math.log(5)),
(7.0, sum(math.log(p) for p in [2, 3, 5, 7])),
(10.0, sum(math.log(p) for p in [2, 3, 5, 7])),
(11.0, sum(math.log(p) for p in [2, 3, 5, 7, 11])),
(20.0, sum(math.log(p) for p in [2, 3, 5, 7, 11, 13, 17, 19])),
]


def test_value_sigma() -> None:
for n, z, expected in SIGMA_VALUES:
Expand Down Expand Up @@ -117,3 +128,57 @@ def test_perfect_power_euler_phi() -> None:
def test_large_input_euler_phi() -> None:
phi = mega.EulerPhi(1_000_000)
assert phi.compute() == 400_000


def test_value_chebyshev() -> None:
for x, expected in CHEBYSHEV_VALUE:
ch = mega.Chebyshev(x)
result = ch.compute()
assert result == pytest.approx(expected, abs=1e-5), f"Chebyshev({x}).compute()"


def test_edge_cases_chebyshev() -> None:
res = mega.Chebyshev(10.9)
assert res.compute() == pytest.approx(
sum(math.log(p) for p in [2, 3, 5, 7]), abs=1e-10
)


def test_setitem_manual_cache_chebyshev() -> None:
res = mega.Chebyshev(10.0)
res[10] = 5.3471
assert res[10] == 5.3471


def test_large_input_chebyshev() -> None:
res = mega.Chebyshev(100.0)
prime = [
2,
3,
5,
7,
11,
13,
17,
19,
23,
29,
31,
37,
41,
43,
47,
53,
59,
61,
67,
71,
73,
79,
83,
89,
97,
]
expected = sum(math.log(p) for p in prime)
res_ch = res.compute()
assert res_ch == pytest.approx(expected, abs=1e-5)