Balken-Element B2

B: Beam, 2: Zwei Knoten

Tags: Web-App

Lineares System

../../../_images/beam-fem1.png

Hinweis

Zu dieser Gleichung:

Andere Notationen

Es gibt verschiedene Notationen für diese Aussagen. Denn man kann in den Aussagen zum Linearen System die Symbole, die Zählrichtungen und die Reihenfolge wählen.

../../../_images/beam-fem_1.png

Grüner Kasten: Reihenfolge wie bei Dr. Bernd Klein. Roter Kasten: Symbole, Zählrichtungen und Reihenfolge wie bei Dr. Manfred Hahn.

Details

Umformung des Gleichungssystems in Matrix-Schreibweise.

  1. Ausgangspunkt ist:

    \[\begin{split}\begin{bmatrix}4 l^{2} & - 6 l & 2 l^{2} & 6 l\\- 6 l & 12 & - 6 l & -12\\2 l^{2} & - 6 l & 4 l^{2} & 6 l\\6 l & -12 & 6 l & 12\end{bmatrix} \begin{bmatrix}\psi_{1}\\w_{1}\\\psi_{2}\\w_{2}\end{bmatrix} = \begin{bmatrix}M_{1}\\F_{1}\\M_{2}\\F_{2}\end{bmatrix}\end{split}\]
  2. Reihenfolge der Aussagen ändern liefert:

    \[\begin{split}\begin{bmatrix}12 & - 6 l & -12 & - 6 l\\- 6 l & 4 l^{2} & 6 l & 2 l^{2}\\-12 & 6 l & 12 & 6 l\\- 6 l & 2 l^{2} & 6 l & 4 l^{2}\end{bmatrix} \begin{bmatrix}w_{1}\\\psi_{1}\\w_{2}\\\psi_{2}\end{bmatrix} = \begin{bmatrix}F_{1}\\M_{1}\\F_{2}\\M_{2}\end{bmatrix}\end{split}\]
  3. Umbenennung der Symbole liefert:

    \[\begin{split}\begin{bmatrix}12 & - 6 l & -12 & - 6 l\\- 6 l & 4 l^{2} & 6 l & 2 l^{2}\\-12 & 6 l & 12 & 6 l\\- 6 l & 2 l^{2} & 6 l & 4 l^{2}\end{bmatrix} \begin{bmatrix}- u_{1}\\\varphi_{1}\\- u_{2}\\\varphi_{2}\end{bmatrix} = \begin{bmatrix}- Q_{1}\\B_{1}\\- Q_{2}\\B_{2}\end{bmatrix}\end{split}\]
  4. Minuszeichen in die Matrix bringen liefert:

    \[\begin{split}\begin{bmatrix}-12 & - 6 l & 12 & - 6 l\\6 l & 4 l^{2} & - 6 l & 2 l^{2}\\12 & 6 l & -12 & 6 l\\6 l & 2 l^{2} & - 6 l & 4 l^{2}\end{bmatrix} \begin{bmatrix}u_{1}\\\varphi_{1}\\u_{2}\\\varphi_{2}\end{bmatrix} = \begin{bmatrix}- Q_{1}\\B_{1}\\- Q_{2}\\B_{2}\end{bmatrix}\end{split}\]
  5. Multiplizieren der 1. und 3. Zeile mit -1 liefert:

    \[\begin{split}\begin{bmatrix}12 & 6 l & -12 & 6 l\\6 l & 4 l^{2} & - 6 l & 2 l^{2}\\-12 & - 6 l & 12 & - 6 l\\6 l & 2 l^{2} & - 6 l & 4 l^{2}\end{bmatrix} \begin{bmatrix}u_{1}\\\phi_{1}\\u_{2}\\\phi_{2}\end{bmatrix} = \begin{bmatrix}Q_{1}\\B_{1}\\Q_{2}\\B_{2}\end{bmatrix}\end{split}\]

SymPy

  1. Copy: Source Code (siehe unten) aufklappen und kopieren.

  2. Paste: Einfügen als Python-Notebook z.B. auf:

  3. Play: Ausführen.

Source Code

from sympy import *

EI, l = var("EI, l")

p1H, w1H, p2H, w2H = var("phi1, u1, phi2, u2")
M1H, F1H, M2H, F2H = var("B1, Q1, B2, Q2")
p1, w1, p2, w2 = var("psi1, w1, psi2, w2")
M1, F1, M2, F2 = var("M1, F1, M2, F2")

sub_list=[
    (p1,  p1H),
    (w1, -w1H),
    (p2,  p2H),
    (w2, -w2H),
    (M1,  M1H),
    (F1, -F1H),
    (M2,  M2H),
    (F2, -F2H),
    ]

u = Matrix([p1, w1, p2, w2])
f = Matrix([M1, F1, M2, F2])

def K(EI, l):
    # l = Element length
    l2, l3 = l*l, l*l*l
    K = EI/l3
    K *= Matrix([
        [  4*l2 ,  -6*l ,  2*l2 ,   6*l  ],
        [ -6*l  ,  12   , -6*l  , -12    ],
        [  2*l2 ,  -6*l ,  4*l2 ,   6*l  ],
        [  6*l  , -12   ,  6*l    , 12   ],
        ])
    return K

print("\nStandard Notation with EI/l³ omitted:")
K = K(EI, l)
fac = EI/(l**3)
pprint(K/fac)
pprint(u)
pprint(f)

pprint("\n1. Change Order:")
def swap(K, u, f, n, m):
    n = n-1
    m = m-1
    K.row_swap(n,m)
    K.col_swap(n,m)
    u.row_swap(n,m)
    f.row_swap(n,m)

# Shift to Hahn Notation:
swap(K, u, f, 1, 2)
swap(K, u, f, 3, 4)
pprint(K/fac)
pprint(u)
pprint(f)

pprint("\n2. Substitute Symbols in u and f:")
u = u.subs(sub_list)
f = f.subs(sub_list)
pprint(u)
pprint(f)

pprint("\n3. Move Negative Signs of u into Matrix:")
K.col_op(0, lambda i, j: i*(-1))
K.col_op(2, lambda i, j: i*(-1))
u.row_op(0, lambda i, j: i*(-1))
u.row_op(2, lambda i, j: i*(-1))
pprint(K/fac)
pprint(u)

pprint("\n4. Multiply Row 1 and 3 by -1:")
K.row_op(0, lambda i, j: i*(-1))
K.row_op(2, lambda i, j: i*(-1))
f.row_op(0, lambda i, j: i*(-1))
f.row_op(2, lambda i, j: i*(-1))
pprint(K/fac)
pprint(f)

# Standard Notation with EI/l³ omitted:
# ⎡   2           2     ⎤
# ⎢4⋅l   -6⋅l  2⋅l   6⋅l⎥
# ⎢                     ⎥
# ⎢-6⋅l   12   -6⋅l  -12⎥
# ⎢                     ⎥
# ⎢   2           2     ⎥
# ⎢2⋅l   -6⋅l  4⋅l   6⋅l⎥
# ⎢                     ⎥
# ⎣6⋅l   -12   6⋅l   12 ⎦
# ⎡ψ₁⎤
# ⎢  ⎥
# ⎢w₁⎥
# ⎢  ⎥
# ⎢ψ₂⎥
# ⎢  ⎥
# ⎣w₂⎦
# ⎡M₁⎤
# ⎢  ⎥
# ⎢F₁⎥
# ⎢  ⎥
# ⎢M₂⎥
# ⎢  ⎥
# ⎣F₂⎦

# 1. Change Order:
# ⎡ 12   -6⋅l  -12  -6⋅l⎤
# ⎢                     ⎥
# ⎢         2          2⎥
# ⎢-6⋅l  4⋅l   6⋅l  2⋅l ⎥
# ⎢                     ⎥
# ⎢-12   6⋅l   12   6⋅l ⎥
# ⎢                     ⎥
# ⎢         2          2⎥
# ⎣-6⋅l  2⋅l   6⋅l  4⋅l ⎦
# ⎡w₁⎤
# ⎢  ⎥
# ⎢ψ₁⎥
# ⎢  ⎥
# ⎢w₂⎥
# ⎢  ⎥
# ⎣ψ₂⎦
# ⎡F₁⎤
# ⎢  ⎥
# ⎢M₁⎥
# ⎢  ⎥
# ⎢F₂⎥
# ⎢  ⎥
# ⎣M₂⎦

# 2. Substitute Symbols in u and f:
# ⎡-u₁⎤
# ⎢   ⎥
# ⎢φ₁ ⎥
# ⎢   ⎥
# ⎢-u₂⎥
# ⎢   ⎥
# ⎣φ₂ ⎦
# ⎡-Q₁⎤
# ⎢   ⎥
# ⎢B₁ ⎥
# ⎢   ⎥
# ⎢-Q₂⎥
# ⎢   ⎥
# ⎣B₂ ⎦

# 3. Move Negative Signs of u into Matrix:
# ⎡-12  -6⋅l   12   -6⋅l⎤
# ⎢                     ⎥
# ⎢        2           2⎥
# ⎢6⋅l  4⋅l   -6⋅l  2⋅l ⎥
# ⎢                     ⎥
# ⎢12   6⋅l   -12   6⋅l ⎥
# ⎢                     ⎥
# ⎢        2           2⎥
# ⎣6⋅l  2⋅l   -6⋅l  4⋅l ⎦
# ⎡u₁⎤
# ⎢  ⎥
# ⎢φ₁⎥
# ⎢  ⎥
# ⎢u₂⎥
# ⎢  ⎥
# ⎣φ₂⎦

# 4. Multiply Row 1 and 3 by -1:
# ⎡12   6⋅l   -12   6⋅l ⎤
# ⎢                     ⎥
# ⎢        2           2⎥
# ⎢6⋅l  4⋅l   -6⋅l  2⋅l ⎥
# ⎢                     ⎥
# ⎢-12  -6⋅l   12   -6⋅l⎥
# ⎢                     ⎥
# ⎢        2           2⎥
# ⎣6⋅l  2⋅l   -6⋅l  4⋅l ⎦
# ⎡Q₁⎤
# ⎢  ⎥
# ⎢B₁⎥
# ⎢  ⎥
# ⎢Q₂⎥
# ⎢  ⎥
# ⎣B₂⎦

  • Lieber mit eigenem SymPy und offline? Es gibt z.B. Anaconda und Miniconda.

  • Statt SymPy lieber anderes CAS? Auswahl hier.

Erweiterte Steifigkeitsmatrix

Die erweiterte Element-Steifigkeitsmatrix bietet den Vorteil, dass die Freiheitsgrade \((\psi_1, w_1, \psi_2, w_2)\) sichtbar sind, auf die sich die Einträge der Steifigkeitsmatrix beziehen.

\begin{align*} k^| = \tfrac{EI}{l^3} \left[ \begin{array}{cccc|c} 4 l^2 & -6 l & 2 l^2 & 6 l & \psi_1 \\ & 12 & -6 l & -12 & w_1 \\ & & 4 l^2 & 6 l & \psi_2 \\ \mathsf{sym} & & & 12 & w_2 \\ \end{array} \right] \end{align*}

Wahl der Knoten-Reihenfolge

Wie beim Stab-Element lässt sich die Knoten-Reihenfolge vertauschen, ohne dass sich die Einträge der Steifigkeitsmatrix ändern:

\begin{align*} \underbrace{ \tfrac{EI}{l^3} \begin{bmatrix} 4 l^2 & -6 l & 2 l^2 & 6 l \\ -6 l & 12 & -6 l & -12 \\ 2 l^2 & -6 l & 4 l^2 & 6 l \\ 6 l & -12 & 6 l & 12 \\ \end{bmatrix} }_k \begin{bmatrix} \psi_1 \\ w_1 \\ \psi_2 \\ w_2 \\ \end{bmatrix} &= \begin{bmatrix} M_1 \\ F_1 \\ M_2 \\ F_2 \\ \end{bmatrix} + q \begin{bmatrix} - \tfrac{1}{12}l^2 \\ \tfrac12 l\\ \tfrac{1}{12}l^2\\ \tfrac12 l \\ \end{bmatrix} \\ \underbrace{ \tfrac{EI}{l^3} \begin{bmatrix} 4 l^2 & -6 l & 2 l^2 & 6 l \\ -6 l & 12 & -6 l & -12 \\ 2 l^2 & -6 l & 4 l^2 & 6 l \\ 6 l & -12 & 6 l & 12 \\ \end{bmatrix} }_k \begin{bmatrix} \psi_2 \\ w_2 \\ \psi_1 \\ w_1 \\ \end{bmatrix} &= \begin{bmatrix} M_2 \\ F_2 \\ M_1 \\ F_1 \\ \end{bmatrix} + q \begin{bmatrix} \tfrac{1}{12}l^2\\ \tfrac12 l \\ - \tfrac{1}{12}l^2 \\ \tfrac12 l\\ \end{bmatrix} \end{align*}

Zugehörige erweiterte Steifigkeitsmatrizen:

\begin{align*} k^| &= \tfrac{EI}{l^3} \left[ \begin{array}{cccc|c} 4 l^2 & -6 l & 2 l^2 & 6 l & \psi_1 \\ -6 l & 12 & -6 l & -12 & w_1 \\ 2 l^2 & -6 l & 4 l^2 & 6 l & \psi_2 \\ 6 l & -12 & 6 l & 12 & w_2 \\ \end{array} \right] \\ k^| &= \tfrac{EI}{l^3} \left[ \begin{array}{cccc|c} 4 l^2 & -6 l & 2 l^2 & 6 l & \psi_2 \\ -6 l & 12 & -6 l & -12 & w_2 \\ 2 l^2 & -6 l & 4 l^2 & 6 l & \psi_1 \\ 6 l & -12 & 6 l & 12 & w_1 \\ \end{array} \right] \end{align*}

Postprocessing / Interpolation

\(N_1, N_2, N_3, N_4, N_5\) sind die dimensionslosen Interpolationsfunktionen:

\[\begin{split}N_1 &= - \xi^{3} + 2 \xi^{2} - \xi\\ N_2 &= 2 \xi^{3} - 3 \xi^{2} + 1 \\ N_3 &= - \xi^{3} + \xi^{2} \\ N_4 &= - 2 \xi^{3} + 3 \xi^{2} \\ N_5 &= \tfrac{\xi^{4}}{24} - \tfrac{\xi^{3}}{12} + \tfrac{\xi^{2}}{24}\\\end{split}\]

Querverschiebung w

Interpolierte Querverschiebung in Matrix-Schreibweise:

\[\begin{split}w &= w_{\psi_1} + w_{w_1} + w_{\psi_2} + w_{w_2} + w_q \\ &= \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} l N_1 \\ N_2 \\ l N_3 \\ N_4 \\ \tfrac{l^4}{B} N_5 \end{bmatrix} \\ &= \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} l \left( - \xi^{3} + 2 \xi^{2} - \xi \right) \\ 2 \xi^{3} - 3 \xi^{2} + 1 \\ l \left( - \xi^{3} + \xi^{2} \right) \\ - 2 \xi^{3} + 3 \xi^{2} \\ \tfrac{l^4}{B} \left( \tfrac{\xi^{4}}{24} - \tfrac{\xi^{3}}{12} + \tfrac{\xi^{2}}{24} \right) \end{bmatrix}\end{split}\]

Bemerkung

  • Eine \(1 \times 5\)-Matrix multipliziert mit einer \(5 \times 1\)-Matrix ergibt eine \(1 \times 1\)-Matrix, also eine Zahl.

  • Man kann das Produkt auch so schreiben:

    \[\begin{split}w &= \underbrace{ \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} }_{1 \times 5} \underbrace{ \begin{bmatrix} l N_1 \\ N_2 \\ l N_3 \\ N_4 \\ \tfrac{l^4}{B} N_5 \end{bmatrix} }_{5 \times 1} \\ &= \underbrace{ \begin{bmatrix} l N_1 & N_2 & l N_3 & N_4 & \tfrac{l^4}{B} N_5 \end{bmatrix} }_{1 \times 5} \underbrace{ \begin{bmatrix} \psi_1 \\ w_1 \\ \psi_2 \\ w_2 \\ q \end{bmatrix} }_{5 \times 1}\end{split}\]
  • \(w\) ist von der Ordnung \(\xi^4\).

Verdrehung ψ, Biegemoment M, Querkraft Q

Mit \(\psi=-w'\) und \(M = - B w''\) und \(Q = - B w'''\) gilt:

\[\begin{split}\psi &= - \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} l N_1' \\ N_2' \\ l N_3' \\ N_4' \\ \tfrac{l^4}{B} N_5' \end{bmatrix} \\ &= \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} 3 \xi^{2} - 4 \xi + 1\\- \tfrac{1}{l} \left(6 \xi^{2} - 6 \xi\right)\\ 3 \xi^{2} - 2 \xi\\ - \tfrac{1}{l} \left(- 6 \xi^{2} + 6 \xi\right)\\ - \tfrac{l^{3}}{B} \left(\tfrac{\xi^{3}}{6} - \tfrac{\xi^{2}}{4} + \tfrac{\xi}{12}\right) \end{bmatrix}\end{split}\]

\[\begin{split}M&= - B \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} l N_1'' \\ N_2'' \\ l N_3'' \\ N_4'' \\ \tfrac{l^4}{B} N_5'' \end{bmatrix} \\ &= \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} - \tfrac{B}{l} \left(- 6 \xi + 4\right) \\ - \tfrac{B}{l^{2}} \left(12 \xi - 6\right) \\ - \tfrac{B}{l} \left(- 6 \xi + 2\right) \\ - \tfrac{B}{l^{2}} \left(- 12 \xi + 6\right) \\ - l^{2} \left(\tfrac{\xi^{2}}{2} - \tfrac{\xi}{2} + \tfrac{1}{12}\right) \end{bmatrix}\end{split}\]

\[\begin{split}Q&= - B \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} l N_1''' \\ N_2''' \\ l N_3''' \\ N_4''' \\ \tfrac{l^4}{B} N_5''' \end{bmatrix} \\ &= \begin{bmatrix} \psi_1 & w_1 & \psi_2 & w_2 & q \end{bmatrix} \begin{bmatrix} \tfrac{6 B}{l^{2}}\\ - \tfrac{12 B}{l^{3}}\\ \tfrac{6 B}{l^{2}}\\ \tfrac{12 B}{l^{3}}\\ - l \left(\xi - \tfrac{1}{2}\right) \end{bmatrix}\end{split}\]

Graphische Darstellung

Hinweise zu den Diagrammen

Zwei Interaktive Diagramme

  • Schieberegler verschieben = Variieren des Knoten-Werts.

  • In Legende Funktion klicken = Aus- und Einblenden des Funktionsgraphen.

  • Horizontale Achse: Werte werden nach rechts größer - wie üblich.

  • Vertikale Achse: Werte werden nach unten größer.

Äquivalente Knotenlasten

Tags: Video

Video

../../../_images/arbitrary-load.png

Knotenlasten aufgrund äußerer Lasten \(q_\xi, M_\xi, F_\xi\) zwischen den Element-Knoten. Oben: Lineares System. Mitte: Generische Lasten. Unten: Äquivalente Knotenlasten,

Details

Generische Lasten -> Generische Rechte Seite

  • Momente \(M_1, M_2\) sowie Kräfte \(F_1, F_2\) an den Knoten.

  • Pro Element konstante verteilte \(q\)-Last.

  • Die FEM-Lösung gleicht der klassischen Lösung, denn das Balken-Element wurde definiert, um diese generischen Lasten abzubilden.

  • Die rechte Seite \(f\) im Gleichungssystem \(ku = f\) ist:

\[\begin{split}f= \begin{bmatrix} M_1 \\ F_1 \\ M_2 \\ F_2 \\ \end{bmatrix} + q \begin{bmatrix} - \tfrac{1}{12}l^2 \\ \tfrac12 l\\ \tfrac{1}{12}l^2\\ \tfrac12 l \\ \end{bmatrix}\end{split}\]

Sonstige Lasten -> Rechte Seite mit Äquivalenten Knotenlasten

  • Momente \(M_\xi\) an beliebiger Position \(\xi_M\) sowie Kräfte \(F_\xi\) an beliebiger Position \(\xi_F\).

  • Beliebig verteilte \(q_\xi\)-Last.

  • Die FEM-Lösung gleicht nicht der klassischen Lösung.

  • Die sonstigen Lasten werden umgerechnet in die Äquivalenten Knotenlasten. 1 Diese Äquivalenten Knotenlasten bilden die rechte Seite \(f\):

\[\begin{split}f = M_\xi \begin{bmatrix} - l N_1' \\ - N_2' \\ - l N_3' \\ - N_4' \\ \end{bmatrix}_{\xi=\xi_M} + F_\xi \begin{bmatrix} l N_1 \\ N_2 \\ l N_3 \\ N_4 \\ \end{bmatrix}_{\xi=\xi_F} + \int_{\xi=0}^1 q_\xi \begin{bmatrix} l N_1 \\ N_2 \\ l N_3 \\ N_4 \\ \end{bmatrix} l \mathsf{d}\xi\end{split}\]

mit Interpolationsfunktionen und deren Ableitungen.

Fußnote:

1

Die bei der Umrechnung verwendete Annahme ist, dass die virtuelle äußere Arbeit der sonstigen Lasten gleich ist der virtuellen äußeren Arbeit der daraus berechneten Äquivalenten Knotenlasten.

Beispiele

../../../_images/ex-1.png

Spaltenmatrizen = Äquivalente Knotenlasten: Berechnung siehe Python-Programm.

SymPy

  1. Copy: Source Code (siehe unten) aufklappen und kopieren.

  2. Paste: Einfügen als Python-Notebook z.B. auf:

  3. Play: Ausführen.

Source Code

from sympy import *

half = S(1)/2

M, F, q, l, xi = var("M, F, q, l, xi")

N1 = -xi**3 + 2*xi**2 - xi
N2 = 2*xi**3 - 3*xi**2 + 1
N3 = -xi**3 + xi**2
N4 = -2*xi**3 + 3*xi**2

N = Matrix([l*N1, N2, l*N3, N4])
# N':
Nx = diff(N, xi) / l

pprint("\nConcentrated Moment Mξ at ξ = 1/2:")
f  = - M * Nx.subs(xi, half)
pprint("\nf / M:")
tmp = f/M
pprint(tmp)

pprint("\nConcentrated Force Fξ at ξ = 1/2:")
f  = F * N.subs(xi, half)
pprint("\nf / F:")
tmp = f/F
pprint(tmp)

pprint("\nPiecewise linear distributed load qξ:")
# Piecewise defined function:
q1, q2 = 2*q*xi, 2*q*(1 - xi)

# Piecewise integration:
f  = integrate(q1*N*l, (xi, 0, half))
f += integrate(q2*N*l, (xi, half, 1))
pprint("\nf / (-ql):")
tmp = f/(-q*l)
pprint(tmp)

pprint("\nSpecial Cases:")
pprint("\nConstant distributed load qξ=q:")
f  = integrate(q*N*l, (xi, 0, 1))
pprint("\nf / (-ql):")
tmp = f/(-q*l)
pprint(tmp)

pprint("\nM1 and F1 at ξ = 0 and M2 and F2 at ξ = 1:")
M1, F1, M2, F2 = var("M1, F1, M2, F2")
f  = - M1 * Nx.subs(xi, 0)
f +=   F1 *  N.subs(xi, 0)
f += - M2 * Nx.subs(xi, 1)
f +=   F2 *  N.subs(xi, 1)
pprint("\nf:")
tmp = f
pprint(tmp)

  • Lieber mit eigenem SymPy und offline? Es gibt z.B. Anaconda und Miniconda.

  • Statt SymPy lieber anderes CAS? Auswahl hier.

Lösung einiger Aufgaben

SymPy

  1. Copy: Source Code (siehe unten) aufklappen und kopieren.

  2. Paste: Einfügen als Python-Notebook z.B. auf:

  3. Play: Ausführen.

Source Code

# -*- coding: utf-8 -*-
from sympy.physics.units import *
from sympy import *

# Rounding:
import decimal
from decimal import Decimal as DX
from copy import deepcopy
def iso_round(obj, pv, rounding=decimal.ROUND_HALF_EVEN):
    import sympy
    """
    Rounding acc. to DIN EN ISO 80000-1:2013-08
    place value = Rundestellenwert
    """
    assert pv in set([
        # place value     #  round to:
        "1",              #  closest (even) integer
        "0.1",            #  1st digit after decimal
        "0.01",           #  2nd
        "0.001",          #  3rd
        "0.0001",         #  4th
        "0.00001",        #  5th
        "0.000001",       #  6th
        "0.0000001",      #  7th
        "0.00000001",     #  8th
        "0.000000001",    #  9th
        "0.0000000001",   # 10th
        ])
    objc = deepcopy(obj)
    try:
        tmp = DX(str(float(objc)))
        objc = tmp.quantize(DX(pv), rounding=rounding)
    except:
        for i in range(len(objc)):
            tmp = DX(str(float(objc[i])))
            objc[i] = tmp.quantize(DX(pv), rounding=rounding)
    return objc

# Units:
(mm, cm)  =  ( m/1000, m/100 )
Newton    =  kg*m/s**2
kN        =  10**3*Newton
Pa        =  Newton/m**2
MPa       =  10**6*Pa
GPa       =  10**9*Pa
deg       =  pi/180
half      =  S(1)/2


kwargs = {}
kwargs["mat_str"] = "bmatrix"
kwargs["mat_delim"] = ""

a, c, l, q, F, M, EI = var("a, c, l, q, F, M, EI")

quants = False
# quants = True
# sublist=[
#     ( M,   10 *Newton*m              ),
#     ( a,    1 *m                     ),
#     ( EI, 200*GPa * 2*mm*6*mm**3/ 12 ),
#     ]

def K(EI, l):
    # l = Element length
    l2, l3 = l*l, l*l*l
    K = EI/l3
    K *= Matrix([
        [  4*l2 ,  -6*l ,  2*l2 ,   6*l  ],
        [ -6*l  ,  12   , -6*l  , -12    ],
        [  2*l2 ,  -6*l ,  4*l2 ,   6*l  ],
        [  6*l  , -12   ,  6*l    , 12   ],
        ])
    return K

def K2(EI, l):
    # l = Element length
    l2, l3 = l*l, l*l*l
    K = EI/l3
    K *= Matrix(
    [
    [  4*l2 ,  -6*l ,  2*l2 ,   6*l ,  0    ,   0    ],
    [ -6*l  ,  12   , -6*l  , -12   ,  0    ,   0    ],
    [  2*l2 ,  -6*l ,  8*l2 ,   0   ,  2*l2 ,   6*l  ],
    [  6*l  , -12   ,  0    ,  24   , -6*l  , -12    ],
    [  0    ,   0   ,  2*l2 ,  -6*l ,  4*l2 ,   6*l  ],
    [  0    ,   0   ,  6*l  , -12   ,  6*l  ,  12    ],
    ]
    )
    return K

def Fq(l):
    return q*l/2
def Mq(l):
    return q*l*l / 12

    

p0,w0,p1,w1,p2,w2 = var("ψ₀,w₀,ψ₁,w₁,ψ₂,w₂")
M0,F0,M1,F1,M2,F2 = var("M₀,F₀,M₁,F₁,M₂,F₂")

# task = "5.7"
# task = "B2.A"
# task = "B2.B"
# task = "B2.C"
# task = "B2.D"
task = "B2.E"
task = "B2.E-MIRROR"
task = "B2.E-POINT"
# task = "B2.F"
# task = "B2.G"
# task = "B2.H"
# task = "B2.I"
# task = "B2.J"

pprint("\nTask: "+task)

if task == "5.7":
    K = K(EI, a)
    #                    1                             2
    #
    #                                                  M
    #                    |------------------------------
    #                                    a             Λ
    M = var("M")
    #
    u = Matrix([       0,0,                          p2,0      ])
    f = Matrix([      M1,F1,                          M,F2     ])
    unks =     [      M1,F1,                         p2,F2     ]

elif task == "B2.A":
    K = K2(EI, a)
    #                    0              1              2
    #                                   F1             F2
    #                    |------------------------------
    u = Matrix([        0,0,          p1,w1,         p2,w2       ])
    f = Matrix([       M0,F0,          0,F1,          0,F2       ])
    unks =     [       M0,F0,         p1,w1,         p2,w2       ]

elif task == "B2.B":
    K = K(EI, a)
    #                    1                             2
    #                     qqqqqqqqqqqqqqqqqqqqqqqqqqqqq
    #                    |------------------------------
    #                                   a              Λ
    #
    u = Matrix([        0,0,                        p2,0         ])
    f = Matrix([ M1-Mq(a),F1+Fq(a),            0+Mq(a),F2+Fq(a)  ])
    unks =     [       M1,F1,                       p2,F2        ]

elif task == "B2.C":
    K = K2(EI, a/4)
    #                    0              1              2
    #                    ◆ qqqqqqqqqqqqqqqqqqqqqqqqqqqq
    #                    -------------------------------
    #                    ◆                             Λ
    #
    u = Matrix([        0,w0,        p1,w1,        p2,0             ])
    f = Matrix([M0-Mq(a/4),Fq(a/4), 0,2*Fq(a/4), Mq(a/4),F2+Fq(a/4) ])
    unks =     [       M0,w0,        p1,w1,        p2,F2            ]

elif task == "B2.D":
    K = K(EI, a)
    #                    1                             2
    #                    F
    #                     qqqqqqqqqqqqqqqqqqqqqqqqqqqqq
    #                    ------------------------------|
    #                    Λ              a
    #                    |
    #                    | cw₁
    #
    u = Matrix([       p1,w1,                        0,0          ])
    f = Matrix([  0-Mq(a),F+Fq(a),           M2+Mq(a),F2+Fq(a)    ])
    # Spring instead of q:
    # f = Matrix([        0,F-c*w1,                       M2,F2   ])
    unks =     [       p1,w1,                       M2,F2         ]

elif task == "B2.E":
    K = K2(EI, a)
    Ml = -q*a*a / 12
    Mr =  q*a*a / 12
    Fq =  q*a/2
    #                    0              1              2
    #                                    qqqqqqqqqqqqqq
    #                    |-----------------------------|
    #
    u = Matrix([        0,0,          p1,w1,             0,0         ])
    f = Matrix([       M0,F0,         Ml,Fq,         M2+Mr,F2+Fq     ])
    unks =     [       M0,F0,         p1,w1,            M2,F2        ]

elif task == "B2.E-MIRROR" or task == "B2.E-POINT":
    K = K(EI, a)

    Ml = -q*a*a / 12 / 2
    Mr =  q*a*a / 12 / 2
    Fq =  q*a/2 / 2

    if task == "B2.E-MIRROR":
        #
        #                    0               1
        #                     q/2   q/2   q/2
        #                    |----------------  Mirror Symm.: ψ1=0, F1=0
        #
        u = Matrix([        0,0,              0,w1                     ])
        f = Matrix([    M0+Ml,F0+Fq,      M1+Mr,0+Fq                   ])
        unks =     [       M0,F0,            w1,M1                     ]
    elif task == "B2.E-POINT":
        #                    0               1
        #                    -q/2  -q/2  -q/2                 
        #                    |----------------  Point Symm.: w1=0, M1=0
        #
        u = Matrix([        0,0,          p1,0                        ])
        f = Matrix([    M0-Ml,F0-Fq,    0-Mr,F1-Fq                    ])
        unks =     [       M0,F0,         p1,F1                       ]

elif task == "B2.F":
    K = K2(EI, a)
    #                    0              1              2
    #                     qqqqqqqqqqqqqq
    #                    ------------------------------|
    #                    Λ
    #
    u = Matrix([       p0,0,          p1,w1,          0,0        ])
    f = Matrix([  0-Mq(a),F0+Fq(a), 0+Mq(a),0+Fq(a), M2,F2       ])
    unks =     [       p0,F0,         p1,w1,         M2,F2       ]

elif task == "B2.G":
    Mq1, Fq1 = -q/12*a*a, q/2*a
    Mq2, Fq2 = -Mq1, Fq1
    #                    0                    1        2
    #                                          qqqqqqqq
    #                    |------------------------------
    #                               2a        Λ    a
    #
    u = Matrix([          0,0,              p1,0,         p2,w2  ])
    f = Matrix([         M0,F0,             Mq1,Fq1+F1,  Mq2,Fq2 ])
    unks =     [         M0,F0,             p1,F1,        p2,w2  ]

    # Elem. 1
    l = 2*a
    l2 = l*l
    l3 = l*l*l
    K = EI/l3 * Matrix(
        [
        [  4*l2 ,   -6*l ,   2*l2 ,    6*l ,   0   ,   0   ],
        [ -6*l  ,   12   ,  -6*l  ,  -12   ,   0   ,   0   ],
        [  2*l2 ,   -6*l ,   4*l2 ,    6*l ,   0   ,   0   ],
        [  6*l  ,  -12   ,   6*l  ,   12   ,   0   ,   0   ],
        [  0    ,    0   ,   0    ,    0   ,   0   ,   0   ],
        [  0    ,    0   ,   0    ,    0   ,   0   ,   0   ],
        ]
        )
    # Elem. 2
    l = a
    l2 = l*l
    l3 = l*l*l
    K += EI/l3 * Matrix(
        [
        [  0   ,   0   ,  0    ,   0   ,  0    ,   0       ],
        [  0   ,   0   ,  0    ,   0   ,  0    ,   0       ],
        [  0   ,   0   ,  4*l2 ,  -6*l ,  2*l2 ,   6*l     ],
        [  0   ,   0   , -6*l  ,  12   , -6*l  , -12       ],
        [  0   ,   0   ,  2*l2 ,  -6*l ,  4*l2 ,   6*l     ],
        [  0   ,   0   ,  6*l  , -12   ,  6*l  ,  12       ],
        ]
        )

elif task == "B2.H":
    M, c = var("M, c")
    Fc = c*w1
    K = K2(EI, a/2)
    #                    0              1              2
    #                                   M              M
    #                    |------------------------------
    #                                   Λ
    #                                   |  Fc = c w₁
    #                                   |
    #
    u = Matrix([        0,0,          p1,w1,          p2,w2      ])
    f = Matrix([       M0,F0,          M,-Fc,          M,0       ])
    unks =     [       M0,F0,         p1,w1,          p2,w2      ]

elif task == "B2.I":
    K = K(EI, a/2)
    #                    1                             2
    #                   F/2
    #
    #                    ◆
    #                    -------------------------------
    #                    ◆             a/2             Λ
    #
    #
    #
    u = Matrix([       0,w1,                       p2,0    ])
    f = Matrix([      M1,F/2,                       0,F2   ])
    unks =     [      M1,w1,                       p2,F2   ]

elif task == "B2.J":
    K = K(EI, a)
    #                    1                             2
    #                                    |
    #                                    | F
    #                                    V            -M
    #                    |------------------------------
    #                                    a             Λ
    #
    #
    #
    # Equiv. Nodal Loads due to F:
    M, F = var("M, F")
    Me1, Fe1, Me2, Fe2 = -F*a/8, F/2, F*a/8, F/2
    #
    u = Matrix([       0,0,                        p2,0      ])
    f = Matrix([  M1+Me1,F1 + Fe1,             -M+Me2,F2+Fe2 ])
    unks =     [      M1,F1,                       p2,F2     ]


eq = Eq(K*u , f)
sol = solve(eq, unks)
pprint("\nSolution")
pprint(sol)

exit()
for s in unks:
    for tmp in [s, sol[s]]:
        tmp = tmp.expand()
        tmp = tmp.simplify()
        pprint(tmp)
        # pprint(latex(tmp,**kwargs))
        pprint("\n")
        if task =="5.7" and quants == True:
            tmp = sol[s]
            tmp = tmp.subs(sub_list)
            if s == M1:
                pprint("\nM1 / Nm:")
                tmp /= Newton*m
            elif s == F1:
                pprint("\nF1 / N:")
                tmp /= Newton
            elif s == p2:
                pprint("\np2 / rad (unrealistic task):")
                tmp = tmp
            elif s == F2:
                pprint("\nF2 / N:")
                tmp /= Newton
            tmp = iso_round(tmp,"1.0")
            pprint(tmp)

if task == "B2.H":
    pprint("\nw2 for ca³ = 10 EI:")
    tmp = sol[w2]
    tmp = tmp.subs(EI, c*a**3/10)
    tmp = N(tmp,3)
    pprint(tmp)

  • Lieber mit eigenem SymPy und offline? Es gibt z.B. Anaconda und Miniconda.

  • Statt SymPy lieber anderes CAS? Auswahl hier.