This post is inspired by minhthien2016's question.

The problem, denoted 2/N/1, for reasons that will appear clearly further on, is to pack N disks into the unit square in such a way that the sum of their radii is maximum.

I replied this problem using Optimization-NLPSolve for N from 1 (obvious solution) to 16, which motivated a few questions, in particular:

• @Carl Love: "Can we confirm that the maxima are global (NLPSolve tends to return local optima)?" 
  Using NLPSolve indeed does not guarantee that the solution found is the (a?) global maximum. In fact packing problems are generaly tackled by using evolutionnary algorithms, greedy algorithms, or specific heuristic strategies.
  Nevertheless, running NLPSolve a large number of times from different initial points may provide several different optima whose the largest one can be reasonably considered as the global maximum you are looking for.
  Of course this may come to a large price in term of computational time.
   
• @acer: "How tight are [some constraints], at different N? Are they always equality?"
  The fact some inequality constraints type always end to equality constraints (which means that in an optimal packing each disk touches at least one other annd, possibly the boundary of the box) seems hard to prove mathematically, but I gave here a sketch of informal proof.

I found 2/N/1 funny enough to spend some time digging into it and looking to some generalizations I will refer to as D/N/M:  How to pack N D-hypersheres into the unit D-hypercube such that the sum of the M-th power of their radii is maximum?
For the sake of simplicity I will say ball instead of disk/sphere/hypersphere and box instead of square/cube/hypercube.

The first point is that problems D/N/1 do not have a unique solution as soon as N > 1 , indeed any solution can be transformed into another one using symmetries with respect to medians and diagonals of the box. Hereafter I use this convention:

Two solutions s and s' are fundamental solutions if:

1. the ordered lists of radii s and s'  contain are identical but there is no composition of symmetries from s to s',
2. or, the ordered lists of radii s and s'  contain are not identical.
    

It is easy to prove that 2/2/1 and 3/2/1, and likely D/2/1, have an infinity of fundamental solutions: see directory FOCUS__2|2|1_and_3|2|1 in the attached zip file..
At the same time 2/N/2, 3/N/3, and likely D/N/D, have only one fundamental solution (see directory FOCUS__2|N|2 for more details and a simple way to characterize these solutions

 (Indeed the strategy ito find the solution of D/N/D  in placing the biggest possible ball in the largest void D/N-1/D contains. Unfortunately this characterization is not algorithmically constructive in the sense that findind this biggest void is a very complex geometrical and combinatorial problem.
 it requires finding the largest void  in a pack of balls)

Let M_{d, 1}(N)  the maximum value of the sum of balls radii for problem d/N/1.
The first question I asked myself is: How does M_{d, 1}(N) grows with N?
 

(M_{d, 1}(N) is obviously a strictly increasing function of N: indeed the solution of problem d/N/1 contains several voids where a ball of strictly positive radius can be placed, then  M_{d+1, 1}(N) > M_{d, 1}(N) )

The answer seems amazing as intensive numerical computations suggest that
                                      

See D|N|M__Growth_law.mw in the attached sip file.
This formula fits very well the set of points  { [n, S_{d, 1}(n) , n=1..48) } for d=2..6.
I have the feeling that this conjecture might be proven (rejected?) by rigourous mathematical arguments.


Fundamental solutions raise several open problems:

• Are D/2/1 problems the only one with more than one fundamental solutions?
  
  The truth is that I have not been capable to find any other example (which does not mean they do not exist).
  A quite strange thing is the behaviour of NLPSolve: as all the solutions of D/2/1 are equally likely, the histogram of the solutions provided by a large number of NLPSolve runs from different initial points is far from being uniform.
  For more detail refer ro directory FOCUS__2|2|1_and_3|2|1 in the attached zip file. 
  I do not understand where this bias comes from: is it due to the implementation of SQP in NLPSolve, or to SQP itself?
   
• For some couples (D, N) the solution of D/N/1 is made of balls of same radius.
  For N from 1 to 48 this is (numerically) the case for 2/1/1 and2/2/1, but the three dimensional case is reacher as it contains  3/1/1, 3/2/1,  3/3/1,  3/4/1 and 3/14/1 (this latter being quite surprising).
  Is there only a finite number of values N such that D/N/1 is made of balls with identical radii?
  If it is so, is this number increasing with D?
  It is worth noting that those values of N mean that the solution of problems D/N/1 are identical to those of a more classic packing problem: "What is the largest radius N identical balls packed in a unit bow may have?".
  For an exhaustive survey of this latter problem see Packomania.
   
• A related question is "How does the number of different radii evolves as N increases dor given values of D and M?
  Displays of 2D and 3D packings show that the set of radii has significantly less elements than N... at least for values of N not too large. So might we expect that solution of, let us say, 2/100/1 can contain 100 balls of 10 different radii, or it is more reasonable to expect it contains 100 balls of 100 different radii?
   
• At the opposite numerical investigations of  2/N/1 and  3/N/1 suggest that the number of different radii a fundamental solution contains increases with N (more a trend than a continuous growth).
  So, is it true that very large values of N correspond to solutions where the number of different radii is also very large?
  Or could it be that the growth of the number of different radii I observed is simply the consequence of partially converged results?
   
• Numerical investigations show that for a given dimension d and a given number of balls n,  solutions of d/n/1 and d/n/M (1 < M < d) problems are rather often the same. Is this a rule with a few exceptions or a false impression due to the fact that I did not pushed the simulations to values of n large enough to draw a more solid picture)?

It is likely that some of the questions above could be adressed by using a more powerful machine than mine.


All the codes and results are gathered in the attached PACKING.zip file a zip file you can download from OneDrive Google  (link at the end of this post, 262 Mb, 570 Mb when unzipped, 1119 files).
Install this zip file in the directory of your choice and unzip it to get a directory named PACKING. 
Within it:

• README.mw contains a description of the different codes and directories
• Repository.rtf must contain a string repesenting the absolute path of directory PACKING. 

Follow this link OneDrive Google