# How to straighten a curve?

What's the easiest way to get the following animation for a given parabola where the arc length is preserved (that is, for two points in the parabola, I straighten it as follows)

My main goal is the following, given a set of random points and a set of points following a parabola, I want to "straighten" the parabola and all the surrounding points according to the straightening of that parabola. Consider the following schematic

For instance, we could consider the points as a starter

f[x_, a_, b_] := a x^2 + b;
pts1 = RandomReal[{-1, 1}, {150, 2}];
ptsf = RandomReal[{-1, 1}, 50];
pts2 = Map[{#, f[#, .5, -.5]} &, ptsf];
Graphics[{Gray, Point /@ pts1, Red, Point /@ pts2}]


Where the red dots follow a parabola. How do change pts1 based on a and b? From such a pattern, I'd expect, after the transformation, to get something like

where the dots are sparser above the line and denser below it. I tried to make it work by using preservation of the arc length somehow, but it's too slow. For more details on the math, here's a more general question.

Edit: I'm pretty satisfied with the answers provided and they do answer part of the question. As for the plane transformation, this seems to be trickier. Consider the following (very messy) image.

What I want, as suggested by @I.M., is to transform the plane such that we get two focusing and defocusing regions. The tricky part is how to "straighten" the surrounding points. If we track three line segments with end points $$\{1,2,3,4,5,6\}$$, the line segment will change size as the curve changes, and I suspect the new points will follow new parabolas (see thiner coloured lines). Not sure if I can prove this. Ultimately, what I'm asking is how to implement the following plane transformation, where the x-axis becomes a parabola

What do you think?

• Just a comment to hopefully point you in the right direction, but possibly you may start with a regularly spaced grid, then figure out how to deform such a grid to be parabolic as you desire. That is, imagine taking the grid from a graph paper and then deforming it as you describe. This would give you an indication of how to properly adjust the points as you have commented you would like to do. It seems that some of your graphical examples don’t necessarily follow this paradigm, which I suspect has encouraged the answers you have gotten so far. Jul 31, 2021 at 17:45
• @CATrevillian I'm not sure I follow, which paradigm are you referring to? The main idea is that the distance and position of the points relative to the parabolic axis is changing. How it should change I don't exactly know, but I suspect a possible way of doing that is by assuming the points follow the level curves of parabolic arclength functions, where its quadratic coefficients are, for example, proportional to the main axis coefficient deformation. If so, straightening the parabola is needed to describe this plane transformation, but I'm not sure which examples not follow this paradigm. Aug 1, 2021 at 12:30
• Sam, as I understand it, pedagogically, it may be better to start with a normal, straightened grid, and map it to the shape you would like. Then the reverse transformation should be trivial. That is, imagine a square section of space, with a regular grid, then what you want to do is essentially the opposite of taking the square space with regular grid & mapping it to the curvature you would like to start with. I don’t myself know how to do this, but I have seen similar things done here with complex mappings. Aug 1, 2021 at 17:36
• – uhoh
Aug 7, 2021 at 0:51

Here is a "homotopy style" transformation between functions $$h(x, q) = (1-q) f(x) + q g(x)$$.

(* set initial data points *)
ClearAll[x] ;
x = Subdivide[-1, 1, 14] ;

(* replace x^2 with q/2 + (1-q) x^2, i.e q = 0 -- initial function, q = 1 -- final function (line) *)
ClearAll[fun];
fun[q_][x_] := 0.5*q + (1.0 - q)*x^2 ;

(* arc length *)
ClearAll[arc];
arc[q_][x_] := 1/4 (2 x Sqrt[1 + 4 (-1 + q)^2 x^2] + ArcSinh[2 (-1 + q) x]*1/(-1 + q)) ;

(* given x and q,find new x with same arc length *)
ClearAll[fuc];
fuc[q_][x_] := Catch[
Block[
{y},
If[q == 1.0, Throw[{arc[0.0][x], fun[q][x]}]] ;
y = y /. FindRoot[arc[0.0][x] == arc[q][y], {y, x}] ;
{y, fun[q][y]}
]
] ;

(* test arc conservation *)
{arc[0.0][1.0], arc[0.25][First[fuc[0.25][1.0]]]}

(* plot *)
Manipulate[
ListLinePlot[
Map[fuc[q], x],
Mesh -> All, AspectRatio -> 1, PlotRange -> {{-2.5, 2.5}, {-2.5, 2.5}}, Background -> Black, PlotStyle -> Cyan, Axes -> False,
PlotLabel -> Style[StringTemplate["q = "][q], White]
],
{q, 0.0, 1.0, 0.01}
]
(* {1.4789428575445975,1.478942857544597} *)


• Brilliant answer! This allows me to deform a curve. I can then extend this idea for all my points pts1 by simply shifting the parabola's height accordingly. However, I cannot reproduce that "coning" effect as sketched by my drawing. That is, points above and below the red parabola get equally sparser or tighter, as opposed to getting a gradual sparsening or tightening depending on how far (and which side) you are from the parabola. There must be a way of getting this effect by tracking the red lines, where points get tighter towards the middle line (y-axis). Any ideas? Jul 29, 2021 at 14:45
• @samwolfe, perhaps the following might be close, for each point find a normal vector to a parabola, track normal during transformation, a point should stay on its normal with preserved distance. This seems to produce a focusing /defocusing effect.
– I.M.
Jul 29, 2021 at 15:17
• If by preserved distance you mean the distance of the point to the parabola (and not the arclength), I suspect this can't happen. The distance must be somehow scaled progressively as the point gets closer to the y-axis. Otherwise the points will cross this axis. Focusing/defocusing is exactly what I mean, but it's happening at different speeds in the plane. See the recent edit for a clarification on this. Hope it's not too confusing. Jul 29, 2021 at 16:04

The arclength is given by this and is the big expression containing ArcSinh later on:

Integrate[Sqrt[1 + (2 a x)^2], {x, -k, k},
Assumptions -> a > 0 && k > 0]


For a given a and len (arc length) find the best k (integration range) that gives that length by doing an NMinimize (for whatever reason, NSolve didn't work that well so I minimized a square error). The error in the NMinimize is extremely small.

Manipulate[
{err, sol} =
NMinimize[(k Sqrt[1 + 4 a^2 k^2] + ArcSinh[2 a k]/(2 a) - len)^2, k];
k1 = k /. sol;
Plot[a t^2, {t, -k1, k1}, Background -> Black,
PlotRange -> {{-1, 1}, {-0.2, 1}}, AspectRatio -> 1, Axes -> False]
, {a, 0.001, 4}, {len, .5, 3}]

f[t_] = t^2;
{l, u} = {-1, 1};


To transition to a line, curves that preserve the arc length and scales the signed curvature by q are used. First such curves are made for a variety of q values:

sols = <|# -> First[{x, y} /. NDSolve[##2]] & @@@ Table[Evaluate[{1 - q,
{D[x'[t]^2 + y'[t]^2 == Total[D[{t, f[t]}, t]^2], t],
(x'[t] y''[t] - y'[t] x''[t])/(x'[t]^2 + y'[t]^2)^(3/2) == (q f''[t])/(1 + f'[t]^2)^(3/2),
x[l] == l, y[l] == f[l], x'[l] == 1, y'[l] == f'[l]}, {x, y}, {t, l, u}}], {q, Subdivide[25]}]|>;


The curves are heuristically adjusted to be similar to the initial curve:

rigids = Last[FindGeometricTransform[#, #2, TransformationClass -> "Rigid"] & @@
Transpose[Table[{{t, f[t]}, Through[#[t]]}, {t, Subdivide[l, u, 15]}]]] & /@ sols;


Linear interpolation between the solved q values is used to allow an arbitrary q:

getKeys = Nearest[Keys[sols]];
straightening[q_, t_] := Interpolation[{#, rigids[#][Through[sols[#][t]]]} & /@ getKeys[q, 2], q, InterpolationOrder -> 1]

Manipulate[ParametricPlot[straightening[q, t], {t, -1, 1},
PlotRange -> {{-2, 2}, {-2, 2}}, AspectRatio -> Automatic], {q, 0, 1}]


Clear[a, x, y, x1, x2]
SeedRandom[2]
npts = 200;
r = 3;
pts1 = RandomReal[{-1.5, 1.5}, {npts, 2}];
parms = {a -> 1/2, b -> -1/2};
f[x_, y_] := y - a x^2 - b /. parms
px[x_] := 1/2 x Sqrt[1 + 4 a^2 x^2] + ArcSinh[2 a x]/(4 a) /. parms
vecn = Grad[f[x, y], {x, y}];
GRS = {};
INTS = {};

For[k = 1, k <= npts, k++,
{x0, y0} = pts1[[k]];
equs = Join[Thread[{x, y} == {x0, y0} + lambda vecn], {f[x, y] == 0}] /. parms;
sols = Quiet@Solve[equs, {x, y, lambda}, Reals];
{x1, y1, l1} = {x, y, lambda} /. sols[[1]];
AppendTo[INTS, {px[x1], Sign[f[x0, y0]] Sqrt[(x1 - x0)^2 + (y1 - y0)^2]}];
vn1 = vecn /. sols[[1]] /. parms;
AppendTo[GRS, ParametricPlot[Thread[{x0, y0} + l vn1], {l, 0, l1}, PlotStyle -> {Thin, Black}]]
];

gr1 = ContourPlot[(f[x, y] /. parms) == 0, {x, -r, r}, {y, -r, r}, ContourStyle -> Red];
gr2 = Graphics[{Black, PointSize[0.01], Point /@ pts1}];
xk = Transpose[INTS][[1]];
yk = Transpose[INTS][[2]] + b /. parms;
INTS0 = Join[{xk}, {yk}] // Transpose;
gr3 = ListPlot[INTS0, PlotStyle -> {Red, PointSize[0.01]}, Filling -> (b /. parms)];
gr4 = Plot[(b /. parms), {x, -r, r}, PlotStyle -> Red];
Show[gr1, gr2, GRS, PlotRange -> {{-r, r}, {-r, r}}, AspectRatio -> 1]
Show[gr4, gr3, PlotRange -> {{-r, r}, {-r, r}}, AspectRatio -> 1]
Show[gr1, gr2, gr3, gr4, GRS, PlotRange -> {{-r, r}, {-r, r}}, AspectRatio -> 1]


EDIT

Included a Manipulateversion. Run at low speed.

Clear["Global*"]
SeedRandom[5]
npts = 10; r = 3; a = 2; b = -1/2;
pts1 = RandomReal[{-2, 2}, {npts, 2}];
L = (x0 - x)^2 + (y0 - y)^2 + lambda (y - a x^2 - b);
sollambda1 = Solve[grad[[1]] == 0, lambda][[1]];
sollambda2 = Solve[grad[[2]] == 0, lambda][[1]];
equ = (lambda /. sollambda1) == (lambda /. sollambda2);
solsdist = Solve[{equ, y - a x^2 - b == 0}, {x, y}];
f[x_, y_, a_] := y - a x^2 - b
px[x_, a_] := 1/2 x Sqrt[1 + 4 a^2 x^2] + ArcSinh[2 a x]/(4 a)
vecn[x_, a_] := {-((2 a x)/Sqrt[1 + 4 (a x)^2]), 1/Sqrt[1 + 4 (a x)^2]}

nearpt[x0_, y0_] := Module[{sols, ptsp, dists, p0, pnear, near}, sols = solsdist // N // Chop;
ptsp = {x, y} /. sols;
p0 = {x0, y0};
dists = Table[{Norm[ptsp[[k]] - p0], k}, {k, 1, 3}];
pnear = Sort[dists][[1]];
near = ptsp[[pnear[[2]]]];
{x1, y1} = near;
Return[{x1, y1, Norm[near - p0] Sign[f[p0[[1]], p0[[2]], a]]}]
]

Manipulate[
GRS1 = {};
INTS1 = {};
For[k = 1, k <= npts, k++,
{x0, y0} = pts1[[k]];
{x1, y1, l1} = nearpt[x0, y0];
dist = px[x1, a];
x2 = x /. Quiet@FindRoot[px[x, a1] - dist == 0, {x, x1}];
y2 = y /. Quiet@Solve[f[x2, y, a1] == 0, y][[1]];
vn2 = vecn[x2, a1];
AppendTo[GRS1, ParametricPlot[{x2, y2} - l vn2, {l, 0, l1}, PlotStyle -> {Thin, Black}]];
AppendTo[INTS1, {x2, y2} - l1 vn2]
];
gr11 = ContourPlot[f[x, y, a1] == 0, {x, -r, r}, {y, -r, r}, PlotPoints -> 30, ContourStyle -> Red];
Show[gr11, ListPlot[INTS1, PlotStyle -> {Black, PointSize[0.02]}], GRS1, PlotRange -> {{-r, r}, {-r, r}}, AspectRatio -> 1],
{{a1, 0.01}, 0.01, a, 0.01}
]
`
• Thanks. This is almost what I want, but I need the distance of the points to the curve to be scaling accordingly. Here, you'd get the x-values of some points changing order and you wouldn't get the focusing/defocusing, if you see what I mean. Jul 30, 2021 at 14:05