考古学家程序员：挖掘现代算法的古代起源

def babylonian_sqrt(n, tolerance=0.0001):
    x = n  # Initial guess
    while True:
        next_x = 0.5 * (x + n / x)
        if abs(next_x - x) < tolerance:
            return next_x
        x = next_x

# Example usage
print(babylonian_sqrt(2))  # Output: ~1.4142135623730951

def euclidean_gcd(a, b):
    while b != 0:
        a, b = b, a % b
    return a

# Example usage
print(euclidean_gcd(48, 18))  # Output: 6

def sieve_of_eratosthenes(n):
    primes = [True] * (n + 1)
    primes[0] = primes[1] = False

    for i in range(2, int(n**0.5) + 1):
        if primes[i]:
            for j in range(i*i, n+1, i):
                primes[j] = False

    return [i for i in range(n+1) if primes[i]]

# Example usage
print(sieve_of_eratosthenes(30))
# Output: [2, 3, 5, 7, 11, 13, 17, 19, 23, 29]

def chinese_remainder_theorem(n, a):
    total = 0
    product = math.prod(n)

    for n_i, a_i in zip(n, a):
        p = product // n_i
        total += a_i * mod_inverse(p, n_i) * p

    return total % product

def mod_inverse(a, m):
    def egcd(a, b):
        if a == 0:
            return b, 0, 1
        else:
            g, y, x = egcd(b % a, a)
            return g, x - (b // a) * y, y

    g, x, _ = egcd(a, m)
    if g != 1:
        raise Exception('Modular inverse does not exist')
    else:
        return x % m

# Example usage
n = [3, 5, 7]
a = [2, 3, 2]
print(chinese_remainder_theorem(n, a))  # Output: 23

import math

def solve_quadratic(a, b, c):
    discriminant = b**2 - 4*a*c
    if discriminant < 0:
        return None
    elif discriminant == 0:
        return -b / (2*a)
    else:
        return ((-b + math.sqrt(discriminant)) / (2*a),
                (-b - math.sqrt(discriminant)) / (2*a))

# Example usage
print(solve_quadratic(1, 5, 6))  # Output: (-2.0, -3.0)

def newton_method(f, df, x0, epsilon=1e-7, max_iter=100):
    x = x0
    for _ in range(max_iter):
        fx = f(x)
        if abs(fx) < epsilon:
            return x
        dfx = df(x)
        if dfx == 0:
            return None  # Derivative is zero, method fails
        x = x - fx / dfx
    return None  # Max iterations reached, method fails

# Example: Finding the square root of 2
f = lambda x: x**2 - 2
df = lambda x: 2*x

result = newton_method(f, df, 1)
print(f"Square root of 2: {result}")
# Output: Square root of 2: 1.4142135623730951

class SimpleTuringMachine:
    def __init__(self, tape, initial_state, final_states, transitions):
        self.tape = list(tape)
        self.head = 0
        self.state = initial_state
        self.final_states = final_states
        self.transitions = transitions

    def step(self):
        if self.state in self.final_states:
            return False

        current_symbol = self.tape[self.head]
        if (self.state, current_symbol) not in self.transitions:
            return False

        new_state, new_symbol, direction = self.transitions[(self.state, current_symbol)]
        self.tape[self.head] = new_symbol
        self.state = new_state
        self.head += 1 if direction == 'R' else -1

        return True

    def run(self):
        while self.step():
            pass
        return ''.join(self.tape)

# Example: A Turing machine that replaces '0's with '1's
tm = SimpleTuringMachine(
    tape='00011010',
    initial_state='q0',
    final_states={'qf'},
    transitions={
        ('q0', '0'): ('q0', '1', 'R'),
        ('q0', '1'): ('q0', '1', 'R'),
        ('q0', ''): ('qf', '', 'R')
    }
)

result = tm.run()
print(f"Result: {result}")
# Output: Result: 11111111

def bubble_sort(arr):
    n = len(arr)
    for i in range(n):
        for j in range(0, n - i - 1):
            if arr[j] > arr[j + 1]:
                arr[j], arr[j + 1] = arr[j + 1], arr[j]
    return arr

# Example usage
print(bubble_sort([64, 34, 25, 12, 22, 11, 90]))
# Output: [11, 12, 22, 25, 34, 64, 90]

def quicksort(arr):
    if len(arr) <= 1:
        return arr
    pivot = arr[len(arr) // 2]
    left = [x for x in arr if x < pivot]
    middle = [x for x in arr if x == pivot]
    right = [x for x in arr if x > pivot]
    return quicksort(left) + middle + quicksort(right)

# Example usage
print(quicksort([64, 34, 25, 12, 22, 11, 90]))
# Output: [11, 12, 22, 25, 34, 64, 90]

import numpy as np

def pagerank(graph, damping=0.85, epsilon=1e-8):
    n = len(graph)

    # Create adjacency matrix
    A = np.array([[1 if j in graph[i] else 0 for j in range(n)] for i in range(n)])

    # Normalize adjacency matrix
    out_degree = np.sum(A, axis=1)
    A_hat = (A / out_degree[:, np.newaxis]).T

    # Initialize PageRank vector
    pr = np.ones(n) / n

    while True:
        pr_next = (1 - damping) / n + damping * A_hat @ pr
        if np.sum(np.abs(pr_next - pr)) < epsilon:
            return pr_next
        pr = pr_next

# Example usage
graph = {0: [1, 2], 1: [2], 2: [0, 1]}
result = pagerank(graph)
print(f"PageRank: {result}")
# Output: PageRank: [0.38461538 0.30769231 0.30769231]

import numpy as np

def linear_regression(X, y, learning_rate=0.01, iterations=1000):
    m, n = X.shape
    theta = np.zeros(n)

    for _ in range(iterations):
        h = np.dot(X, theta)
        gradient = np.dot(X.T, (h - y)) / m
        theta -= learning_rate * gradient

    return theta

# Example usage
X = np.array([[1, 1], [1, 2], [1, 3]])
y = np.array([1, 2, 3])
theta = linear_regression(X, y)
print(f"Theta: {theta}")
# Output: Theta: [0. 1.]

import numpy as np

def quantum_coin_flip():
    # Create a superposition state
    state = np.array([1/np.sqrt(2), 1/np.sqrt(2)])

    # Measure the state
    measurement = np.random.choice(['Heads', 'Tails'], p=state**2)

    return measurement

# Example usage
results = [quantum_coin_flip() for _ in range(1000)]
heads_count = results.count('Heads')
tails_count = results.count('Tails')

print(f"Heads: {heads_count}, Tails: {tails_count}")
# Output will be approximately: Heads: 500, Tails: 500