Things I did in optimization:

1. As we discussed last week, I obtained the initial condition for optimizing the disorder potential from `scipy.linalg.lstsq` in the following manner:

temp_voltages = voltages.copy()
tunable_gates_array = np.hstack(tunable_gates.values())

greens_function_matrix = []

for i in range(len(np.hstack(tunable_gates.values()))):

    temp_voltages = dict.fromkeys(temp_voltages, 0)
    temp_voltages['gate_'+tunable_gates_array[i]] = 1.0
    
    linear_problem, potential = gate_potential(
        poisson_system,
        linear_problem,
        site_coords[:, [0, 1]],
        site_indices,
        temp_voltages,
        charges,
        [],
        mixed_charges,
        random_generator = None,
        general_shape_charges = None,
        update_pos = [],
    )
    total_potential = potential.copy()

    greens_function_matrix.append(np.array(list(total_potential.values())))
    
greens_function_matrix = np.stack(greens_function_matrix).T

rhs = clean_potential - disorder_potential
lhs = greens_function_matrix

from scipy.linalg import lstsq as least_squares_norm

initial_condition = least_squares_norm(lhs, rhs)[0]

initial_condition += list(clean_solution.values())[:len(np.hstack(list(tunable_gates.values())))]

Here I first build the green's function matrix and leave it to the least squares function.

This method seems to work fine for two dots as it moves the initial condition bit closer to ideal solution, but can't find a big marginal improvement for four dots.

2. Next, I split the optimization process into two steps. First bring the onsite energies close to 0meV. Then optimize the tunnel couplings and dot positions. This method seems to perform much better than findng the initial condition. First optimizing the voltages with onsite energies as close as possible to 0eV with any arbitrary tunnel couplings seems to be a good bet.

3. However, during the second step of optimization, when we optimize the onsite energies, tunnel couplings and dot positons all together, optimizer converges to a solution in which the onsite energies are few tens of meV away from 0eV (without two-step method, optimizer is stuck in the region which is 100's meV away from ideal value).

4. To further improve the two-step optimization method, I tried to understand why the onsite energies are moved away from 0eV when optimizing tunnel couplings and dot positions. Instead of choosing a single position for dots, I chose a range along y-axis (other than the gates). This seems to improve but not close to ideal solution yet.

5. Then I changed the weights of the tunnel couplings in the cost function. This brings the onsite energies close around 5meV but leaves interdot hoppings within the range of ~30-50microeV around the ideal value (20microeV).

cost = {
        "onsite": (1, onsite_cost),
        "hopping": (0.01, hopping_cost),
#         "symmetry": (0.01, symmetry_cost )
    }

6. This shows that we need to do some sort of adaptive optimization where we leave the dot positions within certain range as well as the interdot hoppings.

7. I also thought about changing the relative weights of parameters in the cost function such that the optimizer converges somewhat close to the ideal solution within reasonable amount of time.


`Nutshell - I suspect that we are too strict in optimizing the potential that might not even be required.`