from sympy import *
import numpy as np
# Python 3.13.2
u1 = Symbol("u1")
u2 = Symbol("u2")
q3 = Symbol("q3")
u4 = Symbol("u4")
q5 = Symbol("q5")
disp_vector = np.array([u1, u2, q3, u4, q5])
stiffness_matrix = np.array([[1, 0, 0, -1, 0],
[0, 0.12, 0.6, 0, 0.6],
[0, 0.6, 4, 0, 2],
[-1, 0, 0, 1, 0],
[0, 0.6, 2, 0, 4]])
force_vector = np.array([0, 40, -26.7, 0, 0])
I am trying to code static condensation. This allows me to reduce the size of the stiffness matrix above but for that I need to be able to divide the stiffness matrix into following matrices.
krr = np.array([[1, 0, 0],
[0, 0.12, 0.6],
[0, 0.6, 4]])
krc = np.array([[-1, 0],
[0, 0.6],
[0, 2]])
kcr = np.array([[-1, 0, 0],
[0, 0.6, 2]])
kcc = np.array([[1, 0],
[0, 4]])
How would I do this ?
You need to divide your stiffness matrix into four submatrices based on the degrees of freedom you want to r
and c
. So to retain the first 3
rows/columns and condense the last 2
krr: The top left block
krc: The top right block
kcr: The bottom left block
kcc: The bottom right block
The code
n_retain = 3
n_condense = 2
krr = stiffness_matrix[:n_retain, :n_retain]
krc = stiffness_matrix[:n_retain, n_retain:]
kcr = stiffness_matrix[n_retain:, :n_retain]
kcc = stiffness_matrix[n_retain:, n_retain:]
print("krr =\n", krr)
print("\nkrc =\n", krc)
print("\nkcr =\n", kcr)
print("\nkcc =\n", kcc)
submatrices
krr =
[[1. 0. 0. ]
[0. 0.12 0.6 ]
[0. 0.6 4. ]]
krc =
[[-1. 0. ]
[ 0. 0.6]
[ 0. 2. ]]
kcr =
[[-1. 0. 0. ]
[ 0. 0.6 2. ]]
kcc =
[[1. 0.]
[0. 4.]]
Animations Credits - https://ocw.mit.edu/courses/1-571-structural-analysis-and-control-spring-2004/pages/readings/
<!DOCTYPE html>
<html>
<head>
<style>
body {
font-family: Arial, sans-serif;
margin: 20px;
background-color: #f5f5f5;
}
.container {
max-width: 800px;
margin: 0 auto;
background-color: white;
padding: 20px;
border-radius: 8px;
box-shadow: 0 2px 10px rgba(0,0,0,0.1);
}
h1, h2 {
color: #2c3e50;
text-align: center;
}
.matrix-container {
display: flex;
justify-content: center;
align-items: center;
margin: 20px 0;
flex-wrap: wrap;
}
.matrix {
border: 2px solid #3498db;
margin: 10px;
padding: 5px;
background-color: white;
border-radius: 4px;
transition: all 0.5s ease;
}
.matrix-row {
display: flex;
justify-content: center;
}
.matrix-cell {
width: 40px;
height: 40px;
display: flex;
align-items: center;
justify-content: center;
margin: 2px;
font-weight: bold;
transition: all 0.5s ease;
}
.retained {
background-color: #2ecc71;
color: white;
}
.condensed {
background-color: #e74c3c;
color: white;
}
.coupling {
background-color: #f39c12;
color: white;
}
.controls {
display: flex;
justify-content: center;
margin: 20px 0;
}
button {
background-color: #3498db;
color: white;
border: none;
padding: 10px 20px;
margin: 0 10px;
border-radius: 4px;
cursor: pointer;
font-size: 16px;
transition: background-color 0.3s;
}
button:hover {
background-color: #2980b9;
}
button:disabled {
background-color: #bdc3c7;
cursor: not-allowed;
}
.explanation {
margin: 20px 0;
padding: 15px;
background-color: #f8f9fa;
border-left: 4px solid #3498db;
border-radius: 4px;
}
.equation {
font-family: 'Times New Roman', Times, serif;
font-style: italic;
text-align: center;
margin: 15px 0;
font-size: 18px;
}
.force-vector {
display: flex;
flex-direction: column;
align-items: center;
margin: 10px;
}
.force-cell {
width: 40px;
height: 30px;
display: flex;
align-items: center;
justify-content: center;
margin: 2px;
font-weight: bold;
border: 1px solid #3498db;
background-color: #fff;
transition: all 0.5s ease;
}
.progress {
width: 100%;
height: 8px;
background-color: #ecf0f1;
margin-bottom: 20px;
border-radius: 4px;
overflow: hidden;
}
.progress-bar {
height: 100%;
background-color: #3498db;
width: 0%;
transition: width 0.3s ease;
}
</style>
</head>
<body>
<div class="container">
<h1>Static Condensation Animation</h1>
<div class="progress">
<div class="progress-bar" id="progressBar"></div>
</div>
<div class="explanation" id="explanation">
<p>Static condensation is a technique used in structural analysis to reduce the size of the stiffness matrix by eliminating degrees of freedom (DOFs) that are not of primary interest.</p>
<p>Click "Next" to walk through the process step by step.</p>
</div>
<div class="matrix-container" id="matrixDisplay">
<!-- Matrices will be displayed here -->
</div>
<div class="controls">
<button id="prevBtn" disabled>Previous</button>
<button id="nextBtn">Next</button>
<button id="resetBtn">Reset</button>
</div>
</div>
<script>
// Original stiffness matrix and force vector
const stiffnessMatrix = [
[1, 0, 0, -1, 0],
[0, 0.12, 0.6, 0, 0.6],
[0, 0.6, 4, 0, 2],
[-1, 0, 0, 1, 0],
[0, 0.6, 2, 0, 4]
];
const forceVector = [0, 40, -26.7, 0, 0];
// Steps for the animation
const steps = [
{
title: "Original System",
explanation: "We start with the full stiffness matrix K and force vector F. Our goal is to solve Ku = F for the displacement vector u.",
showMatrices: ["K", "F"]
},
{
title: "Partition the System",
explanation: "We partition the matrix into 4 submatrices: Krr (retained DOFs), Kcc (condensed DOFs), and Krc/Kcr (coupling terms). The force vector is similarly split into Fr and Fc.",
showMatrices: ["K_partitioned", "F_partitioned"]
},
{
title: "Extract Submatrices",
explanation: "We extract the four submatrices explicitly: Krr (3×3), Krc (3×2), Kcr (2×3), and Kcc (2×2).",
showMatrices: ["Krr", "Krc", "Kcr", "Kcc", "Fr", "Fc"]
},
{
title: "Condensation Equations",
explanation: "The static condensation process can be represented by these equations:\n\nCondensed stiffness matrix: K* = Krr - Krc × Kcc⁻¹ × Kcr\nCondensed force vector: F* = Fr - Krc × Kcc⁻¹ × Fc",
showMatrices: ["equation"]
},
{
title: "Condensed System",
explanation: "We now have a reduced system K* u_r = F* that only involves the retained DOFs. This smaller system is easier to solve.",
showMatrices: ["K_condensed", "F_condensed"]
},
{
title: "Solve for Retained DOFs",
explanation: "We solve the condensed system to find the values of the retained DOFs (u_r).",
showMatrices: ["u_r"]
},
{
title: "Back-Calculate Condensed DOFs",
explanation: "Once we have u_r, we can back-calculate the condensed DOFs (u_c) using: u_c = Kcc⁻¹ × (Fc - Kcr × u_r)",
showMatrices: ["u_r", "u_c"]
},
{
title: "Complete Solution",
explanation: "Finally, we have the complete solution vector u combining u_r and u_c. We've solved the full system by working with a reduced matrix!",
showMatrices: ["u_full"]
}
];
// Function to create matrix display
function createMatrixDisplay(matrix, title, className = "") {
const matrixDiv = document.createElement("div");
matrixDiv.className = `matrix ${className}`;
const titleDiv = document.createElement("div");
titleDiv.textContent = title;
titleDiv.style.textAlign = "center";
titleDiv.style.marginBottom = "5px";
matrixDiv.appendChild(titleDiv);
for (let i = 0; i < matrix.length; i++) {
const rowDiv = document.createElement("div");
rowDiv.className = "matrix-row";
if (Array.isArray(matrix[i])) {
// It's a 2D matrix
for (let j = 0; j < matrix[i].length; j++) {
const cellDiv = document.createElement("div");
cellDiv.className = "matrix-cell";
cellDiv.textContent = typeof matrix[i][j] === 'number' ?
matrix[i][j].toFixed(2).replace(/\.00$/, "") : matrix[i][j];
// Apply colors based on partitioning for the full matrix
if (title === "K (partitioned)") {
if (i < 3 && j < 3) cellDiv.classList.add("retained");
else if (i >= 3 && j >= 3) cellDiv.classList.add("condensed");
else cellDiv.classList.add("coupling");
}
rowDiv.appendChild(cellDiv);
}
} else {
// It's a vector
const cellDiv = document.createElement("div");
cellDiv.className = "force-cell";
cellDiv.textContent = typeof matrix[i] === 'number' ?
matrix[i].toFixed(2).replace(/\.00$/, "") : matrix[i];
// Apply colors for partitioned force vector
if (title === "F (partitioned)") {
if (i < 3) cellDiv.classList.add("retained");
else cellDiv.classList.add("condensed");
}
rowDiv.appendChild(cellDiv);
}
matrixDiv.appendChild(rowDiv);
}
return matrixDiv;
}
// Function to create force vector display
function createForceVector(vector, title, className = "") {
const vectorDiv = document.createElement("div");
vectorDiv.className = `force-vector ${className}`;
const titleDiv = document.createElement("div");
titleDiv.textContent = title;
titleDiv.style.textAlign = "center";
titleDiv.style.marginBottom = "5px";
vectorDiv.appendChild(titleDiv);
for (let i = 0; i < vector.length; i++) {
const cellDiv = document.createElement("div");
cellDiv.className = "force-cell";
cellDiv.textContent = typeof vector[i] === 'number' ?
vector[i].toFixed(2).replace(/\.00$/, "") : vector[i];
// Apply colors for partitioned force vector
if (title === "F (partitioned)") {
if (i < 3) cellDiv.classList.add("retained");
else cellDiv.classList.add("condensed");
}
vectorDiv.appendChild(cellDiv);
}
return vectorDiv;
}
// Function to display an equation
function createEquation(text) {
const eqDiv = document.createElement("div");
eqDiv.className = "equation";
eqDiv.innerHTML = text;
return eqDiv;
}
// Get DOM elements
const matrixDisplay = document.getElementById("matrixDisplay");
const explanation = document.getElementById("explanation");
const prevBtn = document.getElementById("prevBtn");
const nextBtn = document.getElementById("nextBtn");
const resetBtn = document.getElementById("resetBtn");
const progressBar = document.getElementById("progressBar");
let currentStep = 0;
// Submatrices
const krr = [
[1, 0, 0],
[0, 0.12, 0.6],
[0, 0.6, 4]
];
const krc = [
[-1, 0],
[0, 0.6],
[0, 2]
];
const kcr = [
[-1, 0, 0],
[0, 0.6, 2]
];
const kcc = [
[1, 0],
[0, 4]
];
const fr = [0, 40, -26.7];
const fc = [0, 0];
// Solution (mockup values for demonstration)
const ur = [0.78, 345.58, -9.50];
const uc = [-0.78, 1.49];
const uFull = [0.78, 345.58, -9.50, -0.78, 1.49];
// Condensed matrices (mockup values for demonstration)
const kCondensed = [
[2, 0, 0],
[0, 0.21, 0.6],
[0, 0.6, 3.5]
];
const fCondensed = [0, 40, -26.7];
// Update display based on current step
function updateDisplay() {
// Update progress bar
progressBar.style.width = `${(currentStep / (steps.length - 1)) * 100}%`;
// Clear previous content
matrixDisplay.innerHTML = "";
// Set explanation
explanation.innerHTML = `<h2>${steps[currentStep].title}</h2><p>${steps[currentStep].explanation}</p>`;
// Display relevant matrices
steps[currentStep].showMatrices.forEach(matrixName => {
switch (matrixName) {
case "K":
matrixDisplay.appendChild(createMatrixDisplay(stiffnessMatrix, "K (Stiffness Matrix)"));
break;
case "F":
matrixDisplay.appendChild(createForceVector(forceVector, "F (Force Vector)"));
break;
case "K_partitioned":
matrixDisplay.appendChild(createMatrixDisplay(stiffnessMatrix, "K (partitioned)"));
break;
case "F_partitioned":
matrixDisplay.appendChild(createForceVector(forceVector, "F (partitioned)"));
break;
case "Krr":
matrixDisplay.appendChild(createMatrixDisplay(krr, "Krr (retained)", "retained"));
break;
case "Krc":
matrixDisplay.appendChild(createMatrixDisplay(krc, "Krc (coupling)", "coupling"));
break;
case "Kcr":
matrixDisplay.appendChild(createMatrixDisplay(kcr, "Kcr (coupling)", "coupling"));
break;
case "Kcc":
matrixDisplay.appendChild(createMatrixDisplay(kcc, "Kcc (condensed)", "condensed"));
break;
case "Fr":
matrixDisplay.appendChild(createForceVector(fr, "Fr (retained)", "retained"));
break;
case "Fc":
matrixDisplay.appendChild(createForceVector(fc, "Fc (condensed)", "condensed"));
break;
case "equation":
const eqDiv = document.createElement("div");
eqDiv.style.width = "100%";
eqDiv.style.textAlign = "center";
const eq1 = document.createElement("div");
eq1.className = "equation";
eq1.innerHTML = "K* = K<sub>rr</sub> - K<sub>rc</sub> · K<sub>cc</sub><sup>-1</sup> · K<sub>cr</sub>";
const eq2 = document.createElement("div");
eq2.className = "equation";
eq2.innerHTML = "F* = F<sub>r</sub> - K<sub>rc</sub> · K<sub>cc</sub><sup>-1</sup> · F<sub>c</sub>";
eqDiv.appendChild(eq1);
eqDiv.appendChild(eq2);
matrixDisplay.appendChild(eqDiv);
break;
case "K_condensed":
matrixDisplay.appendChild(createMatrixDisplay(kCondensed, "K* (Condensed Stiffness Matrix)"));
break;
case "F_condensed":
matrixDisplay.appendChild(createForceVector(fCondensed, "F* (Condensed Force Vector)"));
break;
case "u_r":
matrixDisplay.appendChild(createForceVector(ur, "ur (Retained DOFs)", "retained"));
break;
case "u_c":
matrixDisplay.appendChild(createForceVector(uc, "uc (Condensed DOFs)", "condensed"));
break;
case "u_full":
matrixDisplay.appendChild(createForceVector(uFull, "u (Complete Solution)"));
break;
}
});
// Update button states
prevBtn.disabled = currentStep === 0;
nextBtn.disabled = currentStep === steps.length - 1;
}
// Event listeners
prevBtn.addEventListener("click", () => {
if (currentStep > 0) {
currentStep--;
updateDisplay();
}
});
nextBtn.addEventListener("click", () => {
if (currentStep < steps.length - 1) {
currentStep++;
updateDisplay();
}
});
resetBtn.addEventListener("click", () => {
currentStep = 0;
updateDisplay();
});
// Initialize display
updateDisplay();
</script>
</body>
</html>