0
$\begingroup$

I am using Minimize to find a solution for the following set of equations that minimizes my cost function $c + k\cdot r$:

Set of equations

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, r\cdot k), p_2=(2r,r + r\cdot 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:

Result of Minimization

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:

Plot for the polynom of k

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.

$\endgroup$

1 Answer 1

4
$\begingroup$
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.

$\endgroup$
5
  • $\begingroup$ Thank you for your answer! I plotted again with the new working precision and it seems to work just fine. Do you think that there is any hope to find another representation of my number (maybe using a completely different approach) or is the fact that Mathematica gave me this root object the proof that there is no "nicer" representation of my number (which might include sin, cos, tan, pi etc.)? $\endgroup$
    – michip96
    Dec 25, 2020 at 0:35
  • $\begingroup$ @michip96 - the result for c has two Root functions. RootReduce should be able to combine them into a single Root function. I don't know what order the polynomials are, but at least for low-order polynomials, ToRadicals can convert Root functions to radical form. Otherwise, Root functions are pretty compact ways of representing exact results. $\endgroup$
    – Bob Hanlon
    Dec 25, 2020 at 2:03
  • $\begingroup$ ToRadicals did not change the Root objects add all. The problem is that I would like to add the exact values of these objects to my thesis, but the factors of the polynomial are so long that their full expansion uses more than three pages. I found that RootApproximant gives me a good approximation. However, I am curious if I can find a better representation of the exact number. $\endgroup$
    – michip96
    Dec 25, 2020 at 10:39
  • $\begingroup$ @michip96 - Put the long form root expressions in an appendix. For the body use sol//N[#,75]//N. The first N ensures precision in the result, the second N shortens the precise value for display. $\endgroup$
    – Bob Hanlon
    Dec 25, 2020 at 13:18
  • $\begingroup$ Thank you for the advice. I will shorten the polynomials as you suggested. $\endgroup$
    – michip96
    Dec 25, 2020 at 13:39

Your Answer

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

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