# Mathematica can't find Minimum sum under integer constrains

I'm trying to verify a solution to a simple probability problem using Mathematica. Here's the problem:

A drawer contains red socks and black socks. When two socks are drawn at random, the probability that both are red is 1/2. How small can the number of socks in the drawer be if the number of black socks is even?

From manual trial and error, it is trivial to find a solution of red = 15 and black = 6. However, when I try to verify the solution, Minimize and NMinimize fails to find an exact solution.

Minimize simply returns without having found a solution.

Minimize[{r + b, Binomial[r, 2]/Binomial[r + b, 2] == 1/2, r > 0,
b > 0, b/2 == i}, {r \[Element] Integers, b \[Element] Integers,
i \[Element] Integers}]


NMinimize throws an NMinimize:nosat error and proposes a close enough solution of red = 5 and black = 2. The probability works out to be $$\frac{\binom{5}{2}}{\binom{7}{2}}=\frac{10}{21}$$, which is not exactly 1/2.

NMinimize[{r + b, Binomial[r, 2]/Binomial[r + b, 2] == 1/2, r > 0,
b > 0, b/2 == i}, {r \[Element] Integers, b \[Element] Integers,
i \[Element] Integers}]


Is there something wrong with my setup, or is this a bug? I've tried all Methods, but neither made a difference

We solve the constrains at first.

sol = Solve[{Binomial[r, 2]/Binomial[r + b, 2] == 1/2, r > 0, b > 0,
b/2 == i}, {r ∈ Integers, b ∈ Integers,
i ∈ Integers}][[1]]


min=Minimize[{r + b /. sol /. C[1] -> c,
c ∈ PositiveIntegers}, {c}] // Simplify

{r, b} /. sol /. C[1] -> c /. min[[2]] // Simplify


{15, 6}

• That works, thanks! What do you think makes Minimize struggle with this? Am I misusing it? Commented Sep 18, 2022 at 0:03

I suppose that the problem might be difficult for the minimizer because there are many local minima.

A direct search approach is the following:

MinimalBy[
Tuples[Range[20], {2}],
Apply[#1 + #2 + 1000 Boole@OddQ[#2] + 10000 Abs[(Binomial[#1, 2]/Binomial[#1 + #2, 2] - 1/2)]&]
]

(* Out: {{15, 6}} *)


One can get NMinimize to work without previously solving the constraints by modifying the objective function to minimize simultaneously the sum of r and b, and the absolute deviation from the desired probability. In practice, that means minimizing the product of the sum and the absolute value of (the probability - 1/2). Additionally, we provide all possible 2-tuples as starting points for the minimization method. Finally, the constraint that b should be even is treated directly through an auxiliary function rather than introduction another parameter.

Clear[even]
even[x_?NumericQ] := Boole@EvenQ[x]

NMinimize[{
(r + b) Abs[Binomial[r, 2]/Binomial[r + b, 2] - 1/2],
even[b] == 1,
r \[Element] PositiveIntegers,
b \[Element] PositiveIntegers},
{r, b},
Method -> {"NelderMead", "InitialPoints" -> Tuples[Range[20], {2}]}
]

(* Out: {0., {r -> 15, b -> 6}} *)