# -------------------------------------------------------------
# Set up a SiO2-quartz crystal.
# -------------------------------------------------------------
lattice = Hexagonal(4.916*Angstrom, 5.4054*Angstrom)
elements = [Silicon, Silicon, Silicon, Oxygen, Oxygen, Oxygen, Oxygen, Oxygen,
            Oxygen]

fractional_coordinates = [[ 0.4697,  0.0000,  0.00000 ],
                          [ 0.0000,  0.4697,  0.66667 ],
                          [ 0.5303,  0.5303,  0.33333 ],
                          [ 0.4135,  0.2669,  0.11910 ],
                          [ 0.2669,  0.4135,  0.547567],
                          [ 0.7331,  0.1466,  0.785767],
                          [ 0.5865,  0.8534,  0.214233],
                          [ 0.8534,  0.5865,  0.452433],
                          [ 0.1466,  0.7331,  0.8809  ]]

bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates
    )

# Define the TremoloX-calculator with a ReaxFF-potential.
potentialSet = ReaxFF_CHONSiPtZrYBaTi_2013()
calculator = TremoloXCalculator(parameters=potentialSet)

bulk_configuration.setCalculator(calculator)
bulk_configuration.update()

# Set up the PartialCharges analysis object.
partial_charges = PartialCharges(bulk_configuration)

# Print the partial charges to the log file.
nlprint(partial_charges)

# Extract and store the partial charges into an array.
q = partial_charges.evaluate().inUnitsOf(elementary_charge)