Packing squares in a rectangle I need to pack $N$ identical squares arranged in array fashion, $R$ rows of $C$ columns, the last row possibly incomplete, in a $W\times H$ rectangle (reals). The squares must be as large as possible.
With $S$ the side of the squares, this leads to the system of inequations
$$\begin{cases}CS\le W,\\RS\le H,\\CR\ge N.\end{cases}$$
Multiplying the first two and dividing by the third,
$$S^2\le\frac{WH}N.$$
My first naïve approach was to compute $S$ by assuming that the equality holds and deduce
$$C=\left\lfloor\frac WS\right\rfloor,R=\left\lfloor\frac HS\right\rfloor.$$
But this is wrong as the third inequation can be invalidated.
My second approach is by trying every possible value of $C$ from $1$ to $N$, deducing $R$ from $CR\ge N$, then $S$. Hence, finding the $C$ that maximizes
$$\min_C\left(\frac WC,\frac H{\left\lceil\dfrac NC\right\rceil}\right)$$
or minimizes 
$$\max_C\left(\frac HWC,\left\lceil\dfrac NC\right\rceil\right).$$
I suspect that there should be a closed-form solution, probably close to my naïve estimate, but so far I haven't found a line of attack. Can you help ?

Update:
Following @DarthGeek, the solution is found where $\dfrac{\left\lceil\dfrac NC\right\rceil}C$ crosses $\dfrac HW=:a$. This function is piecewise linear and monotonic. As we are only interested in values at integers, a linear or dichotomic search can do.
From
$$\frac NC\le\left\lceil\frac NC\right\rceil<\frac NC+1$$
we deduce that the solution will lie in
$$\left[\sqrt{\frac{N}a},\sqrt{\frac{N}a+\frac1{4a^2}}+\frac1{2a}\right).$$
 A: I don't have a "strong" proof for this solution but I think this is the right one.
For each set $W,H,N$ there are $N$ possible configurations:
$$C = n\quad  R = \left\lceil\dfrac{N}{n}\right\rceil \qquad n = 1,\ldots, N$$
Consider the set of rational numbers $\left\{\dfrac{\lceil N/n\rceil}{n},\  n = 1,\ldots, N\right\}$ and write them from smallest to largest: $\left\{\dfrac{a_1}{b_1},\ldots,\dfrac{a_N}{b_N}\right\}$.
If $\dfrac{H}{W} \leq \dfrac{a_1}{b_1} = \dfrac{1}{N}$ then the optimal solution is:
$$C = N \qquad R = 1 \qquad S = H/N$$
If $\dfrac{a_n}{b_n} \leq \dfrac{H}{W} \leq \dfrac{a_{n+1}}{b_{n+1}}$ then the optimal solution is:
$$C = b_{n+1} \qquad R = a_{n+1} \qquad S = H/R $$
If $\dfrac{H}{W} \geq \dfrac{a_N}{b_N} = N$ then the optimal solution is:
$$ C = 1 \qquad R = N \qquad S = W/N$$
A: Given are $n \in \mathbb{N}$ and $w, h \in \mathbb{R}$, both positive. 
We have variables $s \in \mathbb{R}$ and $c, r \in \mathbb{Z}$.
The problem is:
$$
\begin{array}{lr}
\max  & s \\
\text{w.r.t.} 
& c \, s \le w \\
& r \, s \le h \\
& n \le cr \\
& s > 0 \\
& c \ge 1 \\
& c \le n \\
& r \ge 1 \\
& r \le n
\end{array}
$$
We try to express $r$ in terms of $c$ and $n$, to have just two variables for a graphical solution:
$$
\begin{array}{lr}
\max  & s \\
\text{w.r.t.} 
& s \le w / c \\
& r(c) = \lfloor n / c \rfloor + 
         \Theta\left( n - c \lfloor n / c \rfloor \right) \\ 
& s \le h / r(c) \\
& s > 0 \\
& c \ge 1 \\
& c \le n
\end{array}
$$
where $\Theta$ is the Heaviside function, which can be approximated by $\DeclareMathOperator{sgn}{sgn}\Theta_{1/2}(x) = (\sgn(x) + 1) / 2 $.
A resulting view for $n=9$, $w=18$, $h=10$ is:
 
The $x$-coordinate features the choice of $c$, the $y$-coordinate the choice of $s$. 
The constraint $s \le w/c$ is the light green surface.
The upper border of the constraint $s \le h/r(c)$ is formed by the red lines. The feasible region is the area below them, not drawn alas.
The constraints $c \ge 1$ and $c \le n$ are the dark blue lines.
Finally we have to keep in mind that $c$ is an integer.
The objective is to get a maximal $y$-coordinate within the feasible region.
Here the optimal solution seems be $(c, s) = (5, 18/5)$ (magenta point).
The image suggests to start searching close to the intersection (yellow point) of
$s(c) = w/c$ (green line) and $s(c) = (h/n) c$ (light blue line).
