# -*- coding: utf-8 -*-
# -------------------------------------------------------------
# Bulk Configuration
# -------------------------------------------------------------

# Set up lattice
lattice = FaceCenteredCubic(5.6537*Angstrom)

# Define elements
elements = [Gallium, Arsenic]

# Define coordinates
fractional_coordinates = [[ 0.  ,  0.  ,  0.  ],
                          [ 0.25,  0.25,  0.25]]

# Set up configuration
bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates
    )

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
#----------------------------------------
# Basis Set
#----------------------------------------
basis_set = [
    GGABasis.Gallium_DoubleZetaPolarized,
    GGABasis.Arsenic_DoubleZetaPolarized,
    ]

#----------------------------------------
# Exchange-Correlation
#----------------------------------------
exchange_correlation = GGAHalf.PBE

k_point_sampling = MonkhorstPackGrid(
    na=13,
    nb=13,
    nc=13,
    )
numerical_accuracy_parameters = NumericalAccuracyParameters(
    density_mesh_cutoff=45.0*Hartree,
    k_point_sampling=k_point_sampling,
    occupation_method=FermiDirac(100.0*Kelvin*boltzmann_constant),
    )

calculator = LCAOCalculator(
    basis_set=basis_set,
    exchange_correlation=exchange_correlation,
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    )

bulk_configuration.setCalculator(calculator)
nlprint(bulk_configuration)
bulk_configuration.update()
nlsave('GaAs.hdf5', bulk_configuration)

# -------------------------------------------------------------
# Second Harmonics Generation Susceptibility
# -------------------------------------------------------------
kpoint_grid = MonkhorstPackGrid(
    na=14,
    nb=14,
    nc=14,
    )

second_harmonics_generation_susceptibility = SecondHarmonicsGenerationSusceptibility(
    configuration=bulk_configuration,
    kpoints=kpoint_grid,
    energies=numpy.linspace(0, 5, 101)*eV,
    broadening=0.1*eV,
    bands_below_fermi_level=4,
    bands_above_fermi_level=8,
    use_symmetries=True,
    tensor_indices=None,
    )
nlsave('GaAs.hdf5', second_harmonics_generation_susceptibility)