---
jupyter:
  description: 'The squircle and the superellipse with planex2d: superellipses drawn with px.Curve, their areas against
    4ab Γ(1 + 1/n)² / Γ(1 + 2/n), the squircle''s curvature from its spline against the exact formula, and a rounded
    square of the same area, with arc corners and with curvature-continuous corners, compared by curvature combs.'
---

# The squircle and the superellipse

What makes a squircle look smoother than a rounded square? A superellipse is the curve

$$\left|\frac{x}{a}\right|^n + \left|\frac{y}{b}\right|^n = 1,$$

which Gabriel Lamé studied. With $n$ = 2 it is an ellipse, and as $n$ grows it fills out towards the rectangle around
it. Piet Hein chose $n$ = 2.5 for the Sergels torg roundabout in Stockholm, and the squircle, with $a = b$ and $n$ = 4,
is used for app icons and in product design. Its area is known in closed form,

$$A = 4ab \, \frac{\Gamma(1 + 1/n)^2}{\Gamma(1 + 2/n)}.$$

A rounded square of the same size looks almost the same, but its curvature jumps from zero on the side to $1/r$ on the
corner arc. That jump shows in the reflections on a glossy surface. The squircle has no straight part and no jump, and
its curvature is zero only at the middle of each side and largest at the diagonals.

The example draws superellipses with `px.Curve` and checks their areas, reads the squircle's curvature from its
spline and checks it against the exact formula, and compares it with the rounded square of the same area, with arc
corners and with curvature-continuous corners.

This page is part of a collection of [geometry and design examples](geometry-and-design.md) using planex2d.

<!-- #raw raw_mimetype="text/restructuredtext" -->
.. rst-class:: px-example-links

`Open this example in the playground <https://planex2d.org/playground/>`__, the interactive code editor in the browser, to change the code and run it.
<!-- #endraw -->

![The squircle and a rounded square of the same area, with their curvature combs](../_static/gallery/squircle.png){.result}

```python
import math

import matplotlib.pyplot as plt
import numpy as np
import planex as px
from scipy.interpolate import BSpline
from scipy.optimize import brentq
```

```python nbsphinx="hidden"
from docs_figures import Chart, Row
```

## Superellipses

A superellipse is drawn from its polar form, in which the point at the angle $t$ lies at the distance
$r(t) = \left(|\cos t / a|^n + |\sin t / b|^n\right)^{-1/n}$ from the center. One turn of $t$ gives one closed spline
edge, and the area of its face matches the formula to about 1e-10 (Figure 1):

```python
def superellipse(n, a=1, b=1):
    """the face inside |x/a|^n + |y/b|^n = 1"""
    def point(t):
        r = (abs(math.cos(t) / a) ** n + abs(math.sin(t) / b) ** n) ** (-1 / n)
        return (r * math.cos(t), r * math.sin(t))
    return px.Face([px.Curve(point, (0, 2 * math.pi))])


def area(n, a=1, b=1):
    return 4 * a * b * math.gamma(1 + 1 / n) ** 2 / math.gamma(1 + 2 / n)


shapes = {n: superellipse(n) for n in (2, 2.5, 4, 10)}
for n, face in shapes.items():
    print("n", n, "| edges", len(face.edges), "| area", face.mass, "| formula", area(n))
hein = superellipse(2.5, 1.2, 1)
print("Piet Hein's 6 by 5, n 2.5: area", hein.mass, "| formula", area(2.5, 1.2, 1))
# px.Plot(px.Compound(list(shapes.values())))
```

```python nbsphinx="hide-input"
Row(*[(face, f"n = {n}") for n, face in shapes.items()],
    caption="Superellipses with a = b = 1, from the circle to almost a square")
```

## Their curvature

The curvature of a spline follows from its exact B-spline, as on the page on
[curvature-continuous blends](../core/sketch/blends.md). The function below goes along a face's outline, edge by edge,
and returns points on it, the unit normals pointing out and the curvature there. The curvature is zero on a line and
$1/r$ on an arc, and on a spline it follows from the first and second derivatives:

```python
def along(face, m=100):
    """points along the outline, the unit normals pointing out and the curvature, m per edge"""
    points, tangents, curvatures = [], [], []
    for e in face.edges:
        edge = e.edge
        if isinstance(edge, px.SplineEdge):
            data = edge.bspline
            p, knots = data["degree"], np.asarray(data["knots"])
            curve = BSpline(knots, np.asarray(data["poles"]), p)
            u = np.linspace(knots[p], knots[-p - 1], m)
            xy, d1, d2 = curve(u), curve.derivative(1)(u), curve.derivative(2)(u)
            speed = np.hypot(d1[:, 0], d1[:, 1])
            t, k = d1 / speed[:, None], (d1[:, 0] * d2[:, 1] - d1[:, 1] * d2[:, 0]) / speed**3
        elif isinstance(edge, px.ArcEdge):
            c, r, turn = np.array([edge.center[0], edge.center[1]]), edge.radius, 1 if edge.ccw else -1
            a0, a1 = (math.atan2(q[1] - c[1], q[0] - c[0]) for q in (edge.start, edge.end))
            a = np.linspace(a0, a0 + turn * ((turn * (a1 - a0)) % (2 * math.pi)), m)
            xy = c + r * np.c_[np.cos(a), np.sin(a)]
            t, k = turn * np.c_[-np.sin(a), np.cos(a)], np.full(m, turn / r)
        else:
            s, f = np.array([edge.start[0], edge.start[1]]), np.array([edge.end[0], edge.end[1]])
            xy = s + np.linspace(0, 1, m)[:, None] * (f - s)
            t, k = np.tile((f - s) / np.hypot(*(f - s)), (m, 1)), np.zeros(m)
        points.append(xy)
        tangents.append(t)
        curvatures.append(k)
    t = np.vstack(tangents)
    return np.vstack(points), np.c_[t[:, 1], -t[:, 0]], np.concatenate(curvatures)
```

On the superellipse itself, the curvature at a point $(x, y)$ is

$$\kappa = \frac{(n - 1) \, |x y|^{n - 2}}{\left(|x|^{2n - 2} + |y|^{2n - 2}\right)^{3/2}},$$

which is zero on the axes for $n > 2$ and $(n - 1) \, 2^{1/n - 1/2}$ at the diagonals. The squircle's spline follows it:

```python
squircle = shapes[4]
xy, normals, kappa = along(squircle, 2000)
x, y = np.abs(xy[:, 0]), np.abs(xy[:, 1])
exact = 3 * (x * y) ** 2 / (x**6 + y**6) ** 1.5
print("largest curvature", kappa.max(), "| at the diagonals", 3 * 2 ** (1 / 4 - 1 / 2))
print("smallest", kappa.min(), "| largest difference from the formula", np.abs(kappa - exact).max())
```

## Against a rounded square

A square of side 2 with corner arcs of radius $r$ has the area $4 - (4 - \pi) r^2$, so the rounded square with the
squircle's area has $r$ = 0.583. The same square with curvature-continuous (G2) corners, made with `Filleted` and
`px.Curvature`, has its own area, and SciPy's `brentq` finds the radius that gives it the squircle's area too. The arc
corners stay within 0.017 of the squircle and the G2 corners within 0.010:

```python
square = px.QuickSketch().RectangleC((0, 0), 2, 2).Face()
r_arc = math.sqrt((4 - area(4)) / (4 - math.pi))
r_g2 = brentq(lambda r: square.Filleted(square.vertices, r, px.Curvature).mass - area(4), 0.3, 0.99, xtol=1e-14)
arcs = square.Filleted(square.vertices, r_arc)
g2 = square.Filleted(square.vertices, r_g2, px.Curvature)
for name, face, r in (("arcs", arcs, r_arc), ("curvature-continuous", g2, r_g2)):
    print(name, "| r", r, "| area", face.mass, "| Hausdorff distance from the squircle", face.Hausdorff(squircle))
```

The differences show in the curvature. A curvature comb, as CAD systems draw it, puts a line along the normal at each
point of the outline, as long as the curvature there. The arc corners' comb steps from zero to $1/r$ where the side
ends, the curvature-continuous corners rise from zero, and the squircle's comb has no step and no straight part
(Figure 2):

```python
def comb(face, m, scale=0.15):
    """the curvature comb of the outline: the teeth along the normals and the line through their tips"""
    xy, normals, kappa = along(face, m)
    tips = xy + scale * kappa[:, None] * normals
    teeth = [px.Line(tuple(p), tuple(q)) for p, q in zip(xy, tips) if np.hypot(*(q - p)) > 1e-9]
    keep = np.r_[True, np.hypot(*np.diff(tips, axis=0).T) > 1e-9]
    return px.Compound(teeth + [px.FromPolyline(tips[keep], closed=True)])


combs = {"Squircle": (squircle, 240), "Arc corners": (arcs, 30), "G2 corners": (g2, 30)}
# px.Plot(px.Compound([squircle, comb(squircle, 240)]))
```

```python nbsphinx="hide-input"
Row(*[(px.Compound([face, comb(face, m)]), name, {"lw": 0.8}) for name, (face, m) in combs.items()],
    caption="Curvature combs of the squircle and of the rounded square of the same area")
```

Along a quarter of the outline, the curvature of the arc corners jumps twice, at the directions where the corners start
and end. The G2 corners, of radius 0.724, peak at a curvature of 2.03, and the squircle peaks at 2.52 on its diagonal
(Figure 3):

```python
quarter = {}
for name, (face, m) in combs.items():
    xy, normals, kappa = along(face, 8 * m)
    direction = np.degrees(np.arctan2(xy[:, 1], xy[:, 0])) % 360
    order = np.argsort(direction, kind="stable")
    keep = order[direction[order] <= 90]
    quarter[name] = (direction[keep], kappa[keep])
    print(name, "| largest curvature", kappa.max())
# plt.plot(*quarter["Squircle"]); plt.show()
```

```python nbsphinx="hide-input"
Chart("Direction from the center in degrees", "Curvature",
      lines={"Arc corners": quarter["Arc corners"], "G2 corners": quarter["G2 corners"]},
      curves={"Squircle": quarter["Squircle"]}, title="The curvature along a quarter of the outline")
```

See [Curves from equations](../core/basics/curves.md) for `px.Curve` and
[Curvature-continuous blends](../core/sketch/blends.md) for `Fillet` with `px.Curvature`.

```python nbsphinx="hidden"
import sys
sys.path.insert(0, "..")
from docs_figures import FACE, INK, PLANEX
from gallery import save_figure

fig, ax = plt.subplots(figsize=(6, 4))
for face, m, dx in ((squircle, 240, -1.45), (arcs, 30, 1.45)):
    face.Translated((dx, 0)).Plot(ax=ax, face_col=FACE, line_col=INK, arc_col=INK, spline_col=INK, legend=False)
    comb(face, m).Translated((dx, 0)).Plot(ax=ax, line_col=PLANEX, arc_col=PLANEX, spline_col=PLANEX, legend=False)
for artist in ax.lines + ax.collections + ax.patches:
    artist.set_linewidth(1.0)
ax.set_xlim(-3.05, 3.05)
ax.set_ylim(-2.03, 2.03)
ax.set_aspect("equal", adjustable="box")
ax.axis("off")
save_figure(fig, "squircle")
plt.close(fig)
```