Find tangency of InterpolatingFunction and circle

Question

I am trying to use a system of equations to solve for the common tangent of a (polar) InterpolatingFunction and a circle. Unfortunately, both NSolve and FindRoot don't yield results, even when I "hand-feed" values into the equations into the equations (i.e. in combinations of equations 1,3,4 / 2,3,4 or 1,2,3)

My guess is that this may have to do with "excessive" precision, although I'm not sure.

Question: what function would I best use to get the closest approximations, if it is a precision problem? For my purposes, three decimal places would be sufficient.

(If someone sees a better basic geometric strategy to go at the problem, please do say.)

Please have a look at the graphic and the code:

diagram

(* intf[ang1] is the InterpolatingFunction (green),
rotated about {0,0} until it is tangent to the circle.
Goal: find out ang1, ang3, angrotation *)

intpoints = {{0, 17.9817}, {0.0521424, 18.5701}, {0.105454, 
    19.1892}, {0.15991, 19.8416}, {0.215494, 20.5299}, {0.272199, 
    21.2567}, {0.330029, 22.0243}, {0.388993, 22.8348}, {0.449111, 
    23.6897}, {0.510411, 24.5896}, {0.572931, 25.5338}, {0.636716, 
    26.5201}, {0.701823, 27.5438}, {0.76832, 28.5972}, {0.836286, 
    29.6684}, {0.905813, 30.7409}, {0.977011, 31.7919}, {1.05, 
    32.7921}, {1.12493, 33.7052}, {1.20195, 34.4887}, {1.28123, 
    35.0956}, {1.36295, 35.4775}, {1.44731, 35.5905}};
intf = Interpolation[intpoints]

(* values used in example *)
r1 = 9.5;
x1 = 34.2584;
y1 = 12;

(* 1. parametric equations for a circle of radius R1 with center \
located at {X1,Y1} *)
xcirclefn[ang2_] = r1*Cos[ang2] + x1;
ycirclefn[ang2_] = r1*Sin[ang2] + y1;

(* 2. equations for locating points based on rotation of axes of intf \
*)
xintf[ang1_] = intf[ang1]*Cos[ang1];
yintf[ang1_] = intf[ang1]*Sin[ang1];
xrot[ang1_, angrotation_] = 
  xintf[ang1]*Cos[angrotation] - yintf[ang1]*Sin[angrotation];
yrot[ang1_, angrotation_] = 
  xintf[ang1]*Sin[angrotation] + yintf[ang1]*Cos[angrotation];

(* 3. equations for angles of tangents *)
slopecircle[ang2_] = ycirclefn'[ang2]/xcirclefn'[ang2] ;
angleang3c[ang2_] = ArcTan[-1/slopecircle[ang2]] ;

slopeintf[ang1_] = yintf'[ang1]/xintf'[ang1];
angleang3i[ang1_, angrotation_] = 
  ArcTan[-1/slopeintf[ang1]] + angrotation;
(* angrotation is negative *)

(* Possible identities: *)
equation1 = xcirclefn[ang2] == xrot[ang1, angrotation];
equation2 = ycirclefn[ang2] == yrot[ang1, angrotation];
equation3 = angleang3c[ang2] == angleang3i[ang1, angrotation];
equation4 = intf[ang1] == Sqrt[xcirclefn[ang2]^2 + ycirclefn[ang2]^2];

(* CAD-measured values in radians:
ang1 = 0.714321 (CCW from intf[0.])
 ang2 = 3.09004 (CCW from 0.)
angrotation = -0.247318 (CW from 0.)
ang3 = -(2*Pi)+(Pi/2 + (-0.0515572)) = -4.76395 (CW from 0.) (1.51924 included angle) *)

(* Method 1:
NSolve[equation1&&equation2&&equation3,{ang1,ang2,angrotation}] *)

(* Method 2:
FindRoot[{equation1,equation3,equation4},{{ang1,0},{ang2,0},{\
angrotation,-Pi/2}}] *)

Could you complete your example with the interpolating function and the solving approaches you already tested?
And maybe, if possible, switch to lowercase symbol names throughout? Uppercase sometimes collides with predefined symbols, causing all kinds of "surprising" behavior. — Jinxed, Feb 22 '15 at 23:24
You might also want to define your functions using := (SetDelayed). — Jinxed, Feb 22 '15 at 23:29
Depending on the Interpolation you used, there might not be a usable derivative. See the documentation for examples. — Jinxed, Feb 22 '15 at 23:31
If you could provide me with your intf, I might be able to help. — Jinxed, Feb 23 '15 at 00:42
Welcome to Mathematica.SE! I suggest the following: 0) Browse the common pitfalls question 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Read the [faq]! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! — Dr. belisarius, Feb 23 '15 at 06:58

Dr. belisarius · Accepted Answer · 2015-02-23T16:19:30.100

2

Not sure if I understand, anyway:

intpoints = {{0, 17.9817}, {0.0521424, 18.5701}, {0.105454, 19.1892}, {0.15991, 19.8416}, 
            {0.215494, 20.5299}, {0.272199, 21.2567}, {0.330029, 22.0243}, {0.388993, 22.8348}, 
            {0.449111, 23.6897}, {0.510411, 24.5896}, {0.572931, 25.5338}, {0.636716, 26.5201}, 
            {0.701823, 27.5438}, {0.76832, 28.5972}, {0.836286,  29.6684}, {0.905813, 30.7409}, 
            {0.977011, 31.7919}, {1.05, 32.7921}, {1.12493, 33.7052}, {1.20195, 34.4887}, 
            {1.28123,  35.0956}, {1.36295, 35.4775}, {1.44731, 35.5905}};
intf = Interpolation[intpoints]

(*values used in example*)
r1 = 9.5; x1 = 34.2584; y1 = 12;
exp = RotationMatrix[-ang1].(intf[p] {Cos@p, Sin@p});
sol = NMinimize[{ang1, EuclideanDistance[{x1, y1}, exp] == r1 && 0 < p < Pi/2}, {p,  ang1}];
Show[ContourPlot[(x - x1)^2 + (y - y1)^2 == r1^2, {x, 0, 50}, {y, -10, 40}], 
     ParametricPlot[intf[p] {Cos@p, Sin@p}, {p, 0, Pi/2}], 
     ParametricPlot[exp /. sol[[2]], {p, 0.01, Pi/2}, PlotStyle -> Red]]

Mathematica graphics

edited Feb 23 '15 at 16:19

answered Feb 23 '15 at 06:30

Dr. belisarius

115,881
13
203
453

@Jinxed: info added, thanks! – jharrison Feb 23 '15 at 06:50
belisarius: thanks, will try adapting this method and see if I can pull the angles out. Looks much cleaner than mine. – jharrison Feb 23 '15 at 06:56
@jharrison There you have it done – Dr. belisarius Feb 23 '15 at 06:57
I started fresh and copied your code, and got: NMinimize::nsol: There are no points that satisfy the constraints {}. >> – jharrison Feb 23 '15 at 07:58
@jharrison: Same for me. – Jinxed Feb 23 '15 at 12:30
@Jinxed: if I use the EuclideanDistance[{x1,y1}...] == r1 as a solver equation, and hand-feed it a value for p, NSolve then successfully finds the corresponding RotationMatrix angle in a snap. Hmmm. – jharrison Feb 23 '15 at 12:42
@Jinxed: belisarius is apparently running v9 (see comments in http://mathematica.stackexchange.com/questions/75704/integral-does-not-converge-when-it-should). I was using v10, but just tried the NMinimize in v9 and it works. – jharrison Feb 23 '15 at 13:41
@jharrison Yup.V9 here. Sorry, I can't test it on v10. – Dr. belisarius Feb 23 '15 at 13:58
Hmm - another v10-oddity. :( My approach was based on solving for curve-circle-coincidence, but that did produce weird results, also. I am still investigating why that might be. – Jinxed Feb 23 '15 at 15:39
I'm good with this solution, will stay in v9 for the nonce. – jharrison Feb 23 '15 at 16:03
@jharrison Good to know :). Please re-read the following comment – Dr. belisarius Feb 23 '15 at 16:09
Welcome to Mathematica.SE! I suggest the following: 0) Browse the common pitfalls question 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Read the [faq]! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! – Dr. belisarius Feb 23 '15 at 16:10

Find tangency of InterpolatingFunction and circle

1 Answers1