An Elastic Möbius Strip

webgpuphysicsnumericsmath

An Elastic Möbius Strip

Everyone likes to talk about the topology of a Möbius strip, but you don’t usually hear too much about the geometry, the shape it would actually make if you assembled it out of a strip of elastic material. Steven Strogatz put this thought in my head in the introduction to a 2025 episode of the Quanta Magazine podcast The Joy of Why, Does Form Really Shape Function?. He recalls being charmed that the first work for which he knew his guest L. Mahadevan was a graduate-school analysis of the equilibrium shape of an elastic Möbius strip.

The image of an elastic Möbius strip stuck in my head for a couple reasons. First, I studied Engineering Mechanics in school and spent a lot of time computing shapes, bending modes, and wave propagation in rods, bars, plates, and the like. The details of what I did feel lost to time at this point, but the mechanics of materials will always have a special place in my heart.

Second, we all labored in school for weeks, sometimes months, over problems and MATLAB or Fortran codes we didn’t really understand. Maybe I’m underappreciating the rite of passage that is wrestling with these problems and eventually triumphing, but it never seemed right to me to take a problem that feels—in principle—so viscerally tangible, and reduce it to some opaque Fortran code that spits out NaN for some incomprehensible reason. I wanted to interact with the problem until I felt it enough to trust the answer.

Bret Victor summarizes what gets lost when these feedback loops don’t exist:

If there is any delay in that feedback loop between thinking of something and seeing it and building on it, then there is this whole world of ideas which will never be.

Bret Victor, Inventing on Principle

Fifteen years later, having done little about it myself, I think there’s still a lot of opportunity here.

So, not even finishing the episode, I felt immediately compelled to build something on the web which computes the equilibrium shape and bending modes of a Möbius strip. The equilibrium shape is an energy minimization problem. Bending modes are an eigenvalue problem.

Already though, I knew exactly where I always hit a brick wall when I try these things. Want to put something on the web that requires eigenvalues? ML.js works great. Maybe plug it into an Observable notebook. But do you need generalized eigenvalues? Eigenvalues of a banded matrix? Of a complex Hermitian matrix? Limited eigenvalues of a banded matrix? Well, good luck. If you want to do any serious analysis, these details matter. One option is to rely on a hosted notebook environment like Google Colab, but for a variety of reasons, that doesn’t suit my particular taste.

That’s no jab at the people who have put a lot of work into these problems. It’s genuinely challenging, stubborn work. BLAS (Basic Linear Algebra Subprograms) and LAPACK (Linear Algrebra PACKage) are the bedrock of numerical computing (they still power NumPy, by the way). LAPACK in particular comprises about 2,200 modules, about half of which are unique up to 32-bit vs. 64-bit precision, but they’re written in Fortran which, no matter how you approach it, makes them nontrivial to run on the web. Heavy use of goto-based control flow makes translation delicate and difficult, and compiling Fortran to WebAssembly, opts you into a terrible world of build scripts, opaque errors, and manual memory management.

My sense has always been that there’s no hope but to translate BLAS and LAPACK into idiomatic JavaScript. I’m not alone in thinking this, and there is long-running work under way in the Stdlib project to confront this challenge head-on. I’ve contributed to the project from time to time, but although I’ve contributed translations of a number of modules by hand, I was on track to finish sometime around 2050.

Maybe I lack imagination or determination or just the raw technical ability, but even though the Böhm–Jacopini theorem tells me it’s just an exercise for the reader to untangle LAPACK’s heavy use of goto, it hasn’t actually helped me get the work done. I’ve tried translating modules by hand, I’ve tried to automate the translation, and I’ve tried waiting for someone else to figure it out. Nothing has panned out.

XKCD #292
Obligatory XKCD #292. How bad can thousands of little gotos be?

I might lose some readers here, but earlier this year I threw Claude Code at improving the automated translation routines I’d struggled to build a year prior. My code would parse Fortran modules, run a sequence of carefully constructed transforms, and output JavaScript. I figured AI would have an easy time tightening up the details which had eluded me, but after some confusion on my part, Claude politely informed that my code was only getting in its way. Instead, it was just sidestepping me and performing the translation directly.

I resisted, but darn it if it wasn’t right. All I really had to do was ask for the generalized eigenvalues of a real matrix, let it churn through the LAPACK module dependency tree, and confirm the result actually solved the problem I’d asked for. If there were an infinite amount of work, I’d have held my ground and demanded a process rather than just working output. It’s a finite amount of work though which Claude Code is only becoming increasingly capable of managing.

That work is accumulating in blahpack, which is either a portmanteau or an expression of what it feels like to be a a programmer in 2026, your choice. Progress currently stands at about 878 of some 2,200 modules — the double and complex-double variants — though only about half of the 2,200 are unique anyway, since routines come in single- and double-precision twins.

So new in 2025-2026 is that you can finally just will these capabilities into existence.

Computing the shape of a Möbius strip

The shape of a Möbius strip is a problem with some history. Sadowsky posed its resting shape as an energy-minimization problem in 1930 but didn’t quite get to an answer. Mahadevan and Keller solved it in 1993 in The Shape of a Möbius Band, treating the band as a thin elastic rod and solving the resulting boundary-value problem numerically.

This notebook starts with the elastic rod solution, using a one-shotted AI JavaScript translation of the same COLSYS solver used by Mahadevan and Keller. From there, this notebook cuts up the geometry into a mesh of triangles and models the Möbius strip as a thin elastic sheet. Discrete triangles turn the continuous physics problem into a linear algebra problem. Using the numerical finite element method (FEM), each triangle resists in-plane stretching and bending with respect to its neighbors.

(a) membrane stretch & shear in-plane θ (b) bending curvature away from flat
Triangular shell elements. 3D Vertex positions are the unknowns. (a) Membrane stress resists stretching and shearing within the plane of the triangle; (b) bending stress resists the surface curving away from flat, measured at each vertex from the fold of the surrounding triangles.

Adding up the membrane and bending stresses of each triangle yields an overall elastic energy where is a vector containing the coordinates of all vertex positions. The equilibrium shape the strip is the one that minimizes energy The minimum is found with the L-BFGS method which requires the energy and its gradient The inner loop uses just dot products and vector sums, correponsing to the ddot and daxpy BLAS methods.

Just for fun, I decided to compute the vibrational bending modes as well. Nudge the strip away from its resting shape and Newton’s law governs what happens next, in terms of two matrices: the stiffness , the curvature of the energy at the minimum, and a mass matrix carrying the inertia of the moving vertices. Looking for motions in which every vertex oscillates in step at a single frequency turns that into a generalized eigenvalue problem.

The eigenvalue problem for bending modes

Expanding the energy to second order in a small displacement away from the equilibrium shape gives

(There is no term linear in since the gradient vanishes at equilibrium.) The matrix is the curvature of the energy, representing restoring force per unit of displacement.

That is an energy, though, and a vibration is a motion, so it takes one more step to connect the two. The elastic force is the negative gradient of the energy, and differentiating the expansion above with respect to the displacement drops the constant and leaves the quadratic term,

using for symmetric . Equivalently, expand the gradient itself rather than the energy: , whose first term is zero precisely because the strip sits at equilibrium. Either way the strip pushes back in proportion to how far it has been displaced, which is Hooke’s law with a matrix where the spring constant used to be.

The kinetic energy of the moving vertices, , supplies their inertia through the mass matrix . Newton’s law for the nudged strip is then simply

A vibration mode is a motion in which every vertex oscillates in step at one frequency, the shape itself held fixed. For harmonic motion at frequency the acceleration is just , so the two time derivatives collapse to a factor of , the minus signs cancel from both sides, and what survives is the eigenvalue problem below.

The eigenvalues are the squared frequencies and the eigenvectors are the mode shapes. After some trial and error using different LAPACK solvers which exploit different optimizations like the banded and symmetric nature of the matrices, I settled on ARPACK’s method for computing the lowest-order modes, with the banded and dense from LAPACK kept as cross-checks.

The figure below shows the computation with additional controls over the simulation and visualization.

The simulation is validated in a number of different ways. The elastic rod solution is validated against the paper. Multiple algorithms are compared for the equilibrium and bending mode steps of the computation. The linear algebra routines are compared to the original Fortran throughout the translation process and error residuals checked for converence. The elastic energy is validated for a flat strip and for an untwisted ring as well as across mesh resolutions and triangle aspect ratios.

A few problems were found. A number of linear algebra routines had errors, for example containing a 0- vs 1-indexed integer error in a branch which covered larger matrices. The first attempt at finite element problem setup only worked for equilateral triangular elements and had to be upgraded.

The button below runs a series of checks using the same functions as the above simulation.

Including modules and their dependents for all variants of the algorithms implemented above, this page uses translated routines, though only about ten are called directly. They’ve been validated and optimized. Additionally, about 899 of 1111 BLAS/LAPACK algorithms have been translated in the blahpack repository.

The figure below shows all modules and the dependency chains. Click a node to navigate to the source code.

Conclusion

In the end, it’s exciting and unnerving. I tried for quite a while to translate and enable this sort of computation, failing repeatedly to maintain enough velocity to get very far. On one hand, AI enables so much of what I always wanted to explore and present, on the other hand it devalues it.

I have a hard time understanding the availablity of BLAS and LAPACK on the web as evil, so I probably continue validating and optimizing the routines. They need to be wrapped up in more convenient interfaces. They need human judgement. Or more precisely, I guess, they need a human to unleash a wave of agents to tend to them, as if farming.

The future is weird.

Human batteries in the machine world from The Matrix (1999)
Farm-fresh BLAS/LAPACK modules, in my mental model of how this AI translation works. The Matrix (1999).

See also

All notebooks →