boltz misc support

Boltzmann Misc support stuff PhL 6.6.15 This was taken from the end of the boltz shot 2 document. _______________________________________________________________________________ How do we know that our solution {Ni} corresponds to a maximum of Ω or f ? ************************************************** Consider now the first and second derivatives of f, f = ln Ω = Σi [(1+lngi)Ni - NilnNi ] // (6.2.3) fi = ln(gi/Ni) = lngi - lnNi // (6.2.6) fii = -1/Ni plots ok H = f + λ1a + λ2 b // (6.2.5) Hi = fi + λ1ai + λ2 bi = fi + λ1*1 + λ2*εi // (6.2.6) Hii = fii = -1/Ni Rcurv = 1/Hii = -Ni am happy to here Ω = ef Ωi = Ω fi plot has peak, is OK Ωii = Ωfii+ Ωifi = Ωfii+ Ωfi2 = Ω (fii + fi2) = Ω ( -1/Ni + [ ln(gi/Ni) ]2 ) h = Ω + λ1'a + λ2' b hi = Ωi + λ1'ai + λ2' bi = Ωi + λ1*1 + λ2*εi hii = Ωii Therefore hii = curvature of the Ω curve = Ωii = Ω ( -1/Ni + [ ln(gi/Ni) ]2 ) But this appears to me to be a large positive number so the h curve is cupping up which is wrong! BUT, hii is not really ∂2h/∂2Ni whereas the curvature is given by d2h/d2Ni which is a completely different function! Realization: Consider f(N1,N2,N3) = ln Ω = Σi=13 [(1+lngi)Ni - NilnNi ] // (6.2.3) The object fi is a partial derivative! fi = ∂f/∂Ni = ln(gi/Ni) = lngi - lnNi But the curve I am plotting for f which shows the bell shape is this F(N1) = f(N1,N2(N1),N3(N1)) The slope of the plot is then dF/dN1 = ∂f/∂N1 + ∂f/∂N2 * ∂N2/∂N1 + ∂f/∂N3 * ∂N3/∂N1 = f1 + f2(∂N2/∂N1) + f3(∂N3/∂N1) = slope of the bell curve plot of f Now I solved to get N2 = - (ε2-ε3)-1[ (e1-e3)N1 - U + Me3] N3 = + (ε2-ε3)-1[ (e1-e2)N1 - U + Me2] Therefore (∂N2/∂N1) = - (ε2-ε3)-1[ (e1-e3) = - (-1)-1(-2) = 2 (∂N3/∂N1) = + (ε2-ε3)-1[ (e1-e2) = + (-1)-1(-1) = 1 Therefore the curve slope is dF/dN1 = f1 + 2 f2 + f3 where in all functions we replaced N2 and N3 as shown. So this is a far cry from f1≥ Continue on this line 6/5/15 Let's think about f or Ω by using the fact that the constraints are simple to just remove two variables and then we have something we can plot and deal with. Consider N1 + N2 + Σn=3m Ni = M ε1N1 + ε2N2 + Σn=3m εiNi = U or N1 + N2 = M - Σn=3m Ni ε1N1 + ε2N2 = U - Σn=3m εiNi Solve for N1 and N2 as follows ε1N1 + ε1N2 = ε1M - ε1Σn=3m Ni ε1N1 + ε2N2 = U - Σn=3m εiNi Subtract (ε1- ε2) N2 = ε1M - ε1Σn=3m Ni - U + Σn=3m εiNi = ε1M - U + Σn=3m (εi - ε1)Ni so N2 = (ε1- ε2)-1 [ ε1M - U + Σn=3m (εi - ε1)Ni] = N2(N3, N4....Nm) Now do it the other way ε2N1 + ε2N2 = ε2M - ε2Σn=3m Ni ε1N1 + ε2N2 = U - Σn=3m εiNi Subtract bottom from top (ε2- ε1)N1 = ε2M - ε2Σn=3m Ni - U + Σn=3m εiNi = ε2M - U + Σn=3m(εi- ε2)Ni so then N1 = (ε2- ε1)-1[ ε2M - U + Σn=3m(εi- ε2)Ni] Now assume ε2 > ε1 just by convention. Then write N1 = + (ε2- ε1)-1[ ε2M - U + Σn=3m (εi- ε2)Ni] N2 = - (ε2- ε1)-1[ ε1M - U + Σn=3m (εi - ε1)Ni] Notice that ∂N1/∂N3 = + (ε2- ε1)-1 [ (ε3- ε2) ] ∂N2/∂N3 = - (ε2- ε1)-1 [ (ε3- ε1) ] and in general ∂N1/∂Ni = (ε2- ε1)-1 [ (εi- ε2) ] i = 3,4....m ∂N1/∂Ni = - (ε2- ε1)-1 [ (εi- ε1) ] So this is the first time I have done this. Now let F be some generic function of all the Ni. Then we have F(N1(N3...Nm), N2(N3...Nm), N3, N4....Nm) dF/dN3 = ∂F/∂N1 ∂N1/∂N3 + ∂F/∂N2 ∂N2/∂N3 + ∂F/∂N3 = F1 (ε2- ε1)-1 [ (ε3- ε2) ] + F2 ( - (ε2- ε1)-1 [ (ε3- ε1) ]) + F3 This certainly has an interesting form. Rewrite dF/dN3 = F1 - F2 + F3 More generally, dF/dNi = ∂F/∂N1 ∂N1/∂Ni + ∂F/∂N2 ∂N2/∂Ni + ∂F/∂Ni i = 3,4...m so dF/dNi = F1 - F2 + Fi i = 3,4...m and notice where the i's are located on the right side! So, if you wanted to find the critical points of this unconstrained function F(N3,....Nm) you would set all these derivatives to 0 0 = F1 - F2 + Fi i = 3,4...m Let's try an example. Suppose F =f f = Σi [ Ni ln gi ] - Σi [ ln Ni! ] We then have these equations to solve for the max: fi = ln(gi/Ni) = lngi - lnNi // from above 0 = f1 - f2 + fi or 0 = ln(g1/N1) - ln(g2/N2) + ln(gi/Ni) = 0 i = 3,4...m But N1 = + (ε2- ε1)-1[ ε2M - U + Σn=3m (εi- ε2)Ni] N2 = - (ε2- ε1)-1[ ε1M - U + Σn=3m (εi - ε1)Ni] Suppose we define si ≡ ri ≡ and then our equations are 0 = si ln(g1/N1) - ri ln(g2/N2) + ln(gi/Ni) = 0 i = 3,4...m or 0 = ln(g1/N1)s + ln(g2/N2)-r + ln(gi/Ni) = 0 i = 3,4...m or 0 = ln [ (g1/N1)s (g2/N2)-r(gi/Ni) ] = 0 i = 3,4...m But this says (g1/N1)s (g2/N2)-r(gi/Ni) = 1 i = 3,4...m But these are very complicated equations because we have to use N1 = + (ε2- ε1)-1[ ε2M - U + Σn=3m (εi- ε2)Ni] N2 = - (ε2- ε1)-1[ ε1M - U + Σn=3m (εi - ε1)Ni] So maybe this just shows how difficult this problem is to solve by this brute force method! So this thread leads nowhere useful. Back to the idea of "moving off the peak" Go back to f = Σi=1m [(1+lngi)Ni - NilnNi ] (6.2.3) fi = lngi - lnNi (6.2.6) Suppose we are at some vector {Ni} (not necessarily the peak) and we "move a little". What happens to f? I already know that N1 + N2 + N3 + Σn=4m Ni = M ε1N1 + ε2N2 + ε3N3 + Σn=4m εiNi = U Suppose we concoct a "movement" such that only N1, N2 and N3 change and all other Ni stay constant. Then one surely can write these 2 equations in 3 unknowns, dN1 + dN2 + dN3 = 0 ε1dN1 + ε2dN2 + ε3dN3 = 0 Now do as was done above to solve for dN2 and dN3 both in terms of dN1 ε1dN1 + ε1dN2 + ε1dN3 = 0 ε1dN1 + ε2dN2 + ε3dN3 = 0 Subtract (ε1- ε2)dN2 + (ε1- ε3)dN3 = 0 Then we can replace dN3 = - dN2 and therefore dN1 + dN2 + dN3 = 0 dN1 + dN2 - dN2 = 0 dN1 = dN2 ( - 1 ) = dN2( - ) = dN2 and this shows that dN2 = dN1 Then we know that dN3 = - dN2 = - dN1 = dN1 Now assume that ε1 is the lowest energy, so write with all positives dN2 = - dN1 dN3 = dN1 This says that if you make a free small movement dN1, the other two variables have to change as shown The first ratio is greater than 1, whereas the second could be < or > 1 depending on energy spacing, where again I assume ε1 is lowest. NOW consider: f(N1,N2,N3....Nm) = Σi=1m [(1+lngi)Ni - NilnNi ] (6.2.3) We are sitting at some general vector N. We now make a change {dN1. dN2, dN3} consistent with the rules above. What happens to F? df = f1dN1 + f2dN2 + f3dN3 // true only for differentials = ( lng1 - lnN1) dN1 + ( lng2 - lnN2) [- dN1 ] + ( lng3 - lnN3) [ dN1] = { ( lng1 - lnN1) + ( lng2 - lnN2) [- ] + ( lng3 - lnN3) [] } dN1 = - { ln(N1/g1) - ln(N2/g2) + ln(N3/g3) } dN1 Notice how i = 4,5...m don't enter into this "movement" dN where only 1,2,3 move. If our vector N is in some arbitrary location, there is some df as above. But if N is the solution vector, then we know that ln(Ni/gi) = λ1+ λ2εi = λ1 - βεi (6.2.9) So if we start with the solution vector N we find that df = - { (λ1 - βε1) - (λ1 - βε2) + (λ1 - βε3) } dN1 = - { λ1 [ 1 - + ] - β [ ε1 - ε2 + ε3] } dN1 = - { λ1 [ - + ] - β [ ε1 - ε2 + ε3] } dN1 = - { λ1 [ 0 ] - β [ 0 ] } dN1 = 0 This is the expected answer because we are at a critical point. So this more or less validates what I am doing here. Let's now try for the curvature!! Start as above to get the first derivative Now let F be some generic function of all the Ni. Then we have F(N1(N3...Nm), N2(N3...Nm), N3, N4....Nm) dF/dN3 = ∂F/∂N1 ∂N1/∂N3 + ∂F/∂N2 ∂N2/∂N3 + ∂F/∂N3 = F1 (ε2- ε1)-1 [ (ε3- ε2) ] + F2 ( - (ε2- ε1)-1 [ (ε3- ε1) ]) + F3 This certainly has an interesting form. Rewrite dF/dN3 = F1 - F2 + F3 More generally, dF/dNi = ∂F/∂N1 ∂N1/∂Ni + ∂F/∂N2 ∂N2/∂Ni + ∂F/∂Ni i = 3,4...m so dF/dNi = F1 - F2 + Fi i = 3,4...m Now consider dF1/dNi = F11 - F12 + F1i Fij are all second partials! dF2/dNi = F21 - F22 + F2i Fij are all second partials dFi/dNi = Fi1 - Fi2 + Fii no sum on i of course Then I get d2F/dNi2 = [ F11 - F12 + F1i] - [ F21 - F22 + F2i] + [ Fi1 - Fi2 + Fii] i = 3,4....M Now at once apply this to F = f where we know that fi = lngi - lnNi (6.2.6) fij = ∂fi/∂Nj = -δij (1/Ni) We then write (since i = 3,4....m) d2f/dNi2 = [ f11 ] - [ - f22 ] + [ + fii] = f11 + f22 + fii = (-1/N1) + (-1/N2) + (-1/Ni) = – { [ ]2 + []2 + } i = 3,4....m Amazingly this comes out definitely negative, meaning cupping down, which is correct. This is a first! Here is the conclusion: d2f/dNi2 = – { [ ]2 + []2 + } i = 3,4....m for any vector N, not just the solution vector. In general, I would say d2f/dNi2 ~ - (1/Nk ) where Nk is one of the vectors , just order of Now consider our earlier calculation where NOW I am taking about FULL derivatives, so be very careful with notation!! Ω = ef Ωi = Ω fi Ωii = Ωfii+ Ωifi = Ωfii+ Ωfi2 = Ω (fii + fi2) i = 3,4...m which really means (d2Ω/dNi2) = Ω [ (d2f/dNi2) + (df/dNi)2 ] Now at the peak, we know (df/dNi) = 0 and we even showed it above. so therefore (d2Ω/dNi2) = – Ω { [ ]2 + []2 + } i = 3,4....m and this is a very major result!! (1) it is negative so peak is really a maximum (2) it is a very large number because although 1/Ni is small, Ω is extremely large! In our example, the above would read (N3 is the only other Ni) (d2Ω/dN32) = – Ω { [ ]2 + []2 + } = – Ω { + 4 + } ≈ – Ω { 6/333} = very large so the corresponding curvature radius at the peak would be very small. So exactly how large is Ω at the peak? f = Σi=1m [(1+lngi)Ni - NilnNi ] (6.2.3) Ni = gi eλ eλε . (6.2.9) So write Ni = A gi e-βε Then we have f = Σi=1m [(1+lngi)Ni - NilnNi ] = Σi Ni [ (1+lngi) - lnNi ] = Σi=1m [(1+lngi)( A gi e-βε) - ( A gi e-βε) ln( A gi e-βε) ] ≈ Σi [lngi( A gi e-βε) - ( A gi e-βε) (ln A + lngi - β εi) ] = Σi ( A gi e-βε) [lngi- (ln A + lngi - β εi) ] = Σi ( A gi e-βε) [lngi- ln A - lngi + β εi ] = Σi ( A gi e-βε) [- ln A + β εi ] = Σi Ni [- ln A + β εi ] = -ln(A) M + β U where A = . (6.2.14) But from just the first line I know that f ≈ Σi Ni ln(gi/Ni) + Σi Ni Now suppose we just assume that ln(gi/Ni) = ln(5000/333) = 2.71. Then f1 ≈ 2.709 M = 2.709* 1000 = 2709 f2 ≈ M = 1000 f = f1+ f2 = 2709 + 1000 = 3709 Then Ω = ef = e3709 ≈ (1/2) * 101611 //Maple The Maple plot shows directly that the peak value is 250,000 * 101600 = 2.5 * 101605 So why am I off by something like 100,000 ??? Probably Sterling. Correct formula is f = lnΩ = Σi [ Ni lngi ] - Σi [ ln Ni! ] In my test case we have Ni = 333 and so f = Σi [ Ni ln(5000) ] - Σi [ ln Ni! ] = 1000*ln(5000) - 3*ln[(1000/3)!] = 8517 - 4820 = 3697 Ω = ef = e3697 = 0.4 * 101606 OK, so now I know how large Ω is at least in my example. Now go back to curvature (d2Ω/dNi2) = – Ω { [ ]2 + []2 + } i = 3,4....m What is Ω for a more real-world example? Suppose we just assume that gi = 107Ni as a ballpark number. And maybe Ni = 1020 as a ballpark number for a small amount of gas. Then ≈ (107N1) well not sure this helps anyone. I think best to stick to the example. I notice from my plot that the curvature at the peak made into a semicircle does NOT give a good estimate for the width of the peak. So, again sticking with the same, I have Ω(N1,N2,N3) = N2 = -2N1-U+3M = -2N1+M U = 2M for ∞ temp N3 = N1 - 2M + U = N1 So there you have it Ω(N1) ≡ Ω(N1,N2(N1),N3(N1)) = = I plotted this with g = 5000, M = 1000 as in example and the plot comes out what I expect. Now, where is the half-height of such a plot?? I guess the gM factor cancels and we can write 1/2 = Ωhalf(N1)Ω(peak) = [(N1)!(N1)!(M-2N1)!] / [(333)!(333)!(333)!] So you have to solve [(333)!]3 (1/2) = [(N1)!]2 (1000-2N1)!] What a pain! Let's try it in terms of f instead. f = Σi=1m [(1+lngi)Ni - NilnNi ] (6.2.3) = Σi=1m Ni[ (1+lngi) - lnNi] If the gi are the same this says f = (1+lng) M - Σi=1m Ni lnNi = (1+lng) M - 2 N1lnN1 - (M-2N1)ln (M-2N1) and the plot looks right! So consider Ωh/Ωp = 1/2 => fh - fp = ln(1/2) = - ln(2) = fp - fh = ln(2) = 0.7 So f drops about 1 single unit. Let's write the equation fp = (1+lng) M - 2 N1plnN1p - (M-2N1p)ln (M-2N1p) fh = (1+lng) M - 2 N1hlnN1h - (M-2N1h)ln (M-2N1h) fp - fh = - 2 N1plnN1p - (M-2N1p)ln (M-2N1p) - [ - 2 N1hlnN1h - (M-2N1h)ln (M-2N1h)] = - 2 N1plnN1p - (M-2N1p)ln (M-2N1p) + 2 N1hlnN1h + (M-2N1h)ln (M-2N1h)] = -(667)ln(333) - (333) ln(333) + 2 N1hlnN1h + (M-2N1h)ln (M-2N1h)] = 0.7 Again I see no intuitive solution! Suppose I assume the log terms are the same for half and peak. Then -(667)ln(333) - (333) ln(333) + 2 N1hln(333) + (M-2N1h)ln (333)] = 0.7 [ -667 - 333 + 2N1h + (1000 - 2N1h} ln(333) = 0.7 // answer is lost! OK, write Nh = Np + d d is very small lnN1h = ln(N1p + d) = ln[(N1p)(1 + d/N1p)] = lnN1p ln(1 + d/N1p) ≈ lnN1p * (d/N1p) ln (M-2N1h) = ln (M-2N1p - 2d) = ln(M-2N1p)* [-2d/(M-2N1p) ] Then fp - fh = - 2 N1plnN1p - (M-2N1p)ln (M-2N1p) - [ - 2 N1hlnN1h - (M-2N1h)ln (M-2N1h)] = - 2 N1plnN1p - (M-2N1p)ln (M-2N1p) + 2N1hlnN1h + (M-2N1h)ln (M-2N1h) ≈ - 2 N1plnN1p - (M-2N1p)ln (M-2N1p) STOP (maybe later). Go back to the curvature of f: d2f/dN32 = – { [ ]2 + []2 + } i = 3,4....m = – { + 4 + } = - 6/333. Very good. Now write f(N1p + d) = f(N1p) + (1/2) f" (d)2 or fh = fp + 1/2 * d2 * [-6/333] or fh = fp - (3/333)d2 or fp - fh = 0.7 = d2/111 => d2 = 111*.7 = 78 = expected half distance away from peak d = 9 so the full half width is then about 18. And this looks VERY good! For general M and our simple model problem, we get N = M/3 and so f" = - 6/N = -18/M fh = fp + (1/2) d2 [ -18/M] fp - fh = (1/2) d2 [ 18/M] = 0.7 0.7 = (9/M)d2 d2 = 0.7 M/9 = .078 M d = .28 half width of bell = .56 as a percentage ΔN/N = .56/ [M/3] = .56*3 / = 1.68/ So if M = 1000 this says ΔN/N = .053 check: 18/333 = .054 So if M = 1020 then ΔN/N = 1.68/ 1010 = 1.68 x 10-10 very narrow peak!~