Functional Representations & Smooth Frame Fields in 2D
A frame (cross) is a set of 4 unit vectors with \(\pi/2\) rotational symmetry:
\[f_i = \bigl(\cos(\theta + i\pi/2),\; \sin(\theta + i\pi/2)\bigr), \quad i \in \{0,1,2,3\}\]
A frame field assigns one frame per vertex of a triangle mesh.
Goal: find a smooth frame field aligned with the domain boundary.
{
const w = 340, h = 340, cx = w/2, cy = h/2, r = 120;
const theta = theta_deg_s * Math.PI / 180;
const svg = d3.create("svg").attr("viewBox", `0 0 ${w} ${h}`).style("width", "100%").style("max-width", w + "px");
const circle_pts = d3.range(0, 2*Math.PI + 0.01, 0.02).map(a => [cx + r*Math.cos(a), cy - r*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circle_pts)).attr("fill", "none").attr("stroke", "#ccc").attr("stroke-width", 1);
svg.append("line").attr("x1", cx - r - 15).attr("y1", cy).attr("x2", cx + r + 15).attr("y2", cy).attr("stroke", "#ddd");
svg.append("line").attr("x1", cx).attr("y1", cy - r - 15).attr("x2", cx).attr("y2", cy + r + 15).attr("stroke", "#ddd");
const colors = ["steelblue", "firebrick", "steelblue", "firebrick"];
for (let i = 0; i < 4; i++) {
const angle = theta + i * Math.PI / 2;
const dx = r * Math.cos(angle), dy = -r * Math.sin(angle);
svg.append("line").attr("x1", cx).attr("y1", cy).attr("x2", cx + dx).attr("y2", cy + dy)
.attr("stroke", colors[i]).attr("stroke-width", 3);
const headLen = 10, headAngle = 0.4;
const a1 = Math.atan2(-dy, dx) + Math.PI + headAngle;
const a2 = Math.atan2(-dy, dx) + Math.PI - headAngle;
svg.append("polygon")
.attr("points", `${cx+dx},${cy+dy} ${cx+dx+headLen*Math.cos(a1)},${cy+dy-headLen*Math.sin(a1)} ${cx+dx+headLen*Math.cos(a2)},${cy+dy-headLen*Math.sin(a2)}`)
.attr("fill", colors[i]);
svg.append("text").attr("x", cx + 1.15*dx).attr("y", cy + 1.15*dy)
.attr("text-anchor", "middle").attr("dominant-baseline", "middle")
.attr("font-size", 14).attr("fill", colors[i]).text(`f${i}`);
}
if (theta > 0.01) {
const arc_pts = d3.range(0, theta + 0.01, 0.02).map(a => [cx + 35*Math.cos(a), cy - 35*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(arc_pts)).attr("fill", "none").attr("stroke", "orange").attr("stroke-width", 2);
svg.append("text").attr("x", cx + 50*Math.cos(theta/2)).attr("y", cy - 50*Math.sin(theta/2))
.attr("text-anchor", "middle").attr("font-size", 14).attr("fill", "orange").text("θ");
}
svg.append("text").attr("x", cx).attr("y", 18).attr("text-anchor", "middle").attr("font-size", 14).attr("fill", "#333")
.text(`2D Frame (θ = ${theta_deg_s}°)`);
return svg.node();
}| Term | Definition |
|---|---|
| Frame | 4 unit vectors invariant under \(\pi/2\) rotation, parameterized by angle \(\theta\) |
| Frame field | One frame per mesh vertex |
| Boundary constraint | Boundary frames must align one vector with the boundary normal |
| Rotation angle | \(\Delta\theta \pmod{\pi/2}\) between adjacent frames (minimum absolute value) |
| Field curvature | \(\displaystyle\sum_{\text{edges}} (\Delta\theta_{ij})^2\) |
| Singular triangle | Sum of rotation angles around triangle \(\neq 0\) |
Key insight: represent frames by functions, not angles.
Reference frame \(\tilde{f}\) (aligned with axes) \(\;\longrightarrow\;\) function \(\tilde{F}(\alpha) = \cos(4\alpha)\)
Any rotated frame: \(F(\alpha) = \cos\bigl(4(\alpha - \theta)\bigr)\)
Fourier decomposition via angle addition:
\[F(\alpha) = \underbrace{\cos(4\theta)}_{a_0} \cos(4\alpha) + \underbrace{\sin(4\theta)}_{a_1} \sin(4\alpha) = B(\alpha) \cdot a\]
A frame is fully described by a 2D coefficient vector:
\[a = \begin{pmatrix} \cos 4\theta \\ \sin 4\theta \end{pmatrix} \quad\text{on the unit circle } (a^\top a = 1)\]
Why frequency 4? Frame has 4-fold symmetry (\(\pi/2\) periodicity).
{
const theta = theta2_s * Math.PI / 180;
const W = 900, H = 280;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
// Panel 1: Frame vectors
const p1cx = 110, p1cy = H/2, p1r = 80;
svg.append("text").attr("x", p1cx).attr("y", 18).attr("text-anchor", "middle").attr("font-size", 12).attr("fill", "#333").text("Physical space");
const circ1 = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [p1cx+p1r*Math.cos(a), p1cy-p1r*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ1)).attr("fill", "none").attr("stroke", "#ddd");
const fc = ["steelblue","firebrick","steelblue","firebrick"];
for (let i = 0; i < 4; i++) {
const a = theta + i*Math.PI/2;
svg.append("line").attr("x1",p1cx).attr("y1",p1cy).attr("x2",p1cx+p1r*Math.cos(a)).attr("y2",p1cy-p1r*Math.sin(a))
.attr("stroke",fc[i]).attr("stroke-width",2.5);
}
// Panel 2: F(α)
const p2x = 270, p2w = 300, p2h = 220, p2top = 35;
const plotL = p2x+35, plotR = p2x+p2w-10, plotT = p2top+10, plotB = p2top+p2h-20;
const xScale = d3.scaleLinear().domain([0, 2*Math.PI]).range([plotL, plotR]);
const yScale = d3.scaleLinear().domain([-1.2, 1.2]).range([plotB, plotT]);
svg.append("text").attr("x",(plotL+plotR)/2).attr("y",18).attr("text-anchor","middle").attr("font-size",12).attr("fill","#333")
.text(`F(α) = cos(4(α − ${theta2_s}°))`);
svg.append("line").attr("x1",plotL).attr("y1",yScale(0)).attr("x2",plotR).attr("y2",yScale(0)).attr("stroke","#ccc");
svg.append("line").attr("x1",plotL).attr("y1",plotT).attr("x2",plotL).attr("y2",plotB).attr("stroke","#ccc");
const refPts = d3.range(0, 2*Math.PI, 0.02).map(a => [xScale(a), yScale(Math.cos(4*a))]);
svg.append("path").attr("d", d3.line()(refPts)).attr("fill","none").attr("stroke","#bbb").attr("stroke-width",1.5).attr("stroke-dasharray","4,3");
const rotPts = d3.range(0, 2*Math.PI, 0.02).map(a => [xScale(a), yScale(Math.cos(4*(a - theta)))]);
svg.append("path").attr("d", d3.line()(rotPts)).attr("fill","none").attr("stroke","firebrick").attr("stroke-width",2);
for (let i = 0; i < 4; i++) {
const peakAlpha = theta + i*Math.PI/2;
if (peakAlpha <= 2*Math.PI) {
svg.append("line").attr("x1",xScale(peakAlpha)).attr("y1",plotT).attr("x2",xScale(peakAlpha)).attr("y2",plotB)
.attr("stroke","red").attr("stroke-width",0.8).attr("stroke-dasharray","3,3");
}
}
svg.append("text").attr("x",plotR).attr("y",yScale(0)+14).attr("text-anchor","end").attr("font-size",10).attr("fill","#666").text("α");
// Panel 3: Coefficient space
const p3cx = 740, p3cy = H/2, p3r = 85;
svg.append("text").attr("x",p3cx).attr("y",18).attr("text-anchor","middle").attr("font-size",12).attr("fill","#333").text("Coefficient space");
const circ3 = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [p3cx+p3r*Math.cos(a), p3cy-p3r*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ3)).attr("fill","none").attr("stroke","#ccc").attr("stroke-width",1.5);
svg.append("line").attr("x1",p3cx-p3r-10).attr("y1",p3cy).attr("x2",p3cx+p3r+10).attr("y2",p3cy).attr("stroke","#ddd");
svg.append("line").attr("x1",p3cx).attr("y1",p3cy-p3r-10).attr("x2",p3cx).attr("y2",p3cy+p3r+10).attr("stroke","#ddd");
const a0 = Math.cos(4*theta), a1 = Math.sin(4*theta);
svg.append("line").attr("x1",p3cx).attr("y1",p3cy).attr("x2",p3cx+p3r*a0).attr("y2",p3cy-p3r*a1)
.attr("stroke","firebrick").attr("stroke-width",1.5).attr("stroke-dasharray","3,3");
svg.append("circle").attr("cx",p3cx+p3r*a0).attr("cy",p3cy-p3r*a1).attr("r",6).attr("fill","firebrick");
svg.append("text").attr("x",p3cx+p3r*a0+10).attr("y",p3cy-p3r*a1-8).attr("font-size",11).attr("fill","firebrick")
.text(`a = (${a0.toFixed(2)}, ${a1.toFixed(2)})`);
svg.append("text").attr("x",p3cx+p3r+12).attr("y",p3cy-5).attr("font-size",10).attr("fill","#666").text("a₀");
svg.append("text").attr("x",p3cx+5).attr("y",p3cy-p3r-8).attr("font-size",10).attr("fill","#666").text("a₁");
return svg.node();
}The frame function plotted in polar coordinates: \(r(\alpha) = 1 + F(\alpha)\)
Lobes point where frame vectors point — rotating \(\theta\) rotates the whole clover.
{
const theta = theta_polar_s * Math.PI / 180;
const W = 750, H = 300;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
const panels = [
{cx: 180, cy: H/2, label: "Reference F̃(α) = cos(4α)", th: 0, color: "#888"},
{cx: 570, cy: H/2, label: `Rotated F(α) = cos(4(α−${theta_polar_s}°))`, th: theta, color: "firebrick"}
];
const R = 110;
for (const p of panels) {
svg.append("text").attr("x", p.cx).attr("y", 16).attr("text-anchor", "middle").attr("font-size", 11).attr("fill", "#333").text(p.label);
const unitCirc = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [p.cx + R*0.5*Math.cos(a), p.cy - R*0.5*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(unitCirc)).attr("fill", "none").attr("stroke", "#eee");
svg.append("line").attr("x1", p.cx - R - 10).attr("y1", p.cy).attr("x2", p.cx + R + 10).attr("y2", p.cy).attr("stroke", "#eee");
svg.append("line").attr("x1", p.cx).attr("y1", p.cy - R - 10).attr("x2", p.cx).attr("y2", p.cy + R + 10).attr("stroke", "#eee");
const pts = d3.range(0, 2*Math.PI + 0.01, 0.015).map(a => {
const r = (1 + Math.cos(4*(a - p.th))) * R * 0.5;
return [p.cx + r*Math.cos(a), p.cy - r*Math.sin(a)];
});
svg.append("path").attr("d", d3.line()(pts)).attr("fill", p.color).attr("fill-opacity", 0.1)
.attr("stroke", p.color).attr("stroke-width", 2);
for (let i = 0; i < 4; i++) {
const a = p.th + i * Math.PI/2;
const dx = R * Math.cos(a), dy = -R * Math.sin(a);
svg.append("line").attr("x1", p.cx).attr("y1", p.cy).attr("x2", p.cx + dx).attr("y2", p.cy + dy)
.attr("stroke", "steelblue").attr("stroke-width", 1.5).attr("stroke-dasharray", "4,3");
svg.append("circle").attr("cx", p.cx + dx).attr("cy", p.cy + dy).attr("r", 3).attr("fill", "steelblue");
}
svg.append("circle").attr("cx", p.cx).attr("cy", p.cy).attr("r", 2.5).attr("fill", "#333");
}
svg.append("line").attr("x1", 300).attr("y1", H/2).attr("x2", 370).attr("y2", H/2).attr("stroke", "#aaa").attr("stroke-width", 1.5);
svg.append("polygon").attr("points", `370,${H/2} 362,${H/2-5} 362,${H/2+5}`).attr("fill", "#aaa");
svg.append("text").attr("x", 335).attr("y", H/2 - 10).attr("text-anchor", "middle").attr("font-size", 11).attr("fill", "#888").text("rotate by θ");
return svg.node();
}Difference between frames in function space:
\[E = \sum_{\text{edges}} \int_0^{2\pi} \bigl(F^j(\alpha) - F^i(\alpha)\bigr)^2 \, d\alpha\]
Since the Fourier basis is orthogonal (\(\int B^\top B\,d\alpha = \pi I\)):
\[\boxed{E = \pi \sum_{\text{edges}} \|a^j - a^i\|^2}\]
Tip
The smoothness energy is just Euclidean distance in coefficient space!
{
const tI = thetaI_s * Math.PI / 180, tJ = thetaJ_s * Math.PI / 180;
const W = 450, H = 400, cx = W/2, cy = H/2 + 10, r = 155;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
const circ = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [cx+r*Math.cos(a), cy-r*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ)).attr("fill","none").attr("stroke","#ddd").attr("stroke-width",1.5);
const ai = [Math.cos(4*tI), Math.sin(4*tI)];
const aj = [Math.cos(4*tJ), Math.sin(4*tJ)];
let arcStart = 4*tI, arcEnd = 4*tJ;
if (arcEnd < arcStart) [arcStart, arcEnd] = [arcEnd, arcStart];
const arcPts = d3.range(arcStart, arcEnd+0.01, 0.02).map(a => [cx+r*Math.cos(a), cy-r*Math.sin(a)]);
if (arcPts.length > 1)
svg.append("path").attr("d", d3.line()(arcPts)).attr("fill","none").attr("stroke","orange").attr("stroke-width",4).attr("opacity",0.8);
svg.append("line").attr("x1",cx+r*ai[0]).attr("y1",cy-r*ai[1]).attr("x2",cx+r*aj[0]).attr("y2",cy-r*aj[1])
.attr("stroke","steelblue").attr("stroke-width",3).attr("stroke-dasharray","6,4");
svg.append("circle").attr("cx",cx+r*ai[0]).attr("cy",cy-r*ai[1]).attr("r",7).attr("fill","firebrick");
svg.append("text").attr("x",cx+r*ai[0]-12).attr("y",cy-r*ai[1]-10).attr("font-size",14).attr("fill","firebrick").attr("font-weight","bold").text("aⁱ");
svg.append("circle").attr("cx",cx+r*aj[0]).attr("cy",cy-r*aj[1]).attr("r",7).attr("fill","forestgreen");
svg.append("text").attr("x",cx+r*aj[0]+8).attr("y",cy-r*aj[1]-10).attr("font-size",14).attr("fill","forestgreen").attr("font-weight","bold").text("aʲ");
const arcLen = Math.abs(arcEnd - arcStart);
const chordLen = Math.sqrt((aj[0]-ai[0])**2 + (aj[1]-ai[1])**2);
svg.append("text").attr("x",15).attr("y",H-35).attr("font-size",13).attr("fill","orange").text(`Arc (rotation): ${arcLen.toFixed(3)}`);
svg.append("text").attr("x",15).attr("y",H-15).attr("font-size",13).attr("fill","steelblue").text(`Chord ‖aʲ−aⁱ‖: ${chordLen.toFixed(3)}`);
svg.append("text").attr("x",W-15).attr("y",H-15).attr("text-anchor","end").attr("font-size",12).attr("fill","#666")
.text(`Ratio: ${arcLen > 0.001 ? (chordLen/arcLen).toFixed(3) : "—"}`);
return svg.node();
}Small rotations: arc \(\approx\) chord. Large rotations: chord is shorter \(\Rightarrow\) \(E\) is smoother near singularities.
Each coefficient vector must be a valid frame:
\[a^\top a = 1 \quad \text{(unit circle)}\]
Frames at boundary vertices align with the normal \(\theta^n\):
Eq. 2 (two equations — correct for init):
\[a_0 = \cos 4\theta^n, \quad a_1 = \sin 4\theta^n\]
Eq. 3 (one equation — sufficient with feasibility):
\[a_0 \cos 4\theta^n + a_1 \sin 4\theta^n = 1\]
{
const theta_n = normal_s * Math.PI / 180;
const W = 400, H = 350;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
const rcx = W/2, rcy = H/2 + 10, rr = 130;
svg.append("text").attr("x",rcx).attr("y",18).attr("text-anchor","middle").attr("font-size",13).attr("fill","#333").text("Coefficient space");
const circ = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [rcx+rr*Math.cos(a), rcy-rr*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ)).attr("fill","none").attr("stroke","#ccc").attr("stroke-width",2);
svg.append("line").attr("x1",rcx-rr-15).attr("y1",rcy).attr("x2",rcx+rr+15).attr("y2",rcy).attr("stroke","#eee");
svg.append("line").attr("x1",rcx).attr("y1",rcy-rr-15).attr("x2",rcx).attr("y2",rcy+rr+15).attr("stroke","#eee");
svg.append("text").attr("x",rcx+rr+17).attr("y",rcy-5).attr("font-size",10).attr("fill","#999").text("a₀");
svg.append("text").attr("x",rcx+5).attr("y",rcy-rr-10).attr("font-size",10).attr("fill","#999").text("a₁");
const c0 = Math.cos(4*theta_n), c1 = Math.sin(4*theta_n);
const linePts = [];
for (let t = -2; t <= 2; t += 0.01) {
const px = c0 + t*(-c1), py = c1 + t*c0;
const sx = rcx + rr*px, sy = rcy - rr*py;
linePts.push([sx, sy]);
}
const filtered = linePts.filter(([x,y]) => x > rcx-rr-15 && x < rcx+rr+15 && y > rcy-rr-15 && y < rcy+rr+15);
if (filtered.length > 1)
svg.append("path").attr("d", d3.line()(filtered)).attr("fill","none").attr("stroke","red").attr("stroke-width",1.5).attr("stroke-dasharray","5,4");
svg.append("circle").attr("cx",rcx+rr*c0).attr("cy",rcy-rr*c1).attr("r",7).attr("fill","orange");
svg.append("text").attr("x",rcx+rr*c0+10).attr("y",rcy-rr*c1-10).attr("font-size",11).attr("fill","orange")
.text(`a = (${c0.toFixed(2)}, ${c1.toFixed(2)})`);
svg.append("text").attr("x",rcx).attr("y",H-8).attr("text-anchor","middle").attr("font-size",10).attr("fill","red")
.text("Dashed = Eq. 3 constraint line");
return svg.node();
}Two endpoints locked, two interior frames free. Drag to interpolate:
viewof t0_lock_s = Inputs.range([0, 90], {value: 0, step: 1, label: "θ₀ locked"})
viewof t3_lock_s = Inputs.range([0, 90], {value: 60, step: 1, label: "θ₃ locked"})
viewof t1_free_s = Inputs.range([-45, 90], {value: 20, step: 1, label: "θ₁ free"})
viewof t2_free_s = Inputs.range([-45, 90], {value: 40, step: 1, label: "θ₂ free"}){
const t0 = t0_lock_s*Math.PI/180, t1 = t1_free_s*Math.PI/180;
const t2 = t2_free_s*Math.PI/180, t3 = t3_lock_s*Math.PI/180;
const thetas = [t0, t1, t2, t3];
const W = 800, H = 300;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
// Left: frame chain
const chainY = 120, xs = [40, 120, 200, 280];
const labels = ["P₀","P₁","P₂","P₃"];
const isLocked = [true, false, false, true];
const colrs = ["orange","steelblue","steelblue","orange"];
svg.append("text").attr("x", 160).attr("y", 20).attr("text-anchor","middle").attr("font-size",13).attr("fill","#333").text("Frame chain");
svg.append("line").attr("x1",xs[0]).attr("y1",chainY).attr("x2",xs[3]).attr("y2",chainY).attr("stroke","#bbb").attr("stroke-width",2);
for (let i = 0; i < 4; i++) {
for (let k = 0; k < 4; k++) {
const a = thetas[i] + k*Math.PI/2;
svg.append("line").attr("x1",xs[i]).attr("y1",chainY).attr("x2",xs[i]+25*Math.cos(a)).attr("y2",chainY-25*Math.sin(a))
.attr("stroke",colrs[i]).attr("stroke-width",2);
}
svg.append("circle").attr("cx",xs[i]).attr("cy",chainY).attr("r",4).attr("fill",colrs[i]);
svg.append("text").attr("x",xs[i]).attr("y",chainY+40).attr("text-anchor","middle").attr("font-size",11).attr("fill",colrs[i])
.text(labels[i] + (isLocked[i] ? " 🔒" : ""));
}
// Energy
function rotAngle(a, b) { let d = b - a; d = ((d + Math.PI/4) % (Math.PI/2) + Math.PI/2) % (Math.PI/2) - Math.PI/4; return d; }
const fieldCurv = rotAngle(t0,t1)**2 + rotAngle(t1,t2)**2 + rotAngle(t2,t3)**2;
const coeffs = thetas.map(t => [Math.cos(4*t), Math.sin(4*t)]);
let objE = 0;
for (let e = 0; e < 3; e++) objE += (coeffs[e+1][0]-coeffs[e][0])**2 + (coeffs[e+1][1]-coeffs[e][1])**2;
objE *= Math.PI;
svg.append("text").attr("x",20).attr("y",200).attr("font-size",12).attr("fill","purple").text(`Curvature Σ(Δθ)² = ${fieldCurv.toFixed(4)}`);
svg.append("text").attr("x",20).attr("y",220).attr("font-size",12).attr("fill","teal").text(`E = π Σ‖aʲ−aⁱ‖² = ${objE.toFixed(4)}`);
// Right: coefficient space
const rcx = 600, rcy = H/2, rr = 120;
svg.append("text").attr("x",rcx).attr("y",20).attr("text-anchor","middle").attr("font-size",13).attr("fill","#333").text("Coefficient space");
const circ = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [rcx+rr*Math.cos(a), rcy-rr*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ)).attr("fill","none").attr("stroke","#ddd").attr("stroke-width",1.5);
for (let e = 0; e < 3; e++) {
svg.append("line")
.attr("x1",rcx+rr*coeffs[e][0]).attr("y1",rcy-rr*coeffs[e][1])
.attr("x2",rcx+rr*coeffs[e+1][0]).attr("y2",rcy-rr*coeffs[e+1][1])
.attr("stroke","#aaa").attr("stroke-width",1.5);
}
for (let i = 0; i < 4; i++) {
svg.append("circle").attr("cx",rcx+rr*coeffs[i][0]).attr("cy",rcy-rr*coeffs[i][1]).attr("r",6).attr("fill",colrs[i]);
svg.append("text").attr("x",rcx+rr*coeffs[i][0]+10).attr("y",rcy-rr*coeffs[i][1]-8).attr("font-size",11).attr("fill",colrs[i]).text(labels[i]);
}
return svg.node();
}Both energies share the same minima, but the functional energy \(E\) is smooth (no creases).
{
const t0 = land_t0_s * Math.PI / 180, t3 = land_t3_s * Math.PI / 180;
const N = 70;
const lo = -Math.PI/4, hi = Math.PI/2;
const step = (hi - lo) / (N - 1);
function rotAngle(a, b) { let d = b - a; d = ((d + Math.PI/4) % (Math.PI/2) + Math.PI/2) % (Math.PI/2) - Math.PI/4; return d; }
const curvGrid = [], objGrid = [];
let curvMin = Infinity, objMin = Infinity, curvMax = 0, objMax = 0;
for (let j = 0; j < N; j++) {
curvGrid[j] = []; objGrid[j] = [];
const t2 = lo + j * step;
for (let i = 0; i < N; i++) {
const t1 = lo + i * step;
const curv = rotAngle(t0,t1)**2 + rotAngle(t1,t2)**2 + rotAngle(t2,t3)**2;
const c0 = [Math.cos(4*t0), Math.sin(4*t0)];
const c1 = [Math.cos(4*t1), Math.sin(4*t1)];
const c2 = [Math.cos(4*t2), Math.sin(4*t2)];
const c3 = [Math.cos(4*t3), Math.sin(4*t3)];
const obj = ((c1[0]-c0[0])**2+(c1[1]-c0[1])**2 + (c2[0]-c1[0])**2+(c2[1]-c1[1])**2 + (c3[0]-c2[0])**2+(c3[1]-c2[1])**2) * Math.PI;
curvGrid[j][i] = curv; objGrid[j][i] = obj;
curvMin = Math.min(curvMin, curv); curvMax = Math.max(curvMax, curv);
objMin = Math.min(objMin, obj); objMax = Math.max(objMax, obj);
}
}
function findMinima(grid) {
let mins = [];
for (let j = 1; j < N-1; j++) for (let i = 1; i < N-1; i++) {
const v = grid[j][i];
if (v < grid[j-1][i] && v < grid[j+1][i] && v < grid[j][i-1] && v < grid[j][i+1]
&& v < grid[j-1][i-1] && v < grid[j-1][i+1] && v < grid[j+1][i-1] && v < grid[j+1][i+1])
mins.push({i, j, v, t1: lo + i*step, t2: lo + j*step});
}
return mins.sort((a,b) => a.v - b.v);
}
function viridis(t) {
const r = Math.max(0, Math.min(255, Math.round(68 + t * (253-68))));
const g = Math.max(0, Math.min(255, Math.round(1 + t * (231-1))));
const b = Math.max(0, Math.min(255, Math.round(84 + (t < 0.5 ? t*2*(168-84) : (168 + (t-0.5)*2*(37-168))))));
return `rgb(${r},${g},${b})`;
}
const W = 800, H = 350;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
const cSize = 270, cTop = 30;
const panels = [
{x: 20, grid: curvGrid, min: curvMin, max: curvMax, title: "Field curvature Σ(Δθ)²"},
{x: 420, grid: objGrid, min: objMin, max: objMax, title: "Functional E = π Σ‖aʲ−aⁱ‖²"}
];
for (const p of panels) {
svg.append("text").attr("x", p.x + cSize/2).attr("y", 16).attr("text-anchor","middle").attr("font-size",12).attr("fill","#333").text(p.title);
const cellW = cSize / N, cellH = cSize / N;
const range = p.max - p.min || 1;
for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) {
const t = (p.grid[j][i] - p.min) / range;
svg.append("rect").attr("x", p.x + i * cellW).attr("y", cTop + (N - 1 - j) * cellH)
.attr("width", cellW + 0.5).attr("height", cellH + 0.5).attr("fill", viridis(t));
}
svg.append("text").attr("x", p.x + cSize/2).attr("y", cTop + cSize + 16).attr("text-anchor","middle").attr("font-size",10).attr("fill","#666").text("θ₁");
svg.append("text").attr("x", p.x - 5).attr("y", cTop + cSize/2).attr("text-anchor","middle").attr("font-size",10).attr("fill","#666")
.attr("transform", `rotate(-90, ${p.x - 5}, ${cTop + cSize/2})`).text("θ₂");
const mins = findMinima(p.grid);
const mC = ["teal","purple"];
for (let k = 0; k < Math.min(mins.length, 2); k++) {
const m = mins[k];
const mx = p.x + (m.i / (N-1)) * cSize, my = cTop + (1 - m.j / (N-1)) * cSize;
svg.append("circle").attr("cx", mx).attr("cy", my).attr("r", 6).attr("fill","none").attr("stroke",mC[k]).attr("stroke-width",2.5);
}
}
svg.append("text").attr("x", W/2).attr("y", H - 5).attr("text-anchor","middle").attr("font-size",11).attr("fill","#888")
.text("Left has creases (piecewise quadratic) — right is globally smooth. Same minima.");
return svg.node();
}Step 1 — Initialization (relax feasibility):
Minimize \(\|AX - b\|^2\) where:
Solve \(A^\top A X = A^\top b\) (linear least-squares)
Step 2 — Project onto feasibility:
\[a^i \leftarrow \frac{(X_{2i}, X_{2i+1})^\top}{\|(X_{2i}, X_{2i+1})\|}\]
Step 3 — Smooth (iterate):
Add feasibility penalty rows: \(\lambda(\tilde{a}_0 X_{2i} + \tilde{a}_1 X_{2i+1}) = \lambda\)
Re-solve and re-project each iteration.
{
const t0 = impl_t0_s*Math.PI/180, t3 = impl_t3_s*Math.PI/180;
const nv = 4, lam = 100, lam_f = 100;
const edgeList = [[0,1],[1,2],[2,3]];
const bdry = {0: t0, 3: t3};
function solveLS(n, getRowsFn) {
const AtA = Array.from({length: n}, () => new Float64Array(n));
const Atb = new Float64Array(n);
for (const {coeffs, rhs} of getRowsFn()) {
for (let i = 0; i < n; i++) { for (let j = 0; j < n; j++) AtA[i][j] += coeffs[i] * coeffs[j]; Atb[i] += coeffs[i] * rhs; }
}
const A = AtA.map(r => [...r]); const b = [...Atb];
for (let col = 0; col < n; col++) {
let maxR = col; for (let r = col+1; r < n; r++) if (Math.abs(A[r][col]) > Math.abs(A[maxR][col])) maxR = r;
[A[col], A[maxR]] = [A[maxR], A[col]]; [b[col], b[maxR]] = [b[maxR], b[col]];
if (Math.abs(A[col][col]) < 1e-15) continue;
for (let r = col+1; r < n; r++) { const f = A[r][col]/A[col][col]; for (let c = col; c < n; c++) A[r][c] -= f*A[col][c]; b[r] -= f*b[col]; }
}
const x = new Float64Array(n);
for (let i = n-1; i >= 0; i--) { let s = b[i]; for (let j = i+1; j < n; j++) s -= A[i][j]*x[j]; x[i] = s / A[i][i]; }
return x;
}
function buildRows(a_prev) {
return function*() {
const sqpi = Math.sqrt(Math.PI);
for (const [vi,vj] of edgeList) for (let comp = 0; comp < 2; comp++) {
const c = new Float64Array(2*nv); c[2*vi+comp] = sqpi; c[2*vj+comp] = -sqpi; yield {coeffs: c, rhs: 0};
}
for (const [vi_s, ti] of Object.entries(bdry)) {
const vi = parseInt(vi_s);
const c0 = new Float64Array(2*nv); c0[2*vi] = lam; yield {coeffs: c0, rhs: lam*Math.cos(4*ti)};
const c1 = new Float64Array(2*nv); c1[2*vi+1] = lam; yield {coeffs: c1, rhs: lam*Math.sin(4*ti)};
}
if (a_prev) for (let vi = 0; vi < nv; vi++) {
const c = new Float64Array(2*nv); c[2*vi] = lam_f*a_prev[vi][0]; c[2*vi+1] = lam_f*a_prev[vi][1];
yield {coeffs: c, rhs: lam_f};
}
};
}
const X0 = solveLS(2*nv, buildRows(null));
const a_relaxed = []; for (let i = 0; i < nv; i++) a_relaxed.push([X0[2*i], X0[2*i+1]]);
let a_current = a_relaxed.map(a => { const n = Math.sqrt(a[0]**2+a[1]**2); return n > 1e-10 ? [a[0]/n, a[1]/n] : [1, 0]; });
for (let it = 0; it < n_smooth_s; it++) {
const X = solveLS(2*nv, buildRows(a_current));
a_current = []; for (let i = 0; i < nv; i++) { const v = [X[2*i], X[2*i+1]]; const n = Math.sqrt(v[0]**2+v[1]**2); a_current.push(n > 1e-10 ? [v[0]/n, v[1]/n] : [1,0]); }
}
const W = 380, H = 330;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
const colrs = ["orange","steelblue","forestgreen","orange"];
const lbls = ["P₀","P₁","P₂","P₃"];
const lcx = W/2, lcy = H/2, lr = 130;
svg.append("text").attr("x",lcx).attr("y",16).attr("text-anchor","middle").attr("font-size",12).attr("fill","#333")
.text(n_smooth_s === 0 ? "Init + projection" : `${n_smooth_s} smoothing iter.`);
const circ = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [lcx+lr*Math.cos(a), lcy-lr*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ)).attr("fill","none").attr("stroke","#ddd").attr("stroke-width",2);
if (n_smooth_s === 0) {
for (let i = 0; i < nv; i++) {
const rx = lcx+lr*a_relaxed[i][0], ry = lcy-lr*a_relaxed[i][1];
svg.append("rect").attr("x",rx-4).attr("y",ry-4).attr("width",8).attr("height",8)
.attr("fill","none").attr("stroke",colrs[i]).attr("stroke-width",1.5).attr("transform",`rotate(45,${rx},${ry})`);
}
}
for (const [vi,vj] of edgeList) svg.append("line")
.attr("x1",lcx+lr*a_current[vi][0]).attr("y1",lcy-lr*a_current[vi][1])
.attr("x2",lcx+lr*a_current[vj][0]).attr("y2",lcy-lr*a_current[vj][1]).attr("stroke","#aaa").attr("stroke-width",1.5);
for (let i = 0; i < nv; i++) {
svg.append("circle").attr("cx",lcx+lr*a_current[i][0]).attr("cy",lcy-lr*a_current[i][1]).attr("r",6).attr("fill",colrs[i]);
svg.append("text").attr("x",lcx+lr*a_current[i][0]+10).attr("y",lcy-lr*a_current[i][1]-8).attr("font-size",12).attr("fill",colrs[i]).text(lbls[i]);
}
let objE = 0;
for (const [vi,vj] of edgeList) objE += (a_current[vj][0]-a_current[vi][0])**2 + (a_current[vj][1]-a_current[vi][1])**2;
objE *= Math.PI;
svg.append("text").attr("x",15).attr("y",H-10).attr("font-size",11).attr("fill","teal").text(`E = ${objE.toFixed(4)}`);
if (n_smooth_s === 0) {
svg.append("rect").attr("x",15).attr("y",H-35).attr("width",7).attr("height",7).attr("fill","none").attr("stroke","#888").attr("transform",`rotate(45,18.5,${H-31.5})`);
svg.append("text").attr("x",28).attr("y",H-28).attr("font-size",10).attr("fill","#888").text("relaxed");
svg.append("circle").attr("cx",100).attr("cy",H-31).attr("r",4).attr("fill","#888");
svg.append("text").attr("x",108).attr("y",H-28).attr("font-size",10).attr("fill","#888").text("projected");
}
return svg.node();
}Eq. 3 (single constraint) fails during initialization on non-rectangular domains — it constrains each boundary vertex to a line, not a point. The optimizer slides along that line.
Fix: use Eq. 2 (both components) during initialization.
{
const skew = skew_s * Math.PI / 180;
const nx = 6, ny = 5, nv = nx * ny;
const baseW = 3, baseH = 2;
const verts = [];
for (let j = 0; j < ny; j++) { const ty = j / (ny-1); const shear = ty * baseH * Math.tan(skew);
for (let i = 0; i < nx; i++) { const tx = i / (nx-1); verts.push([tx * baseW + shear, ty * baseH]); }
}
const edges = [];
for (let j = 0; j < ny; j++) for (let i = 0; i < nx-1; i++) edges.push([j*nx+i, j*nx+i+1]);
for (let j = 0; j < ny-1; j++) for (let i = 0; i < nx; i++) edges.push([j*nx+i, (j+1)*nx+i]);
const rightNormal = Math.atan2(-Math.tan(skew), 1), leftNormal = Math.atan2(Math.tan(skew), -1);
const bdryMap = {};
for (let i = 0; i < nx; i++) { bdryMap[i] = 3*Math.PI/2; bdryMap[(ny-1)*nx+i] = Math.PI/2; }
for (let j = 1; j < ny-1; j++) { bdryMap[j*nx] = leftNormal; bdryMap[j*nx+nx-1] = rightNormal; }
const useEq2 = constraint_s.startsWith("Eq. 2");
const lam = 100, lam_f = 100, sqpi = Math.sqrt(Math.PI);
function solveSystem(n, rowGen) {
const AtA = Array.from({length: n}, () => new Float64Array(n)); const Atb = new Float64Array(n);
for (const {coeffs, rhs} of rowGen()) { for (let ii = 0; ii < n; ii++) { for (let jj = 0; jj < n; jj++) AtA[ii][jj] += coeffs[ii]*coeffs[jj]; Atb[ii] += coeffs[ii]*rhs; } }
const A = AtA.map(r => [...r]); const b = [...Atb];
for (let col = 0; col < n; col++) { let maxR = col; for (let r = col+1; r < n; r++) if (Math.abs(A[r][col]) > Math.abs(A[maxR][col])) maxR = r;
[A[col], A[maxR]] = [A[maxR], A[col]]; [b[col], b[maxR]] = [b[maxR], b[col]];
if (Math.abs(A[col][col]) < 1e-15) continue;
for (let r = col+1; r < n; r++) { const f = A[r][col]/A[col][col]; for (let c = col; c < n; c++) A[r][c] -= f*A[col][c]; b[r] -= f*b[col]; }
}
const x = new Float64Array(n);
for (let ii = n-1; ii >= 0; ii--) { let s = b[ii]; for (let jj = ii+1; jj < n; jj++) s -= A[ii][jj]*x[jj]; x[ii] = s / A[ii][ii]; }
return x;
}
function buildRows(a_prev) {
return function*() {
for (const [vi,vj] of edges) for (let comp = 0; comp < 2; comp++) {
const c = new Float64Array(2*nv); c[2*vi+comp] = sqpi; c[2*vj+comp] = -sqpi; yield {coeffs: c, rhs: 0};
}
for (const [vi_s, theta_n] of Object.entries(bdryMap)) {
const vi = parseInt(vi_s);
if (useEq2) {
const c0 = new Float64Array(2*nv); c0[2*vi] = lam; yield {coeffs: c0, rhs: lam*Math.cos(4*theta_n)};
const c1 = new Float64Array(2*nv); c1[2*vi+1] = lam; yield {coeffs: c1, rhs: lam*Math.sin(4*theta_n)};
} else {
const c = new Float64Array(2*nv); c[2*vi] = lam*Math.cos(4*theta_n); c[2*vi+1] = lam*Math.sin(4*theta_n); yield {coeffs: c, rhs: lam};
}
}
if (a_prev) for (let vi = 0; vi < nv; vi++) {
const c = new Float64Array(2*nv); c[2*vi] = lam_f*a_prev[vi][0]; c[2*vi+1] = lam_f*a_prev[vi][1]; yield {coeffs: c, rhs: lam_f};
}
};
}
const X0 = solveSystem(2*nv, buildRows(null));
let a_current = [];
for (let i = 0; i < nv; i++) { const v = [X0[2*i], X0[2*i+1]]; const nm = Math.sqrt(v[0]**2+v[1]**2); a_current.push(nm > 1e-10 ? [v[0]/nm, v[1]/nm] : [1,0]); }
for (let it = 0; it < para_smooth_s; it++) {
const X = solveSystem(2*nv, buildRows(a_current)); a_current = [];
for (let i = 0; i < nv; i++) { const v = [X[2*i], X[2*i+1]]; const nm = Math.sqrt(v[0]**2+v[1]**2); a_current.push(nm > 1e-10 ? [v[0]/nm, v[1]/nm] : [1,0]); }
}
const W = 700, H = 350;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width", "100%").style("max-width", W + "px");
let minX = Infinity, maxX = -Infinity, minY = Infinity, maxY = -Infinity;
for (const v of verts) { minX = Math.min(minX, v[0]); maxX = Math.max(maxX, v[0]); minY = Math.min(minY, v[1]); maxY = Math.max(maxY, v[1]); }
const meshL = 30, meshR = 350, meshT = 40, meshB = 310;
const meshW = meshR - meshL, meshH = meshB - meshT;
const dataW = maxX - minX || 1, dataH = maxY - minY || 1;
const scale = Math.min(meshW / dataW, meshH / dataH) * 0.85;
const offX = meshL + (meshW - dataW * scale) / 2, offY = meshT + (meshH - dataH * scale) / 2;
function tx(v) { return offX + (v[0] - minX) * scale; }
function ty(v) { return meshB - (offY - meshT) - (v[1] - minY) * scale; }
svg.append("text").attr("x", (meshL+meshR)/2).attr("y", 18).attr("text-anchor","middle").attr("font-size",12).attr("fill","#333")
.text(useEq2 ? "Eq. 2: frames follow the boundary ✓" : "Eq. 3: frames ignore the slant ✗");
for (const [vi, vj] of edges) svg.append("line").attr("x1", tx(verts[vi])).attr("y1", ty(verts[vi]))
.attr("x2", tx(verts[vj])).attr("y2", ty(verts[vj])).attr("stroke", "#ddd").attr("stroke-width", 0.8);
let minEdgeLen = Infinity;
for (const [vi, vj] of edges) { const dx = tx(verts[vi]) - tx(verts[vj]); const dy = ty(verts[vi]) - ty(verts[vj]); minEdgeLen = Math.min(minEdgeLen, Math.sqrt(dx*dx + dy*dy)); }
const frameLen = Math.max(4, minEdgeLen * 0.42);
for (let vi = 0; vi < nv; vi++) {
const px = tx(verts[vi]), py = ty(verts[vi]);
const a = a_current[vi]; const th = Math.atan2(a[1], a[0]) / 4;
const isBdry = vi in bdryMap;
for (let k = 0; k < 4; k++) { const ang = th + k*Math.PI/2;
svg.append("line").attr("x1", px).attr("y1", py).attr("x2", px + frameLen*Math.cos(ang)).attr("y2", py - frameLen*Math.sin(ang))
.attr("stroke", isBdry ? "orange" : "steelblue").attr("stroke-width", 1.5);
}
}
// Coefficient space
const rcx = 540, rcy = H/2, rr = 120;
svg.append("text").attr("x", rcx).attr("y", 18).attr("text-anchor","middle").attr("font-size",12).attr("fill","#333").text("Coefficient space");
const circ = d3.range(0, 2*Math.PI+0.01, 0.02).map(a => [rcx+rr*Math.cos(a), rcy-rr*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ)).attr("fill","none").attr("stroke","#ddd").attr("stroke-width",1.5);
for (let vi = 0; vi < nv; vi++) {
const a = a_current[vi]; const isBdry = vi in bdryMap;
svg.append("circle").attr("cx", rcx + rr*a[0]).attr("cy", rcy - rr*a[1]).attr("r", 3)
.attr("fill", isBdry ? "orange" : "steelblue").attr("opacity", 0.8);
}
return svg.node();
}Click any frame to lock/unlock. Drag a locked frame to rotate it — the solver re-runs live.
{
const nInner = sbs_n_inner, nv = nInner + 2, lam = sbs_lambda, nSmooth = sbs_smooth, pi = Math.PI;
const locks = sbs_locks;
const edgeList = []; for (let i = 0; i < nv - 1; i++) edgeList.push([i, i + 1]);
function lockKey(i) { return i === 0 ? "L" : i === nv - 1 ? "R" : i; }
function solveLS(n, rowIter) {
const AtA = Array.from({length: n}, () => new Float64Array(n)), Atb = new Float64Array(n);
for (const {coeffs, rhs} of rowIter) {
const nz = []; for (let i = 0; i < n; i++) if (coeffs[i] !== 0) nz.push(i);
for (const i of nz) { for (const j of nz) AtA[i][j] += coeffs[i] * coeffs[j]; Atb[i] += coeffs[i] * rhs; }
}
const A = AtA.map(r => [...r]), b = [...Atb];
for (let col = 0; col < n; col++) {
let maxR = col; for (let r = col+1; r < n; r++) if (Math.abs(A[r][col]) > Math.abs(A[maxR][col])) maxR = r;
[A[col], A[maxR]] = [A[maxR], A[col]]; [b[col], b[maxR]] = [b[maxR], b[col]];
if (Math.abs(A[col][col]) < 1e-15) continue;
for (let r = col+1; r < n; r++) { const f = A[r][col]/A[col][col]; for (let c = col; c < n; c++) A[r][c] -= f*A[col][c]; b[r] -= f*b[col]; }
}
const x = new Float64Array(n);
for (let i = n-1; i >= 0; i--) { let s = b[i]; for (let j = i+1; j < n; j++) s -= A[i][j]*x[j]; x[i] = s/A[i][i]; }
return x;
}
function* makeRows(bdryMap, a_prev) {
const sqpi = Math.sqrt(pi);
for (const [vi, vj] of edgeList) for (let comp = 0; comp < 2; comp++) {
const c = new Float64Array(2*nv); c[2*vi+comp] = sqpi; c[2*vj+comp] = -sqpi; yield {coeffs: c, rhs: 0};
}
for (const [vi, ti] of bdryMap) {
const c0 = new Float64Array(2*nv); c0[2*vi] = lam; yield {coeffs: c0, rhs: lam*Math.cos(4*ti)};
const c1 = new Float64Array(2*nv); c1[2*vi+1] = lam; yield {coeffs: c1, rhs: lam*Math.sin(4*ti)};
}
if (a_prev) for (let vi = 0; vi < nv; vi++) {
const c = new Float64Array(2*nv); c[2*vi] = lam*a_prev[vi][0]; c[2*vi+1] = lam*a_prev[vi][1]; yield {coeffs: c, rhs: lam};
}
}
function solveField(bdryMap, iters) {
const X0 = solveLS(2*nv, makeRows(bdryMap, null));
let a = []; for (let i = 0; i < nv; i++) { const v0 = X0[2*i], v1 = X0[2*i+1]; const nm = Math.sqrt(v0**2+v1**2); a.push(nm > 1e-10 ? [v0/nm, v1/nm] : [1,0]); }
for (let it = 0; it < iters; it++) {
const X = solveLS(2*nv, makeRows(bdryMap, a)); a = [];
for (let i = 0; i < nv; i++) { const v0 = X[2*i], v1 = X[2*i+1]; const nm = Math.sqrt(v0**2+v1**2); a.push(nm > 1e-10 ? [v0/nm, v1/nm] : [1,0]); }
}
return a;
}
function buildBdry(lo) {
const b = new Map();
if (typeof lo.L === "number") b.set(0, lo.L*pi/180);
if (typeof lo.R === "number") b.set(nv-1, lo.R*pi/180);
for (let i = 1; i <= nInner; i++) if (typeof lo[i] === "number") b.set(i, lo[i]*pi/180);
return b;
}
function computeEnergy(a) { let E = 0; for (const [vi,vj] of edgeList) E += (a[vj][0]-a[vi][0])**2 + (a[vj][1]-a[vi][1])**2; return E*pi; }
const bdry = buildBdry(locks);
let a_current = solveField(bdry, nSmooth);
const W = 900, H = 380, chainW = 520, coeffW = 340;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width","100%").style("max-width",W+"px");
const chainMarginL = 40, chainMarginR = 30;
const chainSpacing = (chainW - chainMarginL - chainMarginR) / Math.max(1, nv-1);
const chainY = H/2, frameR = Math.min(30, chainSpacing * 0.4);
svg.append("text").attr("x", chainW/2).attr("y", 18).attr("text-anchor","middle").attr("font-size",13).attr("fill","#333")
.text("Click to lock/unlock · Drag locked frames to rotate");
svg.append("line").attr("x1",chainMarginL).attr("y1",chainY).attr("x2",chainMarginL+(nv-1)*chainSpacing).attr("y2",chainY).attr("stroke","#bbb").attr("stroke-width",2);
const rcx = chainW + coeffW/2 + 10, rcy = H/2, rr = Math.min(150, coeffW/2 - 20);
svg.append("text").attr("x",rcx).attr("y",18).attr("text-anchor","middle").attr("font-size",13).attr("fill","#333").text("Coefficient space");
const circ = d3.range(0, 2*pi+0.01, 0.02).map(a => [rcx+rr*Math.cos(a), rcy-rr*Math.sin(a)]);
svg.append("path").attr("d", d3.line()(circ)).attr("fill","none").attr("stroke","#ddd").attr("stroke-width",1.5);
svg.append("line").attr("x1",rcx-rr-10).attr("y1",rcy).attr("x2",rcx+rr+10).attr("y2",rcy).attr("stroke","#eee");
svg.append("line").attr("x1",rcx).attr("y1",rcy-rr-10).attr("x2",rcx).attr("y2",rcy+rr+10).attr("stroke","#eee");
const frameGs = [], coeffDots = [], coeffLabels = [], coeffEdgeLines = [];
for (let e = 0; e < edgeList.length; e++) {
const [vi,vj] = edgeList[e];
coeffEdgeLines.push(svg.append("line")
.attr("x1",rcx+rr*a_current[vi][0]).attr("y1",rcy-rr*a_current[vi][1])
.attr("x2",rcx+rr*a_current[vj][0]).attr("y2",rcy-rr*a_current[vj][1]).attr("stroke","#aaa").attr("stroke-width",1.5));
}
for (let i = 0; i < nv; i++) {
const cx = chainMarginL + i*chainSpacing;
const th = Math.atan2(a_current[i][1], a_current[i][0]) / 4;
const isLocked = bdry.has(i);
const color = isLocked ? "#d97706" : "steelblue";
const g = svg.append("g");
for (let k = 0; k < 4; k++) { const ang = th + k*pi/2;
g.append("line").attr("x1",cx).attr("y1",chainY).attr("x2",cx+frameR*Math.cos(ang)).attr("y2",chainY-frameR*Math.sin(ang)).attr("stroke",color).attr("stroke-width",2);
}
frameGs.push(g);
if (isLocked) svg.append("circle").attr("cx",cx).attr("cy",chainY).attr("r",frameR+5).attr("fill","none").attr("stroke","#d97706").attr("stroke-width",1.5).attr("stroke-dasharray","4,3");
svg.append("circle").attr("cx",cx).attr("cy",chainY).attr("r",4).attr("fill",color);
const label = i === 0 ? "P₀" : i === nv-1 ? `P${nv-1}` : `P${i}`;
const angleDeg = isLocked ? locks[lockKey(i)] : null;
svg.append("text").attr("x",cx).attr("y",chainY+frameR+18).attr("text-anchor","middle").attr("font-size",10).attr("fill",color)
.text(isLocked ? `${label} 🔒 ${angleDeg}°` : label);
coeffDots.push(svg.append("circle").attr("cx",rcx+rr*a_current[i][0]).attr("cy",rcy-rr*a_current[i][1]).attr("r",6).attr("fill",color));
coeffLabels.push(svg.append("text").attr("x",rcx+rr*a_current[i][0]+9).attr("y",rcy-rr*a_current[i][1]-9).attr("font-size",11).attr("fill",color).text(label));
const hitTarget = svg.append("circle").attr("cx",cx).attr("cy",chainY).attr("r",frameR+8).attr("fill","transparent").attr("cursor","pointer");
let isDragging = false;
const drag = d3.drag()
.on("start", () => { isDragging = false; })
.on("drag", (event) => {
if (!isLocked) return; isDragging = true;
const dx = event.x - cx, dy = -(event.y - chainY);
const deg90 = ((Math.atan2(dy, dx)*180/pi) % 90 + 90) % 90;
const tl = {...locks}; tl[lockKey(i)] = deg90;
const tb = buildBdry(tl), as = solveField(tb, nSmooth);
for (let v = 0; v < nv; v++) {
const vx = chainMarginL+v*chainSpacing, vth = Math.atan2(as[v][1], as[v][0])/4;
const vc = tb.has(v) ? "#d97706" : "steelblue";
frameGs[v].selectAll("line").remove();
for (let k = 0; k < 4; k++) { const ang = vth+k*pi/2;
frameGs[v].append("line").attr("x1",vx).attr("y1",chainY).attr("x2",vx+frameR*Math.cos(ang)).attr("y2",chainY-frameR*Math.sin(ang)).attr("stroke",vc).attr("stroke-width",2);
}
coeffDots[v].attr("cx",rcx+rr*as[v][0]).attr("cy",rcy-rr*as[v][1]);
coeffLabels[v].attr("x",rcx+rr*as[v][0]+9).attr("y",rcy-rr*as[v][1]-9);
}
for (let e = 0; e < edgeList.length; e++) { const [vi,vj] = edgeList[e];
coeffEdgeLines[e].attr("x1",rcx+rr*as[vi][0]).attr("y1",rcy-rr*as[vi][1]).attr("x2",rcx+rr*as[vj][0]).attr("y2",rcy-rr*as[vj][1]);
}
energyText.text(`E = ${computeEnergy(as).toFixed(4)}`);
})
.on("end", (event) => {
if (isDragging && isLocked) {
const dx = event.x-cx, dy = -(event.y-chainY);
const deg90 = Math.round(((Math.atan2(dy,dx)*180/pi) % 90 + 90) % 90);
const nl = {...mutable sbs_locks}; nl[lockKey(i)] = deg90; mutable sbs_locks = nl;
} else if (!isDragging) {
const nl = {...mutable sbs_locks}, key = lockKey(i);
if (typeof nl[key] === "number") delete nl[key];
else { const ang = Math.atan2(a_current[i][1], a_current[i][0])/4*180/pi; nl[key] = Math.round(((ang%90)+90)%90); }
mutable sbs_locks = nl;
}
});
hitTarget.call(drag);
}
const energyText = svg.append("text").attr("x",chainW+15).attr("y",H-15).attr("font-size",12).attr("fill","teal").text(`E = ${computeEnergy(a_current).toFixed(4)}`);
return svg.node();
}Same interaction on a 2D grid. Boundary frames start locked to outward normals.
viewof rects_nx = Inputs.range([3, 10], {value: 5, step: 1, label: "Columns"})
viewof rects_ny = Inputs.range([3, 8], {value: 4, step: 1, label: "Rows"})
viewof rects_smooth = Inputs.range([0, 10], {value: 3, step: 1, label: "Smoothing iters"})
viewof rects_lambda = Inputs.range([1, 200], {value: 100, step: 1, label: "λ"}){
const nx = rects_nx, ny = rects_ny, nv = nx*ny, lam = rects_lambda, nSmooth = rects_smooth, pi = Math.PI;
const locks = rects_locks;
const baseW = 3, baseH = 2;
const verts = []; for (let j = 0; j < ny; j++) for (let i = 0; i < nx; i++) verts.push([i/(nx-1)*baseW, j/(ny-1)*baseH]);
const edges = [];
for (let j = 0; j < ny; j++) for (let i = 0; i < nx-1; i++) edges.push([j*nx+i, j*nx+i+1]);
for (let j = 0; j < ny-1; j++) for (let i = 0; i < nx; i++) edges.push([j*nx+i, (j+1)*nx+i]);
const defaultBdryAngle = {};
for (let i = 0; i < nx; i++) { defaultBdryAngle[i] = 3*pi/2; defaultBdryAngle[(ny-1)*nx+i] = pi/2; }
for (let j = 1; j < ny-1; j++) { defaultBdryAngle[j*nx] = pi; defaultBdryAngle[j*nx+nx-1] = 0; }
const sqpi = Math.sqrt(pi);
function solveSystem(n, rowIter) {
const AtA = Array.from({length: n}, () => new Float64Array(n)), Atb = new Float64Array(n);
for (const {coeffs, rhs} of rowIter) {
const nz = []; for (let ii = 0; ii < n; ii++) if (coeffs[ii] !== 0) nz.push(ii);
for (const ii of nz) { for (const jj of nz) AtA[ii][jj] += coeffs[ii]*coeffs[jj]; Atb[ii] += coeffs[ii]*rhs; }
}
const A = AtA.map(r => [...r]), b = [...Atb];
for (let col = 0; col < n; col++) {
let maxR = col; for (let r = col+1; r < n; r++) if (Math.abs(A[r][col]) > Math.abs(A[maxR][col])) maxR = r;
[A[col], A[maxR]] = [A[maxR], A[col]]; [b[col], b[maxR]] = [b[maxR], b[col]];
if (Math.abs(A[col][col]) < 1e-15) continue;
for (let r = col+1; r < n; r++) { const f = A[r][col]/A[col][col]; for (let c = col; c < n; c++) A[r][c] -= f*A[col][c]; b[r] -= f*b[col]; }
}
const x = new Float64Array(n);
for (let ii = n-1; ii >= 0; ii--) { let s = b[ii]; for (let jj = ii+1; jj < n; jj++) s -= A[ii][jj]*x[jj]; x[ii] = s/A[ii][ii]; }
return x;
}
function* makeRows(bdryMap, a_prev) {
for (const [vi,vj] of edges) for (let comp = 0; comp < 2; comp++) {
const c = new Float64Array(2*nv); c[2*vi+comp] = sqpi; c[2*vj+comp] = -sqpi; yield {coeffs: c, rhs: 0};
}
for (const [vi, ti] of bdryMap) {
const c0 = new Float64Array(2*nv); c0[2*vi] = lam; yield {coeffs: c0, rhs: lam*Math.cos(4*ti)};
const c1 = new Float64Array(2*nv); c1[2*vi+1] = lam; yield {coeffs: c1, rhs: lam*Math.sin(4*ti)};
}
if (a_prev) for (let vi = 0; vi < nv; vi++) {
const c = new Float64Array(2*nv); c[2*vi] = lam*a_prev[vi][0]; c[2*vi+1] = lam*a_prev[vi][1]; yield {coeffs: c, rhs: lam};
}
}
function solveField(bdryMap, iters) {
const X0 = solveSystem(2*nv, makeRows(bdryMap, null));
let a = []; for (let i = 0; i < nv; i++) { const v0 = X0[2*i], v1 = X0[2*i+1]; const nm = Math.sqrt(v0**2+v1**2); a.push(nm > 1e-10 ? [v0/nm, v1/nm] : [1,0]); }
for (let it = 0; it < iters; it++) {
const X = solveSystem(2*nv, makeRows(bdryMap, a)); a = [];
for (let i = 0; i < nv; i++) { const v0 = X[2*i], v1 = X[2*i+1]; const nm = Math.sqrt(v0**2+v1**2); a.push(nm > 1e-10 ? [v0/nm, v1/nm] : [1,0]); }
}
return a;
}
function buildBdry(lo) {
const b = new Map();
for (let vi = 0; vi < nv; vi++) {
const key = String(vi), isBd = vi in defaultBdryAngle;
if (lo[key] === "free") continue;
else if (typeof lo[key] === "number") b.set(vi, lo[key]*pi/180);
else if (isBd) b.set(vi, defaultBdryAngle[vi]);
}
return b;
}
function computeEnergy(a) { let E = 0; for (const [vi,vj] of edges) E += (a[vj][0]-a[vi][0])**2 + (a[vj][1]-a[vi][1])**2; return E*pi; }
const bdry = buildBdry(locks);
let a_current = solveField(bdry, nSmooth);
const dragSmooth = Math.min(nSmooth, 1);
const W = 900, H = 480;
const svg = d3.create("svg").attr("viewBox", `0 0 ${W} ${H}`).style("width","100%").style("max-width",W+"px");
const meshL = 40, meshR = W-40, meshT = 50, meshB = H-50;
const meshW = meshR-meshL, meshH = meshB-meshT;
const scale = Math.min(meshW/baseW, meshH/baseH) * 0.85;
const offX = meshL + (meshW - baseW*scale)/2, offY = meshT + (meshH + baseH*scale)/2;
function tx(v) { return offX + v[0]*scale; }
function ty(v) { return offY - v[1]*scale; }
svg.append("text").attr("x",W/2).attr("y",20).attr("text-anchor","middle").attr("font-size",13).attr("fill","#333")
.text("Click to lock/unlock · Drag locked frames to rotate");
for (const [vi,vj] of edges) svg.append("line").attr("x1",tx(verts[vi])).attr("y1",ty(verts[vi]))
.attr("x2",tx(verts[vj])).attr("y2",ty(verts[vj])).attr("stroke","#e8e8e8").attr("stroke-width",0.8);
let minEdgeLen = Infinity;
for (const [vi,vj] of edges) { const dx = tx(verts[vi])-tx(verts[vj]), dy = ty(verts[vi])-ty(verts[vj]); minEdgeLen = Math.min(minEdgeLen, Math.sqrt(dx*dx+dy*dy)); }
const frameLen = Math.max(4, minEdgeLen * 0.40);
const frameGs = [], screenX = [], screenY = [];
for (let vi = 0; vi < nv; vi++) { screenX.push(tx(verts[vi])); screenY.push(ty(verts[vi])); }
for (let vi = 0; vi < nv; vi++) {
const px = screenX[vi], py = screenY[vi];
const th = Math.atan2(a_current[vi][1], a_current[vi][0]) / 4;
const isLocked = bdry.has(vi);
const color = isLocked ? "#d97706" : "steelblue";
const fg = svg.append("g");
for (let k = 0; k < 4; k++) { const ang = th+k*pi/2;
fg.append("line").attr("x1",px).attr("y1",py).attr("x2",px+frameLen*Math.cos(ang)).attr("y2",py-frameLen*Math.sin(ang)).attr("stroke",color).attr("stroke-width",1.5);
}
frameGs.push(fg);
if (isLocked) svg.append("circle").attr("cx",px).attr("cy",py).attr("r",frameLen+3).attr("fill","none").attr("stroke","#d97706").attr("stroke-width",1).attr("stroke-dasharray","3,2").attr("opacity",0.6);
svg.append("circle").attr("cx",px).attr("cy",py).attr("r",2).attr("fill",color);
const hitR = Math.max(frameLen+4, 12);
const hitTarget = svg.append("circle").attr("cx",px).attr("cy",py).attr("r",hitR).attr("fill","transparent").attr("cursor","pointer");
let isDragging = false;
const drag = d3.drag()
.on("start", () => { isDragging = false; })
.on("drag", (event) => {
if (!isLocked) return; isDragging = true;
const dx = event.x-px, dy = -(event.y-py);
const deg90 = ((Math.atan2(dy,dx)*180/pi) % 90 + 90) % 90;
const tl = {...locks}; tl[String(vi)] = deg90;
const tb = buildBdry(tl), as = solveField(tb, dragSmooth);
for (let v = 0; v < nv; v++) {
const vth = Math.atan2(as[v][1], as[v][0])/4, vc = tb.has(v) ? "#d97706" : "steelblue";
frameGs[v].selectAll("line").remove();
for (let k = 0; k < 4; k++) { const ang = vth+k*pi/2;
frameGs[v].append("line").attr("x1",screenX[v]).attr("y1",screenY[v])
.attr("x2",screenX[v]+frameLen*Math.cos(ang)).attr("y2",screenY[v]-frameLen*Math.sin(ang)).attr("stroke",vc).attr("stroke-width",1.5);
}
}
energyText.text(`E = ${computeEnergy(as).toFixed(4)} ${nv} vertices ${nSmooth} smoothing iter.`);
})
.on("end", (event) => {
const key = String(vi);
if (isDragging && isLocked) {
const dx = event.x-px, dy = -(event.y-py);
const deg90 = Math.round(((Math.atan2(dy,dx)*180/pi) % 90 + 90) % 90);
const nl = {...mutable rects_locks}; nl[key] = deg90; mutable rects_locks = nl;
} else if (!isDragging) {
const nl = {...mutable rects_locks}, isBd = vi in defaultBdryAngle;
if (isLocked) { if (isBd && typeof nl[key] !== "number") nl[key] = "free"; else { delete nl[key]; if (isBd) nl[key] = "free"; } }
else { const ang = Math.atan2(a_current[vi][1], a_current[vi][0])/4*180/pi; nl[key] = Math.round(((ang%90)+90)%90); }
mutable rects_locks = nl;
}
});
hitTarget.call(drag);
}
const nLocked = bdry.size;
const energyText = svg.append("text").attr("x",40).attr("y",H-10).attr("font-size",12).attr("fill","teal")
.text(`E = ${computeEnergy(a_current).toFixed(4)} ${nv} vertices (${nLocked} locked, ${nv-nLocked} free) ${nSmooth} smoothing iter.`);
return svg.node();
}| Concept | Angle approach | Functional approach |
|---|---|---|
| Frame | Angle \(\theta\) | \(a = (\cos 4\theta, \sin 4\theta)^\top\) |
| Feasibility | Always feasible | \(a^\top a = 1\) |
| Energy | \(\sum (\Delta\theta)^2\) — piecewise quadratic | \(\pi\sum\|a^j - a^i\|^2\) — quadratic |
| Boundary | \(\theta = \theta^n\) | \(a = (\cos 4\theta^n, \sin 4\theta^n)^\top\) |
| Optimization | Nonlinear (angle wrapping) | Linear system + projection |
Key takeaway
The functional approach converts a nonlinear optimization over angles into a sequence of linear least-squares problems followed by projection onto the unit circle.
Thank you!
All visualizations are interactive — drag the sliders to explore.
Built with Quarto + Observable JS + D3.js