Mathematica Asked by michip96 on December 26, 2020
I am using Minimize
to find a solution for the following set of equations that minimizes my cost function $c + kcdot r$:
Essentially, I have two points $p_1,p_2$ with x coordinates $r$ and $2r$ and initial y coordinates $0$ and $r$ that I can move upwards to a new position $p_1=(r, rcdot k), p_2=(2r,r + rcdot k)$ in the y direction. I am looking for a $k$ that minimizes the shortest distance $|p_1q|+|p_2q|$ to a point $q=(x,y)$ on a circle $C$ centered at $(x_0,y_0)=(0,2r)$ with radius $r$. For a given $k$ the point can be found when calculating the intersection point between my fixed circle and the ellipse $E$, centered at $(x_e,y_e)$ induced by $p_1,p_2$:
$$x_e = frac{1}{2}(p^x_1 + p^x_2) = frac{3}{2}r$$
$$y_e = frac{1}{2}(p^y_1 + p^y_2) = frac{1}{2}(r+2kr)$$
Let $f$ be the distance from the center point of $E$ to the foci
$$f=sqrt{(p^x_1-x_e)^2 + (p^y_1-y_e)^2} = sqrt{frac{r^2}{4} + frac{(r+kr)^2}{4}} = frac{1}{2}sqrt{k cdot (k+2) + 2}r$$
The ellipse can now be defined with its center point $(x_e,y_e)$ and the major/minor axis $a,b$.
$$a= frac{1}{2} c$$
$$b= sqrt{a^2 – f^2} = frac{1}{2} sqrt{c^2 – (k cdot (k+2) + 2)r^2}$$
In order to find the minimal solution, I assumed that $r=1$ to simplify the formulas above:
Dist[x_, y_, z_, v_] := Sqrt[(x - z)^2 + (y - v)^2]
r = 1;
xPos[k_] := 3/2* r
yPos[k_] := 1/2 * (r + k*r)
f[k_] := Dist[xPos[k], yPos[k], r, r*k]
a[c_, k_] := 1/2 * c
b[c_, k_] := Sqrt[a[c, k]^2 - f[k]^2]
res = Minimize[{c + k*r, x^2 + (y - 2 *r)^2 == r^2 && (x - xPos[k])^2/a[c, k]^2 + (y -
yPos[k])^2/b[c, k]^2 == 1 && 0 <= k <= 1 && c == Dist[x, y, r, k*r] + Dist[x, y, 2 r, r +
k*r] && 0 <= x <= r && r <= y <= 2 r && c > 0 }, {c, x, y, k }, Reals]
I am sure that this might not be the most elegant solution but after a few hours of computation I get the following result:
As you can see in the image all values for $x,y,k,c,cost$ are Root
objects. I am fine that Mathematica gives me these objects instead of a closed form, but I am worried that my Minimization has an error, because all Root
objects have exponents up to $128$ with very strange factors. Here is for example the beginning of the polynom for $k$:
$$8837684323394741881152434292665551758032896 – 438169602476233057174556091222444606626988032*x + … $$
As all of my variables have these kinds of factors and all roots have exponents up to 128, my guess is that Minimize
uses some numerical approximations that produce these strange numbers and that the precision is somehow limited to an exponent of $128$. I also plotted the Root polynom for $k$ between $0.48$ and $0.49$ which confuses me even more:
Are there different approaches that I can use to verify my result or techniques to simplify my input formulas? I am assuming that this behavior has something to do with the strange form of $c$ that is the sum of two square roots.
Clear["Global`*"]
Dist[x_, y_, z_, v_] := Sqrt[(x - z)^2 + (y - v)^2]
r = 1;
xPos[k_] := 3/2*r
yPos[k_] := 1/2*(r + k*r)
f[k_] := Dist[xPos[k], yPos[k], r, r*k]
a[c_, k_] := 1/2*c
b[c_, k_] := Sqrt[a[c, k]^2 - f[k]^2]
sys = {c + k*r,
x^2 + (y - 2*r)^2 ==
r^2 && (x - xPos[k])^2/a[c, k]^2 + (y - yPos[k])^2/b[c, k]^2 == 1 &&
0 <= k <= 1 && c == Dist[x, y, r, k*r] + Dist[x, y, 2 r, r + k*r] &&
0 <= x <= r && r <= y <= 2 r && c > 0};
Since you said that your results took "a few hours of computation", I have verified with a numeric result.
resN = NMinimize[sys, {c, x, y, k}, WorkingPrecision -> 30] // N
(* {2.6338, {c -> 2.15115, x -> 0.848236, y -> 1.47038, k -> 0.482657}} *)
The numeric results agree with the values shown in the picture of your Root
expressions
Verifying that the result satisfies the system
(sys[[1]] == resN[[1]] && Rest[sys]) /. resN[[2]]
(* {True} *)
You don't show the code that you used for your Plot
; however, I assume that you used machine precision. Since the Root
expressions have coefficients as large as you show, then a large WorkingPrecision
is required to maintain precision. Start with WorkingPrecision -> 75
. Increase if required.
Correct answer by Bob Hanlon on December 26, 2020
Get help from others!
Recent Questions
Recent Answers
© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP