# -------------------------------------------------------------
# Add Doping
# -------------------------------------------------------------
# Define doping density in units of 1/cm**3 and convert to 1/m**3
doping_density = 1e+19 * (1/(0.01*Meter)**3)

# Define the volume of the InAs region and convert in m**3
volume = 1334.15*Ang**3
volume = volume.convertTo(Meter**3)

# Get charge per electrode which is unitless (float)
q = float(doping_density * volume)
# Specify charge so that electrode doping is 10**19 cm**-3
# positive charge in left electrode, p-doping
left_charge = q
# negative charge in right electrode, n-doping
right_charge = -1*q

# Calculate the compensation charge left
# Number of elements in the left electrode
n_left = len(left_electrode.elements())
# Number of Hydrogen atoms in the left electrode
nh_left = sum([i == Hydrogen for i in left_electrode.elements()])
# Average ion charge on each non Hydrogen atom
left_compensation_charge = -left_charge/(n_left-nh_left)
left_electrode.setExternalPotential(
            AtomicCompensationCharge(
                                [('H_As', -0.25),('H_In', 0.25),
                                (Arsenic, left_compensation_charge),
                                (Indium, left_compensation_charge)])
                                    )

# Calculate the compensation charge right
n_right = len(right_electrode.elements())
nh_right = sum([i == Hydrogen for i in right_electrode.elements()])
right_compensation_charge = -right_charge/(n_right-nh_right)


right_electrode.setExternalPotential(
            AtomicCompensationCharge(
                                [('H_As', -0.25),('H_In', 0.25),
                                (Arsenic, right_compensation_charge),
                                (Indium, right_compensation_charge)])
                                    )

# Set doping for the central region
n_central = len(device_configuration.elements())
# Specify the size of the doped left central region
n_central_left = 1718
left_list = [i for i in range(n_central_left) \
                    if  device_configuration.elements()[i] != Hydrogen ]
central_region.addTags('LEFT', left_list)

# Specify the size of the doped right central region
# ( To make a 10 nm drain underlap change this value to n_central_right = 860)
n_central_right = 1718
right_list = [i for i in range(n_central-n_central_right, n_central) \
                    if  device_configuration.elements()[i] != Hydrogen ]
central_region.addTags('RIGHT', right_list)

external_potential_central_region = AtomicCompensationCharge([
                                ('H_As', -0.25),('H_In', 0.25),
                                ('LEFT', left_compensation_charge),
                                ('RIGHT', right_compensation_charge)
                                                             ])
central_region.setExternalPotential(external_potential_central_region)