4
$\begingroup$

"Code used in Mathematica version 13.3.1 on a Microsoft Windows (64-bit) (July 24, 2023) system"

I used @cvgmt's code to find the intersection points from a parametric plot graphically:

$Assumptions = {r \[Element] Reals, r >= 0, rh \[Element] Reals, 
   rh > 0, Q \[Element] Reals, 
   Q >= 0, \[Lambda] \[Element] Reals, \[Lambda] > 1, 
   P \[Element] Reals, P >= 0, L \[Element] Reals, L > 0};
S[rh_] := \[Pi] rh^2;
\[Phi][rh_, Q_] := -((k Q)/rh^(\[Lambda] - 1))
\[CapitalLambda] := -8 \[Pi] P
M[rh_, Q_, P_] := 
 rh/2 - (k Q^2)/(2 rh^(\[Lambda] - 1)) - (rh^3 \[CapitalLambda])/6
T[rh_, Q_, P_] := 
 1/(4 \[Pi]) ((2 M[rh, Q, P])/rh^2 + (k Q^2 \[Lambda])/
    rh^(\[Lambda] + 1) - (2 \[CapitalLambda] rh)/3)
G[rh_, Q_, P_] := M[rh, Q, P] - (T[rh, Q, P] S[rh]) // Simplify
pplt[Q_] := 
 Block[{P = 3/(8 \[Pi]), k = -1, \[Lambda] = 4}, 
  ParametricPlot[Evaluate@{T[rh, Q, P], G[rh, Q, P]}, {rh, 0.001, 10},
    PlotRange -> {{0.265, 0.280}, {0.095, 0.110}}, 
   PlotHighlighting -> None, AspectRatio -> 1/2]]
Clear[intersections];
intersections[Q_] := 
 Region`Mesh`FindSegmentIntersections[
   Cases[Normal@pplt[Q], _Line, -1], "ReturnSegmentIndex" -> True, 
   "Ignore" -> {"EndPointsTouching"}][[1, 1]]
intersections /@ Range[0.0283, 0.0331, 0.0001]

which gives me the answer, but some intersection points are not being evaluated. How do I resolve this issue?

$\endgroup$

1 Answer 1

5
$\begingroup$
  • This method depend on the plotting, it only work for the ideal cases.Here it only work for MaxRecursion -> 2 and Range[0.0283, 0.0331 - 0.0001*4, 0.0001].
Clear["Global`*"];
 $Assumptions = {r ∈ Reals, r >= 0, 
  rh ∈ Reals, rh > 0, Q ∈ Reals, 
  Q >= 0, λ ∈ Reals, λ > 1, 
  P ∈ Reals, P >= 0, L ∈ Reals, L > 0};
S[rh_] := π  rh^2;
ϕ[rh_, Q_] := -((k  Q)/rh^(λ - 1))
Λ := -8  π  P
M[rh_, Q_, P_] := 
 rh/2 - (k  Q^2)/(2  rh^(λ - 1)) - (rh^3  Λ)/6
T[rh_, Q_, P_] := 
 1/(4  π)  ((2  M[rh, Q, P])/rh^2 + (k  Q^2  λ)/
     rh^(λ + 1) - (2  Λ  rh)/3)
G[rh_, Q_, P_] := M[rh, Q, P] - (T[rh, Q, P]  S[rh]) // Simplify
pplt[Q_] := 
 Block[{P = 3/(8  π), k = -1, λ = 4}, 
  ParametricPlot[Evaluate@{T[rh, Q, P], G[rh, Q, P]}, {rh, 0.001, 10},
    PlotRange -> {{0.265, 0.280}, {0.095, 0.110}}, 
   PlotHighlighting -> None, AspectRatio -> 1/2, MaxRecursion -> 2]]
Clear[intersections];
intersections[Q_] := 
 Region`Mesh`FindSegmentIntersections[
   Cases[Normal@pplt[Q], _Line, -1], "ReturnSegmentIndex" -> True, 
   "Ignore" -> {"EndPointsTouching"}][[1, 1]]
intersections /@ Range[0.0283, 0.0331 - 0.0001*4, 0.0001]

{{0.273623, 0.100382}, {0.273536, 0.10049}, {0.273449, 0.100598}, {0.273361, 0.100707}, {0.273274, 0.100816}, {0.273186, 0.100925}, {0.273098, 0.101035}, {0.273009, 0.101145}, {0.27292, 0.101256}, {0.272831, 0.101366}, {0.272742, 0.101478}, {0.272669, 0.101571}, {0.272601, 0.101659}, {0.272534, 0.101746}, {0.272467, 0.101833}, {0.272401, 0.101919}, {0.272336, 0.102004}, {0.27227, 0.102089}, {0.272205, 0.102175}, {0.27214, 0.102259}, {0.272076, 0.102344}, {0.272011, 0.102429}, {0.271946, 0.102514}, {0.271882, 0.102599}, {0.271817, 0.102683}, {0.271753, 0.102768}, {0.271688, 0.102853}, {0.271623, 0.102938}, {0.271552, 0.103026}, {0.271475, 0.103117}, {0.271398, 0.103208}, {0.271321, 0.103299}, {0.271243, 0.10339}, {0.271165, 0.103482}, {0.271086, 0.103574}, {0.271007, 0.103666}, {0.270928, 0.103759}, {0.270849, 0.103852}, {0.27077, 0.103943}, {0.2707, 0.104027}, {0.270635, 0.104107}, {0.270572, 0.104185}, {0.27051, 0.104262}, {0.270449, 0.104338}, {0.270388, 0.104414}}

  • For the original range Range[0.0283, 0.0331, 0.0001],it is recomment to use Solve or FindInstance to solve such problem. ( We Rationalize the equation to infinite precision to overcome s1 < s2.)
f[Q_] := 
 Block[{P = 3/(8   π), k = -1, λ = 4}, 
  Evaluate@{T[rh, Q, P], G[rh, Q, P]}]; 
sol[Q_] := (f[Q] /. rh -> s1) /. 
  FindInstance[{Rationalize[(f[Q] /. rh -> s1) == (f[Q] /. rh -> s2), 
     0], 0 < s1 < s2 < 10}, {s1, s2}][[1]]
sol /@ Range[0.0283, 0.0331, 0.0001]

{{0.273525, 0.100504}, {0.273451, 0.100602}, {0.273376, 0.100699}, {0.273302, 0.100796}, {0.273228, 0.100892}, {0.273155, 0.100988}, {0.273081, 0.101084}, {0.273008, 0.101179}, {0.272934, 0.101273}, {0.272861, 0.101368}, {0.272788, 0.101461}, {0.272715, 0.101555}, {0.272642, 0.101648}, {0.27257, 0.10174}, {0.272497, 0.101833}, {0.272425, 0.101924}, {0.272353, 0.102016}, {0.272281, 0.102107}, {0.272209, 0.102197}, {0.272137, 0.102287}, {0.272066, 0.102377}, {0.271994, 0.102466}, {0.271923, 0.102555}, {0.271852, 0.102644}, {0.27178, 0.102732}, {0.271709, 0.10282}, {0.271639, 0.102907}, {0.271568, 0.102994}, {0.271497, 0.103081}, {0.271427, 0.103167}, {0.271357, 0.103253}, {0.271287, 0.103338}, {0.271217, 0.103423}, {0.271147, 0.103508}, {0.271077, 0.103592}, {0.271007, 0.103676}, {0.270938, 0.10376}, {0.270868, 0.103843}, {0.270799, 0.103926}, {0.27073, 0.104009}, {0.270661, 0.104091}, {0.270592, 0.104173}, {0.270523, 0.104254}, {0.270455, 0.104335}, {0.270386, 0.104416}, {0.270318, 0.104496}, {0.270249, 0.104576}, {0.270181, 0.104656}, {0.270113, 0.104735}}

$\endgroup$
2
  • $\begingroup$ Which procedure gives more accurate values, the first one or the FindInstance one? $\endgroup$
    – codebpr
    Jan 13 at 12:46
  • 1
    $\begingroup$ @codebpr Solve or FindInstance is more accurate than the graphics method although it is slower. ( We can replace FindInstance to Solve in the code, the result are the same.) $\endgroup$
    – cvgmt
    Jan 13 at 12:49

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.