Unable to find root for user-defined function

Question

I created the function getSolution[xVal,yVal] for a geometry problem I'm working on and wish to solve it for getSolution[2,x]==0. I can plot it clearly showing the roots:

However, when I attempt to solve

 FindRoot[getSolution[2, x] == 0, {x, 1.5}]

I receive the errors:

Part::partd: Part specification x$18835[[2]] is longer than depth of object.

Part::partd: Part specification x$18835[[1]] is longer than depth of object.

Which I've tracked down I believe to FindRoot attempting to obtain the derivative of the function for the Newton method. Since if I code:

D[getSolution[2, x], x]

I receive the same message. I was wondering if someone could offer some suggestions to fix this code? Here's my code:

getAngle[a_, b_, c_] := Module[{line1, line2, line3, aVal, bVal, val},
   aVal = a - b;
   bVal = c - b;
   val = (aVal.bVal)/(EuclideanDistance[{0, 0}, 
        aVal] EuclideanDistance[{0, 0}, bVal]);
   ArcCos[val]
   ];

makeEqn[p1_, p2_] := 
  Function[{x}, (p2[[2]] - p1[[2]])/(p2[[1]] - p1[[1]]) (x - 
       p1[[1]]) + p1[[2]]];

getLength[a_, b_] := EuclideanDistance[a, b];

getSolution[cx_, cy_] := 
 Module[{c, m, acLength, macAngle, d, ebLength, bmLength, emLength, 
   acEqn, acAngle, e, deEqn, cbEqn, f, x, a, b},
  a = {0, 0};
  b = {7, 0};
  c = {cx, cy};
  m = {getLength[a, b]/2, 0};
  acLength = getLength[a, c];
  macAngle = getAngle[m, a, c];
  d = {acLength Cos[macAngle], 0};

  bmLength = getLength[b, m];
  acEqn = makeEqn[a, c];

  e = {x, acEqn[x]} /. 
     NSolve[acEqn[x]^2 + (x - m[[1]])^2 == bmLength^2 && x > 0, x] // 
    First;
  deEqn = makeEqn[d, e];
  cbEqn = makeEqn[c, b];
  f = {x, deEqn[x]} /. NSolve[deEqn[x] == cbEqn[x], x] // First;
  getLength[e, b] - getLength[f, b] // N
  ]

Not sure this is needed but here's the problem: Essentially need to find the point C in the diagram such that the triangle FEB is isosceles. When getSolution[x,y]=0, it is. Also, don't wish to just generate a table of values and find it manually for one case as above as I wish to solve the larger problem getSolution[x,y]=0 for all x and y in applicable domain.

Answer 1

In fact, for at least your specific case, NSolve[] is not really necessary, as one can explicitly solve the underlying linear and quadratic equations:

getSolution[cx_, cy_] := Module[{a, b, c, m, acLength, d, bmLength, e, f, x},
   a = {0, 0}; b = {7, 0}; c = {cx, cy}; m = {EuclideanDistance[a, b]/2, 0};
   acLength = EuclideanDistance[a, c];
   d = {acLength Normalize[m - a].Normalize[c - a], 0};
   bmLength = EuclideanDistance[b, m];
   e = (c (cx m[[1]] + Sqrt[(bmLength Norm[c])^2 - cy^2 m[[1]]^2]))/Norm[c]^2;
   f = LinearSolve[{Cross[e - d], Cross[b - c]}, {(d - e).Cross[d], (c - b).Cross[c]}];
   EuclideanDistance[e, b] - EuclideanDistance[f, b]]

ContourPlot[Evaluate[getSolution[cx, cy] == 0 /. Abs -> Identity],
            {cx, 0, 3}, {cy, -4, 4}, FrameLabel -> Automatic, PlotPoints -> 75]

---

One can in fact do some further analysis to separate out the two branches of the solution. Skipping the details,

ContourPlot[{-7 cx + cx^2 + cy^2 == 0,
             -343 cx + 147 cx^2 - 21 cx^3 + cx^4 + 147 cy^2 - 21 cx cy^2 - cy^4 == 0},
            {cx, 0, 3}, {cy, -4, 4}, FrameLabel -> Automatic, PlotPoints -> 75]

Answer 2

Not an answer, but perhaps a small step further on:

You try to solve the equation getSolution[cx, cy]==0 .

First I changed the definition of getSolution[cx, cy] (thanks to @John Doty 's comment)

getSolution[cx_?NumericQ, cy_?NumericQ] := 
Block[{c, m, acLength, macAngle, d, ebLength, bmLength, emLength, acEqn, acAngle, e, deEqn, cbEqn, f, x, a, b}, a = {0, 0};b = {7, 0};c = {cx, cy};m = {getLength[a, b]/2, 0};acLength = getLength[a, c];macAngle = getAngle[m, a, c];d = {acLength Cos[macAngle], 0};bmLength = getLength[b, m];acEqn = makeEqn[a, c];e = {x, acEqn[x]} /. NSolve[acEqn[x]^2 + (x - m[[1]])^2 ==bmLength^2 && x > 0, x] //First;deEqn = makeEqn[d, e];cbEqn = makeEqn[c,b];f = {x, deEqn[x]} /. NSolve[deEqn[x] == cbEqn[x], x] // First;
getLength[e, b] - getLength[f, b]]

If I now try to apply Contourplot to visualize the solution

ContourPlot[getSolution[cx, cy] == 0 , {cx, 0, 3}, {cy, -4, 4}, FrameLabel -> Automatic, Evaluated -> True]

Mathematica v12 gives an error message Part::partd: Part specification x$2040162[[2]] is longer than depth of object. which indicates that's something wrong with the function definition.

Surprisingly

ContourPlot[getSolution[cx, cy] , {cx, 0, 3}, {cy, -4, 4},Contours -> {0}, FrameLabel -> Automatic]

evaluates and shows

four possible solutions cy for given cx.

Answer 3

Thanks guys. The ?NumericQ did it. Here's the plot in 3D with the yellow line representing the points (x,y) which produce the isosceles triangle:

mycontour[t1_] := t2 /. FindRoot[getSolution[t1, t2], {t2, 1.4}];
pp1 = ParametricPlot3D[{t1, mycontour[t1], 0}, {t1, 1, 3}, 
   PlotStyle -> Yellow];
pp2 = Plot3D[getSolution[t, t2], {t, 1, 3}, {t2, 1, 3}, 
   BoxRatios -> {1, 1, 1}, AxesLabel -> {"x", "y", "z"}];
Show[{pp2, pp1}]

And here's a (numerically) exact solution in a Manipulate environement:

Answer 4

If you have triangles inscribed in circles, there is a good chance that synthetic geometry and FindGeometricConjectures will find some simplifying assumptions given your constraints.

The following code suggest that if FEB is isosceles, then triangle CEB is a right triangle.

gsc[cx_, cy_] := Module[{aa, bb, ri},
  aa = {0, 0};
  bb = {7, 0};
  ri = RandomInstance[
    GeometricScene[
     {a -> aa, b -> bb, c -> {cx, cy}, d, e, f, 
      r}, {p1 == Style[Polygon[{b, f, e}], Cyan], 
      EuclideanDistance[b, f] == EuclideanDistance[b, e], 
      EuclideanDistance[r, f] < EuclideanDistance[r, c], 
      GeometricAssertion[{c, f, r, b}, "Collinear"], 
      GeometricAssertion[{d, f, e}, "Collinear"], 
      GeometricAssertion[{a, c, e}, "Collinear"], 
      CircleThrough[{b, e, c, d}, r]}
     ], RandomSeeding -> 4
    ];
  {ri, FindGeometricConjectures[ri]}
  ]
{ri, fc} = gsc[1.75, 1.5]
ri["Points"]

Although it did not fall out from the conjectures, another useful relation is the Intersecting Chords Theorem that can be verified using the points:

{Norm[c - f] Norm[f - b], Norm[d - f] Norm[f - e]} /. 
 ri["Points"]
(* {4.12068, 4.12068} *)

Another constraint to note is that the solution is only valid for $\measuredangle ACB\geq \frac{\pi }{2}$, which you can use Solve to find the relations:

Last[
 cy /. Simplify[
   Solve[{π/2 == 
      Last@PolygonAngle@Triangle[{{0, 0}, {bx, 0}, {cx, cy}}]}, 
    cy], {bx, cx} ∈ Reals]]
First[cx /. 
  Simplify[Solve[{π/2 == 
      Last@PolygonAngle@Triangle[{{0, 0}, {bx, 0}, {cx, cy}}]}, 
    cx], {bx, cy} ∈ Reals]]
(* Sqrt[(bx - cx) cx] *)
(* 1/2 (bx - Sqrt[bx^2 - 4 cy^2]) *)

Knowing the point C defines the center and radius of the circle. You can use RegionIntersection to solve for E. All other points can be determined from E. Below is a simple demonstration.

ep[cx_, cy_, bx_: 7] := 
 Module[{cf, sal, green, pf, a, b, c, d, r, radius, circ, rc, rl, 
   rint, o, e, ebc, f, tri},
  cf = RGBColor[0.39215686274509803`, 0.5843137254901961, 
    0.9294117647058824];
  sal = RGBColor[
   0.9803921568627451, 0.5019607843137255, 0.44705882352941173`];
  green = RGBColor[
   0.19607843137254902`, 0.803921568627451, 0.19607843137254902`];
  pf = {#3, Disk[#1, 0.25], White, Text[#2, #1, Background -> #3]} &;
  a = {0, 0};
  c = {cx, cy};
  b = {bx, 0};
  (* Circle Center and Radius *)
  r = (b + c)/2;
  radius = Norm[b - c]/2;
  circ = Circle[r, radius];
  (* Find point E by Circle Line Intersection *)
  rc = Region[circ];
  rl = InfiniteLine[a, c];
  rint = RegionIntersection[rc, rl][[1]][[1]];
  o = Ordering[Norm[#] & /@ rint];
  e = rint[[Last@o]];
  (* Rotate BE to align with CB to create F *)
  ebc = Pi/2 - PlanarAngle[c -> {b, e}];
  f = RotationTransform[ebc, b][e];
  tri = Triangle[{e, f, b}];
  (* Use Intersecting Chord Theorem to  create D *)
  d = f + Norm[c - f] Norm[f - b] Normalize[f - e]/Norm[f - e];
  Graphics[{Line[{a, b}], Opacity[0.25], 
    EdgeForm[Directive[Thick, Blue]], Opacity[1], Thick, Line[{c, a}],
     Line[{c, e}], Line[{c, f}], Line[{d, f}], Cyan, tri, Black, 
    EdgeForm[Directive[Thin, Black]], circ, pf[a, "A", sal], 
    pf[b, "B", sal], pf[c, "C", green], pf[e, "E", cf], 
    pf[r, "R", cf], pf[f, "F", cf], pf[d, "D", cf]}, 
   PlotRange -> {{-0.5, bx + 0.8}, {-3 bx/7, 2/3 bx}}]
  ]
Manipulate[
 If[cy > 0.99 Sqrt[(bx - cx) cx], cy = 0.99 Sqrt[(bx - cx) cx]]; 
 ep[cx, cy, bx], {{cx, 2}, 1/2 (bx - Sqrt[bx^2 - 4 cy^2]), 4, 
  Appearance -> "Labeled"}, {{cy, 1.5}, 0.001, 
  0.99 Sqrt[(bx - cx) cx], Appearance -> "Labeled"}, {{bx, 7}, 5, 10, 
  Appearance -> "Labeled"}]