COMPUTATIONAL PHYSICS PHZ 5156
University of Central Florida
Popular in Course
Popular in Physics 2
verified elite notetaker
This 65 page Class Notes was uploaded by Humberto Romaguera II on Thursday October 22, 2015. The Class Notes belongs to PHZ 5156 at University of Central Florida taught by Patrick Schelling in Fall. Since its upload, it has received 36 views. For similar materials see /class/227503/phz-5156-university-of-central-florida in Physics 2 at University of Central Florida.
Reviews for COMPUTATIONAL PHYSICS
Report this Material
What is Karma?
Karma is the currency of StudySoup.
Date Created: 10/22/15
Declarations I program to relax a vacancy in a LJ solid IMPLICIT NONE INTEGER PARAMETER nn4 INTEGER PARAMETER Prec14SELECTEDREALKIND14 REALKINDPrec14 PARAMETER L20dO20dO30d0nnftol01d0 INTEGER PARAMETER natoms4nn31ndim3natoms Add a vacancy INTEGER PARAMETER nxnnnynnnznnnrmax50nrmax210nstep20 INTEGER ijixiyizmunumnmatzm1n1ierrnrnr2 REALKINDPrec14 PARAMETER sigma10dOepsilon10dO REAL KINDPrecI4 PARAMETER dt0000002dO REAL KINDPrecI4 PARAMETER rcut30dOsigma REAL KINDPrec14 epotdvdrdvdr2pifaclambdagammadenomnumapha REAL KINDPrec14 sxijsyijszijrijsovrmag REAL KINDPrec14 DIMENSION4 r1r2r3 REAL KINDPrec14 DIMENSION3 5 REAL KINDPrec14 DIMENSIONnatoms rxryrz REAL KINDPrec14 DIMENSIONndimndim hess REAL KINDPrec14 DIMENSIONndim ghbugrad REAL KINDPrecI4 PARAMETER mass39948d0eps1200d0sar34d0 AAAAAAAAAA VVVVVVVVVV Generate the crystal lattice I now construct lattice nO doix1nx doiy1ny doiz1nz doi14 nn1 ifnenatoms then I last atom not placed vacancy site rxnr1idbleix1sigma rynr2idbleiy1sigma rznr3idbleiz1sigma end enddo enddo enddo enddo l scale lengths by box size rxrxnx ryryny rzrznz 9903 General idea Loop on number of relaxation steps Compute Hessian for current configuration Compute and output potential energy Is it going down Given Hessian relax along conjugate gradient directions After several relaxations go back and compute Hessian again Generate the hessian do nr1nrmax Number of times to compute hessian Clear the hessian out hess00dO epot00dO set potential energy to zero g00dO h00dO do i1 natoms1 doji1natoms sxijrxirxj SyiJ39ryiryJ szijrzirzj sxijsxijdbleanintsxij syijsyijdbleanintsyij szijszijdbleanintszij sxijsxijL syijsyijL szijszijL rijdsqrtsxij2syij2szij2 s1sxijrij s2syijrij s3szijrij ifrij trcut then sovrsigmalrij epotepot40dOepsionsovr12sovr6 dvdr240dOepsionrij220dOsovr12sovr6 dvdr2240dOepsionrij2260dOsovr1270d0sovr6 Generate the hessian continued Get initial relaxation vectors gh m1i 13mu gm1gm1 dvdrrijsmu n1 j1 3mu gn1gn1dvdrrijsmu do nu13 mm1 nj13nu ifmueqnu hessmndvdr2smusnudvdrsmusnu10dO ifmunenu hessmndvdr2smusnudvdrsmusnu hessnmhessmn Matrix needs to be symmetric and we are only using each m1i 13mu n1i 13nu hessm1n1hessm1n1hessmn m1 j1 3mu n1 j1 3nu hessm1n1hessm1n1hessmn enddo enddo end enddo enddo End loop over all atoms Recall the CG approach We generate successive g h vectors 3N dimensional gt gt i1 gi AiA hi hi1 i1 Viki Relax along direction given by current hvector Subsequent directions satisfy h A hj llti This way relaxations don t interfere Take initial direction from g0 11O VU This is just the forces on the atoms in their unrelaxed state Parameters for relaxation starting vectors A gt ggtl39 ggtigt yi i1ltgi1 hi A hi gr 3139 get initial direction h hg Save the intial gradient and coordinates bh output current potential energy write6 39Hessian computed potential energy39epot Relax along directions h simple step approach u00dO l Atomic displacements do nr21nrmax2 l Relax along current direction h start with zero displacement do n1nstep l Compute current gradient from b and Hessian gradb mag00dO l magnitude of hvector do i1 ndim magmaghihi do j1 ndim gradigradihessijuj enddo enddo I Find projection of force along the h direction alpha00dO do i1 ndim alphaalphagradihimag enddo Continuned relaxation along directions h uuaphadth enddo end of relaxation steps From stepping along h we get close to minimum in that direction Next need a new conjugate direction h Should not interfere with relaxations we have already done Evolve new 9 h vectors Compute parameter lambda for CG algorithm compute lambda 100 denom00d0 num00d0 do i1ndim numnumgigi do j1 ndim denomdenomhihessijhj enddo enddo lambdalambdaldenom New gvector for CG algorithm Determine new 9 vector do i1ndim do j1 ndim gigiIambdahessijhj enddo enddo New parameter gamma for CG algorithm Compute gamma gamma00d0 do i1ndim gammagammagigi enddo gammagammanum New hvector for CG algorithm Determine new h vector hggammah enddo end of loop for number of conjugate gradient directi Finally move atoms return to compute new Hessian move atoms using u vector do i1natoms J39i13 rxirxiuj1 ryiryiuj2 rzirziuj3 enddo enddo end of loop for number nr of times we compute hessi Vacancy formation energy Remove atom place it into perfect crystal lattice Total energy of N atom perfect crystal EN Total energy of N1 atom crystal with vacancy EN1 1 E EN1EN N VClC Typical value 10eV fcc metal fexpEvackT exp40 103918 Might be as much as 10396 for some systems F Mitetemperature equation of state 2sinh M 2kg Compute vibrational modes frequencies Evaluate at a given volume V Compute F at various temperatures T with above FVT U0 kBTEIn A Need the Hessian Start again with fcc lattice Atoms should be relaxed zero force on atoms Evaluate at a given volume V Compute F at various temperatures T and volumes V Generate nonprimitive basis 4 atoms start with nonprimitive fcc basis r1 100d0 r2100d0 r3100d0 r1 205d0 r2205d0 r3200d0 r1 305d0 r2300d0 r3305d0 r1 400d0 r2405d0 r3405d0 Some parameters for lattice INTEGER PARAMETER nn4 INTEGER PARAMETER Prec14SELECTEDREALKIND14 INTEGER PARAMETER natoms4nn3 INTEGER PARAMETER nxnnnynnnznn REALKINDPrec14 PARAMETER sigma10d0epsilon10d0 4x4x4 lattice of nonprimitive units 64x4256 atoms System has 256x3768 degrees of freedom The Hessian matrix needs to be 768x768 Generate the full lattice n0 doix1nx doiy1ny doiz1nz doiL4 nn1 rxnr1idbeix1sigma rynr2idbeiy1sigma rznr3idbleiz1sigma enddo enddo enddo enddo scale lengths by box size rxrxnx VFW y rzrznz We need more than forces this time Definition of Hessian matrix elements dZU 613M6er iMjv For N atoms this means a 3Nx3N matrix Many elements will be zero due to cutoff First derivative isjust related to force on j E rjvrkv drjv k ark rkj But need higher order derivatives We need more than forces this time 2 0quot U My 613M6er For ij referring to different atoms 072U 072Uriyrjyrtvrjv 2 2 aw 0713 Izj 0712 3 7 zj rm rmxriv rjv Also need terms where ij zU E 2 Urm rkyriv rkv 2 2 2 6r r rik r k arik in iv kzi l rm rkyriv rkv 3 1m rik rik Notice that for each pair of components My we have EjHiMajV 0 or Himiv 2 Himjv j i Can be easily tested in code Also Hessian is symmetric How does this relate to the phonon frequencies Assume harmonic motion so that the displacements are 1 iwkt um T2 aAgiMAe mi A And the potential energy is quadratic in displacements 1 82U U U 0 2 2 6r r imjv I39M JV 0 And the equation of motion for an atom is m d2uiu 2 ZU u l 2 V dt jv 37 er 0 Eigenvalue equation Assuming a single normal mode excited we find from Eq of motio I a 6028 iMjv jvJL A in jv The forceconstant matrix is given by 2 1 a U my 12 mimj driner O Related to Hessian just including the masses of the atoms 9903 General procedure Choose a volume V Compute Hessian gt force constant matrix Diagonalize realsymmetric matrix Compute free energy FVT for various temperatures Return to 1 This should allow us to generate FVT curves Plot FVT vs volume V for various temperatures Find minimum by fitting to Eq of state Thermal expansion Anharmonicity through dependence of modes on V Usually called quasiharmonic approximation Declarations Equation of state code I Program to compute the free energy for a lennard jones solid IMPLICIT NONE INTEGER PARAMETER nn4 INTEGER PARAMETER PrecI4SELECTEDREALKIND14 REALKINDPrec14 PARAMETER LO20dO20dO30d0nntemp1000dO INTEGER PARAMETER nvol25 INTEGER PARAMETER natoms4nn3ndim3natoms INTEGER PARAMETER nxnnnynnnznn INTEGER ijixiyizmunumnmatzm1n1ierr REALKINDPrec14 PARAMETER sigma10dOepsilon10dO REAL KINDPre014 PARAMETER amass10d0 REAL KINDPre014 PARAMETER rcut30d0sigma REAL KINDPre014 epotdvdrdvdr2pifacfreedvovminvmaxL REAL KINDPre014 sxijsyijszijrijsovrvolumenv REAL KINDPre014 DIMENSION4 r1r2r3 REAL KINDPre014 DIMENSION3 5 REAL KINDPre014 DIMENSIONnatoms rxryrz REAL KINDPre014 DIMENSIONndimndim hessz REAL KINDPre014 DIMENSIONndim wfv1fv2 REAL KINDPre014 PARAMETER mass39948d0eps1200d0sar34d0 AAAAAAAAAA VVVVVVVVVV Some preliminaries Need to compute at several volumes one temperature Also need to zero Hessian for calculations pi40d0datan10d0 Compute volume volumeLO3 vminvolume085d0 vmaxvolume105d0 dvolvmaxvminnvol1 Then put the atoms on the fcc lattice as in MD code see several slides previous Do loop on the volumes several volumes do nv1nvo Clear the hessian out hess00d0 epot00d0 set potential energy to zero volumevminnv1dvo Lvoume10dO30d0 Loop on all atom pairs compute for rijltrcut do i1 natoms1 doji1natoms sxijrxi rxj SyiJryiryJ szijrzi rzj sxijsxijdbleanintsxij Added a 1rij here compared to forces sxijsxijL syijsyijL szijszijL rijdsqrtsxij2syij2szij2 s1sxijrij s2syijrij s3szi39 rij ifrijtrcut then sovrsigmalrij epotepot40dOepsilon ovr12 sovr6 dvdr240dOepsionrij220dOsovr12sovr6 dvdr2240dOepsionrij2260dOsovr1270dOsovr6 New lines first and second derivative of V with respect to r Identifying atom numbers with matrix elements 3N degrees of freedom xyz corresponding to p123 Atom 1 123 Atom 2 456 Etc Then for atom number i component mu m3i1mu For atom j component nu n3j1nu Then we can fill matrix element hessij Implement do mu13 do nu13 mi13mu nj13nu ifmueqnu hessmndvdr2smusnudvdrsmusnu10dO ifmunenu hessmndvdr2smusnudvdrsmusnu hessmnhessmnamass hessnmhessmn Matrix needs to be symmetric and we are only using eac Need to add then terms for HWV m1i13mu n1i13nu hessm1n1hessm1n1hessmn m1j13mu n1j13nu hessm1n1hessm1n1hessmn Diagonalize using the subroutine rs once hess is filled What do you get call rsndimndimhesswmatzzfv1fv2ierr Eigenvalues are the m2 for each of the 3N modes First three should be nearly zero might be slightly negative Next ones should all be greater than zero in general Units can be found in THz for example wndsqrtwn l Find frequencies from eigenvalues l Convert frequency to a unit in THz facdsqrt1381d06022d0eps100d0masssar2 wnwnfac20d0pi What do you get Eigenvalues are the m2 for each of the 3N modes First three should be nearly zero might be slightly negative Next ones should all be greater than zero in general Units can be found in THz for example wndsqrtwn Find frequencies from eigenvalues Convert frequency to a unit in THz facdsqrt1381d06022d0eps100dOmasssar2 wnwnfac20d0pi fac1381d41602d0 kB in eVK epotepotepsfac freeepot do n4ndim sum on eigenvalues freefreefactempdog20dOdsinhwn20dOfactemp enddo freefreenatoms write6100 volumesar3natomsfree BirchMurnaghan equation of state Determine equilibrium volume V0 Bulk modulus Murnaghan 07F 10 we av T av T Experimentally it is found that B varies with pressure B B0 B p Combining these and integrating we get W l B V0 p B BirchMurnaghan equation of state A slightly more sophisticated expression is BirchMurnaghan E73 53 V V After integrating we get 23 23 1 1 V V 16 Equilibrium p0 volume is VO B0 is bulk modulus at p0 F0 V0 80 are fitting parameters Monte Carlo code to do fitting exists 2 1Bg 4 will 23 6 45 V 2 3 B0 l Free energy vs volume 70075 a T5K x M T1OK waw ii 0 055 T K 0 095 T1OOK 39 oo 0 Thermal expansion V0 increases with temperature 0 Curvature related to bulk modulus B 0 Vibrational energyentropy After Hessian and intial direction 9 found get initial direction h hg output current potential energy write6 39Hessian computed potential energy39epot compute lambda lambda00dO do n1ndim lambdalambdagngn enddo vec00dO do m1ndim do n1ndim vecmvecmhessmnhn enddo enddo denom00dO do m1ndim denomdenomhmvecm enddo Quantum M0nteCarlo Helium atom H V12V 3 3i rlrz r12 We assume that we can first of all take the trial wave function ITr19r2 ACXp O 7i 73 If we first ignore electronelectron repulsion we get or 2 do Where a0 is the Bohr radius Energy with and without electronelectron repulsion Computing energy without the 1r12 term we get 1088eV Experimental result 7898eV Obviously ee repulsion important Including 1r12 term but not optimizing or we get 748eV Optimizing or169a0 the atom swells a bit With optimized or we get 775eV within about 15eV Can go lower yet Energy with added correlations ITr1r2 AeXp ocr1 r2eXp3r12 or 3 both parameters to be optimized For my code I found an energy within 05eV of experiment Adding correlations accounts for 10eV of difference when i0 Quantum MonteCarlo approach Simply but approximate way to evaluate integral f Jr mm 5 Mm 5 whirlla E E Jr I 31 39iquot1 Tzl39jl pquot39quotquotl1395quotK id iquot EFT This can also be written as Jquot Jquot it H mm jutng arm digit Jr Jr ElsiJIM Elitn a Local energy de nition 1 After some difficult algebra we find F I7 17 05 2 05 2 2 1 EL a232oc3 a3 2 12 l3 r1732 r2732 r1 r2 r12 gt Evaluate for each configuration with electrons 12 at 1 1 7392 9903 CD Outline of code Start with electrons near nucleus Make a random move 01au for electrons 12 Compute ICIT before and after random move If probability density increased after move accept move If probability density decreased after move accept with probability qqTa er before For current configuration compute EL sum for averaging subscripts for beforeafter move After at least 107 moves output the average value of EL Declarations for Helium atom code INTEGER PARAMETER neI2 nstep100000000 INTEGER PARAMETER Prec14SELECTEDREALKIND14 REALKNDPrec14 DIMENSIONneI x y z xn yn zn I coordinates of1 REALKNDPrec14 r1r2 r12 rn1rn2rn12 Fi Fin REALKNDPrec14 P Pn DiffP INTEGER ach REALKNDPrec14 Elocal I Local enegy REALKNDPrec14 PARAMETER alpha180d0 beta025d0 REALKNDPrec14 EIoca SumE ranx rany ranz r REALKNDPrec14 PARAMETER v001d0 v0 is the magnitude of the random moves in au Initial set up x105d0 initial coordinates of the 1st particle y100d0 z100d0 x205d0 initial coordinates ofthe 2nd particle y200d0 z200d0 SumE00dO ach0 call RANDOMSEED Coordinates of atoms in au Bohrs energies will be in Hartrees Initial set up Evaluation of relative probability of con gurations After taking a random step r1 dsqrtx12y12z12 r2dsqrtx22y22z22 r12dSqrtX1X22Y1Y22Z1Z22 rn1dsqrtxn12yn12zn12 rn2dsqrtxn22yn22zn22 rn12dsqrtxn1xn22yn1yn22zn1zn22 Fidexpalphar1dexpalphar2dexpbetar12 amplitude of Phi before mc Findexpalpharn1dexpalpharn2dexpbetarn12 amplitude of Phi after PdabsFiFi probability density before move PndabsFinFin probabality density after move ls Pn greater or smaller than P That determines whether we accept random step In class 4 PHZ 5156 Due Thursday Nov 20 1 The lsing model in the mean eld approximation can also be studied using a Landau free energy density 1 at772 Elm4 7 H77 where H is the applied external eld7 and we can make the connection It Tin In this way7 we will see that t lt 0 corresponds to T lt T07 and likewise t gt 0 corresponds to T gt T0 The 77 is called the order parameter and7 in this case7 can be thought of as the average magnetization per site The a and b are simply phenomenological parameters a Make a plot on the computer of L in the case where H 0 for t lt 07 t 07 and t gt 0 Also7 explain what happens when H lt 0 or H gt 0 Take for simplicity a b 1 in your plots Explain how your pictures agree with the general picture of the phase transition b Assume that the value of 77 is found by taking the partial derivative of and setting the result equal to zero Consider the case where H 0 Show that 77 0 for t gt 07 and 77 iqiaf for t S 0 Hence we see the critical exponent B 12 just as in the Weiss mean eld theory as discussed in lecture c Show that the speci c heat per spin CV 0 for t gt 0 and CV for t lt 0 Again this is the case for H 0 This shows there is a discontinuity7 but no divergence7 in the heat capacity at T Tc d Show that the magnetic equation of state7 found by differentiating with respect to 777 is given by at77 b773 This gives the critical exponent 6 3 we found along the isotherm t 0 from the mean eld theory Hint This is the same thing you did in b above except there we took H 0 e Show that the susceptibility7 de ned as X 3 is given for T gt To by X 2at 1 and for T lt Tc X 740071 This shows the divergence at T T0 with critical exponent y 1 2 This is a short code that explores the Metropolis Monte Carlo approach using the Landau free energy density described above For a spin systern7 take for N spins the free energy L N For sirnplicity7 then we can take L Am2 B774 for H 07 where A 1N and B bN The probability P of a state with a given value of 77 is then P N exp i Ln Write a Metropolis Monte Carlo code that cornputes lt77gt as a function of temperature T You can take A B 1 Start with the system at T gt To and gradually cool to T 0 How do your results compare with the results of problem 1 Explain why you might think lt77gt 0 for all temperatures Think in terms of the free energy diagrarns you made above in 1a Why in practice does this not happen The fact that a nite magnetization is found for T lt To is an example of ergodicz39ty breaking 3 The Hamiltonian HQ in the lsing model used to describer a two dimensional spin system with no applied magnetic eld is given by N N 1 HQ 3 Sij Si1j 5121739 513741 513721 1 As we showed in class the transition temperature of the mean eld theory is given by kBTc 4 for a two dimensional square lattice The transition state is where the system alternates between pararnagnetic high temperatures and ferrornagnetic low tejrnperatures The exact result found by Onsager is kBTc W 22692 Write a code that uses a Metropolis Monte Carlo algorithm to approximately corn pute the average properties of the system for a two dimensional lsing rnodel Write your code in terms of the dimensionless pararneter JkBT Use periodic boundary conditions so that spins on the edge of your square interact with spins on the op posite edge Write a function clusterij which computes the energy in a cluster of spins centered on spin 2397 Write another function which uses clusterij to corn pute the total energy of your square lattice of spins Now here is the Monte Carlo algorithm which sarnples randorn states with probability in accord with the partition function 0 Begin with a random array of spins 1 Compute the energy E1 of the current spin state 2 Choose a spin on the lattice at random and ip its spin state 3 Recompute the energy E2 of the new spin state 4 Determine the energy difference A E1 7 E2 5 If A S 07 accept the current step in the ensemble Compute and average whatever are the quantities of interest7 including the energy and average spin Save the current spins and energy and return to 2 6 If A gt 07 draw a random number between 0 and 1 7 If the random number is less than exp ikBAT7 accept the current step into the ensemble Compute and average whatever are the quantities of interest7 including the energy and average spin Save the current spins and energy Return to step 2 8 If the random number is greater than or equal to exp ikBAT7 then reject the current step Reset the spins to what they were before the spin was ipped in step 2 Compute and average whatever are the quantities of interest7 including the energy and average spin Return to step 2 Evolve a lattice of 32 gtlt 32 spins Make a plot of energy vs iteration You will nd that the system will equilibrate after a certain amount of time After equilibrating7 average properties can be obtained Accumulate the average magnetization Obtain expressions for the spin susceptibility and heat capacity Compute these quantities during your Monte Carlo simulation Establish approximately the transition point kBTcJ For your presentation7 compare to the mean eld and exact results from Onsager7s solution We will discuss and compare critical exponents in each case You should run your code at kBTJ 17 157 1757 207 2257 257 and 30 Finite size effects are important near the transition point Consider also an explanation for why it is much harder to equilibrate near kBTJ 225 For three dimensional systems7 no one has obtained an exact result for the lsing model7 and usually either computational or theoretical techniques using the renor malization group are used A good reference for this is 77 Lectures on Phase Transitions and the Renormalization Group77 by Nigel Goldenfeld Anderson model code H A m 63 Hj jm I 136 1 i lt t W M Wz 136 0 6 gt VVZ Not tridiagonal We will use periodic boundary conditions Anderson model code declarations implicit none INTEGER PARAMETER nn10nn3nn3 I size of squarenumber of site INTEGER PARAMETER nbin1000 I for dospart ratio REAL8 PARAMETER ww50d0 I width of disorder from w2 to w2 REAL8 dimensionnn3 wfv1fv2 REAL8 dimensionnn3nn3 hamiz REAL8 dimensionnbin dospartcount REAL8 rndenran2 REAL8 PARAMETER wminww20d0wmaxww20dO REAL8 PARAMETER delw40d0wwnbin INTEGER matzierrmnnbiiijjjkkkneigh INTEGER PARAMETER mensb50 Initialize diagonal elements from random distribution hamil00d0 Clear everything to zero do n1nn3 set diagonal energies call RANDOMNUMBERrnd haminnrndO5d0ww enddo Offdiagonal elements for nearest neighbors Reference each site coordinates ijk Number of site determined from the ijk trio doi1nn doj1nn do k1nn nnnnG1nnG1k enddo enddo enddo Nearest neighbors The following give nearest neighbors for site ijk Make a loop on the integer variable neigh 16 iii jjj kkk ifneigheq1 iii1 ifneigheq2 iii1 ifneigheq3 jjj1 ifneigheq4 jjj1 ifneigheq5 kkk1 ifneigheq6 kkk1 Need the site number of the neighbor which we get from llnnnnii1 nnjj1kk hamillll10d0 hamillll10d0 Periodic boundary conditions The ii orjj might be less than 1 or bigger than nn ifiit1 iinn ifiigtnn ii1 ifjjgtnn jj1 ifjjt1 jjnn ifkkgtnn kk1 ifkkt1 kknn This get the pbc correct Notice yields a matrix that is not tridiagonal Need pbc in each of the three dimensions Subroutine rs f subroutine rsnmnawmatzzfv1fv2ierr C OOOOOOOOOOOOOOOO integer nnmierrmatz double precision anmnwnznmnfv1nfv2n this subroutine calls the recommended sequence of subroutines from the eigensystem subroutine package eispack to find the eigenvalues and eigenvectors if desired of a real symmetric matrix oninput nm must be set to the row dimension of the twodimensional array parameters as declared in the calling program dimension statement n is the order of the matrix a a contains the real symmetric matrix Subroutine rs f The znm contains the eigenvectors The index m the second index is mode The index n the first index is the amplitude on site n Used to compute participation ratio The wn give the energies of the eigenstates in ascending order Participation ratios measure of localization pgl N ham i1N Delocalized state p1 Localized state p1N goes to zero as Ngt infinity Determine from the eigenvector matrix zmn Density of states number of modes in a range of frequency ga 2500 ooh Measures number of states at a given frequency Must broaden to do numerically As with participation ratio should average over realizations Computing and storing DOS and participation ratios do n1nn3 nbnintwnwmindew ifnbgt0andnbenbin then countnbcountnb10dO dosnbdosnb10d0 do i1nn3 partnbpartnbzin4 enddo endif enddo Will need average participation ratio over all modes found in a given range of frequency Use the count array DOS for W5 8x8x8 lattice 50 realizations of randomness 150 m 120 1 a o o g v a 100 quotw 3 3o 3 o o g m 3 50 v O o 3 x 00 a 3 20 o 39 a Width of the DOS re ects ariiount of disorder hopping Spread below 5 above 5 due to hopping Strong band formation DOS for W20 8x8x8 lattice 50 realizations of randomness 5 10 15 rm 75 Larger spread in DOS due to larger disorder Hopping less important between sites Participation for W5 8x8x8 lattice 50 realizations of randomness 5 re 72 2 o 5 a All states expected to be delocalized forming in nite clusters Participation ratio should increasestay constant as size increases Participation for W20 8x8x8 lattice 50 realizations afrandomne 0015 2 0015 3 0010 0012 Expect all localized states Participation ratio should decrease as system size gets bigger