```python
# https://colab.research.google.com/github/kassbohm/tm-snippets/blob/master/ipynb/TM_A/TM_2/lagrange.ipynb

from sympy.physics.units import *
from sympy import *

a0, a1, a2 = var("a0, a1, a2")
b0, b1, b2 = var("b0, b1, b2")
c0, c1, c2 = var("c0, c1, c2")

xi = var("xi")

L0 = a0 + a1*xi + a2*xi*xi
L1 = b0 + b1*xi + b2*xi*xi
L2 = c0 + c1*xi + c2*xi*xi

e1 = Eq(L0.subs(xi,S(0)/2), 1)
e2 = Eq(L0.subs(xi,S(1)/2), 0)
e3 = Eq(L0.subs(xi,S(2)/2), 0)
e4 = Eq(L1.subs(xi,S(0)/2), 0)
e5 = Eq(L1.subs(xi,S(1)/2), 1)
e6 = Eq(L1.subs(xi,S(2)/2), 0)
e7 = Eq(L2.subs(xi,S(0)/2), 0)
e8 = Eq(L2.subs(xi,S(1)/2), 0)
e9 = Eq(L2.subs(xi,S(2)/2), 1)

eqns = [e1,e2,e3,e4,e5,e6,e7,e8,e9]
unks = [a0,a1,a2,b0,b1,b2,c0,c1,c2]

sol = solve(eqns,unks)
# pprint(sol)
# exit()


l, EA, n = var("l, EA, n")

(x0, x1, x2) = (0, l/2, l)

x = var("x")
u0, u1, u2 = var("u0, u1, u2")
Ax = var("Ax")

# Lagrange-polynomials:
L0  = (x  - x1)*(x  - x2)
L0 /= (x0 - x1)*(x0 - x2)

L1  = (x  - x0)*(x  - x2)
L1 /= (x1 - x0)*(x1 - x2)

L2  = (x  - x0)*(x  - x1)
L2 /= (x2 - x0)*(x2 - x1)

pprint("\nShape-Functions wiht ξ = x/l and Integrals:")
for Li in [L0, L1, L2]:
    pprint("\n")
    # Li = Li.simplify()
    # pprint(Li)
    Li = Li.expand()
    tmp = Li.subs(x/l, xi)
    # pprint(Li)
    pprint(tmp)
    # pprint(latex(Li))
    I = integrate(Li, (x, 0, l))
    pprint(I)
    # pprint("\nChecking x=0, x=l/2, x=l:")
    # pprint(Li.subs(x,0))
    # pprint(Li.subs(x,L/2))
    # pprint(Li.subs(x,L))

L0p = diff(L0,x)
L1p = diff(L1,x)
L2p = diff(L2,x)

pprint("\n\nDerivatives of Shape Functions:")
for Lip in [L0p, L1p, L2p]:
    Lip = Lip.simplify()
    pprint(Lip)

pprint("\n\nStiffness Matrix k / EA:")
k00, k01, k02 = L0p*L0p, L1p*L0p, L2p*L0p
k11, k12 = L1p*L1p, L2p*L1p
k22 = L2p*L2p

k00 = integrate(k00, (x, 0, l))
k01 = integrate(k01, (x, 0, l))
k02 = integrate(k02, (x, 0, l))
k11 = integrate(k11, (x, 0, l))
k12 = integrate(k12, (x, 0, l))
k22 = integrate(k22, (x, 0, l))

k = Matrix([
    [k00, k01, k02],
    [k01, k11, k12],
    [k02, k12, k22]
    ])
pprint(k)
k*=EA

f = Matrix([Ax + n*l/6, 2*n*l/3, n*l/6])
u = Matrix([0, u1, u2])
eq = Eq(k*u,f)

pprint("\n\nSolution:")
sol = solve(eq, [Ax, u1, u2], dict=True)
sol = sol[0]
Ax, u1, u2 = sol[Ax], sol[u1], sol[u2]
pprint("\nAx:")
pprint(Ax)
pprint("\nu1:")
pprint(u1)
pprint("\nu2:")
pprint(u2)

u = Matrix([0, u1, u2])
Lp = Matrix([L0p, L1p, L2p])
eps = u.dot(Lp)

pprint("\nε:")
pprint("\nx=0:")
tmp = eps.subs(x,0)
pprint(tmp)
pprint("\nx=l:")
tmp = eps.subs(x,l)
pprint(tmp)

pprint("\nu(x) / (nl² / EA):")
u0 = 0
ux = L0*u0 + L1*u1 + L2*u2
tmp = ux / ( n*l**2 / EA )
tmp = tmp.simplify()
tmp = tmp.subs(x, xi*l)
tmp = tmp.expand()
pprint(tmp)


# Shape-Functions and Integrals:
# (l - 2⋅x)⋅(l - x)
# ─────────────────
#          2
#         l
# l
# ─
# 6
# 4⋅x⋅(l - x)
# ───────────
#       2
#      l
# 2⋅l
# ───
#  3
# x⋅(-l + 2⋅x)
# ────────────
#       2
#      l
# l
# ─
# 6
#
#
# Derivatives of Shape Functions:
# -3⋅l + 4⋅x
# ──────────
#      2
#     l
# 4⋅(l - 2⋅x)
# ───────────
#       2
#      l
# -l + 4⋅x
# ────────
#     2
#    l
#
#
# Stiffness Matrix k / EA:
# ⎡ 7   -8    1 ⎤
# ⎢───  ───  ───⎥
# ⎢3⋅l  3⋅l  3⋅l⎥
# ⎢             ⎥
# ⎢-8    16  -8 ⎥
# ⎢───  ───  ───⎥
# ⎢3⋅l  3⋅l  3⋅l⎥
# ⎢             ⎥
# ⎢ 1   -8    7 ⎥
# ⎢───  ───  ───⎥
# ⎣3⋅l  3⋅l  3⋅l⎦
#
#
# Solution:
#
# Ax:
# -l⋅n
#
# u1:
#    2
# 3⋅l ⋅n
# ──────
#  8⋅EA
#
# u2:
#  2
# l ⋅n
# ────
# 2⋅EA

```

75c31fb01bc0cf9d15b3614cb0e9e4f4aa45078b
ipynb
Jupyter Notebook
ipynb/TM_A/TM_2/lagrange.ipynb
kassbohm/tm-snippets
5e0621ba2470116e54643b740d1b68b9f28bff12
MIT
ipynb/TM_A/TM_2/lagrange.ipynb
kassbohm/tm-snippets
5e0621ba2470116e54643b740d1b68b9f28bff12
MIT
ipynb/TM_A/TM_2/lagrange.ipynb
kassbohm/tm-snippets
5e0621ba2470116e54643b740d1b68b9f28bff12
MIT
Qwen/Qwen-72B
1. YES
2. YES