import numpy as np
from mpmath import mp
import argparse

# Set very high precision
mp.dps = 300

# Target numbers (ordinal, value)
target_numbers = [
    (1, 1), (2, 3), (3, 7), (4, 8), (5, 21), (6, 49), (7, 76), (8, 224), (9, 467), (10, 514),
    (11, 1155), (12, 2683), (13, 5216), (14, 10544), (15, 26867), (16, 51510),
    (17, 95823), (18, 198669), (19, 357535), (20, 863317), (21, 1811764),
    (22, 3007503), (23, 5598802), (24, 14428676), (25, 33185509),
    (26, 54538862), (27, 111949941), (28, 227634408), (29, 400708894),
    (30, 1033162084), (31, 2102388551), (32, 3093472814), (33, 7137437912),
    (34, 14133072157), (35, 20112871792), (36, 42387769980), (37, 100251560595),
    (38, 146971536592), (39, 323724968937), (40, 1003651412950),
    (41, 1458252205147), (42, 2895374552463), (43, 7409811047825),
    (44, 15404761757071), (45, 19996463086597), (46, 51408670348612),
    (47, 119666659114170), (48, 191206974700443), (49, 409118905032525),
    (50, 611140496167764), (51, 2058769515153876), (52, 4216495639600700),
    (53, 6763683971478124), (54, 9974455244496707), (55, 30045390491869460),
    (56, 44218742292676575), (57, 138245758910846492), (58, 199976667976342049),
    (59, 525070384258266191), (60, 1135041350219496382), (61, 1425787542618654982),
    (62, 3908372542507822062), (63, 8993229949524469768),
    (64, 17799667357578236628), (65, 30568377312064202855),
    (66, 46346217550346335726), (67, 132656943602386256302),
    (68, 219898266213316039825), (70, 970436974005023690481)
]

# Provided alpha values
alpha_values = {
    20: -5.509999999999367, 21: -4.8099999999993575, 22: 6.690000000000806, 23: 4.309999999879222,
    24: -5.109999999999362, 25: -2.4933006701008145, 26: 23.19000000000104, 27: 50.01000000165623,
    28: 50.01000000165623, 29: 4.490000000000775, 30: -4.30999999999935, 31: -4.8036514600379405,
    32: 2.7900000000007505, 33: 10.110000011655663, 34: 8.496609620547714, 35: 1.6905278399975268,
    36: 2.2900000000007434, 37: 6.506918969897979, 38: 1.8912248499932935, 39: 2.9099999998792025,
    40: -2.809999999999329, 41: 5.2945524899731184, 42: 5.690000000000792, 43: -4.691035100114615,
    44: -3.9099999999993447, 45: 2.5900000000007477, 46: 50.01000000165623, 47: -4.50218804004683,
    48: 6.909999999879259, 49: 28.510000011655915, 50: 2.495913989964804, 51: -2.5900000001208756,
    52: -2.5099999999993248, 53: 50.01000000165623, 54: 2.4900000000007463, 55: -5.00999999999936,
    56: 3.8099999998792153, 57: -2.1943843400942242, 58: 8.890000000000837, 59: -3.009999999999332,
    60: -2.2099999999993205, 61: 3.1099999998792054, 62: -6.000280060058447, 63: -2.490000000120874,
    64: -2.809999999999329, 65: -19.689999988344763, 66: 2.7931615399815364, 67: -4.709999999999356,
    68: 8.692996380248756, 70: -24.301635149307526
}

def calculate_prediction(puzzle_number, alpha):
    # Get all previous data points
    ordinals = np.array([x[0] for x in target_numbers if x[0] < puzzle_number], dtype=float)
    values = np.array([x[1] for x in target_numbers if x[0] < puzzle_number], dtype=float)
   
    if len(ordinals) < 2:
        print(f"Not enough data points to make a prediction for puzzle {puzzle_number}")
        return None, None, None
   
    log_values = np.array([float(mp.log(val)) for val in values], dtype=float)
    weights = np.linspace(0.5, 1, len(ordinals))
    coefficients = np.polyfit(ordinals, log_values, 1, w=weights)
   
    a = mp.exp(coefficients[1])
    b = mp.exp(coefficients[0])
   
    correction = (((2 ** puzzle_number) - 1) - (2 ** (puzzle_number - 1))) / alpha
    predicted_value = (a * mp.power(b, puzzle_number)) - correction
   
    return predicted_value, a, b

def interpolate_alpha(puzzle_number):
    """Interpolate alpha for missing puzzle numbers"""
    known_puzzles = sorted(alpha_values.keys())
   
    if puzzle_number < min(known_puzzles) or puzzle_number > max(known_puzzles):
        return 0.0  # Default value if outside known range
   
    # Find closest lower and higher puzzle numbers
    lower = max([p for p in known_puzzles if p < puzzle_number])
    higher = min([p for p in known_puzzles if p > puzzle_number])
   
    # Linear interpolation
    alpha_range = alpha_values[higher] - alpha_values[lower]
    position = (puzzle_number - lower) / (higher - lower)
    return alpha_values[lower] + (alpha_range * position)

def format_mp(value, decimals=None):
    """Format mpmath value with specified decimal places"""
    if decimals is not None:
        return f"{float(value):,.{decimals}f}".replace(',', '')
    return mp.nstr(value, mp.dps)

def main(puzzle_number):
    if puzzle_number == 69:
        # Special case for puzzle 69 - use interpolated alpha
        alpha = interpolate_alpha(69)
        predicted_value, a, b = calculate_prediction(69, alpha)
       
        print("\nPuzzle 69 Prediction:")
        print(f"Using interpolated alpha: {alpha:.6f}")
        print(f"Base exponential model: {format_mp(a)} * {format_mp(b)}^n")
        print(f"Predicted Value: {format_mp(predicted_value, 2)}")
       
        # Compare with n=70
        n70_value = next(x[1] for x in target_numbers if x[0] == 70)
        growth_rate = n70_value / predicted_value
        print(f"Growth to n=70: {float(growth_rate):.6f}x")
        print(f"Actual n=70 value: {format_mp(n70_value)}")
        return
   
    if puzzle_number not in alpha_values:
        print(f"No alpha value available for puzzle {puzzle_number}")
        return
   
    alpha = alpha_values[puzzle_number]
    predicted_value, a, b = calculate_prediction(puzzle_number, alpha)
   
    if predicted_value is None:
        return
   
    real_value = next((x[1] for x in target_numbers if x[0] == puzzle_number), None)
   
    print(f"\nPuzzle {puzzle_number} Calculation:")
    print(f"Using alpha: {alpha}")
    print(f"Base exponential model: {format_mp(a)} * {format_mp(b)}^n")
   
    if real_value is not None:
        print(f"Real Value: {format_mp(real_value)}")
    print(f"Predicted Value: {format_mp(predicted_value, 2)}")
   
    if real_value is not None:
        error_percentage = abs((predicted_value - real_value) / real_value) * 100
        print(f"Percentage error: {format_mp(error_percentage, 8)}%")

if __name__ == "__main__":
    parser = argparse.ArgumentParser(description="Calculate puzzle values using provided alpha values.")
    parser.add_argument("puzzle_number", type=int, help="Puzzle number to calculate")
    args = parser.parse_args()
    main(args.puzzle_number)