The stiffness matrix: from a single spring to the global system.

En el artículo anterior vimos los siete pasos que ejecuta el software en cada cálculo de elementos finitos. Dos de ellos generan la pieza central de todo el proceso: el paso 4 construye la matriz de rigidez de cada elemento por separado; el paso 5 las ensambla en una única matriz de rigidez global. Esa matriz —[K]— es lo que realmente resuelve el ordenador cuando pulsas «Solve». Entramos en ella con el ejemplo más simple posible: un resorte.

La ecuación que sostiene todo

Casi cualquier problema de estática lineal se reduce a una sola relación:

{K} · {u} = {F}

K es la matriz de rigidez: la traducción de la geometría y el material de la pieza al lenguaje del álgebra lineal. u es el vector de desplazamientos —la incógnita que resuelve el software, cuánto se mueve cada nodo. F es el vector de fuerzas aplicadas en esos mismos nodos.

Un término individual k_ij de esa matriz tiene un significado físico preciso: es la fuerza que hay que aplicar en el grado de libertad i para producir un desplazamiento unitario en el grado de libertad j, manteniendo inmóviles todos los demás. La matriz entera no es más que esa misma pregunta, repetida para cada par de grados de libertad del modelo.

K no describe una fuerza ni un desplazamiento. Describe cómo se relacionan todos los desplazamientos posibles con todas las fuerzas posibles, a la vez.

El bloque de construcción: una barra

Antes de una malla con millones de grados de libertad, conviene mirar el elemento más simple que existe: una barra sometida a carga axial. Con área de sección A, longitud L y módulo de Young E, su rigidez axial es:

k = AE / L

Cuanto más gruesa, más corta o de material más rígido, más cuesta estirarla. Para ese único elemento, con dos nodos y un grado de libertad por nodo, la relación entre fuerzas y desplazamientos nodales es una matriz 2×2:

Fig. 1 — Matriz de rigidez local de una barra articulada. Las celdas laten por pares — diagonal y fuera-de-diagonal — mostrando que k₁₂ = k₂₁: la simetria no es casualidad, es reciprocidad elástica.

La matriz local dice algo muy concreto: una fuerza aplicada en un extremo genera una reacción igual y de signo opuesto en el otro. Es, literalmente, la tercera ley de Newton escrita en forma matricial.

Ensamblaje por superposición

Una barra sola no sirve de mucho. El software construye la matriz global uniendo cientos, miles o millones de matrices locales como esta mediante un proceso de ensamblaje por superposición.

Imagina dos barras en serie, compartiendo un nodo central. El software expande cada matriz local de 2×2 al tamaño total del sistema —en este caso, 3×3— y las suma, término a término. En el nodo compartido, la rigidez final es sencillamente la suma de las rigideces de todos los elementos que confluyen en él.

Fig. 2 — Ensamblaje por superposicion: el elemento A aporta su bloque a las filas/columnas 1–2; el B, a las 2–3. En el nodo compartido (2,2) las rigideces se suman; la celda (1,3), sin conexion directa, queda en cero —ahí nace la dispersión de [K].

Ese paso, mecánico y repetitivo, garantiza dos cosas a la vez: que los elementos no se separen entre sí (compatibilidad de desplazamientos) y que las fuerzas internas estén en equilibrio en cada nodo. No es un truco de programación —es la superposición lineal aplicada con total literalidad.

Tres propiedades que delatan la física

La matriz K no es una tabla de números cualquiera. Hereda del problema físico tres propiedades que el software explota activamente.

Simetría. Por reciprocidad elástica, la fuerza que hace falta en A para mover B es idéntica a la que hace falta en B para mover A: k_ij = k_ji. Eso permite al solver almacenar solo la mitad de la matriz, con el consiguiente ahorro de memoria.

Dispersión (sparsity). En un modelo real, la inmensa mayoría de los términos de K son cero. Un nodo en un extremo de un puente no tiene conexión directa con un nodo en el extremo opuesto —solo con sus vecinos inmediatos en la malla. Esa dispersión es la que hace viable resolver sistemas con millones de incógnitas en un ordenador de sobremesa.

Singularidad inicial. Antes de fijar la estructura al suelo, K es inestable: si intentas resolver el sistema tal cual, el software devuelve un error de matriz singular. Tiene sentido físico —sin restricciones, la estructura entera podría moverse como un cuerpo rígido sin deformarse ni un milímetro, y ese movimiento no tiene una solución única.

Reordenar para poder resolver: ancho de banda y skyline

Con millones de grados de libertad, cómo se almacena K importa tanto como lo que contiene. El software reordena la numeración de los nodos —con algoritmos como el de Cuthill-McKee— para agrupar todos los términos distintos de cero lo más cerca posible de la diagonal principal. A esa franja se la llama ancho de banda (bandwidth), y cuanto más estrecha, menos memoria y menos operaciones hacen falta para resolver el sistema.

Fig. 3 — Reordenar la numeración de los nodos agrupa los términos no nulos junto a la diagonal. Menos ancho de banda significa menos memoria y menos operaciones para resolver el mismo sistema.

Un método más avanzado, el almacenamiento skyline, va un paso más allá: por cada columna, solo guarda los términos desde el primer valor no nulo hasta la diagonal, ignorando todo lo demás. Entre reordenación y almacenamiento inteligente, un ordenador de escritorio resuelve hoy problemas que hace treinta años exigían una supercomputadora.

Fijar la estructura y despejar el sistema

Para que K deje de ser singular hay que aplicar condiciones de contorno —apoyos que impiden el movimiento de cuerpo rígido. En términos matriciales, eso significa particionar el sistema: las filas y columnas asociadas a desplazamientos conocidos (típicamente cero, en un apoyo fijo) se separan del resto.

Con el sistema ya reducido, el software despeja las incógnitas con métodos numéricos como la eliminación de Gauss o la descomposición de Cholesky ([L][D][L]^T), aprovechando la simetría y la dispersión de K para hacerlo con la eficiencia que exige un sistema de ese tamaño. Solo entonces —y no antes— el software conoce el desplazamiento de cada nodo. A partir de ahí deriva las deformaciones elementales y, con la ley de comportamiento del material, las tensiones que ves en el mapa de colores.

Una matriz sana es un resultado fiable

La matriz de rigidez no es una abstracción de fondo: es la traducción matemática, literal, de si tu diseño tiene sentido físico. Un mallado con elementos muy distorsionados o unas condiciones de contorno insuficientes no producen un error visible de inmediato —producen una matriz mal condicionada, y una matriz mal condicionada produce resultados que parecen números pero no son física. Entender qué hay dentro de K es lo que separa confiar en un resultado de simplemente creerlo.

FIN

In the previous article we walked through the seven steps the software runs in every finite element calculation. Two of them generate the centrepiece of the whole process: step 4 builds the stiffness matrix of each element on its own; step 5 assembles them into a single global stiffness matrix. That matrix —[K]— is what the computer actually solves when you hit "Solve". Let's get inside it with the simplest possible example: a spring.

The equation that holds it all together

Almost any linear static problem reduces to a single relationship:

{K} · {u} = {F}

K is the stiffness matrix: the translation of the part's geometry and material into the language of linear algebra. u is the displacement vector —the unknown the software solves for, how much each node moves. F is the vector of forces applied at those same nodes.

A single term k_ij in that matrix has a precise physical meaning: it's the force you'd need to apply at degree of freedom i to produce a unit displacement at degree of freedom j, while holding every other node fixed. The whole matrix is nothing more than that same question, asked for every pair of degrees of freedom in the model.

K doesn't describe a force or a displacement. It describes how every possible displacement relates to every possible force, all at once.

The building block: a bar

Before a mesh with millions of degrees of freedom, it's worth looking at the simplest element there is: a bar under axial load. With cross-sectional area A, length L and Young's modulus E, its axial stiffness is:

k = AE / L

The thicker, shorter, or stiffer the material, the harder it is to stretch. For that single element, with two nodes and one degree of freedom per node, the relationship between nodal forces and displacements is a 2×2 matrix:

Fig. 1 — Local stiffness matrix of a bar element. Cells pulse in pairs — diagonal and off-diagonal — showing that k₁₂ = k₂₁: the symmetry isn’t a coincidence, it’s elastic reciprocity.

The local matrix says something very concrete: a force applied at one end produces an equal and opposite reaction at the other. It's, literally, Newton's third law written in matrix form.

Assembly by superposition

A single bar isn't much use on its own. The software builds the global matrix by joining hundreds, thousands, or millions of local matrices like this one through a process of assembly by superposition.

Picture two bars in series, sharing a central node. The software expands each local 2×2 matrix to the full size of the system —3×3, in this case— and adds them, term by term. At the shared node, the final stiffness is simply the sum of the stiffnesses of every element that meets there.

Fig. 2 — Assembly by superposition: element A contributes its block to rows/columns 1–2; element B, to 2–3. At the shared node (2,2) the stiffnesses add up; cell (1,3), with no direct connection, stays at zero — the origin of [K]’s sparsity.

That step, mechanical and repetitive as it is, guarantees two things at once: that elements don't pull apart from each other (displacement compatibility), and that internal forces balance at every node. It's not a programming trick —it's linear superposition, applied completely literally.

Three properties that give away the physics

The K matrix isn't just any table of numbers. It inherits three properties from the underlying physics that the software actively exploits.

Symmetry. By elastic reciprocity, the force needed at A to move B is identical to the force needed at B to move A: k_ij = k_ji. That lets the solver store only half the matrix, saving memory as a result.

Sparsity. In a real model, the vast majority of K's terms are zero. A node at one end of a bridge has no direct connection to a node at the opposite end —only to its immediate neighbours in the mesh. That sparsity is what makes it feasible to solve systems with millions of unknowns on a desktop computer.

Initial singularity. Before the structure is fixed to the ground, K is unstable: try to solve the system as it stands and the software returns a singular matrix error. That makes physical sense —without constraints, the whole structure could move as a rigid body without deforming a single millimetre, and that motion has no unique solution.

Reordering to make it solvable: bandwidth and skyline

With millions of degrees of freedom, how K is stored matters as much as what it contains. The software renumbers the nodes —using algorithms such as Cuthill-McKee— to cluster every non-zero term as close as possible to the main diagonal. That band is called the bandwidth, and the narrower it is, the less memory and fewer operations it takes to solve the system.

Fig. 3 — Renumbering the nodes clusters the non-zero terms next to the diagonal. A smaller bandwidth means less memory and fewer operations to solve the same system.

A more advanced method, skyline storage, goes one step further: for each column, it only keeps the terms from the first non-zero value down to the diagonal, ignoring everything else. Between renumbering and smart storage, a desktop computer today solves problems that would have needed a supercomputer thirty years ago.

Fixing the structure and solving the system

For K to stop being singular, you need to apply boundary conditions —supports that prevent rigid-body motion. In matrix terms, that means partitioning the system: the rows and columns tied to known displacements (typically zero, at a fixed support) are separated from the rest.

With the system reduced, the software solves for the unknowns using numerical methods such as Gaussian elimination or Cholesky decomposition ([L][D][L]^T), exploiting K's symmetry and sparsity to do so with the efficiency a system of that size demands. Only then —not before— does the software know the displacement of every node. From there it derives element strains and, through the material's constitutive law, the stresses you see in the colour map.

A healthy matrix is a trustworthy result

The stiffness matrix isn't background abstraction: it's the literal mathematical translation of whether your design makes physical sense. A mesh with badly distorted elements, or boundary conditions that don't fully constrain the model, won't throw an obvious error right away —they'll produce an ill-conditioned matrix, and an ill-conditioned matrix produces results that look like numbers but aren't physics. Understanding what's inside K is what separates trusting a result from merely believing it.

END