1 00:00:01,787 --> 00:00:04,120 GERBRAND CEDER: --temperature finite temperature methods 2 00:00:04,120 --> 00:00:06,760 of which you did molecular dynamics earlier with Professor 3 00:00:06,760 --> 00:00:07,780 Marzari. 4 00:00:07,780 --> 00:00:09,520 So what we're going to do is sort of 5 00:00:09,520 --> 00:00:15,070 continue that for a bit, showing you Monte Carlo simulations 6 00:00:15,070 --> 00:00:17,080 and how you can essentially use it. 7 00:00:17,080 --> 00:00:20,200 And after that is kind of a series of fairly random topics. 8 00:00:24,460 --> 00:00:26,922 In some cases, clearly we don't even know what they are. 9 00:00:26,922 --> 00:00:28,630 Even when it says case study, that really 10 00:00:28,630 --> 00:00:31,750 means that the instructors don't have a clue yet what they're 11 00:00:31,750 --> 00:00:35,050 going to be teaching about. 12 00:00:35,050 --> 00:00:37,040 Somebody's listening in I think. 13 00:00:37,040 --> 00:00:37,540 OK. 14 00:00:37,540 --> 00:00:45,160 So I want to have a bit of a discussion on how you calculate 15 00:00:45,160 --> 00:00:46,973 certain properties. 16 00:00:46,973 --> 00:00:48,640 This is a good way to see the difference 17 00:00:48,640 --> 00:00:52,720 between molecular dynamic simulation and Monte Carlo. 18 00:00:52,720 --> 00:00:55,810 So in a molecular dynamic simulation, 19 00:00:55,810 --> 00:00:58,450 one way you calculate certain properties 20 00:00:58,450 --> 00:01:00,980 is by tracking them over time. 21 00:01:00,980 --> 00:01:03,880 So if you can calculate a certain property 22 00:01:03,880 --> 00:01:08,530 for a given microscopic state, for a given configuration, 23 00:01:08,530 --> 00:01:10,330 say the energy-- 24 00:01:10,330 --> 00:01:11,740 let's take some simple ones-- 25 00:01:11,740 --> 00:01:15,920 or the volume, you calculate their macroscopic averages, 26 00:01:15,920 --> 00:01:19,840 which are their thermodynamic qualities, really 27 00:01:19,840 --> 00:01:23,360 as an integral over time. 28 00:01:23,360 --> 00:01:26,330 So there's sort of two issues with that, 29 00:01:26,330 --> 00:01:31,900 that obviously in that average you will only include states 30 00:01:31,900 --> 00:01:34,450 that you can reach with the molecular dynamic simulation. 31 00:01:34,450 --> 00:01:37,660 And so since you simulate over a finite time, 32 00:01:37,660 --> 00:01:40,780 if there are excitations, say, which 33 00:01:40,780 --> 00:01:43,600 can determine the average energy in reality that 34 00:01:43,600 --> 00:01:45,340 occur over a longer timescale, you 35 00:01:45,340 --> 00:01:48,140 will not sample them in that average. 36 00:01:48,140 --> 00:01:52,880 The second thing is that if all you care about 37 00:01:52,880 --> 00:01:57,110 are the average properties but not the dynamics of how 38 00:01:57,110 --> 00:01:59,150 you go between them, then sometimes you 39 00:01:59,150 --> 00:02:01,970 can do simpler techniques than molecular dynamics, 40 00:02:01,970 --> 00:02:04,580 or things that are less computing intensive. 41 00:02:04,580 --> 00:02:06,793 You can literally just go and sample those averages. 42 00:02:06,793 --> 00:02:08,210 And so that's basically what we're 43 00:02:08,210 --> 00:02:10,610 going to talk about this lecture and the next lecture, 44 00:02:10,610 --> 00:02:12,570 talk about how you do that sampling. 45 00:02:12,570 --> 00:02:16,580 So if you only care about the averages but not the dynamics 46 00:02:16,580 --> 00:02:17,510 itself. 47 00:02:17,510 --> 00:02:20,090 You'll see that for some models of materials 48 00:02:20,090 --> 00:02:22,670 you'll actually need to do that, because you 49 00:02:22,670 --> 00:02:24,050 won't know the dynamics. 50 00:02:24,050 --> 00:02:25,520 Remember that in molecular dynamics 51 00:02:25,520 --> 00:02:29,240 you basically assume that the atoms move with Newton's 52 00:02:29,240 --> 00:02:30,155 equation of motion. 53 00:02:30,155 --> 00:02:32,060 So you know the dynamic. 54 00:02:32,060 --> 00:02:34,550 There will be certain models for which you have 55 00:02:34,550 --> 00:02:36,140 no idea about the dynamics. 56 00:02:36,140 --> 00:02:38,790 I'll give you an example, the study an Ising model, 57 00:02:38,790 --> 00:02:43,910 which is a spin model of magnetic moments sort of being 58 00:02:43,910 --> 00:02:45,840 oriented in space. 59 00:02:45,840 --> 00:02:47,840 In many case you don't know the spin dynamics. 60 00:02:47,840 --> 00:02:49,770 It's a fairly complicated problem. 61 00:02:49,770 --> 00:02:52,340 So if all you want to know is the average magnetic moment, 62 00:02:52,340 --> 00:02:54,500 there's no way of doing dynamics on the spin. 63 00:02:54,500 --> 00:02:56,270 You'll actually have to sample them. 64 00:03:00,460 --> 00:03:04,360 So I want to give you one quick example of a time scale 65 00:03:04,360 --> 00:03:06,430 problem, because I'm going to use 66 00:03:06,430 --> 00:03:10,510 this example throughout the lecture to illustrate 67 00:03:10,510 --> 00:03:13,090 how we use Monte Carlo. 68 00:03:13,090 --> 00:03:19,240 Here's a binary mixture, sort of green and blue atoms. 69 00:03:19,240 --> 00:03:22,270 You may want to know the average energy of that system. 70 00:03:22,270 --> 00:03:24,940 So in any given configuration, let's 71 00:03:24,940 --> 00:03:27,190 say you have an energy method to calculate 72 00:03:27,190 --> 00:03:28,750 the energy of that system. 73 00:03:28,750 --> 00:03:34,120 So at high temperature the atoms will simply 74 00:03:34,120 --> 00:03:36,770 sort of hop around and give you some average energy. 75 00:03:36,770 --> 00:03:39,340 So I want to give you an idea of what you need-- 76 00:03:39,340 --> 00:03:40,910 do you mind passing around-- 77 00:03:40,910 --> 00:03:42,260 Thank you. 78 00:03:42,260 --> 00:03:46,630 Thank you, [INAUDIBLE]---- of what kind of timescale you need 79 00:03:46,630 --> 00:03:48,110 to get that average right. 80 00:03:48,110 --> 00:03:48,610 Thanks. 81 00:03:52,880 --> 00:03:56,060 And you saw, from Professor Marzari's last lecture, 82 00:03:56,060 --> 00:03:57,800 the relation between the diffusion 83 00:03:57,800 --> 00:04:00,380 constant and the root mean square displacement. 84 00:04:00,380 --> 00:04:06,350 Essentially, the root mean square displacement scales 85 00:04:06,350 --> 00:04:08,180 linear with the diffusion constant, 86 00:04:08,180 --> 00:04:10,950 where d is the dimensionality of the system-- 87 00:04:10,950 --> 00:04:12,230 so in this case 2-- 88 00:04:12,230 --> 00:04:14,045 and so D is the diffusion constant. 89 00:04:21,220 --> 00:04:22,510 OK. 90 00:04:22,510 --> 00:04:27,610 To get an idea of the kind of hopping rate you need, 91 00:04:27,610 --> 00:04:31,810 you can assume for a second that the atoms do a random walk. 92 00:04:31,810 --> 00:04:33,550 If you do a random walk, then the root 93 00:04:33,550 --> 00:04:35,320 mean square displacement you take 94 00:04:35,320 --> 00:04:36,850 is the number of jumps you take-- 95 00:04:36,850 --> 00:04:41,260 N times the jump distance squared. 96 00:04:41,260 --> 00:04:43,490 You may remember that for a random walk. 97 00:04:43,490 --> 00:04:47,200 So that means that now we know we 98 00:04:47,200 --> 00:04:53,110 can relate the jump rate, the N dt to the diffusion constant. 99 00:04:53,110 --> 00:04:57,530 And I'll call the N dt, which is the jump rate, gamma. 100 00:04:57,530 --> 00:05:03,900 So if you know how long your simulation can last, 101 00:05:03,900 --> 00:05:06,030 you sort of know what jump you need to see 102 00:05:06,030 --> 00:05:08,230 atoms jump in that simulation. 103 00:05:08,230 --> 00:05:10,050 So you can relate that to-- 104 00:05:10,050 --> 00:05:12,105 excuse me-- the diffusivity you need. 105 00:05:12,105 --> 00:05:14,685 And let me sort of do that on the next page. 106 00:05:22,260 --> 00:05:22,760 OK. 107 00:05:25,310 --> 00:05:32,160 We know that the jump rate is some vibrational frequency 108 00:05:32,160 --> 00:05:34,320 times a success factor. 109 00:05:34,320 --> 00:05:38,106 When atoms have to jump over a barrier, 110 00:05:38,106 --> 00:05:42,060 if Ea is the activation barrier, we essentially 111 00:05:42,060 --> 00:05:45,660 assume that they attempt with some vibrational frequency, 112 00:05:45,660 --> 00:05:47,910 which is, you could say, the frequency with which they 113 00:05:47,910 --> 00:05:50,680 sort of try to go up the hill times the success rate. 114 00:05:50,680 --> 00:05:52,620 And the success rate is that Boltzmann factor, 115 00:05:52,620 --> 00:05:56,130 the exponential minus Ea over kT. 116 00:05:56,130 --> 00:05:59,520 If you assume that the vibrational frequencies 117 00:05:59,520 --> 00:06:02,190 of the order of 10 to the 13 Hertz, which typically we 118 00:06:02,190 --> 00:06:05,310 are, anywhere from 10 to 12, 10 to the 14, 119 00:06:05,310 --> 00:06:08,340 let's say you want to get a jump rate of 10 to the 10 hertz. 120 00:06:08,340 --> 00:06:10,320 Why would we want 10 to the 10 12? 121 00:06:10,320 --> 00:06:14,280 So 10 to the 10, that's one jump or 100 picoseconds. 122 00:06:14,280 --> 00:06:14,875 Is that right? 123 00:06:14,875 --> 00:06:15,540 Yeah. 124 00:06:15,540 --> 00:06:18,210 0.1 jump for 0.1 nanoseconds. 125 00:06:18,210 --> 00:06:21,000 So if you did a 10-nanosecond simulation, 126 00:06:21,000 --> 00:06:23,320 you'd have about 100 jumps per atom. 127 00:06:23,320 --> 00:06:26,310 Which would be reasonable for some equilibration, 128 00:06:26,310 --> 00:06:28,870 but in many cases still not enough. 129 00:06:28,870 --> 00:06:29,370 OK. 130 00:06:29,370 --> 00:06:32,340 This is sort of a reasonable jump rate to expect. 131 00:06:32,340 --> 00:06:35,220 That means that your success rate, your Boltzmann factor, 132 00:06:35,220 --> 00:06:37,150 has to be greater than 10 to the minus 3. 133 00:06:37,150 --> 00:06:41,070 So 1 out of 1,000 times the atom has to jump across the barrier. 134 00:06:41,070 --> 00:06:43,860 So you can deduce from that what the activation energy should 135 00:06:43,860 --> 00:06:44,453 be. 136 00:06:44,453 --> 00:06:45,870 Essentially, the activation energy 137 00:06:45,870 --> 00:06:48,810 should be about less than 7 kT. 138 00:06:48,810 --> 00:06:50,820 And then you'll see jumps of the order 139 00:06:50,820 --> 00:06:53,550 of once per 100 picoseconds. 140 00:06:53,550 --> 00:06:56,830 And just to put these numbers in your head, 141 00:06:56,830 --> 00:06:59,610 because they're sort of important to remember, 142 00:06:59,610 --> 00:07:04,230 at 300 k, 7 kT is 180 millielectronvolt-- 143 00:07:04,230 --> 00:07:06,270 that's a low barrier. 144 00:07:06,270 --> 00:07:08,080 If you simulate at 1,000 degrees, 145 00:07:08,080 --> 00:07:10,290 you can do much higher barriers. 146 00:07:10,290 --> 00:07:13,890 Now, 7 kT is about 600 millielectronvolt. 147 00:07:13,890 --> 00:07:18,270 Now you're starting to be in the range of typical to low 148 00:07:18,270 --> 00:07:20,170 barriers in solids. 149 00:07:20,170 --> 00:07:21,600 Of course, in liquids you'll have 150 00:07:21,600 --> 00:07:23,940 much, much faster transport. 151 00:07:23,940 --> 00:07:25,560 So in liquids you'll have no problem. 152 00:07:25,560 --> 00:07:27,420 If you do the liquid at 1,000 degrees, 153 00:07:27,420 --> 00:07:30,570 you'll have a tremendous amount of hopping. 154 00:07:30,570 --> 00:07:33,630 Solid barriers are anywhere typically from 0.5 to 2 155 00:07:33,630 --> 00:07:37,020 electronvolt. So to put-- 156 00:07:37,020 --> 00:07:40,200 it's kind of important, I think, that you kind of get an idea 157 00:07:40,200 --> 00:07:43,150 of these scales that you need. 158 00:07:43,150 --> 00:07:47,140 And so, like I said in the beginning, 159 00:07:47,140 --> 00:07:51,090 if the average is all you care about, 160 00:07:51,090 --> 00:07:53,890 it may be better to simply do a sampling technique. 161 00:07:53,890 --> 00:07:56,220 And that's essentially what Monte Carlo is. 162 00:07:56,220 --> 00:08:01,350 Monte Carlo is often mistaken for another form of dynamics, 163 00:08:01,350 --> 00:08:06,310 but it is really nothing more than a sampling method. 164 00:08:06,310 --> 00:08:08,130 So now remember what the issue is. 165 00:08:08,130 --> 00:08:11,320 So if you want to get a sample-- 166 00:08:11,320 --> 00:08:13,110 this is in any statistical methods-- 167 00:08:13,110 --> 00:08:15,240 there are two things you have to decide, 168 00:08:15,240 --> 00:08:17,370 which population you sample from-- 169 00:08:17,370 --> 00:08:21,060 how do you set constraints on the population of states 170 00:08:21,060 --> 00:08:22,600 that you sample from-- 171 00:08:22,600 --> 00:08:26,790 and with what probability do you sample them, 172 00:08:26,790 --> 00:08:29,108 or should you sample them unbiased for example. 173 00:08:29,108 --> 00:08:30,900 So there are two issues that always come up 174 00:08:30,900 --> 00:08:33,690 with sampling, the probability or the rate 175 00:08:33,690 --> 00:08:36,419 with which you sample a particular state and the group 176 00:08:36,419 --> 00:08:39,630 of states, the constraints on the states, over which you 177 00:08:39,630 --> 00:08:42,450 sample. 178 00:08:42,450 --> 00:08:47,520 So before I get into that, I want 179 00:08:47,520 --> 00:08:50,370 to briefly review the statistical mechanics 180 00:08:50,370 --> 00:08:52,285 and thermodynamics altogether in 15 minutes. 181 00:08:52,285 --> 00:08:53,910 And I know a lot of you have seen that. 182 00:08:53,910 --> 00:08:56,050 Some of you have not seen that so much. 183 00:08:56,050 --> 00:09:00,000 I sort of want to try to put us all on the same page. 184 00:09:00,000 --> 00:09:02,220 You know, I find when I see Monte Carlo simulation 185 00:09:02,220 --> 00:09:05,160 that people do, the biggest mistake people always 186 00:09:05,160 --> 00:09:07,012 make is at the outset it has nothing 187 00:09:07,012 --> 00:09:09,220 to do with the technical implications of Monte Carlo, 188 00:09:09,220 --> 00:09:11,290 it's defining the boundary conditions, 189 00:09:11,290 --> 00:09:13,890 which is essentially defining the constraints on the states 190 00:09:13,890 --> 00:09:15,240 you sample. 191 00:09:15,240 --> 00:09:17,820 Those constraints come from the thermodynamic boundary 192 00:09:17,820 --> 00:09:21,930 conditions you impose on your system. 193 00:09:21,930 --> 00:09:24,870 I have sometimes really funny stories. 194 00:09:24,870 --> 00:09:29,190 I've seen discussions at meetings where people argue 195 00:09:29,190 --> 00:09:31,920 strongly back and forth about problems, 196 00:09:31,920 --> 00:09:34,890 and all it is a misdefining of their thermodynamic boundary 197 00:09:34,890 --> 00:09:36,205 conditions. 198 00:09:36,205 --> 00:09:37,830 About two years ago, I was at a meeting 199 00:09:37,830 --> 00:09:40,290 where people had enormous discussions about where 200 00:09:40,290 --> 00:09:42,370 hydrogen sits on palladium. 201 00:09:42,370 --> 00:09:44,700 Palladium is a pretty good hydrogen absorber, 202 00:09:44,700 --> 00:09:47,560 and people would argue back and forth. 203 00:09:47,560 --> 00:09:48,570 It sits on the surface-- 204 00:09:48,570 --> 00:09:51,660 and somebody says, no, it's at subsurface, layer subsurface. 205 00:09:51,660 --> 00:09:53,860 All that is-- it actually can sit both, 206 00:09:53,860 --> 00:09:55,590 it's just a matter of thermodynamics. 207 00:09:55,590 --> 00:09:58,650 If your hydrogen chemical potential, your external field, 208 00:09:58,650 --> 00:10:00,510 is high enough, it sits everywhere. 209 00:10:00,510 --> 00:10:02,040 Period. 210 00:10:02,040 --> 00:10:04,470 And so thinking about your thermodynamic boundary 211 00:10:04,470 --> 00:10:09,370 conditions before you make any conclusions is quite important. 212 00:10:09,370 --> 00:10:14,340 So for those of you-- 213 00:10:14,340 --> 00:10:16,370 let me just remind you of the relation 214 00:10:16,370 --> 00:10:18,770 between the thermodynamic boundary conditions 215 00:10:18,770 --> 00:10:20,270 and the statistical ensemble. 216 00:10:20,270 --> 00:10:21,895 The statistical ensemble is essentially 217 00:10:21,895 --> 00:10:24,630 the group of states over which you sample. 218 00:10:24,630 --> 00:10:27,020 So we're always going to make the relation 219 00:10:27,020 --> 00:10:30,050 between a microscopic ensemble-- 220 00:10:30,050 --> 00:10:32,840 so this is the group of states for which you sample. 221 00:10:35,780 --> 00:10:38,600 You could think of it, it's the boundary conditions 222 00:10:38,600 --> 00:10:41,450 on the states over which you sample-- 223 00:10:41,450 --> 00:10:43,460 should you sample any configuration, should you 224 00:10:43,460 --> 00:10:45,560 constrain them-- 225 00:10:45,560 --> 00:10:48,380 and the macroscopic thermodynamics. 226 00:10:48,380 --> 00:10:52,160 And the relation here is made between statistical mechanics. 227 00:10:52,160 --> 00:10:53,810 In the end, it's very easy. 228 00:10:53,810 --> 00:10:56,220 Any fixed thermodynamic variable-- 229 00:10:56,220 --> 00:10:58,220 I'll come to that in the second, what they are-- 230 00:10:58,220 --> 00:11:01,730 sets a constraint on how you can sample. 231 00:11:01,730 --> 00:11:03,470 And more particularly, the constraints 232 00:11:03,470 --> 00:11:05,930 come from what are called the extensive thermodynamic 233 00:11:05,930 --> 00:11:06,503 variables. 234 00:11:06,503 --> 00:11:08,420 So remember, extensive thermodynamic variables 235 00:11:08,420 --> 00:11:13,100 are things at scale with size, things like volume, 236 00:11:13,100 --> 00:11:15,200 number of atoms of a given type. 237 00:11:15,200 --> 00:11:17,390 And intensive variables are the conjugates. 238 00:11:17,390 --> 00:11:19,910 Those are the ones that don't scale that size-- temperature, 239 00:11:19,910 --> 00:11:23,840 pressure, chemical potential. 240 00:11:23,840 --> 00:11:26,630 You'll see in a second that the probability with which you 241 00:11:26,630 --> 00:11:29,360 sample depends on the relevant Hamiltonian, 242 00:11:29,360 --> 00:11:33,050 and the macroscopic pound counterpart of that 243 00:11:33,050 --> 00:11:34,470 is a free energy function. 244 00:11:34,470 --> 00:11:36,860 So a Hamiltonian in microscopic space 245 00:11:36,860 --> 00:11:41,390 will correspond to a free energy function in macroscopic space. 246 00:11:44,110 --> 00:11:45,220 OK. 247 00:11:45,220 --> 00:11:48,550 So here's all I'm going to say about thermal is sort of one 248 00:11:48,550 --> 00:11:50,830 slide that I think represents pretty much all you need 249 00:11:50,830 --> 00:11:55,570 to know for doing sampling. 250 00:11:55,570 --> 00:11:58,270 In thermodynamics variables appear in groups, 251 00:11:58,270 --> 00:12:00,020 what we call conjugate pairs. 252 00:12:00,020 --> 00:12:01,880 And one way you can see how they belong, 253 00:12:01,880 --> 00:12:08,000 it's the couples that appear together in the first law. 254 00:12:08,000 --> 00:12:09,400 First law tells you, essentially, 255 00:12:09,400 --> 00:12:11,470 that the change of internal energy of a system 256 00:12:11,470 --> 00:12:15,970 is the net energy flow in or out of that system. 257 00:12:15,970 --> 00:12:18,550 And you can write that from different contributions. 258 00:12:18,550 --> 00:12:20,500 TdS you could think of as the heat flow, 259 00:12:20,500 --> 00:12:24,040 pdV is the mechanical or volume work done, 260 00:12:24,040 --> 00:12:25,983 mu dN as the chemical work. 261 00:12:25,983 --> 00:12:27,400 So if u is the chemical potential, 262 00:12:27,400 --> 00:12:30,460 N is the number of atoms of a given type. 263 00:12:30,460 --> 00:12:32,800 So remember, these are conjugate pairs. 264 00:12:32,800 --> 00:12:37,390 S is an extensive, T is an intensive, V is an extensive, 265 00:12:37,390 --> 00:12:43,318 p is an intensive, N is an extensive, mu is an intensive. 266 00:12:43,318 --> 00:12:44,860 Essentially, the rule is very simple. 267 00:12:44,860 --> 00:12:46,360 To define your boundary conditions, 268 00:12:46,360 --> 00:12:49,465 you always need to specify one out of each pair. 269 00:12:49,465 --> 00:12:50,840 So you need to say whether you're 270 00:12:50,840 --> 00:12:52,820 at constant entropy or constant temperature, 271 00:12:52,820 --> 00:12:55,295 whether you're at constant pressure or constant volume, 272 00:12:55,295 --> 00:12:57,170 whether you're at constant chemical potential 273 00:12:57,170 --> 00:12:59,780 or constant number of particles. 274 00:12:59,780 --> 00:13:04,205 You could not leave one unspecified. 275 00:13:04,205 --> 00:13:06,080 If you do that, you're actually always making 276 00:13:06,080 --> 00:13:07,730 an implicit assumption. 277 00:13:07,730 --> 00:13:08,600 OK. 278 00:13:08,600 --> 00:13:12,230 If you have a closed box and nothing can go in and out, 279 00:13:12,230 --> 00:13:16,880 you're obviously under constant number of particles. 280 00:13:16,880 --> 00:13:19,100 A really classic mistake is simulating 281 00:13:19,100 --> 00:13:20,810 on the constant number of particles, 282 00:13:20,810 --> 00:13:24,470 when you're really, in reality, under constant chemical field. 283 00:13:24,470 --> 00:13:27,450 Surfaces are a really good example of that. 284 00:13:27,450 --> 00:13:33,302 If you study a surface, you will often study it in simulation 285 00:13:33,302 --> 00:13:34,760 under constant number of particles. 286 00:13:34,760 --> 00:13:36,540 It's sort of the obvious thing to do. 287 00:13:36,540 --> 00:13:38,870 But real surfaces are actually equilibrated 288 00:13:38,870 --> 00:13:40,440 with their environment. 289 00:13:40,440 --> 00:13:44,060 So they're under constant chemical potential field 290 00:13:44,060 --> 00:13:46,940 of the environment or under constant chemical potential 291 00:13:46,940 --> 00:13:49,550 from their bulk. 292 00:13:49,550 --> 00:13:51,900 And in some cases that doesn't matter. 293 00:13:51,900 --> 00:13:54,290 But in other cases that actually will. 294 00:13:54,290 --> 00:13:56,660 So the number of atoms on a surface can change. 295 00:13:59,930 --> 00:14:02,660 In statistical mechanics, we tend to not use this 296 00:14:02,660 --> 00:14:06,950 formulation of the first law, but we tend interchange 297 00:14:06,950 --> 00:14:09,680 the role of the internal energy and the entropy literally 298 00:14:09,680 --> 00:14:13,250 by moving this to the left-hand side and dividing by T-- 299 00:14:13,250 --> 00:14:16,077 so riding the differential of the entropy instead. 300 00:14:16,077 --> 00:14:17,660 And you'll see in a second that that's 301 00:14:17,660 --> 00:14:19,790 a much more useful thermodynamic notation 302 00:14:19,790 --> 00:14:23,130 to correlate to statistical mechanics. 303 00:14:23,130 --> 00:14:25,230 So this is called the entropy formulation. 304 00:14:25,230 --> 00:14:28,220 And in the entropy formulation, U and 1 over T 305 00:14:28,220 --> 00:14:29,330 are conjugate variables. 306 00:14:29,330 --> 00:14:32,210 U is the extensive, 1 over T is the intensive, 307 00:14:32,210 --> 00:14:36,560 V is the the extensive, p over T is the extensive. 308 00:14:36,560 --> 00:14:39,050 I hope you can see this is just a rewrite 309 00:14:39,050 --> 00:14:40,490 of the standard first law. 310 00:14:47,150 --> 00:14:47,650 OK. 311 00:14:51,192 --> 00:14:53,150 So how do the thermodynamic boundary conditions 312 00:14:53,150 --> 00:14:56,660 relate to the ensemble's over which you average? 313 00:14:56,660 --> 00:15:01,610 The ensembles are determined by the extensive variables 314 00:15:01,610 --> 00:15:03,240 you keep constant. 315 00:15:03,240 --> 00:15:06,170 So if you think about, if you go back for a second 316 00:15:06,170 --> 00:15:08,540 the previous slide, the simplest is 317 00:15:08,540 --> 00:15:12,142 that you keep the energy, the volume, the number 318 00:15:12,142 --> 00:15:13,100 of particles, constant. 319 00:15:13,100 --> 00:15:15,050 That's all the extensive variables. 320 00:15:17,800 --> 00:15:20,110 All the extensive variables, that's 321 00:15:20,110 --> 00:15:23,560 the ensemble at constant energy, volume, number of particles. 322 00:15:23,560 --> 00:15:26,470 That's what's called the microcanonical ensemble. 323 00:15:26,470 --> 00:15:31,450 A realization of that would be if you did a system 324 00:15:31,450 --> 00:15:35,330 in a box, Newtonian dynamics. 325 00:15:35,330 --> 00:15:37,330 If it's in a closed box, the number of particles 326 00:15:37,330 --> 00:15:39,640 can't change, the volume is fixed. 327 00:15:39,640 --> 00:15:43,660 And if you do Newtonian dynamics, the energy is fixed. 328 00:15:43,660 --> 00:15:47,890 So the basic form of molecular dynamics, you could say, 329 00:15:47,890 --> 00:15:51,050 goes through a microcanonical ensemble. 330 00:15:51,050 --> 00:15:55,400 If you start thermostating your energy to fix the temperature, 331 00:15:55,400 --> 00:15:59,450 then you go out of the microcanonical ensemble. 332 00:15:59,450 --> 00:16:05,540 So if you release the constraint of fixed energy, 333 00:16:05,540 --> 00:16:07,610 then you have to control the average energy. 334 00:16:07,610 --> 00:16:09,230 And you do that with the temperature. 335 00:16:09,230 --> 00:16:14,140 Then you end up in what's called a canonical ensemble. 336 00:16:14,140 --> 00:16:17,860 If you release the constraint of fixed number of particles, 337 00:16:17,860 --> 00:16:20,650 again, then you have to take into account 338 00:16:20,650 --> 00:16:23,870 the field that sets your average number of particles. 339 00:16:23,870 --> 00:16:27,940 So then you have an ensemble with T, V, and mu. 340 00:16:27,940 --> 00:16:30,830 That's called a grand canonical ensemble. 341 00:16:30,830 --> 00:16:34,540 You can keep on making other ensembles by switching 342 00:16:34,540 --> 00:16:36,812 in intensive variables for extensive, 343 00:16:36,812 --> 00:16:38,270 it's just that we run out of names. 344 00:16:38,270 --> 00:16:41,230 These are three classic names-- microcanonical, canonical, 345 00:16:41,230 --> 00:16:42,220 and grand canonical. 346 00:16:42,220 --> 00:16:45,370 After that, we sort of have no names for them anymore, 347 00:16:45,370 --> 00:16:46,270 but they do exist. 348 00:16:52,970 --> 00:16:54,830 OK. 349 00:16:54,830 --> 00:16:58,760 So that defines the states from which you should sample. 350 00:16:58,760 --> 00:17:01,857 The question is, now, how do you sample? 351 00:17:01,857 --> 00:17:04,190 Well, first I'm going to tell you how you should sample, 352 00:17:04,190 --> 00:17:08,150 and then I'll tell you how you really sample in practice. 353 00:17:08,150 --> 00:17:10,380 You should sample with the correct probability. 354 00:17:10,380 --> 00:17:13,849 And the probability of a state in an ensemble 355 00:17:13,849 --> 00:17:16,460 is proportional to the exponential of minus beta, 356 00:17:16,460 --> 00:17:18,589 the Hamiltonian-- 357 00:17:18,589 --> 00:17:21,780 again, remember, beta is 1 over kT-- 358 00:17:21,780 --> 00:17:23,869 and then normalized by the sum of that 359 00:17:23,869 --> 00:17:26,369 over all the states, which is called the partition function, 360 00:17:26,369 --> 00:17:26,869 PrZ. 361 00:17:31,550 --> 00:17:34,100 And the free energy is essentially the log 362 00:17:34,100 --> 00:17:35,630 of that partition function. 363 00:17:35,630 --> 00:17:37,520 Now, what should that Hamiltonian be? 364 00:17:37,520 --> 00:17:40,850 You've probably seen this kind of probability function a lot 365 00:17:40,850 --> 00:17:43,760 already when you substitute in the energy. 366 00:17:43,760 --> 00:17:45,200 In more general terms, that should 367 00:17:45,200 --> 00:17:49,460 be the relevant microscopic Hamiltonian, 368 00:17:49,460 --> 00:17:51,920 and it should include all the things that 369 00:17:51,920 --> 00:17:55,700 can fluctuate in your system. 370 00:17:55,700 --> 00:17:58,490 Typically, the way I sort of remember 371 00:17:58,490 --> 00:18:01,460 what goes into Hamiltonian, you think 372 00:18:01,460 --> 00:18:03,770 of your thermodynamic boundary conditions, 373 00:18:03,770 --> 00:18:07,040 you take the relevant Legendre transform of your entropy, 374 00:18:07,040 --> 00:18:08,630 and it's essentially those pieces 375 00:18:08,630 --> 00:18:10,280 that go in your Hamiltonian. 376 00:18:10,280 --> 00:18:15,620 For example, if you work at fixed N, V, and T, 377 00:18:15,620 --> 00:18:19,400 the relevant free energy is the Helmholtz free energy, 378 00:18:19,400 --> 00:18:23,650 which is S minus E over T. This is the part that 379 00:18:23,650 --> 00:18:25,405 goes into your Hamiltonian. 380 00:18:25,405 --> 00:18:29,740 It's whatever is the Legendre transform part of the entropy. 381 00:18:29,740 --> 00:18:31,630 And that's what you see. 382 00:18:31,630 --> 00:18:35,680 In a canonical ensemble, you weigh states 383 00:18:35,680 --> 00:18:39,070 with the exponential of minus beta the energy. 384 00:18:39,070 --> 00:18:40,790 That's the one you're used to. 385 00:18:40,790 --> 00:18:44,705 But let's say you go to a grand canonical ensemble. 386 00:18:47,630 --> 00:18:50,700 So now you're at fixed chemical potential 387 00:18:50,700 --> 00:18:52,907 so the number of states can fluctuate. 388 00:18:52,907 --> 00:18:54,740 So now the relevant thermodynamic potential, 389 00:18:54,740 --> 00:18:57,020 the relevant Legendre transform of the entropy 390 00:18:57,020 --> 00:19:00,730 has a minus E over T and a plus mu over T in it. 391 00:19:00,730 --> 00:19:03,140 So that's essentially the part that 392 00:19:03,140 --> 00:19:05,490 appears in your Hamiltonian. 393 00:19:05,490 --> 00:19:09,290 So now you weigh with minus beta times 394 00:19:09,290 --> 00:19:14,090 the energy minus mu N. Just sort of think of that as a Legendre 395 00:19:14,090 --> 00:19:17,580 transform Hamiltonian almost. 396 00:19:17,580 --> 00:19:20,960 And if you were to simulate that constant pressure , 397 00:19:20,960 --> 00:19:27,420 you would have a P over T times V term in your Hamiltonian. 398 00:19:27,420 --> 00:19:29,170 If you sort think about it, it makes sense 399 00:19:29,170 --> 00:19:31,958 that you need these terms in your Hamiltonian. 400 00:19:31,958 --> 00:19:34,000 Because if they don't appear in your Hamiltonian, 401 00:19:34,000 --> 00:19:38,350 you have no way of controlling the number of particles. 402 00:19:38,350 --> 00:19:38,850 OK. 403 00:19:44,960 --> 00:19:48,940 It's sort of going backwards. 404 00:19:48,940 --> 00:19:50,440 So this is here sort of the summary. 405 00:19:53,140 --> 00:19:55,540 The macroscopic boundary conditions 406 00:19:55,540 --> 00:19:59,500 are set by the constant extensive thermodynamic 407 00:19:59,500 --> 00:20:01,390 variables. 408 00:20:01,390 --> 00:20:05,350 That in turn defines the relevant ensemble 409 00:20:05,350 --> 00:20:09,382 and the probability which you sample over an ensemble. 410 00:20:09,382 --> 00:20:10,840 And when you do that all correctly, 411 00:20:10,840 --> 00:20:13,810 you can get the averages of properties, 412 00:20:13,810 --> 00:20:17,930 such as the energy and the volume 413 00:20:17,930 --> 00:20:20,280 with those probabilities. 414 00:20:20,280 --> 00:20:20,780 OK. 415 00:20:24,350 --> 00:20:26,570 Just one brief aside-- 416 00:20:26,570 --> 00:20:28,100 there's an implicit assumption you 417 00:20:28,100 --> 00:20:33,290 make in here, which is common to statistical mechanics. 418 00:20:33,290 --> 00:20:37,310 We've essentially said that following a system in time 419 00:20:37,310 --> 00:20:43,250 will give the same averages then sampling over its ensemble. 420 00:20:43,250 --> 00:20:45,890 The assumption in there is that when 421 00:20:45,890 --> 00:20:49,490 you follow a system in time, it will reach essentially 422 00:20:49,490 --> 00:20:52,250 all the states in its ensemble. 423 00:20:52,250 --> 00:20:56,770 Or, maybe simpler said, it will reach all the states it can-- 424 00:20:56,770 --> 00:20:59,870 it will reach all the states that are within the boundary 425 00:20:59,870 --> 00:21:01,760 conditions. 426 00:21:01,760 --> 00:21:03,440 So if you put-- 427 00:21:03,440 --> 00:21:04,442 think about this. 428 00:21:04,442 --> 00:21:05,900 What you're really saying is if you 429 00:21:05,900 --> 00:21:10,250 put a system in a box, fixed number of particle, 430 00:21:10,250 --> 00:21:14,720 say, at constant temperature, not constant energy, 431 00:21:14,720 --> 00:21:16,880 you're essentially saying that that system will, 432 00:21:16,880 --> 00:21:20,240 if you wait long enough, reach any possible configuration 433 00:21:20,240 --> 00:21:21,770 of atoms. 434 00:21:21,770 --> 00:21:25,870 Because that's the one that will be allowed in your ensemble. 435 00:21:25,870 --> 00:21:27,993 And that's called the ergodicity principle. 436 00:21:27,993 --> 00:21:30,410 The ergodicity principle says that if you wait long enough 437 00:21:30,410 --> 00:21:33,612 a system will reach, finally, all its states . 438 00:21:33,612 --> 00:21:36,070 Because if you know that it's going to reach all its states 439 00:21:36,070 --> 00:21:37,670 and you know the probability of them, 440 00:21:37,670 --> 00:21:38,920 then you can just sum over it. 441 00:21:38,920 --> 00:21:41,200 And that's what you do in statistical mechanics, 442 00:21:41,200 --> 00:21:44,710 rather than doing the full dynamics. 443 00:21:44,710 --> 00:21:47,500 This is where we usually have to give you the caveat-- 444 00:21:47,500 --> 00:21:51,880 the lawyer's caveat-- that not all systems are ergodic. 445 00:21:51,880 --> 00:21:53,800 Most practical ones are. 446 00:21:53,800 --> 00:21:56,230 The examples that are often given as non-ergodic system 447 00:21:56,230 --> 00:21:57,310 are a bit fake. 448 00:21:57,310 --> 00:22:01,840 And one that I like to bring up is the harmonic oscillator. 449 00:22:01,840 --> 00:22:04,600 If you think of a harmonic oscillator, 450 00:22:04,600 --> 00:22:11,100 think of something, a mass, that is 451 00:22:11,100 --> 00:22:15,000 in contact with a rigid object to a spring. 452 00:22:15,000 --> 00:22:17,460 If you think that the two coordinates-- 453 00:22:17,460 --> 00:22:20,830 what are the two coordinates of this system? 454 00:22:20,830 --> 00:22:24,480 Well, there's a position coordinate, where this mass is, 455 00:22:24,480 --> 00:22:27,000 and there's a velocity coordinate. 456 00:22:27,000 --> 00:22:29,070 Well, because it's a harmonic oscillator, 457 00:22:29,070 --> 00:22:31,570 the two are actually correlated. 458 00:22:31,570 --> 00:22:34,260 So if you think of the phase pace of the system as-- 459 00:22:34,260 --> 00:22:39,970 here's the spatial coordinate, here's the velocity coordinate. 460 00:22:39,970 --> 00:22:42,160 The system actually goes through an ellipse. 461 00:22:42,160 --> 00:22:46,330 It actually never goes to states like this and this and this. 462 00:22:46,330 --> 00:22:49,870 So people say the harmonic oscillator is not ergodic. 463 00:22:49,870 --> 00:22:54,500 I'd say it's not ergodic because you took a stupid phase space. 464 00:22:54,500 --> 00:22:58,850 It's because you essentially said that the system has 465 00:22:58,850 --> 00:23:01,415 to coordinates. 466 00:23:01,415 --> 00:23:03,415 So in a two-dimensional space, it's not ergodic. 467 00:23:03,415 --> 00:23:06,540 In a one-dimensional space, it is ergodic. 468 00:23:06,540 --> 00:23:08,490 This system has only one coordinate. 469 00:23:08,490 --> 00:23:11,570 That's its normal mode. 470 00:23:11,570 --> 00:23:14,690 So if you write it in the normal mode, it is ergodic. 471 00:23:14,690 --> 00:23:17,190 In one dimension, it is ergodic. 472 00:23:17,190 --> 00:23:19,880 So along that trajectory, if you say that's your phase space, 473 00:23:19,880 --> 00:23:21,990 the system is ergodic. 474 00:23:21,990 --> 00:23:23,960 So it's a bit of a matter of definition. 475 00:23:23,960 --> 00:23:27,260 In the obvious coordinate space, this one is not ergodic. 476 00:23:31,110 --> 00:23:38,070 And then there are systems that for practical reasons 477 00:23:38,070 --> 00:23:39,720 are essentially non-ergodic. 478 00:23:39,720 --> 00:23:44,520 They just never get there going through all their phase spaces. 479 00:23:44,520 --> 00:23:51,090 Glasses are really difficult to do thermodynamics on, 480 00:23:51,090 --> 00:23:53,400 to integrate phase space over, because you really 481 00:23:53,400 --> 00:23:56,460 don't know over what region to integrate. 482 00:23:56,460 --> 00:23:59,800 If I give you atoms-- 483 00:23:59,800 --> 00:24:02,820 I say, give me the free energy of a glass. 484 00:24:02,820 --> 00:24:04,320 Really, when you think about it, you 485 00:24:04,320 --> 00:24:06,750 don't know what ensemble to integrate over. 486 00:24:06,750 --> 00:24:08,550 Because if you say, well, what is a glass? 487 00:24:08,550 --> 00:24:10,170 Well, it's sort of random. 488 00:24:10,170 --> 00:24:16,350 So can I just average over all the possible positions 489 00:24:16,350 --> 00:24:17,920 in the box? 490 00:24:17,920 --> 00:24:21,680 Then I'm actually getting the free energy of a liquid. 491 00:24:21,680 --> 00:24:24,620 The glass actually never samples all the possible positions. 492 00:24:24,620 --> 00:24:27,760 So you could argue that it's a liquid that's not ergodic. 493 00:24:27,760 --> 00:24:29,260 And if you sort of think of it, it's 494 00:24:29,260 --> 00:24:35,360 a liquid that's frozen [? in. ?] 495 00:24:35,360 --> 00:24:35,930 OK. 496 00:24:35,930 --> 00:24:38,430 So finally we're going to get to the Monte Carlo simulation. 497 00:24:38,430 --> 00:24:41,090 Before that, I want to introduce one toy model. 498 00:24:41,090 --> 00:24:43,940 And you'll see the relation with the mixing problem. 499 00:24:43,940 --> 00:24:46,680 And that's the Ising model, which 500 00:24:46,680 --> 00:24:49,510 is essentially a spin model. 501 00:24:49,510 --> 00:24:51,660 So we have lots of sites. 502 00:24:51,660 --> 00:24:55,380 At every site we can have an up spin or a down spin. 503 00:24:55,380 --> 00:24:59,190 So there's only a binary variable at every site. 504 00:24:59,190 --> 00:25:03,660 We'll mark that with a plus 1 or minus 1-- say plus 1 505 00:25:03,660 --> 00:25:06,270 for up and minus 1 for down. 506 00:25:06,270 --> 00:25:09,090 And so the simplest Hamiltonian or energy function 507 00:25:09,090 --> 00:25:11,550 you can write for this system is one 508 00:25:11,550 --> 00:25:15,430 that counts the nearest neighbor bonds and associates 509 00:25:15,430 --> 00:25:17,920 some interaction with them. 510 00:25:17,920 --> 00:25:23,000 So I and J are only nearest neighbors in this model. 511 00:25:23,000 --> 00:25:29,080 So let's say if J is positive, then, 512 00:25:29,080 --> 00:25:30,580 since there's a minus sign here, you 513 00:25:30,580 --> 00:25:33,390 will want to make this positive so you'll get a ferromagnet. 514 00:25:33,390 --> 00:25:35,140 That means if you have a plus 1 somewhere, 515 00:25:35,140 --> 00:25:39,310 you'll want a plus 1 next to it to get flow energy. 516 00:25:39,310 --> 00:25:44,490 If J is negative, you will want the spin product sigma I sigma 517 00:25:44,490 --> 00:25:46,680 J to be negative. 518 00:25:46,680 --> 00:25:48,540 That means if you you have a plus spin, 519 00:25:48,540 --> 00:25:50,680 you'll want a negative spin next to it. 520 00:25:50,680 --> 00:25:51,990 So you have antiferromagnet. 521 00:25:51,990 --> 00:25:55,650 So it's the simplest possible magnetic model you can make. 522 00:25:55,650 --> 00:25:57,900 And that was Ising's PhD thesis. 523 00:25:57,900 --> 00:26:01,020 And his advisors really hated it. 524 00:26:01,020 --> 00:26:03,300 His advisor was so mad at him for that PhD thesis 525 00:26:03,300 --> 00:26:06,750 that Ising left science and he became a businessman. 526 00:26:06,750 --> 00:26:09,700 So he never knew he became famous later on. 527 00:26:12,480 --> 00:26:14,760 But what you're going to see is that the Ising model 528 00:26:14,760 --> 00:26:16,020 is sort of a conduit. 529 00:26:16,020 --> 00:26:19,060 It's a model I will use for any two-state system. 530 00:26:19,060 --> 00:26:21,330 So any time you have a model where 531 00:26:21,330 --> 00:26:23,640 there's sort of two possibilities on the lattice 532 00:26:23,640 --> 00:26:25,020 side, we'll use this model. 533 00:26:25,020 --> 00:26:26,850 Because I call it a magnetic spin. 534 00:26:26,850 --> 00:26:30,850 I mean, you could call it low fat, 535 00:26:30,850 --> 00:26:33,880 high fat cheese, whatever-- any time you have two states. 536 00:26:33,880 --> 00:26:36,550 So we're going to use it for binary alloys, for example. 537 00:26:36,550 --> 00:26:38,880 If we have a mixture of A and B on a lattice, 538 00:26:38,880 --> 00:26:43,800 you could say the plus 1 is A, the minus 1 is B. 539 00:26:43,800 --> 00:26:46,020 So that's why it is such a useful model. 540 00:26:49,900 --> 00:26:51,940 OK. 541 00:26:51,940 --> 00:26:53,693 So how would you sample a model like that? 542 00:26:53,693 --> 00:26:55,360 Well, let me first tell you a little bit 543 00:26:55,360 --> 00:26:59,170 about where the ideas of sampling and Monte Carlo 544 00:26:59,170 --> 00:27:02,320 sampling came from. 545 00:27:02,320 --> 00:27:04,960 It's typically credited to the people 546 00:27:04,960 --> 00:27:09,460 at Los Alamos, Monte Carlo samples, who built the ENIAC 547 00:27:09,460 --> 00:27:11,830 computer, which is essentially, I think, 548 00:27:11,830 --> 00:27:14,540 the first computer or one of the first computers that was built. 549 00:27:14,540 --> 00:27:16,123 One of these things with vacuum tubes, 550 00:27:16,123 --> 00:27:18,250 and you could calculate for 10 seconds 551 00:27:18,250 --> 00:27:20,020 and then one of the tubes would break. 552 00:27:20,020 --> 00:27:22,353 And then you'd replace the tube and calculate another 10 553 00:27:22,353 --> 00:27:23,770 seconds. 554 00:27:23,770 --> 00:27:27,940 But I believe idea of Monte Carlo sampling really 555 00:27:27,940 --> 00:27:29,920 goes back to Fermi. 556 00:27:29,920 --> 00:27:32,890 People before that sort of had sampling ideas already. 557 00:27:32,890 --> 00:27:35,140 They sort of knew how to use sampling for things that 558 00:27:35,140 --> 00:27:37,210 were difficult to integrate. 559 00:27:37,210 --> 00:27:41,290 And the first example that I thought is known is 560 00:27:41,290 --> 00:27:44,080 by the Comte de Buffon-- which I think is a beautiful name., 561 00:27:44,080 --> 00:27:46,600 buffon is like-- 562 00:27:46,600 --> 00:27:50,830 who thought of this as a method for integrating functions. 563 00:27:50,830 --> 00:27:54,160 And the idea is very simple, that if you 564 00:27:54,160 --> 00:27:57,280 have this function that you can't integrate, 565 00:27:57,280 --> 00:27:59,380 or maybe it's numerically defined, 566 00:27:59,380 --> 00:28:05,030 a really simple way to integrate is to draw a box 567 00:28:05,030 --> 00:28:07,250 and throw darts in that box. 568 00:28:07,250 --> 00:28:09,417 Randomly take points in that box, 569 00:28:09,417 --> 00:28:11,750 and you'll hit here once in a while, and here, and here, 570 00:28:11,750 --> 00:28:12,570 and here. 571 00:28:12,570 --> 00:28:17,180 And if you think about it, if you track the number of times 572 00:28:17,180 --> 00:28:19,490 that you're below the curve, that's 573 00:28:19,490 --> 00:28:23,330 telling you something about that integral under the curve. 574 00:28:23,330 --> 00:28:36,980 Because the number of times below the curve is the integral 575 00:28:36,980 --> 00:28:38,390 divided by the total area. 576 00:28:41,120 --> 00:28:45,970 So the area of the square is the total area, 577 00:28:45,970 --> 00:28:49,190 and the integral is a fraction of that. 578 00:28:49,190 --> 00:28:51,430 And so that fraction is the probability 579 00:28:51,430 --> 00:28:52,820 that you hit below the curve. 580 00:28:52,820 --> 00:28:55,720 So by just knowing how many times you're below the curve, 581 00:28:55,720 --> 00:28:58,343 you have a reasonable estimate of the integral. 582 00:28:58,343 --> 00:29:00,010 And if you do this with a lot of points, 583 00:29:00,010 --> 00:29:04,522 you'll get a good approximation to the integral. 584 00:29:04,522 --> 00:29:07,210 This is how kids in school learn to integrate. 585 00:29:07,210 --> 00:29:09,130 My daughter's in third grade, and the way 586 00:29:09,130 --> 00:29:12,580 they learn to integrate is by putting grids on this 587 00:29:12,580 --> 00:29:15,770 and counting how many grid points are below the curve. 588 00:29:15,770 --> 00:29:19,060 Which you could say is a sampling method-- it's just 589 00:29:19,060 --> 00:29:21,220 not a random sampling methods. 590 00:29:21,220 --> 00:29:24,640 You use sample with a fixed grid rather than with a random grid. 591 00:29:27,280 --> 00:29:29,340 So I was going to give you an assignment. 592 00:29:33,120 --> 00:29:35,930 This is one of these assignments for when you're 593 00:29:35,930 --> 00:29:37,260 sitting in a boring lecture-- 594 00:29:37,260 --> 00:29:40,100 which is not this one, but another boring lecture-- 595 00:29:40,100 --> 00:29:42,800 or you're sitting on an airplane and they hold you 596 00:29:42,800 --> 00:29:46,760 for two hours on the tarmac with nowhere to go and no food. 597 00:29:46,760 --> 00:29:49,880 You should try to determine pi. 598 00:29:49,880 --> 00:29:52,700 And here's a way of determining pi. 599 00:29:52,700 --> 00:29:55,520 You know that pi essentially relates 600 00:29:55,520 --> 00:30:00,800 to the ratio of the area of a circle divided 601 00:30:00,800 --> 00:30:06,200 by the area of the square that's circumscribed. 602 00:30:06,200 --> 00:30:08,360 Because this area here, if you just 603 00:30:08,360 --> 00:30:16,060 take, say, a quadrant of this, is 1 over 4 pi r squared. 604 00:30:16,060 --> 00:30:19,040 And the area of the square is r squared. 605 00:30:19,040 --> 00:30:24,750 So the ratio of those two is 1/4 times pi. 606 00:30:24,750 --> 00:30:30,650 So if you throw darts randomly at this thing, if you count 607 00:30:30,650 --> 00:30:33,500 the fraction of the times you're within the circle divided 608 00:30:33,500 --> 00:30:38,070 by the total number of points, you'll get an estimate for pi. 609 00:30:38,070 --> 00:30:41,042 So you should try this, like throw darts. 610 00:30:41,042 --> 00:30:43,500 You can't bring darts anymore in an airplane unfortunately, 611 00:30:43,500 --> 00:30:46,500 so you'll have to kind of like hold your pen randomly. 612 00:30:46,500 --> 00:30:49,950 And you should see once how quickly or slowly you 613 00:30:49,950 --> 00:30:51,450 actually approach pi. 614 00:30:51,450 --> 00:30:53,490 You get the first few digits pretty quickly 615 00:30:53,490 --> 00:30:54,640 actually-- not so bad. 616 00:30:54,640 --> 00:30:57,120 But then, after that, of course, it goes really slowly. 617 00:30:57,120 --> 00:30:59,970 Because, essentially, to get the third digit, 618 00:30:59,970 --> 00:31:03,930 you're going to have to almost do 1,000 points. 619 00:31:03,930 --> 00:31:06,710 But when you're bored, it's not a bad thing to do. 620 00:31:11,830 --> 00:31:16,420 So let's talk more about sampling. 621 00:31:16,420 --> 00:31:19,090 There are essentially many ways of sampling, 622 00:31:19,090 --> 00:31:21,350 but I'm going to give you two extremes. 623 00:31:21,350 --> 00:31:22,610 One is called simple sampling. 624 00:31:22,610 --> 00:31:24,940 One is called importance sampling. 625 00:31:24,940 --> 00:31:26,620 Later, in one of the later lectures, 626 00:31:26,620 --> 00:31:28,953 I'll sort of show you that you could kind of do anything 627 00:31:28,953 --> 00:31:31,490 in between as well. 628 00:31:31,490 --> 00:31:34,480 So what is simple sampling? 629 00:31:34,480 --> 00:31:36,940 It's kind of what these previous methods were. 630 00:31:36,940 --> 00:31:42,160 You would randomly pick points in your phase space. 631 00:31:42,160 --> 00:31:47,040 So if you randomly pick points, that means now states of atoms 632 00:31:47,040 --> 00:31:48,790 or whatever the states of your system are, 633 00:31:48,790 --> 00:31:52,843 you have to weigh them with the correct probability. 634 00:31:52,843 --> 00:31:53,760 But you could do that. 635 00:31:53,760 --> 00:31:56,580 So you could randomly pick states. 636 00:31:56,580 --> 00:31:59,370 Let's say the state is new here. 637 00:31:59,370 --> 00:32:02,370 The property in the state is A. So the way 638 00:32:02,370 --> 00:32:04,260 you get the average of A is by weighing 639 00:32:04,260 --> 00:32:07,230 that A with the probability of that state. 640 00:32:07,230 --> 00:32:09,432 And we know how to calculate the probability. 641 00:32:09,432 --> 00:32:10,890 Probability we know is proportional 642 00:32:10,890 --> 00:32:14,220 to that exponential of minus beta of the Hamiltonian. 643 00:32:14,220 --> 00:32:16,950 What I don't know is the partition function 644 00:32:16,950 --> 00:32:18,880 that I have to normalize with. 645 00:32:18,880 --> 00:32:20,610 But what I can do is after I sampled 646 00:32:20,610 --> 00:32:22,410 a bunch of states I can literally 647 00:32:22,410 --> 00:32:25,020 define an approximation to my partition function 648 00:32:25,020 --> 00:32:27,990 as the sum of that exponential over all the states that I've 649 00:32:27,990 --> 00:32:30,120 sampled [INAUDIBLE]. 650 00:32:30,120 --> 00:32:32,130 And I need that, because I need to normalize 651 00:32:32,130 --> 00:32:33,390 my probabilities to 1. 652 00:32:36,210 --> 00:32:41,550 Simple sampling, it's called simple because, as I said, 653 00:32:41,550 --> 00:32:42,630 economists use it. 654 00:32:45,480 --> 00:32:49,170 That was a joke, but, you know. 655 00:32:49,170 --> 00:32:51,780 It works for small spaces-- 656 00:32:51,780 --> 00:32:52,690 small phase spaces. 657 00:32:52,690 --> 00:32:54,810 It does not work at all for materials for a reason 658 00:32:54,810 --> 00:32:57,540 I'll tell you in a second. 659 00:32:57,540 --> 00:32:59,850 And you think that's obvious. 660 00:32:59,850 --> 00:33:02,430 Did you think that was obvious that simple sampling does not 661 00:33:02,430 --> 00:33:05,010 work for materials? 662 00:33:05,010 --> 00:33:07,540 Who thought that was obvious? 663 00:33:07,540 --> 00:33:08,800 Nobody. 664 00:33:08,800 --> 00:33:09,940 So it's not obvious. 665 00:33:09,940 --> 00:33:10,990 Good. 666 00:33:10,990 --> 00:33:13,420 Because otherwise, it was a trick question. 667 00:33:13,420 --> 00:33:16,630 I'll show you in a second that it's not obvious. 668 00:33:16,630 --> 00:33:17,730 But why does it not work? 669 00:33:22,710 --> 00:33:25,560 You essentially, with simple sampling, 670 00:33:25,560 --> 00:33:32,150 pick states in proportion to their degeneracy. 671 00:33:32,150 --> 00:33:36,040 So if you think of your spin model, 672 00:33:36,040 --> 00:33:38,830 let's say you're spinning at a low temperature 673 00:33:38,830 --> 00:33:41,470 where this thing wants to be a ferromagnet with all 674 00:33:41,470 --> 00:33:43,060 the spins aligned. 675 00:33:43,060 --> 00:33:46,030 How would you pick spins? 676 00:33:46,030 --> 00:33:47,770 Well, if you pick them randomly, you 677 00:33:47,770 --> 00:33:50,590 would never end up with a ferromagnetic configuration. 678 00:33:50,590 --> 00:33:53,320 Because you picked the spin here up, here down, here up, 679 00:33:53,320 --> 00:33:54,070 here down. 680 00:33:54,070 --> 00:33:57,130 So you would sample a lot of states in the phase space 681 00:33:57,130 --> 00:33:59,370 but never the relevant one. 682 00:33:59,370 --> 00:34:00,540 OK. 683 00:34:00,540 --> 00:34:06,180 And typically you pick states with higher entropy. 684 00:34:06,180 --> 00:34:09,719 And you know that, because the entropy is essentially 685 00:34:09,719 --> 00:34:12,719 the log of the multiplicity. 686 00:34:12,719 --> 00:34:14,610 The entropy at a given energy is the log 687 00:34:14,610 --> 00:34:17,940 of the number of states of that energy. 688 00:34:17,940 --> 00:34:20,060 So remember, if you do simple sampling, 689 00:34:20,060 --> 00:34:24,389 you're going to pick states proportional to omega. 690 00:34:24,389 --> 00:34:26,699 So you're always going to pick states with high energy, 691 00:34:26,699 --> 00:34:31,239 because we know that d log dE is 1 over kT is a positive number. 692 00:34:31,239 --> 00:34:36,219 So that means the number of states with a given energy 693 00:34:36,219 --> 00:34:37,900 is an increasing function of energy. 694 00:34:40,679 --> 00:34:44,427 So that means that if you simple sample 695 00:34:44,427 --> 00:34:46,260 you're almost never going to pick down here, 696 00:34:46,260 --> 00:34:48,000 because there are almost no states. 697 00:34:48,000 --> 00:34:49,590 You're always going to pick up here, 698 00:34:49,590 --> 00:34:51,770 because there are a lot of states there. 699 00:34:51,770 --> 00:34:54,980 See the parallel with your spin model is-- 700 00:34:54,980 --> 00:34:57,020 see, these are the ferromagnetic states. 701 00:34:57,020 --> 00:34:59,540 You have low energy, everything's nicely lined up. 702 00:34:59,540 --> 00:35:02,970 These are all the random states, the paramagnetic states. 703 00:35:02,970 --> 00:35:05,490 In reality, your system lives somewhere in between. 704 00:35:08,550 --> 00:35:11,100 At some temperature it has some average energy, 705 00:35:11,100 --> 00:35:15,090 and then it has a fluctuation of energy around it. 706 00:35:15,090 --> 00:35:17,490 So your real system kind of samples these 707 00:35:17,490 --> 00:35:20,775 states, and you're always hanging out here-- 708 00:35:20,775 --> 00:35:22,400 you're always sampling the wrong thing. 709 00:35:28,490 --> 00:35:30,990 If you had said that this was obvious that it wouldn't work, 710 00:35:30,990 --> 00:35:32,850 I would have thrown this at you. 711 00:35:32,850 --> 00:35:35,280 Here's a paper from a famous physicist 712 00:35:35,280 --> 00:35:37,530 that I respect very well, but he was so wrong. 713 00:35:40,930 --> 00:35:42,472 This was somebody who's been great 714 00:35:42,472 --> 00:35:43,680 in density functional theory. 715 00:35:43,680 --> 00:35:46,013 He has done beautiful work in it and at some point tried 716 00:35:46,013 --> 00:35:47,460 to do a phase diagram. 717 00:35:47,460 --> 00:35:49,860 And realize, of course, I blanked out 718 00:35:49,860 --> 00:35:51,810 the name out of respect. 719 00:35:51,810 --> 00:35:53,350 But of course you could find out. 720 00:35:58,480 --> 00:36:00,940 They were really good at doing quantum mechanics 721 00:36:00,940 --> 00:36:02,350 and getting energies of states. 722 00:36:02,350 --> 00:36:04,660 And they wanted to get free energies now 723 00:36:04,660 --> 00:36:05,710 at finite temperature. 724 00:36:05,710 --> 00:36:07,495 This was on the aluminum lithium system. 725 00:36:07,495 --> 00:36:09,070 I mean this is 1988. 726 00:36:09,070 --> 00:36:11,590 This was like if there was any heydays ever of aluminum 727 00:36:11,590 --> 00:36:12,700 lithium, this was it. 728 00:36:12,700 --> 00:36:14,860 Everybody was working on aluminum lithium, 729 00:36:14,860 --> 00:36:17,560 because everybody thought that was the alloy of the future. 730 00:36:17,560 --> 00:36:19,990 Because lithium is the lightest metal. 731 00:36:19,990 --> 00:36:22,520 It's the lightest solid in the periodic table. 732 00:36:22,520 --> 00:36:25,720 So any atom you can replace by lithium 733 00:36:25,720 --> 00:36:27,800 will make your alloys lighter. 734 00:36:27,800 --> 00:36:29,080 So that was a big dream. 735 00:36:29,080 --> 00:36:31,880 You could make really light airplanes and whatever. 736 00:36:31,880 --> 00:36:34,250 And [INAUDIBLE] built anomalous capacity 737 00:36:34,250 --> 00:36:36,510 for making aluminum lithium alloys, et cetera, 738 00:36:36,510 --> 00:36:38,180 and it never panned out. 739 00:36:38,180 --> 00:36:40,050 Turned out you couldn't weld them very well, 740 00:36:40,050 --> 00:36:42,045 and they have really bad fatigue resistance. 741 00:36:42,045 --> 00:36:44,420 So you really don't want to make an airplane out of them. 742 00:36:44,420 --> 00:36:46,850 Ironically, they make the doors of airplanes out of them. 743 00:36:46,850 --> 00:36:48,960 I never know how to figure that one out. 744 00:36:48,960 --> 00:36:53,000 And the seat frames, there's aluminum. 745 00:36:53,000 --> 00:36:56,210 But the seat is supported on is out of aluminum lithium. 746 00:36:56,210 --> 00:36:58,520 But the doors always scares me, when I sit by the door 747 00:36:58,520 --> 00:37:00,440 and I know this is aluminum lithium, which 748 00:37:00,440 --> 00:37:02,840 has poor fatigue resistance and very poor 749 00:37:02,840 --> 00:37:03,920 welding characteristics. 750 00:37:03,920 --> 00:37:08,900 But anyway, everybody worked on aluminum lithium. 751 00:37:08,900 --> 00:37:11,310 And so these guys were good at doing quantum mechanics. 752 00:37:11,310 --> 00:37:12,920 So how did they get the phase diagram? 753 00:37:12,920 --> 00:37:19,440 That they actually took a box of atoms and randomly made up 754 00:37:19,440 --> 00:37:22,812 configurations and then calculated the energy along 755 00:37:22,812 --> 00:37:23,770 with quantum mechanics. 756 00:37:23,770 --> 00:37:26,520 So this was essentially simple sampling. 757 00:37:26,520 --> 00:37:28,897 And then they summed that and got the partition function. 758 00:37:28,897 --> 00:37:30,480 If you have the partition functioning, 759 00:37:30,480 --> 00:37:31,522 you have the free energy. 760 00:37:31,522 --> 00:37:34,920 And this was the phase diagram that came out to the right. 761 00:37:34,920 --> 00:37:37,110 It sort of violates all kinds of phase rules, 762 00:37:37,110 --> 00:37:39,980 and it's completely wrong. 763 00:37:39,980 --> 00:37:43,330 And so it's not that obvious that simple sampling 764 00:37:43,330 --> 00:37:45,720 doesn't work, because even the best people have done it. 765 00:37:54,680 --> 00:37:55,740 I like this example. 766 00:37:55,740 --> 00:37:57,440 This comes out of a great book. 767 00:37:57,440 --> 00:37:58,910 This is the book by Franklin Smith. 768 00:37:58,910 --> 00:38:01,022 And I'll show you a list of references at the end. 769 00:38:01,022 --> 00:38:02,480 This is a great book on Monte Carlo 770 00:38:02,480 --> 00:38:04,650 if you want to delve deeper into it. 771 00:38:04,650 --> 00:38:07,040 They have this great analogy of simple sample 772 00:38:07,040 --> 00:38:09,860 versus importance sampling. 773 00:38:09,860 --> 00:38:13,880 The question is, how do you measure the average depth 774 00:38:13,880 --> 00:38:16,070 of the Nile River? 775 00:38:16,070 --> 00:38:18,020 And simple sampling is essentially, well, 776 00:38:18,020 --> 00:38:22,070 you know the Nile is in Africa, so you just walk around Africa 777 00:38:22,070 --> 00:38:26,690 with a dipstick and you sort of measure how deep the Nile is. 778 00:38:26,690 --> 00:38:29,820 So of course, most of the time you're not in the Nile, 779 00:38:29,820 --> 00:38:32,570 so you get a lot of zeros. 780 00:38:32,570 --> 00:38:34,520 So what's a smarter way to do it? 781 00:38:34,520 --> 00:38:38,340 Well, a smarter ways to first find the Nile. 782 00:38:38,340 --> 00:38:40,410 That's putting a bias on your sample. 783 00:38:40,410 --> 00:38:42,395 You really want to get to the Nile. 784 00:38:42,395 --> 00:38:43,770 And then once you're in the Nile, 785 00:38:43,770 --> 00:38:45,870 you want to stay in the Nile. 786 00:38:45,870 --> 00:38:48,657 So that means you still want to have some amount of randomness 787 00:38:48,657 --> 00:38:50,490 to walk around in the Nile, otherwise you're 788 00:38:50,490 --> 00:38:52,860 not getting the average depth of the Nile. 789 00:38:52,860 --> 00:38:56,810 But you want to be biased, so you stay in it. 790 00:38:56,810 --> 00:38:59,350 And that's the idea of importance sampling. 791 00:38:59,350 --> 00:39:01,320 Remember, what is importance sampling? 792 00:39:01,320 --> 00:39:03,560 You want to sample the important states. 793 00:39:03,560 --> 00:39:06,020 Those are the ones with low energy. 794 00:39:06,020 --> 00:39:10,430 If you want to set up a thermodynamically-relevant 795 00:39:10,430 --> 00:39:12,800 ensemble, obviously you know that your system 796 00:39:12,800 --> 00:39:15,310 spends most of its time in the states with low energy 797 00:39:15,310 --> 00:39:16,610 and then fluctuates around it. 798 00:39:16,610 --> 00:39:20,070 So you want to get to the low energy and stay around it. 799 00:39:20,070 --> 00:39:21,830 So you want a sample with bias and then 800 00:39:21,830 --> 00:39:23,883 correct your probabilities for that bias. 801 00:39:23,883 --> 00:39:24,800 That's the whole idea. 802 00:39:31,950 --> 00:39:33,150 OK. 803 00:39:33,150 --> 00:39:38,790 So here's the idea of importance sampling. 804 00:39:38,790 --> 00:39:42,450 Remember, in a random sample, you pick randomly. 805 00:39:42,450 --> 00:39:46,650 And to average a property, you weigh 806 00:39:46,650 --> 00:39:51,020 the property in each state with the correct probability. 807 00:39:51,020 --> 00:39:53,480 Now, if you think about it, maybe the best 808 00:39:53,480 --> 00:39:55,880 sample we could make-- 809 00:39:55,880 --> 00:39:59,520 couldn't we sample right away with that probability? 810 00:39:59,520 --> 00:40:02,330 So could we sample with a probability that's 811 00:40:02,330 --> 00:40:05,270 directly proportional to the real probability the state 812 00:40:05,270 --> 00:40:06,800 should have in the ensemble? 813 00:40:06,800 --> 00:40:10,430 So if you do that, if you could sample these states 814 00:40:10,430 --> 00:40:13,460 with the propr probability, then their average is just 815 00:40:13,460 --> 00:40:18,318 determined by summing them, by a trivial average. 816 00:40:18,318 --> 00:40:19,860 So does everybody see the difference? 817 00:40:19,860 --> 00:40:21,750 In simple sampling, you randomly pick things. 818 00:40:21,750 --> 00:40:23,380 In importance sampling, you right away 819 00:40:23,380 --> 00:40:26,130 try to shoot for getting them with the proper probability 820 00:40:26,130 --> 00:40:28,560 so the low energy state's a lot more than the high energy 821 00:40:28,560 --> 00:40:29,520 state. 822 00:40:29,520 --> 00:40:30,790 Keep that in mind. 823 00:40:30,790 --> 00:40:33,510 There will be properties where you need the high energy 824 00:40:33,510 --> 00:40:34,260 states. 825 00:40:34,260 --> 00:40:36,780 And then you're going to have to adjust your Monte Carlo 826 00:40:36,780 --> 00:40:39,360 sampling algorithm to go there. 827 00:40:39,360 --> 00:40:43,080 Because in some cases maybe some properties 828 00:40:43,080 --> 00:40:45,750 are determined by extreme in your distribution. 829 00:40:45,750 --> 00:40:47,590 Maybe when it gets to a certain energy 830 00:40:47,590 --> 00:40:49,530 your system does something weird. 831 00:40:49,530 --> 00:40:51,720 If you want to see that, then you 832 00:40:51,720 --> 00:40:54,138 should bias your system towards there. 833 00:40:54,138 --> 00:40:56,430 But if we're interested in just thermodynamic behavior, 834 00:40:56,430 --> 00:41:00,000 we want to sample around the energy minimum. 835 00:41:00,000 --> 00:41:03,990 So there are various ways that you 836 00:41:03,990 --> 00:41:07,740 can construct the probability-weighted sample. 837 00:41:07,740 --> 00:41:10,300 A common way is what's called a metropolis algorithm. 838 00:41:10,300 --> 00:41:13,500 And this is the one we'll use pretty much 839 00:41:13,500 --> 00:41:15,660 exclusively in the [INAUDIBLE]. 840 00:41:15,660 --> 00:41:18,420 And the idea is that you construct 841 00:41:18,420 --> 00:41:20,310 what's called a Markov chain, which is really 842 00:41:20,310 --> 00:41:26,100 just a sequence of states, where, over a long time, if you 843 00:41:26,100 --> 00:41:28,140 do it long enough, you visit each state 844 00:41:28,140 --> 00:41:30,727 with the proper probability. 845 00:41:30,727 --> 00:41:32,310 So I'm going to define two quantities, 846 00:41:32,310 --> 00:41:36,232 and then I'm going to tell you what rules they should obey. 847 00:41:36,232 --> 00:41:38,190 So obviously you have to start with some state. 848 00:41:38,190 --> 00:41:39,648 Well, in a metropolis algorithm you 849 00:41:39,648 --> 00:41:41,020 can pick that state randomly. 850 00:41:41,020 --> 00:41:43,500 So I call that state i. 851 00:41:43,500 --> 00:41:49,050 Then you have to pick a new state called j from i. 852 00:41:49,050 --> 00:41:52,660 And the rate at which we pick that state 853 00:41:52,660 --> 00:41:56,880 I call omega 0 i to j. 854 00:41:56,880 --> 00:41:58,710 But you don't always accept that state. 855 00:41:58,710 --> 00:42:02,940 You accept that with some probability p i going to j. 856 00:42:02,940 --> 00:42:08,010 So the total transfer rate, how many times you pick j from i, 857 00:42:08,010 --> 00:42:09,990 and put it in the Markov chain is determined 858 00:42:09,990 --> 00:42:12,540 by the product of these two. 859 00:42:12,540 --> 00:42:15,780 It's times you actually pick it as a possible transfer 860 00:42:15,780 --> 00:42:19,800 state times the success rate in getting there. 861 00:42:19,800 --> 00:42:20,554 That's p i j. 862 00:42:35,080 --> 00:42:37,180 The two conditions that a metropolis algorithm 863 00:42:37,180 --> 00:42:40,363 has to satisfy are typically these two. 864 00:42:40,363 --> 00:42:42,280 The first one isn't actually a real condition, 865 00:42:42,280 --> 00:42:45,790 it's just sort of a practical one to work with. 866 00:42:45,790 --> 00:42:48,280 Usually people look for what's called equal 867 00:42:48,280 --> 00:42:51,220 a priori probabilities. 868 00:42:51,220 --> 00:42:54,340 What that means is that the rate at which you pick j froim i 869 00:42:54,340 --> 00:42:58,150 should be the same as the rate at which you pick i from j. 870 00:42:58,150 --> 00:43:00,260 That's actually not absolutely necessary. 871 00:43:00,260 --> 00:43:03,820 But if you don't satisfy this, you have to correct later on. 872 00:43:03,820 --> 00:43:06,490 And we'll do that two lectures from now-- we'll 873 00:43:06,490 --> 00:43:07,910 give you an example of that. 874 00:43:07,910 --> 00:43:11,650 But it's easiest to work with equal a priori probabilities. 875 00:43:11,650 --> 00:43:14,260 In most this is a condition that's 876 00:43:14,260 --> 00:43:16,180 amazingly easily satisfied. 877 00:43:16,180 --> 00:43:17,230 You'll see that. 878 00:43:17,230 --> 00:43:20,500 Next lecture, we'll do some examples of complicated Monte 879 00:43:20,500 --> 00:43:23,530 Carlo moves, and it's hard to violate this one. 880 00:43:23,530 --> 00:43:24,670 You have to work on it. 881 00:43:24,670 --> 00:43:26,320 But it does happen. 882 00:43:26,320 --> 00:43:29,620 The more important one is what's called detailed balance. 883 00:43:29,620 --> 00:43:33,100 Detailed balance is that the rate at which you actually 884 00:43:33,100 --> 00:43:37,210 transfer from i to jay and the rate at which you transfer 885 00:43:37,210 --> 00:43:41,590 from j to i has to be inversely proportional to the probability 886 00:43:41,590 --> 00:43:43,280 that you're in the initial state. 887 00:43:43,280 --> 00:43:44,980 So essentially the rate at which you 888 00:43:44,980 --> 00:43:49,570 go from i to j over the rate at which you go from j to i 889 00:43:49,570 --> 00:43:54,820 has to equal p j over p i, or these are the probabilities 890 00:43:54,820 --> 00:43:57,280 that you are in that state. 891 00:43:57,280 --> 00:43:59,650 It's actually fairly easy to convince yourself 892 00:43:59,650 --> 00:44:03,910 that if that's true you're going to end up 893 00:44:03,910 --> 00:44:06,160 with a correct probability distribution. 894 00:44:06,160 --> 00:44:07,910 I mean, think about it. 895 00:44:07,910 --> 00:44:11,710 If you look at this carefully, what is this really saying? 896 00:44:11,710 --> 00:44:17,440 That if you were doing a lot of samplings in parallel, 897 00:44:17,440 --> 00:44:20,380 this gives you a steady state distribution of states. 898 00:44:20,380 --> 00:44:22,660 Because if you look at the left-hand side 899 00:44:22,660 --> 00:44:25,780 of the green box, that's the number of times 900 00:44:25,780 --> 00:44:29,480 you're going out if i into j. 901 00:44:29,480 --> 00:44:32,600 Because it's the probability that you're in i times 902 00:44:32,600 --> 00:44:35,580 the rate at which you transfer to j. 903 00:44:35,580 --> 00:44:38,820 And the right-hand side is the probability 904 00:44:38,820 --> 00:44:41,740 at the rate at which you go from j to i. 905 00:44:41,740 --> 00:44:43,980 Which is, first of all, you have to be in j-- 906 00:44:43,980 --> 00:44:46,050 otherwise you can't go from j to i-- 907 00:44:46,050 --> 00:44:52,640 times the rate W j i, the rate at which you go from j to i. 908 00:44:52,640 --> 00:44:54,980 So when these two are equal, you sort of 909 00:44:54,980 --> 00:44:58,890 have a net zero flow of probability density. 910 00:44:58,890 --> 00:45:00,920 So you have a steady state distribution. 911 00:45:00,920 --> 00:45:04,182 If this holds for all pair of states i j, 912 00:45:04,182 --> 00:45:06,390 your probability distribution isn't changing anymore. 913 00:45:06,390 --> 00:45:09,600 So you end up with a steady state distribution. 914 00:45:25,460 --> 00:45:27,380 You can be fairly creative about how 915 00:45:27,380 --> 00:45:31,910 you make metropolis algorithms. 916 00:45:31,910 --> 00:45:35,090 But here's a very simple one, which you'll often see 917 00:45:35,090 --> 00:45:37,760 a lot in a Monte Carlo codes. 918 00:45:37,760 --> 00:45:45,320 So you first have to design some pick rate, some rate at which 919 00:45:45,320 --> 00:45:47,315 you choose new states for i. 920 00:45:47,315 --> 00:45:50,390 So that's omega 0. 921 00:45:50,390 --> 00:45:52,190 So you have to have some mechanism 922 00:45:52,190 --> 00:45:56,110 for picking potential moves. 923 00:45:56,110 --> 00:45:59,320 So omega 0 i to j. 924 00:45:59,320 --> 00:46:02,140 And then you have to have an acceptance criteria. 925 00:46:02,140 --> 00:46:03,940 And a typical acceptance criteria 926 00:46:03,940 --> 00:46:08,140 is that you accept the probability from i to j is 1. 927 00:46:08,140 --> 00:46:11,270 That means you definitely accept when the energy goes down-- 928 00:46:11,270 --> 00:46:14,020 so when the new state has lower energy than i, 929 00:46:14,020 --> 00:46:15,550 than the previous state. 930 00:46:15,550 --> 00:46:20,260 And if the new state has higher energy than the initial state, 931 00:46:20,260 --> 00:46:23,500 then you accept with a Boltzmann factor. 932 00:46:23,500 --> 00:46:27,340 So then you accept with a certain probability. 933 00:46:27,340 --> 00:46:30,530 You can easily show that this satisfies detailed balance. 934 00:46:30,530 --> 00:46:33,220 Because if you-- let's see if I can squeeze this in here. 935 00:46:33,220 --> 00:46:38,830 Detailed balances that essentially p i p i to j 936 00:46:38,830 --> 00:46:50,710 has to equal p j p j to I. Well p i 937 00:46:50,710 --> 00:47:01,670 is the exponential of minus beta i over the partition function. 938 00:47:01,670 --> 00:47:06,180 p i to J, let's assume that i is the state with low energy is 1. 939 00:47:06,180 --> 00:47:10,950 That has to equal p j, which is the exponential. 940 00:47:10,950 --> 00:47:13,784 Sorry for my poor writing here. 941 00:47:13,784 --> 00:47:18,850 This pen is displaced from where I see it. 942 00:47:18,850 --> 00:47:23,920 This is exponential minus beta E j over z, 943 00:47:23,920 --> 00:47:28,390 times p i j, which is the exponential. 944 00:47:28,390 --> 00:47:29,950 I'm going to write it below here. 945 00:47:29,950 --> 00:47:36,280 The exponential of minus beta Ei minus Ej. 946 00:47:39,710 --> 00:47:43,430 This Ej cancels with this one, and so you 947 00:47:43,430 --> 00:47:45,890 get that the exponential is minus beta Ei equals 948 00:47:45,890 --> 00:47:47,480 the exponential minus beta Ei. 949 00:47:47,480 --> 00:47:54,150 So this set of transfer rates satisfies detailed balance 950 00:47:54,150 --> 00:47:55,430 in an almost trivial way. 951 00:47:58,480 --> 00:48:00,370 So how would you practically implement 952 00:48:00,370 --> 00:48:04,110 a metropolis algorithm for sampling phase space? 953 00:48:04,110 --> 00:48:07,300 You would start with some configuration you would choose 954 00:48:07,300 --> 00:48:08,930 a perturbation of the system. 955 00:48:08,930 --> 00:48:10,900 So that's essentially determine your mechanism 956 00:48:10,900 --> 00:48:13,660 for picking omega 0 i j. 957 00:48:13,660 --> 00:48:16,630 You compute the energy for that perturbation. 958 00:48:16,630 --> 00:48:20,880 If the change is negative, you accept the perturbation. 959 00:48:20,880 --> 00:48:24,870 If it's positive, you accept it with some probability. 960 00:48:24,870 --> 00:48:27,330 And then you go on and you choose the next perturbation. 961 00:48:27,330 --> 00:48:30,270 And so you essentially cycle this. 962 00:48:30,270 --> 00:48:32,047 So you build up a sequence of states. 963 00:48:32,047 --> 00:48:32,880 This is your sample. 964 00:48:57,170 --> 00:49:05,430 This must be a Pentium 286. 965 00:49:05,430 --> 00:49:06,850 It's checking the equations. 966 00:49:06,850 --> 00:49:11,190 That's-- OK. 967 00:49:11,190 --> 00:49:15,300 So what you will get then is some property 968 00:49:15,300 --> 00:49:17,980 as a function of sample size or sampling, 969 00:49:17,980 --> 00:49:20,830 which is often called Monte Carlo time. 970 00:49:20,830 --> 00:49:23,340 Which is very confusing, because this is not 971 00:49:23,340 --> 00:49:26,400 a dynamical trajectory. 972 00:49:26,400 --> 00:49:29,700 You could think of it as a trajectory through phase space, 973 00:49:29,700 --> 00:49:31,560 but there's nothing here that says that it's 974 00:49:31,560 --> 00:49:33,523 a dynamical trajectory. 975 00:49:33,523 --> 00:49:35,190 I'm going to later confuse you even more 976 00:49:35,190 --> 00:49:37,830 and talk about what's called kinetic Monte Carlo, which 977 00:49:37,830 --> 00:49:41,520 is a way of mimicking dynamics with Monte Carlo. 978 00:49:41,520 --> 00:49:43,530 But standard Monte Carlo simply has 979 00:49:43,530 --> 00:49:45,850 nothing to do with kinetics. 980 00:49:45,850 --> 00:49:48,160 It's an efficient way to sample phase space. 981 00:49:48,160 --> 00:49:51,870 So you'll go through some property, 982 00:49:51,870 --> 00:49:55,740 and then you just essentially track the average. 983 00:49:55,740 --> 00:49:58,230 Typically what people do is they often, 984 00:49:58,230 --> 00:49:59,820 just for practical consideration, 985 00:49:59,820 --> 00:50:02,280 cut off the first part and remove it. 986 00:50:02,280 --> 00:50:06,630 Because, remember, you started with a random state. 987 00:50:06,630 --> 00:50:09,870 You may be very far away from the low energy state. 988 00:50:09,870 --> 00:50:13,800 And so by cutting that off, the first part, you often 989 00:50:13,800 --> 00:50:17,580 get a much faster relaxation of the average. 990 00:50:17,580 --> 00:50:19,800 Because the initial state may have 991 00:50:19,800 --> 00:50:24,330 really very high energy and thereby some weird property. 992 00:50:24,330 --> 00:50:27,570 And so you only start sampling after the system has 993 00:50:27,570 --> 00:50:30,360 relaxed towards the minimum. 994 00:50:30,360 --> 00:50:32,490 So you often talk about, in Monte Carlo, people 995 00:50:32,490 --> 00:50:35,550 using a certain amount of equilibration passes 996 00:50:35,550 --> 00:50:37,470 and after a certain amount of sampling passes. 997 00:50:37,470 --> 00:50:39,417 And you only average over the sampling pass. 998 00:50:44,520 --> 00:50:48,210 So let's go through this for the very simple magnetic model. 999 00:50:48,210 --> 00:50:52,280 How will you do it in a magnetic model, in a simple Ising model? 1000 00:50:52,280 --> 00:50:54,627 You could randomly initialize the spin. 1001 00:50:54,627 --> 00:50:55,710 You could start with that. 1002 00:50:55,710 --> 00:50:58,170 You could also start with a ferromagnetic configuration. 1003 00:50:58,170 --> 00:51:00,740 The closer you start to your equilibrium state, 1004 00:51:00,740 --> 00:51:03,180 the faster this is going to go. 1005 00:51:03,180 --> 00:51:06,360 An obvious perturbation is you pick a spin 1006 00:51:06,360 --> 00:51:08,170 and you flip it over. 1007 00:51:08,170 --> 00:51:09,600 So if it was up, you make it down. 1008 00:51:09,600 --> 00:51:11,220 If it was down, you make it up. 1009 00:51:11,220 --> 00:51:13,440 Calculate the energy for perturbation. 1010 00:51:13,440 --> 00:51:15,370 If it goes down, you accept it. 1011 00:51:15,370 --> 00:51:17,790 If it doesn't, you accept it with a certain Boltzman 1012 00:51:17,790 --> 00:51:19,870 probability and you go back. 1013 00:51:19,870 --> 00:51:22,610 And as you can see, you can do this very fast. 1014 00:51:22,610 --> 00:51:26,130 If you think of computer code now, how you write this, 1015 00:51:26,130 --> 00:51:28,100 these are trivial operations. 1016 00:51:28,100 --> 00:51:29,850 You have to calculate a random number 1017 00:51:29,850 --> 00:51:32,490 to randomly pick a spin in the lattice, 1018 00:51:32,490 --> 00:51:34,298 or two if it's a two dimensional. 1019 00:51:34,298 --> 00:51:36,465 You have to compute the energy, but your Hamiltonian 1020 00:51:36,465 --> 00:51:38,790 is just like j sigma i sigma j. 1021 00:51:38,790 --> 00:51:39,720 That goes really fast. 1022 00:51:39,720 --> 00:51:42,570 It's like multiplying a few numbers. 1023 00:51:42,570 --> 00:51:45,090 And then you have to compute that exponential also. 1024 00:51:45,090 --> 00:51:47,880 So you can cycle this through literally millions of times 1025 00:51:47,880 --> 00:51:52,110 per second on a modern computer. 1026 00:51:52,110 --> 00:51:54,600 Even on a very basic computer, if you 1027 00:51:54,600 --> 00:51:56,430 did just this two-dimensional Ising model, 1028 00:51:56,430 --> 00:51:59,400 if you coded it well, you could do probably 50 million 1029 00:51:59,400 --> 00:52:01,440 of these per second. 1030 00:52:01,440 --> 00:52:04,290 So that's how fast you can sample through your ensemble. 1031 00:52:07,980 --> 00:52:09,120 So here's what comes out. 1032 00:52:09,120 --> 00:52:12,400 Here's real output. 1033 00:52:12,400 --> 00:52:16,590 So this is the 2D square lattice model I showed you. 1034 00:52:16,590 --> 00:52:19,500 The temperature is in units of j over k-- 1035 00:52:19,500 --> 00:52:22,270 so the interaction constant or the Boltzmann constant. 1036 00:52:22,270 --> 00:52:26,310 So kT over j equals 2, that's actually below the period 1037 00:52:26,310 --> 00:52:27,670 temperature of this model. 1038 00:52:27,670 --> 00:52:30,280 So this is a ferromagnetic model. 1039 00:52:30,280 --> 00:52:33,420 And so if you heat it up, it would go paramagnetic 1040 00:52:33,420 --> 00:52:34,260 at some point. 1041 00:52:34,260 --> 00:52:37,650 This is still below that temperature, but it's not far. 1042 00:52:37,650 --> 00:52:41,010 The transition temperature would be about 2.23 1043 00:52:41,010 --> 00:52:43,360 or something like that, 2.24 in this model. 1044 00:52:43,360 --> 00:52:46,180 So you're below the paramagnetic transition temperature, 1045 00:52:46,180 --> 00:52:47,410 but not much. 1046 00:52:47,410 --> 00:52:51,900 So I started it up randomly I think. 1047 00:52:51,900 --> 00:52:53,760 And so the energy is sort of high. 1048 00:52:53,760 --> 00:52:55,230 It kind of fluctuates for a while 1049 00:52:55,230 --> 00:52:57,700 and then starts going down. 1050 00:52:57,700 --> 00:53:00,060 And here it seems to have relaxed fairly well, 1051 00:53:00,060 --> 00:53:04,920 and then it just kind of fluctuates around the average. 1052 00:53:04,920 --> 00:53:09,450 And you average this and you'd get the energy. 1053 00:53:09,450 --> 00:53:12,090 If you want to know the magnetization, 1054 00:53:12,090 --> 00:53:14,400 well, you would keep track of the magnetization, which 1055 00:53:14,400 --> 00:53:16,483 is sort of the difference between how much up spin 1056 00:53:16,483 --> 00:53:18,390 do I have and how much down spin do I have. 1057 00:53:18,390 --> 00:53:20,910 You track that, and this is that quantity 1058 00:53:20,910 --> 00:53:23,160 as a function of temperature, and see 1059 00:53:23,160 --> 00:53:24,300 it's clearly ferromagnetic. 1060 00:53:24,300 --> 00:53:27,176 You have preferred up spin. 1061 00:53:27,176 --> 00:53:29,860 It slowly goes down with temperature. 1062 00:53:29,860 --> 00:53:32,170 And then this is the transition temperature. 1063 00:53:32,170 --> 00:53:34,180 It's the paramagnetic transition temperature, 1064 00:53:34,180 --> 00:53:37,990 beautiful second-order transition. 1065 00:53:37,990 --> 00:53:40,720 So just to give an idea, doing this with Monte Carlo 1066 00:53:40,720 --> 00:53:44,830 is a one-page code without the bells and whistles. 1067 00:53:44,830 --> 00:53:47,530 And it runs before you blink. 1068 00:53:47,530 --> 00:53:52,730 I mean, it's so fast that you would sample this. 1069 00:53:52,730 --> 00:53:55,530 OK. 1070 00:53:55,530 --> 00:53:58,170 How do you detect phase transitions? 1071 00:53:58,170 --> 00:54:02,100 This is sort of a critical issue in any simulation method 1072 00:54:02,100 --> 00:54:03,720 to be honest. 1073 00:54:03,720 --> 00:54:06,300 Whether you do Monte Carlo or molecular dynamics, 1074 00:54:06,300 --> 00:54:08,460 it's one thing to look at the atoms, 1075 00:54:08,460 --> 00:54:12,270 and it's another thing to figure out exactly what they're doing. 1076 00:54:12,270 --> 00:54:15,090 You may have seen that if you try to do things 1077 00:54:15,090 --> 00:54:16,650 like melting with molecular dynamics, 1078 00:54:16,650 --> 00:54:19,622 it's not trivial to exactly find where the atoms melt. 1079 00:54:19,622 --> 00:54:20,580 Because what do you do? 1080 00:54:20,580 --> 00:54:23,070 You just stare at them. 1081 00:54:23,070 --> 00:54:26,130 But it seems kind of a fairly random decision about-- 1082 00:54:26,130 --> 00:54:30,060 let's say I gave you this lattice model, 1083 00:54:30,060 --> 00:54:33,050 and I asked you, where does it become disordered? 1084 00:54:33,050 --> 00:54:36,060 Where is the paramagnetic transition temperature? 1085 00:54:36,060 --> 00:54:38,070 You could just stare at the spins. 1086 00:54:38,070 --> 00:54:40,440 But I'll tell you, when you're close to that transition, 1087 00:54:40,440 --> 00:54:42,932 just below and just above, it looks pretty similar. 1088 00:54:42,932 --> 00:54:44,640 And one is still the ferromagnetic phase. 1089 00:54:44,640 --> 00:54:46,180 The other one's the antiferromagnetic phase. 1090 00:54:46,180 --> 00:54:47,760 So you need a more systematic way 1091 00:54:47,760 --> 00:54:49,930 of looking at phase transitions. 1092 00:54:49,930 --> 00:54:51,420 And the beauty is you can do that 1093 00:54:51,420 --> 00:54:54,720 through the thermodynamic quantities you sampled. 1094 00:54:54,720 --> 00:54:57,620 If you sample the energy, you typically 1095 00:54:57,620 --> 00:54:59,510 use the fact that the energy is discontinuous 1096 00:54:59,510 --> 00:55:02,260 at first-order transitions. 1097 00:55:02,260 --> 00:55:05,050 It's not so at second-order transitions. 1098 00:55:05,050 --> 00:55:07,300 If you work at constant chemical potential-- and I'll 1099 00:55:07,300 --> 00:55:08,482 come to that later. 1100 00:55:08,482 --> 00:55:10,440 It's often a useful thing to do in Monte Carlo. 1101 00:55:10,440 --> 00:55:12,610 it tends to equilibrate faster-- 1102 00:55:12,610 --> 00:55:14,830 you will have concentration discontinuities 1103 00:55:14,830 --> 00:55:16,625 at first-order transition. 1104 00:55:19,740 --> 00:55:22,460 Another thing to usually track is the heat capacity 1105 00:55:22,460 --> 00:55:23,810 of the system. 1106 00:55:23,810 --> 00:55:27,800 The heat capacity is the fluctuation of the energy. 1107 00:55:27,800 --> 00:55:29,210 If you go back-- 1108 00:55:29,210 --> 00:55:29,900 can I go back? 1109 00:55:33,730 --> 00:55:35,770 The heat capacity is essentially a measure 1110 00:55:35,770 --> 00:55:39,280 of how much you fluctuate around the average. 1111 00:55:39,280 --> 00:55:40,480 It's your sigma squared. 1112 00:55:44,457 --> 00:55:46,540 So the heat capacity is defined as the fluctuation 1113 00:55:46,540 --> 00:55:49,030 of the energy-- so E squared average minus E 1114 00:55:49,030 --> 00:55:51,640 average squared. 1115 00:55:51,640 --> 00:55:55,090 That tends to show singularities at second-order transitions. 1116 00:55:55,090 --> 00:55:58,300 But it's even useful often to spot first-order transitions, 1117 00:55:58,300 --> 00:55:59,950 even though it's a little harder to. 1118 00:55:59,950 --> 00:56:03,130 I'll show you a lot of sort of case studies and examples 1119 00:56:03,130 --> 00:56:05,412 on more complicated systems later on. 1120 00:56:05,412 --> 00:56:07,120 But for example, here's the heat capacity 1121 00:56:07,120 --> 00:56:12,360 of that 2D model as a function of normalized temperature. 1122 00:56:12,360 --> 00:56:15,360 So notice how the heat capacity the transition 1123 00:56:15,360 --> 00:56:18,180 is here of course, about 2.24 or so. 1124 00:56:18,180 --> 00:56:20,760 [INAUDIBLE] has a beautiful singularity. 1125 00:56:20,760 --> 00:56:23,800 This is actually a logarithmic divergence to infinity. 1126 00:56:23,800 --> 00:56:27,540 So you can't miss this transition. 1127 00:56:27,540 --> 00:56:29,800 Even if you tried, you wouldn't miss it. 1128 00:56:34,400 --> 00:56:37,070 There's a great Java program on the web that 1129 00:56:37,070 --> 00:56:38,663 does two-dimensional lattice files 1130 00:56:38,663 --> 00:56:39,830 if you want to play with it. 1131 00:56:39,830 --> 00:56:43,130 And you can literally type in temperatures and interaction 1132 00:56:43,130 --> 00:56:45,840 strength, and it just simulates on the screen for you 1133 00:56:45,840 --> 00:56:48,960 and plots like thermodynamic quantities. 1134 00:56:48,960 --> 00:56:49,460 Pretty cool. 1135 00:56:53,430 --> 00:56:53,930 OK. 1136 00:56:56,560 --> 00:56:57,940 So in the last bit of time I want 1137 00:56:57,940 --> 00:57:03,220 to go through some variants of this lattice model, how 1138 00:57:03,220 --> 00:57:04,540 we use it for other fields. 1139 00:57:04,540 --> 00:57:07,060 And then, in the next lecture, I'll pick up on that. 1140 00:57:10,150 --> 00:57:14,110 Like I said before, you can use this for anything 1141 00:57:14,110 --> 00:57:18,920 that has essentially two states on a fixed topology. 1142 00:57:18,920 --> 00:57:23,200 So like I said, you could use it for A and B atoms on a lattice. 1143 00:57:23,200 --> 00:57:26,420 It's very often used for surface absorption as well. 1144 00:57:26,420 --> 00:57:29,620 If you say you have a surface-- 1145 00:57:29,620 --> 00:57:31,888 let's say this is a triangular lattice, 1146 00:57:31,888 --> 00:57:33,680 so that could be the surface sites on a 1 1 1147 00:57:33,680 --> 00:57:38,710 1 FCC metal, like platinum or palladium. 1148 00:57:38,710 --> 00:57:43,420 If stuff can absorb on the lattice side-- 1149 00:57:43,420 --> 00:57:47,050 now you can really see that my pen is miscalibrated. 1150 00:57:47,050 --> 00:57:51,070 I'm drawing on the latest points. 1151 00:57:51,070 --> 00:57:53,250 So if you want to know if stuff absorbs, 1152 00:57:53,250 --> 00:57:55,320 well, you see that becomes a binary problem. 1153 00:57:55,320 --> 00:57:57,630 At each site you want to know, is there 1154 00:57:57,630 --> 00:57:59,440 an absorbent there or not? 1155 00:57:59,440 --> 00:58:01,740 So you can define a binary variable and whether it's 1156 00:58:01,740 --> 00:58:03,870 plus or minus 1 doesn't matter. 1157 00:58:03,870 --> 00:58:08,160 People often write things in not spin variables 1158 00:58:08,160 --> 00:58:11,250 but occupation variables, which are maybe 1159 00:58:11,250 --> 00:58:12,330 a little more natural. 1160 00:58:12,330 --> 00:58:14,460 They're mathematically less elegant, 1161 00:58:14,460 --> 00:58:16,350 but you could say the occupation variable is 1162 00:58:16,350 --> 00:58:20,490 1 when the side is occupied and is 0 when the side is not. 1163 00:58:20,490 --> 00:58:23,880 The advantage then is that you only get interaction 1164 00:58:23,880 --> 00:58:27,930 when two things are occupied. 1165 00:58:27,930 --> 00:58:31,440 So one occupied side with one occupied side 1166 00:58:31,440 --> 00:58:33,780 gives you no contribution to the Hamiltonian, 1167 00:58:33,780 --> 00:58:36,370 because one of the p i's is 0. 1168 00:58:36,370 --> 00:58:40,060 Whereas if you use spins, then this is a plus 1, 1169 00:58:40,060 --> 00:58:41,260 and this is minus 1. 1170 00:58:41,260 --> 00:58:43,812 So it gives you a contribution to the Hamiltonian. 1171 00:58:43,812 --> 00:58:45,520 Mathematically it all works out, but it's 1172 00:58:45,520 --> 00:58:47,020 a little more confusing. 1173 00:58:50,460 --> 00:58:53,570 So you can even write that out for a binary alloy. 1174 00:58:53,570 --> 00:58:58,160 You can sum over the probability that you have A B bonds on i j 1175 00:58:58,160 --> 00:59:00,530 and multiply with the A B interaction. 1176 00:59:00,530 --> 00:59:04,490 You can get the probability that it's A A 1177 00:59:04,490 --> 00:59:06,380 and multiply it with the AA interaction, 1178 00:59:06,380 --> 00:59:08,630 and you get the probability that it's B B and multiply 1179 00:59:08,630 --> 00:59:10,486 with the B B interaction. 1180 00:59:16,630 --> 00:59:17,130 OK. 1181 00:59:21,000 --> 00:59:22,890 I think the coolest thing about Monte Carlo 1182 00:59:22,890 --> 00:59:26,910 is that you have anomalous degrees of freedom 1183 00:59:26,910 --> 00:59:28,705 about the moves you do. 1184 00:59:28,705 --> 00:59:30,330 These are essentially the perturbations 1185 00:59:30,330 --> 00:59:34,870 you attempt to get to the next state in the Markov chain. 1186 00:59:34,870 --> 00:59:37,210 And the one thing you should never do 1187 00:59:37,210 --> 00:59:40,450 is let yourself be guided there or constrained 1188 00:59:40,450 --> 00:59:42,490 by physical principles. 1189 00:59:42,490 --> 00:59:46,990 People often do moves that they think are physical. 1190 00:59:46,990 --> 00:59:49,540 But remember, it's not the point to do a physical move. 1191 00:59:49,540 --> 00:59:53,270 The point is to efficiently sample phase space. 1192 00:59:53,270 --> 00:59:54,760 And so on the spin model, well, you 1193 00:59:54,760 --> 00:59:57,850 think, a spin going from up to down 1194 00:59:57,850 --> 00:59:59,710 or vice versa is a physical move. 1195 00:59:59,710 --> 01:00:02,260 Spins flip up and down all the time. 1196 01:00:02,260 --> 01:00:03,760 But in some cases, for example, it's 1197 01:00:03,760 --> 01:00:07,300 much more efficient to flip over whole patches of spin. 1198 01:00:07,300 --> 01:00:10,750 When you're near a second-order transition in a material, 1199 01:00:10,750 --> 01:00:14,260 the fluctuations in the system are very long range. 1200 01:00:14,260 --> 01:00:16,000 So it's sort of like the wind blowing 1201 01:00:16,000 --> 01:00:18,550 to a wheat field or something. 1202 01:00:18,550 --> 01:00:21,230 There are whole patches of spin that move together. 1203 01:00:21,230 --> 01:00:23,200 So you'll get much better convergence there 1204 01:00:23,200 --> 01:00:25,300 if you flip whole patches of spin. 1205 01:00:25,300 --> 01:00:26,980 So rather than attempting one at a time, 1206 01:00:26,980 --> 01:00:30,100 you could say I'm going to attempt 10 at a time close 1207 01:00:30,100 --> 01:00:32,900 to each other. 1208 01:00:32,900 --> 01:00:34,720 So you can always pick your dynamic. 1209 01:00:34,720 --> 01:00:36,750 So I've given you an example. 1210 01:00:36,750 --> 01:00:38,560 Here's an example of two types of dynamics. 1211 01:00:38,560 --> 01:00:40,477 Let's say you have A and B atoms on a lattice. 1212 01:00:44,170 --> 01:00:46,240 So if you want to equilibrate that, 1213 01:00:46,240 --> 01:00:47,830 you could pick the obvious dynamics. 1214 01:00:47,830 --> 01:00:52,750 You could interchange A and B. So if you have a lot of us 1215 01:00:52,750 --> 01:00:57,240 and you have A here and B here, you could interchange them. 1216 01:00:57,240 --> 01:01:00,010 That could be a possible move. 1217 01:01:00,010 --> 01:01:01,540 And it's sort of an obvious move, 1218 01:01:01,540 --> 01:01:03,520 because it sort of mimics diffusion. 1219 01:01:03,520 --> 01:01:06,700 Atoms interchange and equilibrate like that. 1220 01:01:06,700 --> 01:01:11,170 That's called Glauber dynamics on the lattice model. 1221 01:01:11,170 --> 01:01:16,120 Turns out that that's useful, but a much faster move 1222 01:01:16,120 --> 01:01:18,160 to equilibrate is actually exchanging 1223 01:01:18,160 --> 01:01:20,210 the identity of the atoms. 1224 01:01:20,210 --> 01:01:23,380 So rather than interchanging A and B at different positions, 1225 01:01:23,380 --> 01:01:26,920 you take an atom, and if it's A, you try to turn it into B. 1226 01:01:26,920 --> 01:01:28,300 And if it's B, you try to turn it 1227 01:01:28,300 --> 01:01:33,520 into A. That's called Kawasaki dynamics. 1228 01:01:33,520 --> 01:01:36,990 Now, you see that if you do Kawasaki dynamics, 1229 01:01:36,990 --> 01:01:39,660 you will need some kind of chemical potential. 1230 01:01:39,660 --> 01:01:41,820 Because in Glauber dynamics, you constrain-- 1231 01:01:41,820 --> 01:01:44,270 the composition is fixed. 1232 01:01:44,270 --> 01:01:46,700 You just interchange the positions of the A and B, 1233 01:01:46,700 --> 01:01:50,720 but you're not changing the number of A or B. 1234 01:01:50,720 --> 01:01:52,580 You conserve the concentration. 1235 01:01:52,580 --> 01:01:57,120 In Kawasaki dynamics, you don't change the number of particles, 1236 01:01:57,120 --> 01:01:59,210 but you do change the concentration of the A and B 1237 01:01:59,210 --> 01:01:59,710 species. 1238 01:01:59,710 --> 01:02:03,090 Because whenever I flip A to B, I'm changing the concentration. 1239 01:02:03,090 --> 01:02:05,480 So you will need some kind of chemical potential term 1240 01:02:05,480 --> 01:02:08,090 in your Hamiltonian, and that chemical potential term 1241 01:02:08,090 --> 01:02:11,120 is the difference in chemical potentials of A and B. 1242 01:02:11,120 --> 01:02:14,580 So that difference will drive the concentration. 1243 01:02:14,580 --> 01:02:16,220 So it's essentially multiplied with-- 1244 01:02:16,220 --> 01:02:17,803 this is essentially the concentration. 1245 01:02:17,803 --> 01:02:21,800 It's the sum of the probability that you have A on each side. 1246 01:02:21,800 --> 01:02:24,220 That's essentially the concentration. 1247 01:02:24,220 --> 01:02:25,720 OK. 1248 01:02:25,720 --> 01:02:30,190 You'll find that Kawasaki dynamics equilibrates 100 times 1249 01:02:30,190 --> 01:02:32,830 faster easily in many cases. 1250 01:02:32,830 --> 01:02:36,850 The reason is that, if you think about it, if you started 1251 01:02:36,850 --> 01:02:42,640 with a situation where you had too much A in one area and too 1252 01:02:42,640 --> 01:02:45,768 much B in another, if you do Glauber dynamics, 1253 01:02:45,768 --> 01:02:47,560 you can really only equilibrate that system 1254 01:02:47,560 --> 01:02:49,710 by slowly diffusing these. 1255 01:02:49,710 --> 01:02:52,150 A has to go this way, B has to go that way. 1256 01:02:52,150 --> 01:02:56,902 In Kawasaki dynamics, the A here feels right away 1257 01:02:56,902 --> 01:02:58,360 that there's too much of it, and it 1258 01:02:58,360 --> 01:03:01,060 starts turning to B. The B here starts turning to A, 1259 01:03:01,060 --> 01:03:02,860 and you're done. 1260 01:03:02,860 --> 01:03:05,690 So it tends to be a much more efficient way. 1261 01:03:05,690 --> 01:03:07,150 This is sort of a generic thing. 1262 01:03:07,150 --> 01:03:13,220 Open systems tend to equilibrate faster than closed systems. 1263 01:03:13,220 --> 01:03:16,250 So if you can at all make your system open 1264 01:03:16,250 --> 01:03:19,970 and end up with the same end result, it'll go much faster. 1265 01:03:19,970 --> 01:03:22,400 And it's essentially because you give open systems 1266 01:03:22,400 --> 01:03:25,220 more degrees of freedom. 1267 01:03:25,220 --> 01:03:27,200 You can change the chemistry of them. 1268 01:03:30,250 --> 01:03:31,013 Yes. 1269 01:03:31,013 --> 01:03:33,180 AUDIENCE: I mean, principally [INAUDIBLE] complexity 1270 01:03:33,180 --> 01:03:34,390 of the system. 1271 01:03:34,390 --> 01:03:37,702 First of all, you're getting [INAUDIBLE].. 1272 01:03:37,702 --> 01:03:39,660 GERBRAND CEDER: Not in the thermodynamic limit. 1273 01:03:39,660 --> 01:03:42,180 That's sort of an idea of statistical mechanics 1274 01:03:42,180 --> 01:03:44,754 is that you get the same averages. 1275 01:03:44,754 --> 01:03:45,770 AUDIENCE: [INAUDIBLE]. 1276 01:03:51,250 --> 01:03:52,250 GERBRAND CEDER: You can. 1277 01:03:52,250 --> 01:03:53,570 But but you still-- 1278 01:03:53,570 --> 01:03:55,430 so that'll definitely make it go faster 1279 01:03:55,430 --> 01:04:01,790 if you do random exchanges, but you still lose efficiency. 1280 01:04:01,790 --> 01:04:02,540 I'll tell you why. 1281 01:04:02,540 --> 01:04:08,430 First of all, you have to pick a dissimilar pair. 1282 01:04:08,430 --> 01:04:11,240 So if you think about it, half of the time 1283 01:04:11,240 --> 01:04:12,650 you'll pick a similar pair. 1284 01:04:12,650 --> 01:04:14,463 So you will pick A here and A here, 1285 01:04:14,463 --> 01:04:16,130 and so you have to give up on that move, 1286 01:04:16,130 --> 01:04:17,610 because it's not going to do anything. 1287 01:04:17,610 --> 01:04:18,443 So that's one thing. 1288 01:04:18,443 --> 01:04:21,900 That's not a major thing. 1289 01:04:21,900 --> 01:04:25,020 But the fact that you can fluctuate the concentration 1290 01:04:25,020 --> 01:04:27,090 of your system will in many cases 1291 01:04:27,090 --> 01:04:29,550 still allow you to equilibrate faster. 1292 01:04:29,550 --> 01:04:31,800 OK. 1293 01:04:31,800 --> 01:04:36,348 And usually this is not-- 1294 01:04:36,348 --> 01:04:37,890 it's a good thing you pointed it out, 1295 01:04:37,890 --> 01:04:39,307 because I haven't talked about it. 1296 01:04:41,780 --> 01:04:44,570 The choice of the ensemble should not 1297 01:04:44,570 --> 01:04:49,220 matter for thermodynamic properties. 1298 01:04:49,220 --> 01:04:51,410 So any time you're doing something 1299 01:04:51,410 --> 01:04:52,850 with periodic boundary conditions 1300 01:04:52,850 --> 01:04:55,970 and your objective is to study the bulk of something, 1301 01:04:55,970 --> 01:04:57,760 you're OK. 1302 01:04:57,760 --> 01:04:59,558 When you go to finite systems-- 1303 01:04:59,558 --> 01:05:01,600 you're studying, I don't know, the thermodynamics 1304 01:05:01,600 --> 01:05:04,290 of nanoparticles-- then the ensembles do matter. 1305 01:05:07,010 --> 01:05:08,930 The problem that you run into there 1306 01:05:08,930 --> 01:05:11,990 is that it's not obvious what the right ensemble is. 1307 01:05:11,990 --> 01:05:15,730 Like if you study a nanoparticle, now you have to-- 1308 01:05:15,730 --> 01:05:18,020 there is essentially no thermodynamic limit 1309 01:05:18,020 --> 01:05:18,770 of a nanoparticle. 1310 01:05:18,770 --> 01:05:21,063 It's a finite system by definition. 1311 01:05:21,063 --> 01:05:22,730 And you just have to sort of figure out, 1312 01:05:22,730 --> 01:05:26,400 can it exchange things within the environment or not. 1313 01:05:26,400 --> 01:05:29,930 But in general, the choice of the ensemble doesn't matter. 1314 01:05:29,930 --> 01:05:31,010 Let me sort of skip this. 1315 01:05:35,830 --> 01:05:38,830 I've given you examples on a lattice model right now. 1316 01:05:38,830 --> 01:05:43,360 But you can, of course, do off-lattice Monte Carlo. 1317 01:05:43,360 --> 01:05:45,130 For example, this the example of a liquid. 1318 01:05:45,130 --> 01:05:48,580 Let's say a liquid is a bunch of atoms in a box. 1319 01:05:48,580 --> 01:05:51,310 And again you could study the electrodynamics and track them. 1320 01:05:51,310 --> 01:05:53,500 You could also get their thermodynamic averages 1321 01:05:53,500 --> 01:05:55,480 just by sampling all the possible positions 1322 01:05:55,480 --> 01:05:56,110 of the atoms. 1323 01:05:58,660 --> 01:06:01,290 So which perturbations would you pick? 1324 01:06:01,290 --> 01:06:03,790 Again, anything that's consistent with the degrees 1325 01:06:03,790 --> 01:06:07,100 of freedom of your ensemble is fine. 1326 01:06:07,100 --> 01:06:11,650 So here an obvious way to do this would be pick an atom 1327 01:06:11,650 --> 01:06:15,580 and you then pick a random displacement vector. 1328 01:06:15,580 --> 01:06:19,030 So you say, I'm going to pick a random displacement vector, 1329 01:06:19,030 --> 01:06:20,230 and that's my perturbation. 1330 01:06:20,230 --> 01:06:22,900 I calculate the energy, it goes down, I accept. 1331 01:06:22,900 --> 01:06:25,480 If it doesn't, I don't accept it. 1332 01:06:25,480 --> 01:06:27,700 OK. 1333 01:06:27,700 --> 01:06:33,180 The question is, do you truly pick this displacement vector 1334 01:06:33,180 --> 01:06:33,870 random? 1335 01:06:33,870 --> 01:06:36,900 But what about its magnitude? 1336 01:06:36,900 --> 01:06:38,890 You probably pick its orientation random, 1337 01:06:38,890 --> 01:06:40,832 but what about its magnitude? 1338 01:06:40,832 --> 01:06:44,070 Do you pick it out of the whole box length? 1339 01:06:44,070 --> 01:06:46,180 Do you pick it out of a very small length? 1340 01:06:46,180 --> 01:06:48,020 So typically you'll do it somehow, 1341 01:06:48,020 --> 01:06:53,550 like you pick a delta x-coordinate out of range-- 1342 01:06:53,550 --> 01:07:00,560 call it minus a over 2 to plus a over 2, 1343 01:07:00,560 --> 01:07:03,890 and then the same for a delta y and z-- 1344 01:07:03,890 --> 01:07:06,410 and then you randomly pick a coordinate 1345 01:07:06,410 --> 01:07:10,400 in that range for x, y, and z. 1346 01:07:10,400 --> 01:07:13,410 But the question is, what do you take A? 1347 01:07:13,410 --> 01:07:16,200 If you take A really big, you're not 1348 01:07:16,200 --> 01:07:18,870 going to accept a lot of your perturbations. 1349 01:07:18,870 --> 01:07:21,840 Because if you take A of the order of the box length, 1350 01:07:21,840 --> 01:07:24,510 then many times if this is a dense liquid 1351 01:07:24,510 --> 01:07:26,640 your atoms is going to end up very much 1352 01:07:26,640 --> 01:07:28,860 on top of something else. 1353 01:07:28,860 --> 01:07:31,080 The density of liquids is what? 1354 01:07:31,080 --> 01:07:33,390 In many cases, 10% less than a solid. 1355 01:07:33,390 --> 01:07:35,010 At best they're still dense. 1356 01:07:35,010 --> 01:07:39,830 So most space is filled by atoms. 1357 01:07:39,830 --> 01:07:43,540 So if you take A very big, you're 1358 01:07:43,540 --> 01:07:45,850 going to sort of equilibrate in principle 1359 01:07:45,850 --> 01:07:47,090 the composition easily-- 1360 01:07:47,090 --> 01:07:48,460 I mean the density-- 1361 01:07:48,460 --> 01:07:52,180 but you're not going to accept a lot of your moves. 1362 01:07:52,180 --> 01:07:56,050 So you say, well then, the solution is I pick A small. 1363 01:07:56,050 --> 01:08:00,120 But if you pick A too small, then you're probably going-- 1364 01:08:00,120 --> 01:08:01,870 you have a high acceptance rate. 1365 01:08:01,870 --> 01:08:04,105 So if you just displace A-- 1366 01:08:04,105 --> 01:08:08,470 if A is like 100 of an Angstrom, then you're 1367 01:08:08,470 --> 01:08:10,070 atom is just moving a very little bit. 1368 01:08:10,070 --> 01:08:12,130 So most time, you're going to accept it. 1369 01:08:12,130 --> 01:08:14,350 But what's the problem now is you have high degree 1370 01:08:14,350 --> 01:08:16,018 of correlation in your average. 1371 01:08:16,018 --> 01:08:18,310 Essentially all the states you put in your Markov chain 1372 01:08:18,310 --> 01:08:20,790 look all the same. 1373 01:08:20,790 --> 01:08:23,290 If you barely this place your atom, they all look the same. 1374 01:08:23,290 --> 01:08:26,880 So it's going to take a really long time to get good averages. 1375 01:08:26,880 --> 01:08:30,450 Actually, you could almost think of Monte Carlo with very 1376 01:08:30,450 --> 01:08:32,819 small perturbations on the system, 1377 01:08:32,819 --> 01:08:35,310 it's almost like molecular dynamics. 1378 01:08:35,310 --> 01:08:37,950 It's like the atoms sort of move slowly 1379 01:08:37,950 --> 01:08:39,930 with a kind of random force, but they still 1380 01:08:39,930 --> 01:08:42,689 move slowly along trajectories and are real trajectories, 1381 01:08:42,689 --> 01:08:46,368 but they don't wildly go through phase space. 1382 01:08:46,368 --> 01:08:47,910 The nice thing about models like this 1383 01:08:47,910 --> 01:08:51,930 is that you can actually, during the simulation, adjust A. 1384 01:08:51,930 --> 01:08:53,310 And this is often what's done. 1385 01:08:53,310 --> 01:08:55,649 Often people will adjust it so that they get 1386 01:08:55,649 --> 01:08:58,050 a reasonable acceptance rate. 1387 01:08:58,050 --> 01:09:02,427 If you set A too big, then you get low acceptance rate. 1388 01:09:02,427 --> 01:09:05,010 If you set A too small, you get a high acceptance rate people. 1389 01:09:05,010 --> 01:09:06,745 Typically shoot for anywhere from 0.3 1390 01:09:06,745 --> 01:09:09,399 to 0.5 acceptance rate. 1391 01:09:09,399 --> 01:09:12,810 So you want to be somewhere where the order of a 1/3 to 1/2 1392 01:09:12,810 --> 01:09:15,312 of the moves you try are accepted. 1393 01:09:15,312 --> 01:09:16,770 Because if they don't get accepted, 1394 01:09:16,770 --> 01:09:20,200 you're wasting all your effort calculating energies. 1395 01:09:20,200 --> 01:09:22,200 And this depends a little bit on the cost 1396 01:09:22,200 --> 01:09:24,569 of your energy function. 1397 01:09:24,569 --> 01:09:27,337 If you're doing this with quantum mechanics, 1398 01:09:27,337 --> 01:09:29,670 then you really want to optimize your Monte Carlo moves, 1399 01:09:29,670 --> 01:09:33,729 because any energy evaluation costs you a lot of time. 1400 01:09:33,729 --> 01:09:36,580 If you're doing this with maybe some simple hard sphere 1401 01:09:36,580 --> 01:09:40,720 model or something, you can afford to have a lot of moves-- 1402 01:09:40,720 --> 01:09:43,450 calculate a lot of energies on which you don't execute. 1403 01:09:55,770 --> 01:09:57,670 So the last thing I want to say is, 1404 01:09:57,670 --> 01:10:00,540 how do you pick random numbers? 1405 01:10:00,540 --> 01:10:03,120 How do you implement probabilities? 1406 01:10:03,120 --> 01:10:05,730 Because it's sort of not obvious in a computer. 1407 01:10:05,730 --> 01:10:08,640 Computers are not probabilistic. 1408 01:10:08,640 --> 01:10:12,150 And that's done by random number generation. 1409 01:10:12,150 --> 01:10:16,470 Sort of maybe obvious, but I want to show it anyway. 1410 01:10:16,470 --> 01:10:18,390 So remember, what you want to do is 1411 01:10:18,390 --> 01:10:22,410 you want to implement a move with a given probability 1412 01:10:22,410 --> 01:10:24,180 with this Boltzmann factor. 1413 01:10:24,180 --> 01:10:25,830 Remember, if delta E is positive, 1414 01:10:25,830 --> 01:10:29,160 that Boltzmann factor is between 0 and 1. 1415 01:10:29,160 --> 01:10:31,570 So let's say this is the value. 1416 01:10:31,570 --> 01:10:34,770 Well, what you do is you pick random numbers between 0 and 1. 1417 01:10:37,490 --> 01:10:40,490 And what's the probability that the random number 1418 01:10:40,490 --> 01:10:43,150 is less than this exponential? 1419 01:10:43,150 --> 01:10:47,540 Well, that probability is that exponential of minus beta E i. 1420 01:10:47,540 --> 01:10:49,390 We'll call this p i. 1421 01:10:49,390 --> 01:10:53,020 If you take a random number and the probability 1422 01:10:53,020 --> 01:10:58,140 that that random number is below p i is p i. 1423 01:10:58,140 --> 01:11:03,880 So you pick the random number if it's below that exponential. 1424 01:11:03,880 --> 01:11:04,780 You execute. 1425 01:11:04,780 --> 01:11:06,910 If it's above that exponential, you don't execute. 1426 01:11:06,910 --> 01:11:09,430 So you execute now with the correct probability. 1427 01:11:09,430 --> 01:11:15,570 So random numbers give you a way of implementing probabilities. 1428 01:11:15,570 --> 01:11:18,210 There's a few caveats. 1429 01:11:18,210 --> 01:11:20,160 A lot of cheap random number generators 1430 01:11:20,160 --> 01:11:22,050 are not at all random. 1431 01:11:22,050 --> 01:11:25,165 I used to not believe this until I tried it. 1432 01:11:25,165 --> 01:11:27,540 I always thought that was kind of like people write books 1433 01:11:27,540 --> 01:11:29,940 and they just talk about it because it's kind 1434 01:11:29,940 --> 01:11:31,410 of cool to say how non-random. 1435 01:11:31,410 --> 01:11:34,560 They actually are very non-random. 1436 01:11:34,560 --> 01:11:39,030 You will often find that if you plot 1437 01:11:39,030 --> 01:11:41,250 the density of random numbers you generate, 1438 01:11:41,250 --> 01:11:44,220 that there are discrete holes where 1439 01:11:44,220 --> 01:11:45,470 you have very little coverage. 1440 01:11:45,470 --> 01:11:47,678 And it's because of the sort of fractional algorithms 1441 01:11:47,678 --> 01:11:49,760 by which they work. 1442 01:11:49,760 --> 01:11:52,550 In many cases, that really doesn't matter. 1443 01:11:52,550 --> 01:11:54,650 For a lot of practical things, it doesn't matter. 1444 01:11:54,650 --> 01:11:58,530 I'll tell you where I've seen it matter. 1445 01:11:58,530 --> 01:12:03,230 It matters when you're in this regime here, really close to 0. 1446 01:12:03,230 --> 01:12:05,600 Because essentially, a random number 1447 01:12:05,600 --> 01:12:08,000 generator, the way it generates random numbers between 0 1448 01:12:08,000 --> 01:12:08,420 and 1-- 1449 01:12:08,420 --> 01:12:09,795 it doesn't do it between 0 and 1. 1450 01:12:09,795 --> 01:12:12,440 It generates between 0-- 1451 01:12:12,440 --> 01:12:14,943 or is it-- between 0 or 1. 1452 01:12:14,943 --> 01:12:16,235 No, between 1 and some integer. 1453 01:12:16,235 --> 01:12:19,970 No, sorry, between 0 and some maximum integer, like 2 1454 01:12:19,970 --> 01:12:22,440 to the 16 or something like that. 1455 01:12:22,440 --> 01:12:24,980 And then it divides by 2 to 16. 1456 01:12:24,980 --> 01:12:27,140 And that's how you get it between 0 and 1. 1457 01:12:27,140 --> 01:12:31,130 But that means your set is discretized. 1458 01:12:31,130 --> 01:12:36,410 So often what's important is the spacing between 0 1459 01:12:36,410 --> 01:12:40,190 and the smallest random number. 1460 01:12:40,190 --> 01:12:43,520 Because if you work with very high energy 1461 01:12:43,520 --> 01:12:46,430 barriers, very high energy excitations, 1462 01:12:46,430 --> 01:12:48,590 they fall in that first gap. 1463 01:12:48,590 --> 01:12:54,230 And so you have a completely discretized probability there. 1464 01:12:54,230 --> 01:12:56,570 And there are problems. 1465 01:12:56,570 --> 01:12:59,060 It's the kind of thing I thought I'd never run into, 1466 01:12:59,060 --> 01:13:00,830 and I've actually run into several times 1467 01:13:00,830 --> 01:13:02,210 within that [INAUDIBLE] research. 1468 01:13:05,420 --> 01:13:07,760 If your system is stuck and it doesn't 1469 01:13:07,760 --> 01:13:10,250 have a lot of low energy excitations left, 1470 01:13:10,250 --> 01:13:13,430 it tends to execute high energy excitations with a way 1471 01:13:13,430 --> 01:13:15,330 too higher probability. 1472 01:13:15,330 --> 01:13:16,330 Because think about it-- 1473 01:13:16,330 --> 01:13:19,330 if the biggest number is 2 to the 16th, 1474 01:13:19,330 --> 01:13:24,630 your random number hits 0 with a probability of 1 over 2 1475 01:13:24,630 --> 01:13:26,780 to the 16th. 1476 01:13:26,780 --> 01:13:27,890 That's way too high. 1477 01:13:27,890 --> 01:13:30,190 It should be 1 over infinity. 1478 01:13:30,190 --> 01:13:34,000 You should hit 0 with a rate of 1 over infinity 1479 01:13:34,000 --> 01:13:36,340 if you had a perfect continue of random numbers, 1480 01:13:36,340 --> 01:13:40,250 but you hit 0 with a rate 1 over 2 to the 16. 1481 01:13:40,250 --> 01:13:43,120 That means if you hit that random number you 1482 01:13:43,120 --> 01:13:47,010 execute any move, even if it's 100 electronvolt 1483 01:13:47,010 --> 01:13:49,140 above the ground state. 1484 01:13:49,140 --> 01:13:51,060 Because that exponential is not 0. 1485 01:13:51,060 --> 01:13:54,540 It's small, but it's not 0, So you will execute. 1486 01:13:54,540 --> 01:13:57,730 And you'd be surprised how often this is a problem, 1487 01:13:57,730 --> 01:13:59,250 especially if you have a wild phase 1488 01:13:59,250 --> 01:14:01,470 space in which you can go. 1489 01:14:01,470 --> 01:14:03,360 So how do you solve this? 1490 01:14:03,360 --> 01:14:04,770 We solve it in a trivial way. 1491 01:14:04,770 --> 01:14:09,050 We add a small epsilon to a random number. 1492 01:14:09,050 --> 01:14:11,290 Which is so, trivial but you have to know it I. 1493 01:14:11,290 --> 01:14:13,930 I've never seen this in books actually. 1494 01:14:13,930 --> 01:14:17,130 But see you should not hit 0. 1495 01:14:17,130 --> 01:14:20,410 Because when you hit 0, you execute everything. 1496 01:14:20,410 --> 01:14:22,190 OK. 1497 01:14:22,190 --> 01:14:22,690 OK. 1498 01:14:22,690 --> 01:14:24,130 So I'm going to stop here. 1499 01:14:24,130 --> 01:14:26,170 And let's see, today's Thursday. 1500 01:14:26,170 --> 01:14:30,060 So yeah, I'll see you on Tuesday for the rest of Monte Carlo.