1 #!/usr/bin/env python
   2 
   3 from Dacapo import Dacapo
   4 from ASE import Atom, ListOfAtoms
   5 
   6 from Numeric import *
   7 import sys
   8 
   9 if len(sys.argv)<2:
  10   print """
  11   Makes an estimation of the magnetic moment of the individual atoms
  12   based on the atom projected local density of states. 
  13   
  14   USAGE: python LMagMom.py [nc.file]
  15   
  16   NOTE: Works only if 
  17    !calc.CalculateAtomicDOS()!
  18   is specified for the dacapo calculation
  19   """
  20   sys.exit()  
  21 atoms=Dacapo.ReadAtoms(sys.argv[1])
  22 calc=atoms.GetCalculator()
  23 
  24 
  25 spin0=calc.GetDensityArray(spin=0)
  26 spin1=calc.GetDensityArray(spin=1)
  27 ngx,ngy,ngz=shape(spin0)
  28 gridpoints=ngx*ngy*ngz
  29 VolumeElement=atoms.GetUnitCellVolume()/float(gridpoints)
  30 Nup=sum(sum(sum(spin0)))*VolumeElement
  31 Ndown=sum(sum(sum(spin1)))*VolumeElement
  32 print "Nup", Nup,"Ndown", Ndown,"Diff", Nup-Ndown
  33 
  34 for i in range(len(atoms)):
  35   dosu=calc.GetLDOS(atoms=[i+1],angularchannels=["s","p","d"],spin=[0])
  36   dosus=calc.GetLDOS(atoms=[i+1],angularchannels=["s","p","d"],spin=[0],cutoffradius = 'short')
  37   dosd=calc.GetLDOS(atoms=[i+1],angularchannels=["s","p","d"],spin=[1])
  38   dosds=calc.GetLDOS(atoms=[i+1],angularchannels=["s","p","d"],spin=[1],cutoffradius = 'short')
  39   print "Atom_"+str(i)+" "+atoms[i].GetChemicalSymbol()
  40   print "N_inf=",dosd.GetIntegratedDOS()+dosu.GetIntegratedDOS(),"N_short=", dosus.GetIntegratedDOS()+dosds.GetIntegratedDOS()  
  41   print "M_inf=",dosu.GetIntegratedDOS()-dosd.GetIntegratedDOS(), "M_short=",dosus.GetIntegratedDOS()-dosds.GetIntegratedDOS()