import pylab
# Try k-point sampling grids in the range 1x1x1 to 19x19x19
k_values = range(2,20) 
# Create a list to store the computed total energies
energies = []
# The configuration to compute is bulk gold
lattice = FaceCenteredCubic(4.07825*Angstrom) 
elements = [Gold] 
cartesian_coordinates = [[ 0.,  0.,  0.]]*Angstrom
bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice, 
    elements=elements, 
    cartesian_coordinates=cartesian_coordinates
)
# For the first step there is no initial state
initial_state = None
# Loop over the different k-point sampling grids
for Nk in k_values:     
    numerical_accuracy_parameters = NumericalAccuracyParameters(
        k_point_sampling=(Nk, Nk, Nk)
    ) 
    iteration_control_parameters = IterationControlParameters(
        tolerance=4e-07,
    )
    calculator = LCAOCalculator(
    numerical_accuracy_parameters=numerical_accuracy_parameters,
        iteration_control_parameters=iteration_control_parameters,
    )
    # This is where we use the previous calculation as starting guess!
    # Note that None is a valid value for the initial_state
    bulk_configuration.setCalculator(         
        calculator,         
        initial_state=initial_state
    )
    # Compute and store the total energy
    energies.append(TotalEnergy(bulk_configuration).evaluate())
    # Save the converged result to use as initial state for the next value of Nk
    initial_state = bulk_configuration
# Print results (on the master node only)
if processIsMaster():
    for Nk,energy in zip(k_values,energies):
        print(Nk, energy)

pylab.plot(k_values,energies,'bo-')
pylab.xlabel('K points',fontsize=16)
pylab.ylabel('Energy (eV)',fontsize=16 )
pylab.title('K-point Scanning')

# Show the plot.
pylab.show()