# Nudged Elastic Band configuration

# Setup the configuration list
configuration_list = []

# Define elements
elements = [Carbon,]*2 + [Hydrogen,]*6

# Define initial configuration
cartesian_coordinates = [[  0.0000000,  -1.2257e-05,   0.7716335],
                         [  0.0000000,   1.2257e-05,  -0.7716335],
                         [  0.0000000,   1.02162746,   1.1513662],
                         [ -0.8847539,  -0.51081476,   1.1513615],
                         [  0.8847539,  -0.51081476,   1.1513615],
                         [  0.0000000,  -1.02162746,  -1.1513662],
                         [ -0.8847539,   0.51081476,  -1.1513615],
                         [  0.8847539,   0.51081476,  -1.1513615]]*Angstrom

initial_configuration = MoleculeConfiguration(
    elements=elements,
    cartesian_coordinates=cartesian_coordinates
    )

# Append to the list of configurations
configuration_list.append(initial_configuration)

# Define final configuration
cartesian_coordinates = [[ 0.07475139, -0.04996446,  0.65332658],
                         [-0.07215805,  0.03485788, -0.6533275 ],
                         [ 0.66099367,  3.98448079, -0.09787249],
                         [-0.62932721,  0.35668865,  1.38033322],
                         [ 0.91714293, -0.53575571,  1.1465225 ],
                         [ 0.63207421, -0.37155137, -1.38032121],
                         [-0.91406328,  0.52147801, -1.1465376 ],
                         [ 1.05253733,  3.75843821, -0.68559951]]*Angstrom

final_configuration = MoleculeConfiguration(
    elements=elements,
    cartesian_coordinates=cartesian_coordinates
    )

# Append to the list of configurations
configuration_list.append(final_configuration)

# Construct the Nudged Elastic Band object
neb_configuration = NudgedElasticBand(
    configuration_list,
    image_distance=0.2*Ang)

# Define a calculator and set it on the NEB
neb_configuration.setCalculator(LCAOCalculator())

# Find the reaction pathway
optimized_neb = OptimizeNudgedElasticBand(
        neb_configuration,
        max_forces=0.05*eV/Ang)