def gamma(n):
    import math
    import Science
    def gamma_ln(z):
        if z <= 0:
            raise ValueError("Invalid input: z must be positive.")
        elif z < 0.5:
            # Use the reflection formula for better accuracy
            # ln(Gamma(z)) = ln(pi) - ln(sin(pi * z)) - ln(Gamma(1 - z))
            ln_pi = Science.ln(Science.pi)
            sin_pi_z = Science.sin(Science.pi * z)
            ln_sin_pi_z = Science.ln(sin_pi_z)
            ln_gamma_1_minus_z = gamma_ln_lanczos(1 - z)
            ln_gamma_z = ln_pi - ln_sin_pi_z - ln_gamma_1_minus_z
            return ln_gamma_z
        elif z >= 100:
            # Use Stirling's approximation for large z
            ln_gamma_z = gamma_ln_stirling(z)
            return ln_gamma_z
        else:
            # Use Lanczos approximation directly
            ln_gamma_z = gamma_ln_lanczos(z)
            return ln_gamma_z

    def gamma_ln_lanczos(z):
        # Lanczos approximation coefficients
        p = [
            0.99999999999980993,
            676.5203681218851,
            -1259.1392167224028,
            771.32342877765313,
            -176.61502916214059,
            12.507343278686905,
            -0.13857109526572012,
            9.9843695780195716e-6,
            1.5056327351493116e-7,
        ]

        # Adjust z for the approximation
        z_minus_one = z - 1.0
        x = p[0]  # Initial value for the sum

        # Manually compute the sum without using loops
        x += p[1] / (z_minus_one + 1.0)
        x += p[2] / (z_minus_one + 2.0)
        x += p[3] / (z_minus_one + 3.0)
        x += p[4] / (z_minus_one + 4.0)
        x += p[5] / (z_minus_one + 5.0)
        x += p[6] / (z_minus_one + 6.0)
        x += p[7] / (z_minus_one + 7.0)
        x += p[8] / (z_minus_one + 8.0)

        t = z_minus_one + 7.5  # g + 0.5 where g = 7

        # Compute ln(Gamma(z)) using math functions
        ln_sqrt_two_pi = 0.5 * Science.ln(2.0 * Science.pi)
        ln_gamma = (
            (z_minus_one + 0.5) * Science.ln(t)
            - t
            + Science.ln(x)
            + ln_sqrt_two_pi
        )
        return ln_gamma

    def gamma_ln_stirling(z):
        # Stirling's series approximation for ln(Gamma(z)) for large z
        ln_z = Science.ln(z)
        ln_sqrt_2_pi = 0.5 * Science.ln(2.0 * Science.pi)
        z_inv = 1.0 / z

        # Compute inverse powers of z manually
        z_inv1 = z_inv
        z_inv2 = z_inv1 * z_inv1
        z_inv3 = z_inv2 * z_inv1
        z_inv5 = z_inv3 * z_inv2
        z_inv7 = z_inv5 * z_inv2
        z_inv9 = z_inv7 * z_inv2
        z_inv11 = z_inv9 * z_inv2

        # Coefficients for the Stirling series expansion
        c1 = 1.0 / 12.0
        c2 = -1.0 / 360.0
        c3 = 1.0 / 1260.0
        c4 = -1.0 / 1680.0
        c5 = 1.0 / 1188.0
        c6 = -691.0 / 360360.0

        # Compute ln(Gamma(z)) using the Stirling series expansion
        ln_gamma = (
            (z - 0.5) * ln_z
            - z
            + ln_sqrt_2_pi
            + c1 * z_inv1
            + c2 * z_inv3
            + c3 * z_inv5
            + c4 * z_inv7
            + c5 * z_inv9
            + c6 * z_inv11
        )
        return ln_gamma

    # Compute the actual gamma function value
    ln_gamma_n = gamma_ln(n)
    # Handle large values to prevent overflow
    try:
        result = Science.exp(ln_gamma_n)
    except OverflowError:
        result = float('inf')  # Return infinity if the result overflows
    return result

import math
n = 171
print("math.gamma({}) = {}".format(n, math.gamma(n)))
print("gamma({}) = {}".format(n, gamma(n)))

math.gamma(171) = 7.257415615307998e+306
gamma(171) = 6.993130794311394e+306

попробуй исправить ошибку не используя встроенные библиотеки кроме math , нельзя использовать циклы и рекурсии, нельзя использовать встроенные ** и ^ функции , также нельзя использовать float.hex подобные функции