1
$\begingroup$

I begin by defining an interpolation of some data I have stored at this link for brevity. The interpolating function is a function of two variables, x and y

S = (*data*)
f = Interpolation[S]
Plot3D[f[x, y], {x, -5, 4}, {y, -5, 4}]

I then want to define a one-dimensional surface by when this interpolated function intersects zero

R = ImplicitRegion[f[x, y] == 0.0]
RegionPlot[R, PlotRange -> {{-5, 4}, {-5, 4}}

This looks like it has worked quite well, but now when I try to calculate ArcLength[R], I don't get an answer, even when I try restricting the domain of the implicit region

R = ImplicitRegion[f[x, y] == 0.0, {{x, -5, 4}, {y, -5, 4}}]

The ellipse looks pretty well-defined, is there a way of getting its arc Length?

$\endgroup$
1
  • 2
    $\begingroup$ Here's one way folks can import the data: S = Import["https://gist.github.com/jhwilson/10467024a973e3b8cf0acf4990a1a21d", "Data"][[1, 2]] // ToExpression $\endgroup$
    – Michael E2
    Feb 11, 2021 at 16:22

2 Answers 2

2
$\begingroup$

I don't know how accurately you want the arc length, but going along the lines you started, one can do this:

R = ImplicitRegion[f[x, y] == 0.0, {{x, -5, 4}, {y, -5, 4}}];
DiscretizeRegion[R, MaxCellMeasure -> 0.01] // ArcLength

(*  19.7927  *)

A more accurate determination:

Module[{tf},
 NDSolveValue[{{x'[t], y'[t]} == 
    Cross@D[f[x[t], y[t]], {{x[t], y[t]}}],
   x[0] == (x /. FindRoot[f[x, 0], {x, 2}]), y[0] == 0,
   s'[t] == Sqrt[x'[t]^2 + y'[t]^2], s[0] == 0,
   WhenEvent[y[t] == 0 && x[t] > 0, {tf = t; "StopIntegration"}]},
  s[tf], {t, 0, 10000}]
 ]

(*  19.7955  *)

It takes a bit of work to get there with DiscretizeRegion:

DiscretizeRegion[R, MaxCellMeasure -> 0.0001] // ArcLength

(*  19.7954  *)
$\endgroup$
2
$\begingroup$

This is a bit of a hack. Assume we have your plot: pl from the data:

R = ImplicitRegion[f[x, y] == 0.0, {x, y}]
pl = RegionPlot[R, PlotRange -> {{-5, 4}, {-5, 4}}]

We may no dissect pl to get the relevant information. Unfortunately, the line is not stored directly, but first the coordinates of the points are stored and then the line is specified by the indices of the coordinates. To get the coordinates of the line and to check if it is o.k. we write:

coor = g[[1, 1, 1]];
ind = g[[1, 1, 2, 1, 1, 3, 1]];
line = coor[[ind]];
ListLinePlot[line]

enter image description here

Now the length of the line is easily calculated by taking the successive differences and the norm therefrom and finally adding up:

Total[Norm /@ Differences[line]]
(* 19.7268 *)
$\endgroup$
8
  • 1
    $\begingroup$ Slightly less hacky: Total@Cases[Normal@pl, l_Line :> ArcLength[l], Infinity]. $\endgroup$
    – Michael E2
    Feb 11, 2021 at 18:17
  • $\begingroup$ Michael E2 Hi, Line gives indices, not coordinates. Reason: GraphicsComplex $\endgroup$ Feb 11, 2021 at 19:30
  • $\begingroup$ Try my code. I got 19.7268 when I ran it. $\endgroup$
    – Michael E2
    Feb 11, 2021 at 19:32
  • $\begingroup$ Hmmm, I get 306.008. Look at pl[[1, 1, 2]] then you see that Line[..] is a list of indices. Only to make sure, pl is the output of Plot3D. $\endgroup$ Feb 11, 2021 at 19:47
  • 1
    $\begingroup$ It's the last point under "Details and Options" of the doc page for GraphicsComplex. $\endgroup$
    – Michael E2
    Feb 11, 2021 at 21:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.