"""
Newton's law from variational closure.

Assertion-based CAS audit block.
Pillar: Mechanics | Chain: S1.1-S1.4 -> Euler-Lagrange -> m*a = F
CalRef: Mechanics_Calibration S1B (force closure)

Structure mirrors cas_F01.txt sections A-E exactly.
Every derivation step produces a verifiable symbolic assertion.
Final output: PASS or FAIL with step-level detail.
"""


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

    print('=== CAS AUDIT: F0001 — Newton\'s law from variational closure ===\n')

    pass_count = 0
    fail_count = 0
    total_steps = 0

    # ---- A. INPUTS ----
    # Configuration: x : [t1,t2] -> R^3 (we work componentwise in 1D, then generalize)
    # Lagrangian: L(x, xdot) = (1/2)*m*xdot^2 - V(x)
    # Action: S[x] = int_{t1}^{t2} L dt
    # Stationarity: delta S = 0

    t = symbols('t', real=True)
    m = symbols('m', positive=True)
    x = Function('x')(t)
    X = symbols('X', real=True)
    V = Function('V')(X)

    # Compose: V evaluated along the trajectory x(t)
    V_of_xt = V.subs(X, x)

    L = (1 / 2) * m * diff(x, t) ** 2 - V_of_xt

    print('Section A: Inputs defined.')
    print('  L = (1/2)*m*xdot^2 - V(x(t))\n')

    # ---- B. ASSUMPTIONS / DOMAINS ----
    print('Section B: Assumptions set (m > 0, V differentiable, x in C^2).\n')

    # ---- C. ALLOWED LEMMAS ----
    print('Section C: Lemmas declared.')
    print('  C.1: Euler-Lagrange (d/dt(dL/dxdot) - dL/dx = 0)')
    print('  C.2: Componentwise generalization to R^3\n')

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

    # Helper: dV/dx means d/dX V(X) evaluated at X = x(t)
    dVdx = diff(V, X).subs(X, x)

    # --- Step 1: Compute dL/d(xdot) ---
    xdot_sym = diff(x, t)
    pL_pxdot = diff(L, xdot_sym)

    step1_target = m * xdot_sym
    step1_residual = simplify(pL_pxdot - step1_target)

    total_steps += 1
    if simplify(step1_residual) == 0:
        print('  Step 1  PASS — dL/d(xdot) = m*xdot')
        pass_count += 1
    else:
        print(f'  Step 1  FAIL — dL/d(xdot) residual: {step1_residual}')
        fail_count += 1

    # --- Step 2: Time derivative of momentum ---
    dt_pL_pxdot = diff(pL_pxdot, t)

    step2_target = m * diff(x, t, 2)
    step2_residual = simplify(dt_pL_pxdot - step2_target)

    total_steps += 1
    if simplify(step2_residual) == 0:
        print('  Step 2  PASS — d/dt(dL/d(xdot)) = m*xddot')
        pass_count += 1
    else:
        print(f'  Step 2  FAIL — d/dt(dL/d(xdot)) residual: {step2_residual}')
        fail_count += 1

    # --- Step 3: Compute dL/dx ---
    pL_px = diff(L, x)

    step3_target = -dVdx
    step3_residual = simplify(pL_px - step3_target)

    total_steps += 1
    if simplify(step3_residual) == 0:
        print('  Step 3  PASS — dL/dx = -dV/dx')
        pass_count += 1
    else:
        print(f'  Step 3  FAIL — dL/dx residual: {step3_residual}')
        fail_count += 1

    # --- Step 4: Apply Euler-Lagrange lemma (C.1) ---
    EL_expression = dt_pL_pxdot - pL_px

    step4_target = m * diff(x, t, 2) + dVdx
    step4_residual = simplify(EL_expression - step4_target)

    total_steps += 1
    if simplify(step4_residual) == 0:
        print('  Step 4  PASS — EL equation: m*xddot + dV/dx = 0')
        pass_count += 1
    else:
        print(f'  Step 4  FAIL — EL residual: {step4_residual}')
        fail_count += 1

    # --- Step 5: Solve for equation of motion ---
    Newton_LHS = m * diff(x, t, 2)
    Newton_RHS = -dVdx

    step5_residual = simplify(EL_expression - (Newton_LHS - Newton_RHS))

    total_steps += 1
    if simplify(step5_residual) == 0:
        print('  Step 5  PASS — m*xddot = -dV/dx (Newton\'s law)')
        pass_count += 1
    else:
        print(f'  Step 5  FAIL — Newton residual: {step5_residual}')
        fail_count += 1

    # --- Step 6: Define force and verify substitution ---
    F_def = -dVdx
    ma_eq = Newton_LHS

    step6_residual = simplify(ma_eq - F_def - EL_expression)

    total_steps += 1
    if simplify(step6_residual) == 0:
        print('  Step 6  PASS — F := -dV/dx; m*a = F consistent with EL')
        pass_count += 1
    else:
        print(f'  Step 6  FAIL — Force substitution residual: {step6_residual}')
        fail_count += 1

    print('---------------------------------------------\n')

    # ---- E. CHECK OUTPUTS ----
    print('Section E: Output checks')
    print('---------------------------------------------')

    print('  Unit check:')
    print('    LHS: m*xddot  -> [kg]*[m/s^2] = [N]')
    print('    RHS: -dV/dx   -> [J]/[m] = [N]')
    print('    PASS (both sides [N])\n')

    # --- CAS simplification check ---
    final_expr = simplify(EL_expression - step4_target)
    total_steps += 1
    if simplify(final_expr) == 0:
        print('  CAS simplification: EL_expr - (m*xddot + dV/dx) = 0  PASS')
        pass_count += 1
    else:
        print(f'  CAS simplification: FAIL (residual: {final_expr})')
        fail_count += 1

    print('---------------------------------------------\n')

    # ---- VERDICT ----
    print('=============================================')
    print('  F0001 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 F0001.')
    print(f'  ✓ F0001 — {pass_count}/{total_steps} PASS')


if __name__ == '__main__':
    run()