# Solving of KU=F leads to numpy.linalg.LinAlgError: Singular matrix

Solver.py:

import numpy as np
import math

E=3000000
A=2.0
truss = [
np.array([[0, 0], [0, 8]]),
np.array([[0, 0], [8, 8]]),
np.array([[0, 8], [8, 8]])
]
lengths=[]
for i, t in enumerate(truss):
lengths.append(math.hypot(t[0][0] - t[1][0], t[0][1] - t[1][1]))

n_joints = len(truss) + 1  # assuming the first point is also a joint
n_reactions = 0

k_w=[]
k_w.append((E * A / lengths[0]) * np.array([
[1, 0, -1, 0],
[0, 0, 0, 0],
[-1, 0, 1, 0],
[0, 0, 0, 0]
]))

k_w.append((E * A / lengths[1]) * np.array([
[1, 0, -1, 0],
[0, 0, 0, 0],
[-1, 0, 1, 0],
[0, 0, 0, 0]
]))

k_w.append((E * A / lengths[2]) * np.array([
[1, 0, -1, 0],
[0, 0, 0, 0],
[-1, 0, 1, 0],
[0, 0, 0, 0]
]))

# Initialize the global stiffness matrix
K = np.zeros((4*len(truss), 4*len(truss)))
theta = []

for i in range(len(truss)):
dy = truss[i][0][0] - truss[i][1][0]
dx = truss[i][0][1] - truss[i][1][1]
theta.append(math.atan2(dy, dx))

# Assemble the global stiffness matrix
for i, t in enumerate(truss):
L = np.array([
[1, 0, 0, 0],
[0, 1, 0, 0],
[0, 0, math.cos(theta[i]), math.sin(theta[i])],
[0, 0, -math.sin(theta[i]), math.cos(theta[i])]
])

# Apply the transformation to the local stiffness matrix
k_local = np.dot(np.dot(L.T, k_w[i]), L)

# Add the transformed local stiffness matrix to the global stiffness matrix
K[i * 4:(i + 1) * 4, i * 4:(i + 1) * 4] += k_local

# Check if the current truss member has a support reaction
if i == 0 or i == 4:
n_reactions += 2

P = np.zeros(4*len(truss))
P[0]=1
P[1]=1
# Set boundary conditions
K[2, :] = 0
K[3, :] = 0
K[:, 2] = 0
K[:, 3] = 0
K[2, 2] = 1
K[3, 3] = 1
P[2] = 0
P[3] = 0

# Add rows and columns for reaction forces
#K = np.vstack((np.hstack((K, np.zeros((n_reactions, 4*len(truss))))), np.zeros((4*len(truss), n_reactions))))
#P = np.hstack((P, np.zeros(n_reactions)))

# Impose equilibrium
#for i in range(n_reactions):
#    row = 4*len(truss) + i
#    K[row, row] = 1
#    if i < 2:
#        K[row, i*4+1] = -1
#    else:
#        K[row, i*4+3] = -1

# Solve for the displacements
#displacements = np.linalg.solve(K, P)

# Add rows and columns for reaction forces
#n_reactions = 2  # assuming there are two reaction forces
#K = np.vstack((np.hstack((K, np.zeros((n_reactions, 4*len(truss))))), np.zeros((4*len(truss), n_reactions))))
#P = np.hstack((P, np.zeros(n_reactions)))

# Impose equilibrium
#for i in range(n_reactions):
#    row = 4*len(truss) + i
#    K[row, row] = 1
#    if i < 2:
#        K[row, i*4+1] = -1
#    else:
#        K[row, i*4+3] = -1

# Solve for the displacements
displacements = np.linalg.solve(K, P)

if(len(truss)+n_reactions<(n_joints*2)): print('statically determinant...')


So I am attempting so far to assign Dirichlet boundary conditions to nodes at 0,0 and 0,8 as fixed or zero reaction. As far as assembling K I think it's generally on the right track but I don't feel sure about it yet. Mostly at this point what I am trying to solve is the displacements whatever they start out as so far. So far the use of solve leads to numpy.linalg.LinAlgError: Singular matrix. and I am really not quite sure yet why that is....As far as the commented out bit about reactions and equilibrium I could maybe use a bit of help understading that too. Mostly I don'y understand yet why the singular matrix?