# --- Create the initial configuration ---
# Set up lattice
vector_a = [8.1714, 0.0, 0.0]*Angstrom
vector_b = [0.0, 8.1714, 0.0]*Angstrom
vector_c = [0.0, 0.0, 8.1714]*Angstrom
lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
elements = [Silver, Silver, Silver, Silver, Silver, Silver, Silver, Silver,
            Silver, Silver, Silver, Silver, Silver, Silver, Silver, Silver,
            Silver, Silver, Silver, Silver, Silver, Silver, Silver, Silver,
            Silver, Silver, Silver, Silver, Silver, Silver, Silver, Silver]

# Define coordinates
fractional_coordinates = [[ 0.  ,  0.  ,  0.  ],
                          [ 0.  ,  0.  ,  0.5 ],
                          [ 0.  ,  0.5 ,  0.  ],
                          [ 0.  ,  0.5 ,  0.5 ],
                          [ 0.5 ,  0.  ,  0.  ],
                          [ 0.5 ,  0.  ,  0.5 ],
                          [ 0.5 ,  0.5 ,  0.  ],
                          [ 0.5 ,  0.5 ,  0.5 ],
                          [ 0.25,  0.25,  0.  ],
                          [ 0.25,  0.25,  0.5 ],
                          [ 0.25,  0.75,  0.  ],
                          [ 0.25,  0.75,  0.5 ],
                          [ 0.75,  0.25,  0.  ],
                          [ 0.75,  0.25,  0.5 ],
                          [ 0.75,  0.75,  0.  ],
                          [ 0.75,  0.75,  0.5 ],
                          [ 0.25,  0.  ,  0.25],
                          [ 0.25,  0.  ,  0.75],
                          [ 0.25,  0.5 ,  0.25],
                          [ 0.25,  0.5 ,  0.75],
                          [ 0.75,  0.  ,  0.25],
                          [ 0.75,  0.  ,  0.75],
                          [ 0.75,  0.5 ,  0.25],
                          [ 0.75,  0.5 ,  0.75],
                          [ 0.  ,  0.25,  0.25],
                          [ 0.  ,  0.25,  0.75],
                          [ 0.  ,  0.75,  0.25],
                          [ 0.  ,  0.75,  0.75],
                          [ 0.5 ,  0.25,  0.25],
                          [ 0.5 ,  0.25,  0.75],
                          [ 0.5 ,  0.75,  0.25],
                          [ 0.5 ,  0.75,  0.75]]

# Set up configuration
bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates
    )

# --- Generate a substitutional alloy ---
# Create a list of selected atoms. A percentage of these will be
# randomly substituted by the new element.
selection = [0, 1, 2, 3, 4, 5,
             6, 7, 8, 9, 10, 11,
             12, 13, 14, 15, 16,
             17, 18, 19, 20, 21,
             22, 23, 24, 25, 26,
             27, 28, 29, 30, 31]

configuration = substitutionalAlloy(
    configuration=bulk_configuration,
    indices=numpy.array(selection),
    algorithm=NormalDistribution,
    percentage=50.0,
    new_element=Gold
)