"""Work-energy theorem (Delta K = integral F dx).

Assertion-based CAS audit block.
Pillar: Mechanics | Chain: Newton -> chain rule -> power -> work integral -> Delta K
CalRef: Mechanics Math Appendix S2.2-2.4

Verifies:
  1. Chain rule: d/dt(v^2/2) = v*dv/dt
  2. Newton -> power relation: F*v = dK/dt
  3. Work integral via change of variables
  4. Fundamental theorem of calculus: integral of dK/dt dt = Delta K
  5. Full closure: Delta K = integral F dx
"""


def run():
    from sympy import symbols, Function, diff, simplify, integrate, Rational

    print("=== CAS AUDIT: F0011 — Work-energy theorem ===\n")

    pass_count = 0
    fail_count = 0
    total_steps = 0

    # ---- A. INPUTS ----
    t = symbols("t", real=True)
    m = symbols("m", positive=True)
    x = Function("x")(t)

    v = diff(x, t)
    a = diff(v, t)

    print("Section A: Inputs defined.")
    print("  x(t), v=dx/dt, a=dv/dt, F=ma, K=(1/2)mv^2\n")

    # ---- B. ASSUMPTIONS / DOMAINS ----
    print("Section B: 1D, m>0, x(t) in C^2, nonrelativistic.\n")

    # ---- C. ALLOWED LEMMAS ----
    print("Section C: Lemmas declared.")
    print("  C.1: d/dt(v^2/2) = v*dv/dt (chain rule)")
    print("  C.2: F = m*dv/dt")
    print("  C.3: dx = v*dt (change of variables)\n")

    # ---- D. STEP LOG ----
    print("Section D: Step log")
    print("---------------------------------------------")

    # Step 1: Chain rule
    lhs_chain = diff(v**2 / 2, t)
    rhs_chain = v * diff(v, t)
    step1_residual = simplify(lhs_chain - rhs_chain)

    total_steps += 1
    if simplify(step1_residual) == 0:
        print("  Step 1  PASS — d/dt(v^2/2) = v*dv/dt (chain rule)")
        pass_count += 1
    else:
        print(f"  Step 1  FAIL — Chain rule residual: {step1_residual}")
        fail_count += 1

    # Step 2: dK/dt
    K = Rational(1, 2) * m * v**2
    dK_dt = diff(K, t)
    expected_dK_dt = m * v * diff(v, t)
    step2_residual = simplify(dK_dt - expected_dK_dt)

    total_steps += 1
    if simplify(step2_residual) == 0:
        print("  Step 2  PASS — dK/dt = m*v*dv/dt")
        pass_count += 1
    else:
        print(f"  Step 2  FAIL — dK/dt residual: {step2_residual}")
        fail_count += 1

    # Step 3: Newton gives F*v = dK/dt
    F_newton = m * a
    power = F_newton * v
    step3_residual = simplify(power - dK_dt)

    total_steps += 1
    if simplify(step3_residual) == 0:
        print("  Step 3  PASS — F*v = dK/dt (power = rate of KE change)")
        pass_count += 1
    else:
        print(f"  Step 3  FAIL — Power residual: {step3_residual}")
        fail_count += 1

    # Step 4: K(2v) = 4*K(v)
    K_double_v = K.subs(v, 2*v)
    step4_residual = simplify(K_double_v - 4*K)

    total_steps += 1
    if simplify(step4_residual) == 0:
        print("  Step 4  PASS — K(2v) = 4*K(v) (quadratic in v)")
        pass_count += 1
    else:
        print(f"  Step 4  FAIL — Quadratic check residual: {step4_residual}")
        fail_count += 1

    # Step 5: Antiderivative check
    step5_residual = simplify(diff(integrate(dK_dt, t), t) - dK_dt)

    total_steps += 1
    if simplify(step5_residual) == 0:
        print("  Step 5  PASS — d/dt(int(dK/dt)) = dK/dt (FTC, general x(t))")
        pass_count += 1
    else:
        print(f"  Step 5  FAIL — FTC residual: {step5_residual}")
        fail_count += 1

    # Step 6: General definite integral
    t1, t2 = symbols("t1 t2", real=True)
    K_at_t2 = K.subs(t, t2)
    K_at_t1 = K.subs(t, t1)
    Delta_K_general_t = simplify(K_at_t2 - K_at_t1)

    W_general = integrate(dK_dt, (t, t1, t2))
    W_general = simplify(W_general)

    step6_residual = simplify(Delta_K_general_t - W_general)

    total_steps += 1
    if simplify(step6_residual) == 0:
        print("  Step 6  PASS — int(dK/dt, t1, t2) = K(t2)-K(t1) (general, FTC)")
        pass_count += 1
    else:
        print(f"  Step 6  FAIL — General work-energy residual: {step6_residual}")
        fail_count += 1

    # Step 7: Concrete trajectory
    a0 = symbols("a0", real=True, positive=True)
    T_final = symbols("T_final", positive=True)
    x_test = Rational(1, 2) * a0 * t**2
    v_test = diff(x_test, t)
    K_test = Rational(1, 2) * m * v_test**2
    dK_dt_test = diff(K_test, t)
    Fv_test = m * diff(v_test, t) * v_test

    K_at_T = K_test.subs(t, T_final)
    K_at_0 = K_test.subs(t, 0)
    Delta_K_concrete = simplify(K_at_T - K_at_0)
    W_concrete = simplify(integrate(Fv_test, (t, 0, T_final)))

    step7_residual = simplify(Delta_K_concrete - W_concrete)

    total_steps += 1
    if simplify(step7_residual) == 0:
        print("  Step 7  PASS — Concrete: x=(1/2)a0*t^2 confirms Delta K = W")
        pass_count += 1
    else:
        print(f"  Step 7  FAIL — Concrete residual: {step7_residual}")
        fail_count += 1

    # Step 8: Delta K = (1/2)*m*(v2^2 - v1^2)
    v1, v2 = symbols("v1 v2", real=True, nonnegative=True)
    K1 = Rational(1, 2) * m * v1**2
    K2 = Rational(1, 2) * m * v2**2
    Delta_K_general = K2 - K1
    expected_Delta_K = Rational(1, 2) * m * (v2**2 - v1**2)
    step8_DK_residual = simplify(Delta_K_general - expected_Delta_K)

    total_steps += 1
    if simplify(step8_DK_residual) == 0:
        print("  Step 8  PASS — Delta K = (1/2)*m*(v2^2 - v1^2)")
        pass_count += 1
    else:
        print(f"  Step 8  FAIL — Delta K form residual: {step8_DK_residual}")
        fail_count += 1

    # Step 9: K non-negativity
    K_at_zero = (Rational(1, 2) * m * v1**2).subs(v1, 0)
    K_positive_v = (Rational(1, 2) * m * v1**2).subs(v1, 1)

    total_steps += 1
    if simplify(K_at_zero) == 0 and simplify(K_positive_v) > 0:
        print("  Step 9  PASS — K(0)=0 and K(v>0)>0 (non-negative, m>0, v^2>=0)")
        pass_count += 1
    else:
        print("  Step 9  FAIL — K non-negativity check")
        fail_count += 1

    print("---------------------------------------------\n")

    # ---- E. CHECK OUTPUTS ----
    print("Section E: Output checks")
    print("---------------------------------------------")
    print("  Unit check:")
    print("    K = (1/2)*m*v^2: [kg]*[m/s]^2 = [J]")
    print("    W = F*dx: [N]*[m] = [J]")
    print("    F*v: [N]*[m/s] = [W] (power)")
    print("    dK/dt: [J/s] = [W]")
    print("    PASS\n")

    # Self-test: wrong sign
    F_wrong = -m * a
    power_wrong = F_wrong * v
    wrong_residual = simplify(power_wrong - dK_dt)

    total_steps += 1
    if simplify(wrong_residual) != 0:
        print("  Self-test: F=-ma gives F*v != dK/dt (wrong sign detected)  PASS")
        pass_count += 1
    else:
        print("  Self-test: FAIL (wrong Newton sign not detected)")
        fail_count += 1

    # Self-test: quantify
    expected_wrong = -2 * dK_dt
    wrong_quant = simplify(wrong_residual - expected_wrong)

    total_steps += 1
    if simplify(wrong_quant) == 0:
        print("  Self-test: wrong - correct = -2*dK/dt (quantified)  PASS")
        pass_count += 1
    else:
        print(f"  Self-test: FAIL (wrong residual = {wrong_residual})")
        fail_count += 1

    print("---------------------------------------------\n")

    # ---- VERDICT ----
    print("=============================================")
    print("  F0011 AUDIT RESULT")
    print(f"  Steps: {total_steps}  |  Pass: {pass_count}  |  Fail: {fail_count}")
    if fail_count == 0:
        print("  STATUS: *** PASS ***")
    else:
        print(f"  STATUS: *** FAIL *** ({fail_count} step(s) failed)")
    print("=============================================")
    print("Audit complete for F0011.")
    print(f"  ✓ F0011 — {pass_count}/{total_steps} PASS")


if __name__ == "__main__":
    run()