# balk.py - reacties, dwarskracht, moment en doorbuiging van een balk op twee steunpunten (numeriek)
import numpy as np

L, E, b, h = 2.40, 11000e6, 0.045, 0.145           # vloerbalk (m), vuren C24 (Pa), doorsnede (m)
I, W = b * h**3 / 12, b * h**2 / 6
x = np.linspace(0, L, 2401)                         # 1 mm per stap


def integreer(y):
    """Trapeziumregel van 0 tot elke x (hoofdstuk 3)."""
    return np.concatenate([[0.0], np.cumsum((y[1:] + y[:-1]) / 2 * np.diff(x))])


def balk(q, puntlasten):
    """q: verdeelde last (N/m) in elk punt; puntlasten: [(F, a)]. Geeft RA, RB, V, M, w."""
    F_tot = integreer(q)[-1] + sum(F for F, a in puntlasten)
    M_A = integreer(q * x)[-1] + sum(F * a for F, a in puntlasten)
    RB = M_A / L                                    # momenten om A
    RA = F_tot - RB
    V = RA - integreer(q) - sum(F * (x >= a) for F, a in puntlasten)
    M = integreer(V)                                # dM/dx = V
    hoek = integreer(M / (E * I))                   # w'' = M/(EI): twee keer integreren
    w = integreer(hoek)
    w = w - w[-1] * x / L                           # een rechte optellen zodat w(L) = 0 (steunpunt B)
    return RA, RB, V, M, -w                         # doorbuiging naar beneden positief


reservoir = np.where(abs(x - 1.00) <= 0.32, 65 * 9.81 / 0.64, 0.0)       # 65 kg over 0,64 m
eigen = np.full_like(x, 75.0)
persoon = [(80 * 9.81, 1.40)]

for naam, gG, gQ in (("gebruik (karakteristiek)", 1.0, 1.0), ("rekenwaarde 1,35 G + 1,5 Q", 1.35, 1.5)):
    RA, RB, V, M, w = balk(gG * (eigen + reservoir), [(gQ * F, a) for F, a in persoon])
    i = np.argmax(M)
    print(f"{naam}: RA = {RA:6.1f} N, RB = {RB:6.1f} N, Mmax = {M[i]:6.1f} N m op {x[i]:.2f} m, "
          f"sigma = {M[i] / W / 1e6:4.2f} MPa, doorbuiging {1000 * w.max():4.2f} mm")