Hello, bonjour, mingalarbar, kem che, Welcome to the 9th tutorial of Statistical Mechanics: Algorithms and Computations from the Physics Department of Ecole Normale Supérieure. The Monte Carlo method that we have mastered in the last two months and a half in classical and quantum systems has far outgrown its origin. In this week's lecture already we studied time-dependent phenomena and faster-than-time algorithms. Today, Monte Carlo methods are literally everywhere in science and engineering and time has come to present at least one application in a wider context of these methods. The application that we have chosen is concerned with optimization methods. For concreteness we consider a problem that is familiar to all of you and to all of us you know what it is? Well, it is the packing problem that builds the logo of our course. What you see here is a red unit sphere and a bunch of blue unit spheres around it. All of the blue ones touch the red one but none of the blue ones overlap with each other. So the question is how many of the blue unit spheres we can place around the red unit sphere? So, let's find the first solution to this packing problem. What you see here are six blue spheres grouped around a central red sphere in its equatorial plane. Now we add two more layers one on the top and one on the bottom, in fact we put three more spheres on the top and three more spheres on the bottom all of these spheres touch each other and touch the red sphere, and we have 3+6+3: we have found a packing with 12 spheres As you already saw in the logo we can distort this packing a little bit and create space around the blue spheres Let us find a second solution now not with three layers, but with four layers. We start with one sphere at the bottom we add five spheres on the second level then five spheres on the third level and then a final sphere on top. Now, we have 1+5+5+1 spheres all of them touch the red sphere but in this icosahedral packing none of them touches each other. Again, we make a little space and the question arises whether we have enough space to add a 13th blue sphere to this ensemble of 12 spheres. This is the famous 13-spheres problem and it was first discussed by Newton and Gregory in 1694. So clearly there are different possible packings, and you understand that a sampling problem is just around the corner. So, from a wider perspective the subject of this tutorial is the relationship between the sampling approach and for once not the integration but searching, and orientation in complicated configuration spaces. Next, in this tutorial Michael.. next in this tutorial, Michael will show us how to formulate the packing problem not in terms of spheres on the surface of a sphere, but in terms of disks on the surface of the sphere. Then we will go on to discuss strategies to find solutions that were once so well-hidden that before the advent of Monte Carlo it took years and decades to find them. Let us transform the problem a little bit and search for the maximal radius of 13 blue outer spheres of radius R around a central red sphere. If this maximal radius is larger or equal to 1 then there is enough space for Newton's 13th sphere otherwise not. As Werner just explained, we can simplify the problem a bit by considering the projections of the 13 blue outer spheres onto disks inside the red sphere. Let us start by considering the simple case with only two blue outer spheres touching each other and the central red sphere, as is illustrated here. The centers of the two blue disks are a distance of 2R apart while the centers of the projection disks are a distance of 2r apart. You can see that 2R relates to 2r like (1+R) relates to 1 you might remember this relation as Thales' theorem. It gives us a simple relation between the radii of disks touching on a sphere and spheres touching on a sphere, it is given by R = 1 / (1/r - 1) likewise we can see that the disks have to be at least a distance 2r apart, lest they overlap. So let us now follow our winning strategy of the last two months we randomly sample the 13 positions of the blue disks on the surface of the unit sphere. In a first naive approach, we uniformly sample these positions and check for overlaps. Of course, we are only interested in legal configurations that is configurations of non-overlapping disks, and we reject in case of an overlap. In Python, this direct-sampling approach can be implements as is shown here Of course this program works and will eventually give us legal configurations of 13 disks on a sphere. This works very well for small radii r but it suffers from a huge rejection rate for larger disks. Already for disk radius r=0.32 the program will generate 10000 illegal configurations before finding a legal one for disks of radius r=0.34 the program will generate hundreds of thousands of illegal configurations and for r=0.36 it will generate more than 10000000 illegal configurations. This value, r=0.36, corresponds to a sphere radius R=0.5625 and you will see in a few minutes that this is still far away from the close-packing configuration. The problem of finding legal configurations of disks on a sphere is related to the problem of finding legal configurations of hard disks in a box, as we considered in week 2 of this course. You will remember that in tutorial 2 we discussed the program direct_disks_any.py which allowed us to first place the points in the plane and only then compute the minimal distance and density. We can adapt this program to our needs by placing the points on the sphere using the Gaussian random numbers we learned about in week 4 of this course and then computing the minimal distance. This gives the python program shown here. So, is there a smarter way to generate large-radius configurations? Of course there is.. Alberto will in a few moments explain a method to use a Markov chain algorithm for this purpose, but before moving on, please take a moment to download, run and modify the programs direct_sphere_disks.py and its alternative version direct_sphere_disks_any.py that we discussed in this section. Now, after the direct sampling let's go ahead with the Markov chain sampling, following our nine week old strategy for tackling difficult problems. In Python, this gives markov_sphere_disks.py where here we start from a legal configuration of points on the sphere, and propose little moves which are accepted if there is no overlap, so that during our Markov chain simulation the disks never touch. This suggests that we should slightly swell the disks at certain time during the simulation, by a percentage gamma of the maximal increase that would keep the configuration legal. In Python this is done in the function resize_disks.py, that is shown here where we can set the percentage gamma of the maximal increase that we wish to use. Alternating long runs of markov_sphere_disks.py with resize_disks.py allows to construct a very efficient program which is actually a very general optimization strategy. In our example of 13 spheres we can see that we can reach densities that we never dreamt of using direct sampling. The idea of combining long Markov chains with careful resizing is reminiscent of what in metallurgy we call annealing which means that we cool down from high temperature a metal very slowly in order to drive out the grain boundaries and the other imperfections and this makes the metal less brittle. On the computer, the name simulated annealing was coined by Kirkpatrick, Gelatt and Vecchi more than 30 years ago in a very famous paper. In our problem of 13 disks on a sphere each time step is composed by a resizing and 10000 iterations of the Markov chain After each resizing the step of the Markov chain which is set by the standard deviation sigma of the Gaussians is automatically adjusted in order to keep the acceptance probability in accordance with the 1/2 rule. We see that the packing density eta slowly increases during a run until the disks settle in a jammed configuration and the parameter sigma goes to zero. Naturally there are many inequivalent jammed configurations and some of them are non-optimal but they can trap the simulating annealing routine forever and this is a clear limitation of the method. However if you run the simulated annealing routine many times with different random numbers and keeping gamma small you will observe that the solution is in most of the cases the highest density solution this means that the optimal packing configuration has a large basin of  attraction. On our problem of 13 disks on a unit sphere the maximal value of R that we can obtain and actually you will obtain in this week's homework, is R=0.916467.. with our method we cannot actually exclude the value R=1 in fact this is exclude by a proof given by the mathematicians Schutte and van der Waerden in a famous paper of 1951. Before listening to Vivien, (welcome back Vivien) please download, run and modify the program simulated_annealing.py that we discussed in this section. Please make sure to understand how the Markov chain is implemented and the function resize_disks. Let us now generalize what Alberto has explained to us, and this will allow us to illustrate the successes of the simulated annealing method, but also its limitations. As a first generalization, let us consider the case not of 13 disks on the sphere but the case of any number N of spheres around the inner sphere. The program simulated_annealing.py can be easily modified, and you obtain the following result for the maximal density eta_max as a function of N. You observe that eta_max grows slowly with N but it seems that it does not go to the value eta_close_packing = pi / 2 sqrt(3) which is the maximal density for the close packing in two dimensions However we would expect to reach this asymptotic since a lot of small disks on a sphere would feel like living in a two-dimensional plane and thus adopt the hexagonal ordering. Yet, in fact, this method will still reach this asymptotic and it can really be proven that the convergence is very slow The ratio between the maximal density eta_max and eta_close_packing behaves as 1 - const / N^(1/3) Let us now make an analogy we can imagine the jammed configuration as local minima in a vast and complex landscape the optimal configuration will correspond to the global minimum, that we can call the ground state. In fact, many physical problems can be recast in this description and in general the method of the simulated annealing will work because if one takes an equilibrium system and send the temperature to zero, it will adopt the ground state. However, this argument is purely formal because in practice simulated annealing cannot be done at a very very small velocity which we would require for this method to work always. The nice behavior we have discussed is in fact limited to simple cases: disks, spheres, .. For instance in the case of polydisperse disks the method of simulated annealing does not find the optimal configuration, but still finds good configurations. In real life, simulated annealing must be practiced with care especially in cases where there are imperfections and disorder because then there are a lot of metastable states and the system may be stuck in them. If the imperfections are too pronounced the result of the algorithm is very often a glass rather than a crystalline matter. For the case of the polydisperse disks on a sphere, that is to say disks with varied radius, this is what happens. The disorder prevents the system from reaching its optimal configuration and leads the system into one of the jammed configurations whose basins of attraction is much larger than in the case of monodisperse disks Before continuing, let me grab the object of the week the flip die, and let's flip faster than the clock. So, in conclusion, we have studied in this Tutorial 9 one of the many connection between statistical mechanics and the rest of life. You will continue on the subject of packing in this week's homework but don't forget that there is a connection between sampling not only with integration as we have studied throughout this course but also with searching and orientation in complicated spaces So finally good luck and have a lot of fun with Homework session 9, the final homework session of this course. Next week, we will have a lecture on the alpha and the omega of Monte Carlo and the final tutorial will be devoted to a review and to a big party. And thanks again for your interest, and for the very active part that you're taking in our course.