Symbolic derivation of 1-D FEM stiffness matrix

First import everything. Note how we now import sympy as sy to make it clear which functions we are using

import numpy as np
import matplotlib.pyplot as plt
import sympy as sy

Now we define the symbols, compute the shape function, and build the stiffness matrix of a single elelemt.

# symbols
#dx, x, k, xi = symbols('dx x k xi')

#N1  = 1-(x-xi)/(dx)
#N2  = (x-xi)/(dx)

#E11 = ???
#E12 = ???
#E21 = ???
#E22 = ???

#Ael = Matrix([[?, ?], [?, ?]])
#Ael

# hide
# symbols
dx, x, k, xi = sy.symbols('dx x k xi')

N1  = 1-(x-xi)/(dx)
N2  = (x-xi)/(dx)

E11 = sy.integrate(sy.diff(N1,x)*k*sy.diff(N1,x), ( x, 0, dx))
E12 = sy.integrate(sy.diff(N1,x)*k*sy.diff(N2,x), ( x, 0, dx))
E21 = sy.integrate(sy.diff(N2,x)*k*sy.diff(N1,x), ( x, 0, dx))
E22 = sy.integrate(sy.diff(N2,x)*k*sy.diff(N2,x), ( x, 0, dx))

Ael = sy.Matrix([[E11, E12], [E21, E22]])
Ael

$\displaystyle \left[\begin{matrix}\frac{k}{dx} & - \frac{k}{dx}\\- \frac{k}{dx} & \frac{k}{dx}\end{matrix}\right]$

Note how we made a little simplification here and integrated from 0 to dx and not from x(i) to x(i+1); the result is the same.

Now we build the global stiffness matrix by adding all element stiffness matrices to the big global stiffness matrix. Unfortunately, the indexing in sympy matrices is a bit different to normal numpy indexing.

We there have to use a little trick that is described here: https://www.reddit.com/r/learnpython/comments/f7uxn6/question_how_to_add_matrix_to_submatrix_in_sympy/

We simply identify the upper left index on the left-hand side and use list indexing on the right-hand side when adding the small matrices to the big one.

# numerical paramters
nel     = 4       # number of elements
nnodel  = 2       # nodes per element
nnod    = nel+1   # number of nodes

EL2NOD  = np.array([np.arange(1,nnod), np.arange(2,nnod+1)])  # connectivity matrix
print(EL2NOD)

# global matrix assemly
A = sy.zeros(nnod, nnod)
for iel in range(0,4):
    A[EL2NOD[0,iel]-1,EL2NOD[0,iel]-1] = A[list(EL2NOD[:,iel]-1), list(EL2NOD[:,iel]-1)] + Ael

A

[[1 2 3 4]
 [2 3 4 5]]

$\displaystyle \left[\begin{matrix}\frac{k}{dx} & - \frac{k}{dx} & 0 & 0 & 0\\- \frac{k}{dx} & \frac{2 k}{dx} & - \frac{k}{dx} & 0 & 0\\0 & - \frac{k}{dx} & \frac{2 k}{dx} & - \frac{k}{dx} & 0\\0 & 0 & - \frac{k}{dx} & \frac{2 k}{dx} & - \frac{k}{dx}\\0 & 0 & 0 & - \frac{k}{dx} & \frac{k}{dx}\end{matrix}\right]$

Does this matrix look familiar, or at least has a familiar structure?