
# Create a 3x3x3 supercell of alpha-silicon
lattice = SimpleCubic(5.4306*Angstrom)
elements = [Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon]
fractional_coordinates = [[ 0.  ,  0.  ,  0.  ],
                          [ 0.25,  0.25,  0.25],
                          [ 0.5 ,  0.5 ,  0.  ],
                          [ 0.75,  0.75,  0.25],
                          [ 0.5 ,  0.  ,  0.5 ],
                          [ 0.75,  0.25,  0.75],
                          [ 0.  ,  0.5 ,  0.5 ],
                          [ 0.25,  0.75,  0.75]]

alpha_silicon = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates
)
alpha_silicon = alpha_silicon.repeat(3, 3, 3)

# Use a Stillinger-Weber potential for silicon to calculate the changes in energy.
potentialSet = StillingerWeber_Si_1985()
calculator = TremoloXCalculator(parameters=potentialSet)
alpha_silicon.setCalculator(calculator)

# Define the Monte Carlo method. Note here the temperature changes from the initial 50500 to
# 500K
monte_carlo_method = ContinuousRandomNetwork(
    temperature=50500.0*Kelvin,
    heating_rate=-500.0*Kelvin,
    chain_length=1
)

# Perform the Monte Carlo simulation
monte_carlo_trajectory = MetropolisMonteCarlo(
    configuration=alpha_silicon,
    trajectory_filename='Amorphous_Silicon.hdf5',
    steps=100,
    max_trial_moves=5000,
    log_interval=1,
    method=monte_carlo_method,
    trajectory_object_id=None
)
last_image = monte_carlo_trajectory.lastImage()
nlsave('Amorphous_Silicon.hdf5', last_image)
