3
$\begingroup$

For the toroidal helix defined by

$Assumptions = {r > 0, R > 0, n \[Element] Integers, n > 0, t >= 0, t <= 2*Pi}
x[t_] := (R + r Cos[n t])*Cos[t]
y[t_] := (R + r Cos[n t])*Sin[t]
z[t_] := r Sin[n t]

the differential arc-length is:

l[t_] := FullSimplify[Sqrt[D[x[t], t]^2 + D[y[t], t]^2 + D[z[t], t]^2]]

and the arc-length becomes:

L = Integrate[l[t], t]

Unfortunately, when calculating the length over 1 period:

L /. t -> 0
FullSimplify[L /. t -> 2 Pi]

We get the output 0 for both evaluations. This is because Mathematica introduces some kind of modulo equivalence when integrating. This can be seen by plotting an example arc-length curve:

rule = {R -> 6, r -> 2, n -> 5}
Plot[L /. rule, {t, 0, 2 Pi}]

enter image description here

Still, by noticing that the integrand is $2\pi/ n$ periodic, we can try to extract the solution by finding the difference between, say:

FullSimplify[Limit[L, t -> Pi/n, Direction -> "FromBelow"]]
FullSimplify[Limit[L, t -> Pi/n, Direction -> "FromAbove"]]

and multiplying it by $n$. This is taking too long for my computer and I can't seem to make any progress. Is there a way to solve this problem in the form presented? Would analytical work help to get rid of this modulo equivalence (for example the Weierstrass substitution and further manipulations)?

$\endgroup$
6
  • $\begingroup$ I can't get any results from your Integrate expression in a reasonable time. How long did you have to wait? What does the result look like? Finally, the numerical calculation seems to work fine: Plot[NIntegrate[l[t] /. rule, {t, 0, tmax}], {tmax, 0, 2 Pi}] giving a positive, almost linearly increasing value throughout. $\endgroup$
    – MarcoB
    Commented Jul 12, 2022 at 7:24
  • $\begingroup$ x,y,z are periodic functions of t. Therefore, your definition for the arc length will also be periodic. You have to introduce aperiodicity "by hand" $\endgroup$ Commented Jul 12, 2022 at 8:27
  • $\begingroup$ @MarcoB I waited 30 seconds. Checked with Timing[] $\endgroup$ Commented Jul 12, 2022 at 10:29
  • $\begingroup$ I get the integral in terms of a complex expression containing arcsinh[g(t)] and complex square roots which because they're multivalued, Matheamtica cannot (continuously) integrate across the branch cuts producing the jump-discontinuity exhibited in your plot. I suspect if we replaced these with analytically-continuous versions, we could then use the resulting antiderivative to evaluate the definite integral. Might be an interesting project. Would be tough though. $\endgroup$
    – josh
    Commented Jul 12, 2022 at 13:08
  • $\begingroup$ @josh This makes sense. I will try to manually manipulate the integral to hopefully get some more insight. In the meanwhile, I still do not know why Mathematica is not able to evaluate the limits at the end of the question. $\endgroup$ Commented Jul 12, 2022 at 14:59

3 Answers 3

7
$\begingroup$

Doing integration with the rule-based integrator 'Rubi' https://rulebasedintegration.org/, you get an analytical solution.

'Rubi' gives an antiderivative with ArcTan[ ... Tan[...]], where it is easier to find a continuation over the branch cuts.

Doing it here for given paramters.

$Assumptions = {r > 0, R > 0, n \[Element] Integers, n > 0, t >= 0};
x[t_] = (R + r Cos[n t])*Cos[t];
y[t_] = (R + r Cos[n t])*Sin[t];
z[t_] = r Sin[n t];

rule0 = {R -> 6, r -> 2, n -> 5};

igd0 = FullSimplify[Sqrt[D[x[t], t]^2 + D[y[t], t]^2 + D[z[t], t]^2]]

igd = igd0 /. rule0

(*   Sqrt[276 + 48 Cos[5 t] + 4 Cos[10 t]]/Sqrt[2]   *)

(* rint[t_] = Int[igd, t]  Rubi integration: Get   *)

rint[t_] = 
  6/5 ArcTan[(2 Tan[(5 t)/2])/Sqrt[
     41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4]] + 
   2/5 Cos[(5 t)/2] Sin[(5 t)/2] Sqrt[
    41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4] - (
   2 Sqrt[1189] Tan[(5 t)/2] Sqrt[
    41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4])/(
   5 (41 + Sqrt[1189] Tan[(5 t)/2]^2)) + (2 1189^(1/4)
       EllipticE[2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], (
       1189 - 33 Sqrt[1189])/
       2378] (41 + Sqrt[1189] Tan[(5 t)/2]^2) Sqrt[(
      41 + 66 Tan[(5 t)/2]^2 + 
       29 Tan[(5 t)/2]^4)/(41 + 
        Sqrt[1189] Tan[(5 t)/2]^2)^2])/(5 Sqrt[
      41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4]) - (12 1189^(1/4)
       EllipticF[2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], (
       1189 - 33 Sqrt[1189])/
       2378] (41 + Sqrt[1189] Tan[(5 t)/2]^2) Sqrt[(
      41 + 66 Tan[(5 t)/2]^2 + 
       29 Tan[(5 t)/2]^4)/(41 + 
        Sqrt[1189] Tan[(5 t)/2]^2)^2])/(5 (41 - Sqrt[1189]) Sqrt[
      41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4]) - (29^(
      1/4) (41 - Sqrt[1189]) EllipticF[
       2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], (1189 - 33 Sqrt[1189])/
       2378] (41 + Sqrt[1189] Tan[(5 t)/2]^2) Sqrt[(
      41 + 66 Tan[(5 t)/2]^2 + 
       29 Tan[(5 t)/2]^4)/(41 + Sqrt[1189] Tan[(5 t)/2]^2)^2])/(5 41^(
      3/4) Sqrt[41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4]) + (6 41^(
      1/4) (29 + Sqrt[1189]) EllipticPi[(1189 - 35 Sqrt[1189])/2378, 
       2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], (1189 - 33 Sqrt[1189])/
       2378] (41 + Sqrt[1189] Tan[(5 t)/2]^2) Sqrt[(
      41 + 66 Tan[(5 t)/2]^2 + 
       29 Tan[(5 t)/2]^4)/(41 + Sqrt[1189] Tan[(5 t)/2]^2)^2])/(5 29^(
      3/4) (41 - Sqrt[1189]) Sqrt[
      41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4]);

Plot[rint[t], {t, 0, 5}]

srint[t_] = 
 rint[t] // PowerExpand[#, Assumptions -> Element[t, Reals]] & // 
  Simplify[#, Assumptions -> Element[t, Reals]] &

(-1)^Floor[...] terms are ==1 except some singular points, which is not relevant. Set them to 1 later.

Plot[(-1)^
 Floor[(\[Pi] + 2 Arg[41 + Sqrt[1189] Tan[(5 t)/2]^2] - 
   Arg[41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4])/(2 \[Pi])], {t, 0,
   5}]

Reduce[(-1)^
   Floor[(\[Pi] + 2 Arg[41 + Sqrt[1189] Tan[(5 t)/2]^2] - 
     Arg[41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4])/(2 \[Pi])] == 
   1 && 0 < t < 5, t]

(*   0 < t < \[Pi]/5 || \[Pi]/5 < t < (3 \[Pi])/5 || (3 \[Pi])/5 < 
  t < \[Pi] || \[Pi] < t < (7 \[Pi])/5 || (7 \[Pi])/5 < t < 5   *)

Continuation rule for ArcTan[aa_ Tan[bb_]] terms is

ruleArcTan = 
 ArcTan[aa_ Tan[bb_]] -> 
  ArcTan[aa Tan[bb]] + bb - 2 ArcTan[Cot[bb] (-1 + Sqrt[Sec[bb]^2])]

Plot[Evaluate[
  bb - 2 ArcTan[Cot[bb] (-1 + Sqrt[Sec[bb]^2])] /. bb -> 5/2 t], {t, 
  0, 5}]

Plot[ArcTan[2 Tan[5/2 t]], {t, 0, 5}]

Plot[Evaluate[ArcTan[2 Tan[5/2 t]] /. ruleArcTan], {t, 0, 5}]

First term of srint[t] has to be excluded form ArcTan transformation, would give wrong result.

(sr3 = List @@ (srint[t] /. (-1)^
         Floor[(\[Pi] + 2 Arg[41 + Sqrt[1189] Tan[(5 t)/2]^2] - 
           Arg[41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4])/(
          2 \[Pi])] -> 1 // Expand)) // TableForm;

Plot[#, {t, 0, 5}] & /@ sr3

Get analytical integral int[t] and compare with numerical integration

nint[tmax_?NumericQ] := NIntegrate[igd, {t, 0, tmax}]

int[t_] = 
 sr3[[1]] + (Total@sr3[[2 ;; Length[sr3]]] /. ruleArcTan) // 
  Simplify[#, Assumptions -> Element[t, Reals]] &

(*   1/290 (348 ArcTan[(2 Tan[(5 t)/2])/Sqrt[
     41 + 66 Tan[(5 t)/2]^2 + 29 Tan[(5 t)/2]^4]] + 
   116 1189^(1/4)
     EllipticE[
     5 t - 4 ArcTan[Cot[(5 t)/2] (-1 + Sqrt[Sec[(5 t)/2]^2])] + 
      2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], 
     1/2 - 33/(2 Sqrt[1189])] - 
   116 1189^(1/4)
     EllipticF[
     5 t - 4 ArcTan[Cot[(5 t)/2] (-1 + Sqrt[Sec[(5 t)/2]^2])] + 
      2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], 
     1/2 - 33/(2 Sqrt[1189])] - (1/(-41 + Sqrt[1189]))
   348 1189^(1/4)
     EllipticPi[1/2 - 35/(2 Sqrt[1189]), 
     5 t - 4 ArcTan[Cot[(5 t)/2] (-1 + Sqrt[Sec[(5 t)/2]^2])] + 
      2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], 
     1/2 - 33/(2 Sqrt[1189])] - (1/(-41 + Sqrt[1189]))
   12 1189^(3/4)
     EllipticPi[1/2 - 35/(2 Sqrt[1189]), 
     5 t - 4 ArcTan[Cot[(5 t)/2] (-1 + Sqrt[Sec[(5 t)/2]^2])] + 
      2 ArcTan[(29/41)^(1/4) Tan[(5 t)/2]], 
     1/2 - 33/(2 Sqrt[1189])] + 
   29 Sqrt[2] Sqrt[(69 + 12 Cos[5 t] + Cos[10 t]) Sec[(5 t)/2]^4]
     Sin[5 t] - (
   58 Sqrt[2378] Sqrt[(69 + 12 Cos[5 t] + Cos[10 t]) Sec[(5 t)/2]^4]
     Tan[(5 t)/2])/(41 + Sqrt[1189] Tan[(5 t)/2]^2))   *)

Limit[int[t], t -> 0, Direction -> -1]   (*   0   *)

Plot[{nint[t], int[t]}, {t, 0, 5}, 
 PlotStyle -> {{Opacity[.3], Thickness[.02], Red}, Black}]

enter image description here

Plot[nint[t] - int[t], {t, 0, 5}, PlotStyle -> {Black}, 
 PlotRange -> 10^-8, PlotPoints -> 200]
$\endgroup$
2
  • $\begingroup$ This is a masterpiece. However, NIntegrate is simpler and more convenient. $\endgroup$
    – user64494
    Commented Jul 13, 2022 at 11:45
  • $\begingroup$ @Akku14: Outstanding! Beautiful example extending the Fundamendal Theorem of Calculus to multi-valued antiderivatives via analytic continuation. The concept should become a standard addendum to this important theorem. $\endgroup$
    – josh
    Commented Jul 13, 2022 at 12:33
3
$\begingroup$

Since Mathematica will ignore some constant for everywhere when we integral a function, so Newton-Leibniz formula not always work.

Clear[n, R, r, x, y, z];
n = 5;
R = 2;
r = 1;
x[t_] = (R + r Cos[n t])*Cos[t];
y[t_] = (R + r Cos[n t])*Sin[t];
z[t_] = r Sin[n t];
ArcLength[{x[t], y[t], z[t]}, {t, N[0], N[2 π]}]
plot = ParametricPlot3D[{x[t], y[t], z[t]}, {t, 0, 2 π}];
DiscretizeGraphics[plot] // ArcLength
NIntegrate[D[{x[t], y[t], z[t]}, t] // Norm, {t, 0, 2 π}]

34.0869 34.0784 34.0869 enter image description here

$\endgroup$
4
  • $\begingroup$ It looks like the OP needs to confirm if an analytic form for the arclength function is really necessary, or if a numerical evaluation suffices. $\endgroup$ Commented Jul 12, 2022 at 10:27
  • $\begingroup$ An exact formula 2 *Pi*n *Sqrt[R^2/n^2 + r^2] is claimed here. Unfortunately, the one is in discordance with numerical results. $\endgroup$
    – user64494
    Commented Jul 12, 2022 at 11:38
  • $\begingroup$ @J. M. a numerical evaluation is unfortunately not sufficient. Explicit results would be of importance in the field of plasma confinement. $\endgroup$ Commented Jul 12, 2022 at 14:52
  • $\begingroup$ @user64494 Thank you, I have seen this and was puzzled. It must be an approximation. $\endgroup$ Commented Jul 12, 2022 at 14:54
2
$\begingroup$

This is the result obtained from Mathematica with general parameters following tedious simplification and taking care of the ArcTan jumps with the Floor function (still quite messy...):

r = 2; R = 6; n = 5; Plot[{NIntegrate[
Sqrt[r^2 + 2*n^2*r^2 + 2*R^2 + 4*r*R*Cos[n*tt] + r^2*Cos[2*n*tt]]/
Sqrt[2], {tt, 0, t}], 
  (1/
 n)*((Sqrt[I*r + n*r + I*R]*Sqrt[I*r + n*r + I*R]*
    Sqrt[(-I)*r + n*r + I*R]*
    EllipticF[
     ArcTan[(Sqrt[(-I)*r + n*r + I*R]*Tan[(n*t)/2])/
        Sqrt[I*r + n*r + I*R]] + 
      Pi*(1 + Floor[(n*t)/(2*Pi) - 1/2]), 
              -((4*I*n*r^2)/((-I + n)^2*r^2 + R^2))])/
  Sqrt[(-I)*r + n*r - I*R] - (2*I*
    R*(Sqrt[I*r + n*r + I*R]*Sqrt[I*r + n*r + I*R])*
    EllipticPi[(2*r)/((1 + I*n)*r - R), 
              
     ArcTan[(Sqrt[(-I)*r + n*r + I*R]*Tan[(n*t)/2])/
        Sqrt[I*r + n*r + I*R]] + 
      Pi*(1 + Floor[(n*t)/(2*Pi) - 1/2]), -((4*I*n*
          r^2)/((-I + n)^2*r^2 + R^2))])/(Sqrt[(-I)*r + n*r + 
      I*R]*Sqrt[(-I)*r + n*r - I*R]) + 
       
 Sqrt[(-I)*r + n*r - I*R]*Sqrt[(-I)*r + n*r + I*R]*
  EllipticE[
   ArcTan[(Sqrt[(-I)*r + n*r + I*R]*Tan[(n*t)/2])/
      Sqrt[I*r + n*r + I*R]] + Pi*(1 + Floor[(n*t)/(2*Pi) - 1/2]), 
           -((4*I*n*r^2)/((-I + n)^2*r^2 + R^2))] + (Sqrt[
      r^2 + 2*n^2*r^2 + 2*R^2 + 4*r*R*Cos[n*t] + r^2*Cos[2*n*t]]/
     Sqrt[2] - (Sqrt[
        n*r - I*R - I*r*Cos[n*t]]*((-I)*r + n*r + I*R))/
     Sqrt[n*r + I*R + I*r*Cos[n*t]])*
         Tan[(n*t)/2])}, {t, 0, 5}, PlotPoints -> 300, 
  PlotStyle -> {Blue, Dashed}]

enter image description here

Edit: Compare this with the arc length of a linear helix (green dots)

Sqrt[1 + n^2*r^2]*t

for n=5: enter image description here

and n=15: enter image description here

The simple linear curve approaches the complex elliptic expression for higher n.

$\endgroup$

Your Answer

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

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