# -------------------------------------------------------------
# Bulk configuration
# -------------------------------------------------------------
# Set up lattice
lattice = Rhombohedral(5.06138065551*Angstrom, 33.5573097619*Degrees)
# Define elements
elements = [Cobalt, Cobalt, Oxygen, Oxygen]
# Define coordinates
fractional_coordinates = [[ 0.002, 0.002, -0.006],
[ 0.502, 0.502, 0.494],
[ 0.25 , 0.25 , 0.25 ],
[ 0.75 , 0.75 , 0.75 ]]
# Set up configuration
bulk_configuration = BulkConfiguration(
bravais_lattice=lattice,
elements=elements,
fractional_coordinates=fractional_coordinates
)
# Add tags
bulk_configuration.addTags('cobalt1', [1])
# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
#----------------------------------------
# Basis Set
#----------------------------------------
cobalt_3d = ConfinedOrbital(
principal_quantum_number=3,
angular_momentum=2,
radial_cutoff_radius=4.07498342334*Bohr,
confinement_start_radius=3.25998673867*Bohr,
additional_charge=0,
confinement_strength=24.5399771266*Hartree,
confinement_power=1,
radial_step_size=0.001*Bohr,
)
cobalt_3d_split = AnalyticalSplit(cobalt_3d, split_norm=0.15)
cobalt_4s = ConfinedOrbital(
principal_quantum_number=4,
angular_momentum=0,
radial_cutoff_radius=6.93481368132*Bohr,
confinement_start_radius=5.54785094506*Bohr,
additional_charge=0,
confinement_strength=14.4199980844*Hartree,
confinement_power=1,
radial_step_size=0.001*Bohr,
)
cobalt_4s_polarization = PolarizationOrbital(cobalt_4s)
cobalt_4s_split = AnalyticalSplit(cobalt_4s, split_norm=0.15)
CobaltBasis = BasisSet(
element=PeriodicTable.Cobalt,
orbitals=[ cobalt_3d , cobalt_4s , cobalt_3d_split , cobalt_4s_split , cobalt_4s_polarization ],
occupations=[ 7.0 , 2.0 , 0.0 , 0.0 , 0.0],
hubbard_u=[ 4.1 , 0.0 , 4.1 , 0.0 , 0.0]*eV,
filling_method=SphericalSymmetric,
pseudopotential=NormConservingPseudoPotential("normconserving/CO.GGAPBE.zip"),
)
cobalt_tag1_3d = ConfinedOrbital(
principal_quantum_number=3,
angular_momentum=2,
radial_cutoff_radius=4.07498342334*Bohr,
confinement_start_radius=3.25998673867*Bohr,
additional_charge=0,
confinement_strength=24.5399771266*Hartree,
confinement_power=1,
radial_step_size=0.001*Bohr,
)
cobalt_tag1_3d_split = AnalyticalSplit(cobalt_tag1_3d, split_norm=0.15)
cobalt_tag1_4s = ConfinedOrbital(
principal_quantum_number=4,
angular_momentum=0,
radial_cutoff_radius=6.93481368132*Bohr,
confinement_start_radius=5.54785094506*Bohr,
additional_charge=0,
confinement_strength=14.4199980844*Hartree,
confinement_power=1,
radial_step_size=0.001*Bohr,
)
cobalt_tag1_4s_polarization = PolarizationOrbital(cobalt_tag1_4s)
cobalt_tag1_4s_split = AnalyticalSplit(cobalt_tag1_4s, split_norm=0.15)
Tag1_CobaltBasis = BasisSet(
element=PeriodicTable.Cobalt,
orbitals=[ cobalt_tag1_3d , cobalt_tag1_4s , cobalt_tag1_3d_split , cobalt_tag1_4s_split , cobalt_tag1_4s_polarization ],
occupations=[ 7.0 , 2.0 , 0.0 , 0.0 , 0.0],
hubbard_u=[ 3.1 , 0.0 , 3.1 , 0.0 , 0.0]*eV,
filling_method=SphericalSymmetric,
pseudopotential=NormConservingPseudoPotential("normconserving/CO.GGAPBE.zip"),
)
basis_set = [
GGABasis.Oxygen_DoubleZetaPolarized,
CobaltBasis,
('cobalt1', Tag1_CobaltBasis ),
]
#----------------------------------------
# Exchange-Correlation
#----------------------------------------
exchange_correlation = SGGAU.PBE
numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(3, 3, 3),
density_mesh_cutoff=80.0*Hartree,
)
iteration_control_parameters = IterationControlParameters(
damping_factor=0.3,
max_steps=200,
tolerance=0.0002,
number_of_history_steps=22,
)
calculator = LCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=numerical_accuracy_parameters,
iteration_control_parameters=iteration_control_parameters,
)
bulk_configuration.setCalculator(calculator)
nlprint(bulk_configuration)
bulk_configuration.update()
nlsave('differentu.nc', bulk_configuration)
# -------------------------------------------------------------
# Mulliken population
# -------------------------------------------------------------
mulliken_population = MullikenPopulation(bulk_configuration)
nlsave('differentu.nc', mulliken_population)
nlprint(mulliken_population)
# -------------------------------------------------------------
# Density of states
# -------------------------------------------------------------
density_of_states = DensityOfStates(
configuration=bulk_configuration,
kpoints=MonkhorstPackGrid(5,5,5),
energy_zero_parameter=FermiLevel,
bands_above_fermi_level=None,
)
nlsave('differentu.nc', density_of_states)
nlprint(density_of_states)