seminar:cu_bulk

# Band計算

## スクラッチから計算する

""" Make a bandstructure plot using the Harris functional.
"""
import sys
from Dacapo import Dacapo
from ASE import Atom,ListOfAtoms
from ASE.Utilities.List import List
import Numeric as num
import os

bulk = ListOfAtoms([Atom('Cu', (0,    0,     0))] )
a0_fcc = 3.61
b = a0_fcc/2.
bulk.SetUnitCell([(0, b, b), (b, 0, b), (b, b, 0)])

calc = Dacapo()
calc.SetBZKPoints((10,10,10))        # set the k-points Monkhorst-Pack
calc.SetPlaneWaveCutoff(340)         # planewavecutoff in eV
calc.SetDensityCutoff(400)
calc.SetNumberOfBands(8)             # set the number of electronic bands
calc.SetSymmetryOn()                 # use symmetry to reduce the k-point set
calc.SetNetCDFFile('make_density_converge.nc')
calc.SetTxtFile('make_density_converge.text')

bulk.SetCalculator(calc)

calc.StayAliveOff()

energy = calc.GetPotentialEnergy()
fermilevel = calc.GetFermiLevel()

# setup the Harris calculation

# setup the line from the Gamme point (0,0,0) to the X point (0.5,0.5,0)
N = 80
Nf = float(N)
kpts = [(x/Nf,x/Nf,0) for x in range(N/2)]
harris = Dacapo()
harris.SetBZKPoints(kpts)          # set the k-points along one direction
harris.SetPlaneWaveCutoff(340)     # planewavecutoff in eV
harris.SetNumberOfBands(6)         # set the number of electronic bands
harris.SetSymmetryOn()             # use symmetry to reduce the k-point set
harris.SetNetCDFFile('calc_eigenvalues.nc')
harris.SetTxtFile('calc_eigenvalues.text')

bulk.SetCalculator(harris)

harris.StayAliveOff()
harris.SetChargeMixing(False)

# set the density found using the uniform k-point sampling
harris.SetDensityArray(calc.GetDensityArray())

energy1 = harris.GetPotentialEnergy()
kpoints     = harris.GetIBZKPoints()
eigenvalues = [harris.GetEigenvalues(kpt=kpt) for kpt in range(len(kpoints))]

# using the fermilevel from uniform k-point sampling
print fermilevel
print kpoints
print eigenvalues

## 一回収斂させた電荷密度を使って計算する

""" Make a bandstructure plot using the Harris functional.
"""
import sys
from Dacapo import Dacapo
from ASE import Atom,ListOfAtoms
from ASE.Utilities.List import List
import Numeric as num
import os

calc = bulk.GetCalculator()

energy = calc.GetPotentialEnergy()
fermilevel = calc.GetFermiLevel()

# setup the Harris calculation

# setup the line from the Gamme point (0,0,0) to the X point (0.5,0.5,0)
N = 80
Nf = float(N)
kpts = [(x/Nf,x/Nf,0) for x in range(N/2)]
harris = Dacapo()
harris.SetBZKPoints(kpts)          # set the k-points along one direction
harris.SetPlaneWaveCutoff(340)     # planewavecutoff in eV
harris.SetNumberOfBands(6)         # set the number of electronic bands
harris.SetSymmetryOn()             # use symmetry to reduce the k-point set
harris.SetNetCDFFile('calc_eigenvalues_only.nc')
harris.SetTxtFile('calc_eigenvalues_only.text')

bulk.SetCalculator(harris)

harris.StayAliveOff()
harris.SetChargeMixing(False)

# set the density found using the uniform k-point sampling
harris.SetDensityArray(calc.GetDensityArray())

energy1 = harris.GetPotentialEnergy()
kpoints     = harris.GetIBZKPoints()
eigenvalues = [harris.GetEigenvalues(kpt=kpt) for kpt in range(len(kpoints))]

# using the fermilevel from uniform k-point sampling
print fermilevel
print kpoints
print eigenvalues
# -*- coding: iso-8859-15 -*-
# 日本語でコメントを書くために上の行が必要
from math import sqrt  # sqrt(x)で√xが計算できる
##################################################
ncfile = "a2.nc"  # このファイルは上書きされてしまうので
# 必ずコピーを取っておくこと

M = 10
Mf = float(M)
wk1 = [(2*n/(3*Mf),n/(3*Mf),0) for n in range(M)]

from Dacapo import Dacapo