3
$\begingroup$

Let $m$ be a large integer and $z\in (0,1)$. I want to find $0 < k < m$ such that

$\frac{\sqrt{5-\frac{4k}{m}}-1}{2} + \frac{m-\frac{k}{2}\sqrt{5-\frac{4k}{m}}}{m^2\sqrt{5-\frac{4k}{m}}}$ is closest to $z$. So I set up an equation by letting the preceding expression equal $z$ and solved it for $k$. That is,

Solve[(Sqrt[5 - 4 k/m] - 1)/2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) == z, k]

and it gave some complicated symbolic solutions, as the equation can be easily reduced to a cubic equation. Then I plugged in the actual values $m=20000$ and $z=0.1$,

N[% /. {m -> 20000, z -> 0.1}]

which evaluates as

{{k -> 143373. + 0.0527344 I}, 
 {k -> -3.2001*10^13 + 0. I}, 
 {k -> -100572. - 0.0546875 I}}

which is clearly wrong. There should be a $k$ that is less than $m$, because the initial expression is at least $1/(2m)$ (plugging in $k = m$) and my choice of $z$ is greater than $1/(2m)$.

However, if I plug in the numbers directly into the equation,

Solve[(Sqrt[5 - 4 k/m] - 1)/2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) 
    == z /. {m -> 20000, z -> 0.1}, k]

it gives the correct answer

{{k -> 17800.5}, {k -> 25000.}}

So I am wondering whether this is a bug.

$\endgroup$
  • $\begingroup$ Block[{z = 1/10, m = 20000}, Solve[eq == z, k, Reals]] or in longer form Block[{m = 20000, z = 1/10}, Solve[(Sqrt[5 - 4 k/m] - 1)/ 2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) == z, k, Reals]] // N $\endgroup$ – ciao Sep 29 '15 at 23:35
5
$\begingroup$

You have not encountered a bug but a weakness in machine precision arithmetic. If you switch to Mathematica's more accurate arbitrary precision arithmetic, you can get away with an unconstrained (sloppy) problem statement.

sol = 
  Solve[
    (Sqrt[5 - 4 k/m] - 1)/2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) == z, 
    k];
N[sol /. {m -> 20000, z -> 1/10}, 25] // Chop
{{k -> 24999.99996527657202992317}, 
 {k -> -3.200096001780046597611895*10^13}, 
 {k -> 17800.46601084237600423031}}

But putting the problem to Solve with proper constraints (as others have observed) is a much better way to proceed because Solve, working under the constraints, give a solution that evaluates well even at machine precison.

sol = 
  Solve[
      (Sqrt[5 - 4 k/m] - 1)/2 + (m - k/2 Sqrt[5 - 4 k/m]) /
        (Sqrt[5 - 4 k/m] m^2) == z 
    && 
      0 < k < m, 
    k, Reals]
 N[sol /. {m -> 20000, z -> .1}]
{{k -> 17800.5}, {k -> Undefined}}
$\endgroup$
5
$\begingroup$

You need to build your assumptions into the problem. Reduce appears to work better than Solve for your problem. Let's look at both.

Solve

If one limits the solution to real numbers and wraps Solve with Assumptions which apply your inequalities a solution is found

sol = Assuming[
  m > k > 0 && z > 0,
  Solve[(Sqrt[5 - 4 k/m] - 1)/
      2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) == z, k, 
   Reals
   ]
  ]

Substituting {m, z} = {20000., 0.1} results in

sol /. {m -> 20000., z -> 0.1}

{{k -> Undefined}, {k -> 17800.5}, {k -> 25000.}}

The middle rule matches your expected result.

Reduce

It appears to go smoother using Reduce.

sol = Reduce[(Sqrt[5 - 4 k/m] - 1)/
      2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) == z && 
   m > k > 0 && z > 0, k, Reals]

Reduce produces a single answer

sol /. {m -> 20000., z -> 0.1}

k == 17800.5
$\endgroup$
3
$\begingroup$

Use Minimize. Since you want k such that the expression is "closest" you need to allow for the difference being positive or negative. I chose to square the difference to handle this.

diff[z_, m_] = 
  z - (Sqrt[5 - 4 k/m] - 1)/
    2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2);

soln = With[{m = 20000, z = 0.1}, Minimize[{diff[z, m]^2, 0 < k < m}, k]]

(*  {6.83353*10^-19, {k -> 17799.5}}  *)

The difference is

diff[z, m] /. {z -> 0.1, m -> 20000, soln[[-1, 1]]}

(*  -8.26651*10^-10  *)

Compared with the original result

diff[z, m] /. {z -> 0.1, m -> 20000, 
  Solve[(Sqrt[5 - 4 k/m] - 1)/
        2 + (m - k/2 Sqrt[5 - 4 k/m])/(Sqrt[5 - 4 k/m] m^2) == 
      z /. {m -> 20000, z -> 0.1}, k][[1, 1]]}

(*  0.0000388349  *)

Minimize gives the "closest"

$\endgroup$
  • $\begingroup$ Your answer gets my vote. Clearly answers question about closest. $\endgroup$ – Jack LaVigne Sep 30 '15 at 14:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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