# -*- coding: utf-8 -*-
# -------------------------------------------------------------
# Two-probe Configuration
# -------------------------------------------------------------

# -------------------------------------------------------------
# Left Electrode
# -------------------------------------------------------------

# Set up lattice
vector_a = [9.0, 0.0, 0.0]*Angstrom
vector_b = [0.0, 9.0, 0.0]*Angstrom
vector_c = [0.0, 0.0, 10.16]*Angstrom
left_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
left_electrode_elements = [Gold, Gold, Gold, Gold]

# Define coordinates
left_electrode_coordinates = [[ 5.  ,  5.  ,  1.27],
                              [ 5.  ,  5.  ,  3.81],
                              [ 5.  ,  5.  ,  6.35],
                              [ 5.  ,  5.  ,  8.89]]*Angstrom

# Set up configuration
left_electrode = BulkConfiguration(
    bravais_lattice=left_electrode_lattice,
    elements=left_electrode_elements,
    cartesian_coordinates=left_electrode_coordinates
    )

# -------------------------------------------------------------
# Right Electrode
# -------------------------------------------------------------

# Set up lattice
vector_a = [9.0, 0.0, 0.0]*Angstrom
vector_b = [0.0, 9.0, 0.0]*Angstrom
vector_c = [0.0, 0.0, 10.16]*Angstrom
right_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
right_electrode_elements = [Gold, Gold, Gold, Gold]

# Define coordinates
right_electrode_coordinates = [[ 5.  ,  5.  ,  1.27],
                               [ 5.  ,  5.  ,  3.81],
                               [ 5.  ,  5.  ,  6.35],
                               [ 5.  ,  5.  ,  8.89]]*Angstrom

# Set up configuration
right_electrode = BulkConfiguration(
    bravais_lattice=right_electrode_lattice,
    elements=right_electrode_elements,
    cartesian_coordinates=right_electrode_coordinates
    )

# -------------------------------------------------------------
# Central Region
# -------------------------------------------------------------

# Set up lattice
vector_a = [9.0, 0.0, 0.0]*Angstrom
vector_b = [0.0, 9.0, 0.0]*Angstrom
vector_c = [0.0, 0.0, 32.86]*Angstrom
central_region_lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
central_region_elements = [Gold, Gold, Gold, Gold, Gold, Gold, Hydrogen, Hydrogen, Gold, Gold,
                           Gold, Gold, Gold, Gold]

# Define coordinates
central_region_coordinates = [[  5.   ,   5.   ,   1.27 ],
                              [  5.   ,   5.   ,   3.81 ],
                              [  5.   ,   5.   ,   6.35 ],
                              [  5.   ,   5.   ,   8.89 ],
                              [  5.   ,   5.   ,  11.43 ],
                              [  5.005,   5.005,  14.102],
                              [  5.   ,   5.   ,  15.972],
                              [  5.   ,   5.   ,  16.888],
                              [  5.005,   5.005,  18.758],
                              [  5.   ,   5.   ,  21.43 ],
                              [  5.   ,   5.   ,  23.97 ],
                              [  5.   ,   5.   ,  26.51 ],
                              [  5.   ,   5.   ,  29.05 ],
                              [  5.   ,   5.   ,  31.59 ]]*Angstrom

# Set up configuration
central_region = BulkConfiguration(
    bravais_lattice=central_region_lattice,
    elements=central_region_elements,
    cartesian_coordinates=central_region_coordinates
    )

device_configuration = DeviceConfiguration(
    central_region,
    [left_electrode, right_electrode]
    )

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
#----------------------------------------
# Basis Set
#----------------------------------------
basis_set = [
    GGABasis.Hydrogen_DoubleZetaPolarized,
    GGABasis.Gold_DoubleZetaPolarized,
    ]

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

#----------------------------------------
# Numerical Accuracy Settings
#----------------------------------------
left_electrode_k_point_sampling = MonkhorstPackGrid(
    nc=93,
    )
left_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
    k_point_sampling=left_electrode_k_point_sampling,
    )

right_electrode_k_point_sampling = MonkhorstPackGrid(
    nc=93,
    )
right_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
    k_point_sampling=right_electrode_k_point_sampling,
    )

device_k_point_sampling = MonkhorstPackGrid(
    nc=93,
    )
device_numerical_accuracy_parameters = NumericalAccuracyParameters(
    k_point_sampling=device_k_point_sampling,
    )

#----------------------------------------
# Poisson Solver Settings
#----------------------------------------
left_electrode_poisson_solver = FastFourier2DSolver(
    boundary_conditions=[[PeriodicBoundaryCondition(),PeriodicBoundaryCondition()],
                         [PeriodicBoundaryCondition(),PeriodicBoundaryCondition()],
                         [PeriodicBoundaryCondition(),PeriodicBoundaryCondition()]]
    )

right_electrode_poisson_solver = FastFourier2DSolver(
    boundary_conditions=[[PeriodicBoundaryCondition(),PeriodicBoundaryCondition()],
                         [PeriodicBoundaryCondition(),PeriodicBoundaryCondition()],
                         [PeriodicBoundaryCondition(),PeriodicBoundaryCondition()]]
    )

#----------------------------------------
# Electrode Calculators
#----------------------------------------
left_electrode_calculator = LCAOCalculator(
    basis_set=basis_set,
    exchange_correlation=exchange_correlation,
    numerical_accuracy_parameters=left_electrode_numerical_accuracy_parameters,
    poisson_solver=left_electrode_poisson_solver,
    )

right_electrode_calculator = LCAOCalculator(
    basis_set=basis_set,
    exchange_correlation=exchange_correlation,
    numerical_accuracy_parameters=right_electrode_numerical_accuracy_parameters,
    poisson_solver=right_electrode_poisson_solver,
    )

#----------------------------------------
# Device Calculator
#----------------------------------------
calculator = DeviceLCAOCalculator(
    basis_set=basis_set,
    exchange_correlation=exchange_correlation,
    numerical_accuracy_parameters=device_numerical_accuracy_parameters,
    electrode_calculators=
        [left_electrode_calculator, right_electrode_calculator],
    )

device_configuration.setCalculator(calculator)
nlprint(device_configuration)
device_configuration.update()
nlsave('au-h2-au.nc', device_configuration)