How to use the approximate model¶
You have already seen a first example of the setup_sum_of_squares_problem function in the Getting started tutorial. In this guide, we are going to revisit that function in more detail.
Setting up a simple model problem¶
For this guide, we will reuse the simple model problem from the Getting started tutorial; the Ising model on a line of 10 sites:
\[\hat{\mathcal{H}}_{\text{Ising}} = \sum_{i=1}^{9} J_{i,(i+1)} Z_i Z_{(i+1)} + \sum_{i=1}^{10} h_i X_i \, ,\]
where \(J\) is the coupling strength between two sites and \(h\) is the external magnetic field.
[1]:
from qiskit.transpiler import CouplingMap
from qiskit_addon_utils.problem_generators import generate_xyz_hamiltonian
coupling_map = CouplingMap.from_line(10, bidirectional=False)
hamiltonian = generate_xyz_hamiltonian(
coupling_map,
coupling_constants=(0.0, 0.0, 1.0),
ext_magnetic_field=(0.4, 0.0, 0.0),
)
print(hamiltonian)
SparsePauliOp(['IIIIIIIZZI', 'IIIIIZZIII', 'IIIZZIIIII', 'IZZIIIIIII', 'IIIIIIIIZZ', 'IIIIIIZZII', 'IIIIZZIIII', 'IIZZIIIIII', 'ZZIIIIIIII', 'IIIIIIIIIX', 'IIIIIIIIXI', 'IIIIIIIXII', 'IIIIIIXIII', 'IIIIIXIIII', 'IIIIXIIIII', 'IIIXIIIIII', 'IIXIIIIIII', 'IXIIIIIIII', 'XIIIIIIIII'],
coeffs=[1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j,
1. +0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j,
0.4+0.j, 0.4+0.j, 0.4+0.j])
Choosing \(k_j\)¶
In this guide, we are going to focus on the tweaks that we can apply to the cvxpy.Problem and its solver that are returned by the setup_sum_of_squares_problem. Therefore, we are not concerned with choosing any specific \(k_j\) values and can simply reuse the ones from the Getting started tutorial.
[2]:
trotter_steps = (8, 12, 19)
Computing \(x\) analytically¶
First, we have to construct our LSE, \(Ax=b\). Once again, we reuse the settings from the Getting started tutorial.
[3]:
from qiskit_addon_mpf.static import setup_static_lse
lse = setup_static_lse(trotter_steps, order=2, symmetric=True)
Next, we compute the analytical solution of \(x\) and its L1-norm as our reference values.
[4]:
coeff_analytical = lse.solve()
print(coeff_analytical)
[ 0.17239057 -1.19447005 2.02207947]
[5]:
import numpy as np
print(np.linalg.norm(coeff_analytical, ord=1))
3.388940092165899
Setting up the approximate model¶
Now, we will construct the approximate model. The idea of using approximate solutions \(\tilde{x}\) in order to ensure a smaller L1-norm was proposed by Zhuk et al., 2023. The model is already explained in detail in the API documentation of setup_sum_of_squares_problem so we refrain from repeating it
here. Suffice to say, that the model consists of two constraints which ensure that the coefficients, \(\tilde{x}\), sum to \(1\) and it enforces their L1-norm to be smaller than the user-provided max_l1_norm value. The optimization objective is to minimize the deviation of \(A\tilde{x}=b\).
We can see all of this in the output below. Here, we are setting max_l1_norm=3 because we already know that the L1-norm of the exact solution is approximately \(3.4\) and we want to find a \(\tilde{x}\) that has a smaller L1-norm.
[6]:
from qiskit_addon_mpf.costs import setup_sum_of_squares_problem
model, coeff_var = setup_sum_of_squares_problem(lse, max_l1_norm=3)
print(model)
minimize quad_over_lin(Vstack([1. 1. 1.] @ x + -1.0, [0.015625 0.00694444 0.00277008] @ x + -0.0, [2.44140625e-04 4.82253086e-05 7.67336039e-06] @ x + -0.0), 1.0)
subject to Sum(x, None, False) == 1.0
norm1(x) <= 3.0
Solving the model with default settings¶
First, let us see what happens when we solve the model with only default settings.
[7]:
model.solve(verbose=True)
===============================================================================
CVXPY
v1.5.3
===============================================================================
(CVXPY) Sep 05 09:28:06 AM: Your problem has 3 variables, 2 constraints, and 0 parameters.
(CVXPY) Sep 05 09:28:06 AM: It is compliant with the following grammars: DCP, DQCP
(CVXPY) Sep 05 09:28:06 AM: (If you need to solve this problem multiple times, but with different data, consider using parameters.)
(CVXPY) Sep 05 09:28:06 AM: CVXPY will first compile your problem; then, it will invoke a numerical solver to obtain a solution.
(CVXPY) Sep 05 09:28:06 AM: Your problem is compiled with the CPP canonicalization backend.
-------------------------------------------------------------------------------
Compilation
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:06 AM: Compiling problem (target solver=OSQP).
(CVXPY) Sep 05 09:28:06 AM: Reduction chain: CvxAttr2Constr -> Qp2SymbolicQp -> QpMatrixStuffing -> OSQP
(CVXPY) Sep 05 09:28:06 AM: Applying reduction CvxAttr2Constr
(CVXPY) Sep 05 09:28:06 AM: Applying reduction Qp2SymbolicQp
(CVXPY) Sep 05 09:28:06 AM: Applying reduction QpMatrixStuffing
(CVXPY) Sep 05 09:28:06 AM: Applying reduction OSQP
(CVXPY) Sep 05 09:28:06 AM: Finished problem compilation (took 1.053e-02 seconds).
-------------------------------------------------------------------------------
Numerical solver
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:06 AM: Invoking solver OSQP to obtain a solution.
-----------------------------------------------------------------
OSQP v0.6.3 - Operator Splitting QP Solver
(c) Bartolomeo Stellato, Goran Banjac
University of Oxford - Stanford University 2021
-----------------------------------------------------------------
problem: variables n = 9, constraints m = 11
nnz(P) + nnz(A) = 33
settings: linear system solver = qdldl,
eps_abs = 1.0e-05, eps_rel = 1.0e-05,
eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
rho = 1.00e-01 (adaptive),
sigma = 1.00e-06, alpha = 1.60, max_iter = 10000
check_termination: on (interval 25),
scaling: on, scaled_termination: off
warm start: on, polish: on, time_limit: off
iter objective pri res dua res rho time
1 0.0000e+00 1.00e+00 2.00e+02 1.00e-01 6.26e-05s
200 4.0661e-09 7.52e-06 1.44e-08 1.22e-05 1.28e-04s
status: solved
solution polish: unsuccessful
number of iterations: 200
optimal objective: 0.0000
run time: 1.53e-04s
optimal rho estimate: 2.09e-06
-------------------------------------------------------------------------------
Summary
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Problem status: optimal
(CVXPY) Sep 05 09:28:07 AM: Optimal value: 4.092e-09
(CVXPY) Sep 05 09:28:07 AM: Compilation took 1.053e-02 seconds
(CVXPY) Sep 05 09:28:07 AM: Solver (including time spent in interface) took 1.028e-03 seconds
[7]:
4.092215822096244e-09
We can see from the final value of the cost function that the sum of squares of \(A\tilde{x}-b\) is almost zero.
In the output above, we also get a lot more information from cvxpy.Problem.solve, some of which we will get back to later.
For now, let us inspect the final values of \(\tilde{x}\) and its L1-norm:
[8]:
coeff_approx_plain = coeff_var.value
print(coeff_approx_plain)
print(np.linalg.norm(coeff_approx_plain, ord=1))
[-0.39646076 0.55556835 0.84089365]
1.7929227539822414
We can clearly see that the L1-norm has decreased compared to the analytical solution that we have computed before.
Nonetheless, let us see if we can tweak our solution a bit. In the end, we will compare the results obtained for this set of coefficients with all the other ones.
Tweaking the approximate model¶
Rather than tweaking the input parameters like the chosen trotter_steps and max_l1_norm values, we will focus on the solver settings in this guide.
Using a different solver¶
Notice, that we used the default solver method: OSQP. We can change that solver to see whether this has any effect.
[9]:
model.solve(verbose=True, solver="CLARABEL")
===============================================================================
CVXPY
v1.5.3
===============================================================================
(CVXPY) Sep 05 09:28:07 AM: Your problem has 3 variables, 2 constraints, and 0 parameters.
(CVXPY) Sep 05 09:28:07 AM: It is compliant with the following grammars: DCP, DQCP
(CVXPY) Sep 05 09:28:07 AM: (If you need to solve this problem multiple times, but with different data, consider using parameters.)
(CVXPY) Sep 05 09:28:07 AM: CVXPY will first compile your problem; then, it will invoke a numerical solver to obtain a solution.
(CVXPY) Sep 05 09:28:07 AM: Your problem is compiled with the CPP canonicalization backend.
-------------------------------------------------------------------------------
Compilation
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Compiling problem (target solver=CLARABEL).
(CVXPY) Sep 05 09:28:07 AM: Reduction chain: Dcp2Cone -> CvxAttr2Constr -> ConeMatrixStuffing -> CLARABEL
(CVXPY) Sep 05 09:28:07 AM: Applying reduction Dcp2Cone
(CVXPY) Sep 05 09:28:07 AM: Applying reduction CvxAttr2Constr
(CVXPY) Sep 05 09:28:07 AM: Applying reduction ConeMatrixStuffing
(CVXPY) Sep 05 09:28:07 AM: Applying reduction CLARABEL
(CVXPY) Sep 05 09:28:07 AM: Finished problem compilation (took 8.657e-03 seconds).
-------------------------------------------------------------------------------
Numerical solver
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Invoking solver CLARABEL to obtain a solution.
-------------------------------------------------------------
Clarabel.rs v0.8.1 - Clever Acronym
(c) Paul Goulart
University of Oxford, 2022
-------------------------------------------------------------
problem:
variables = 9
constraints = 11
nnz(P) = 3
nnz(A) = 30
cones (total) = 2
: Zero = 1, numel = 4
: Nonnegative = 1, numel = 7
settings:
linear algebra: direct / qdldl, precision: 64 bit
max iter = 200, time limit = Inf, max step = 0.990
tol_feas = 1.0e-8, tol_gap_abs = 1.0e-8, tol_gap_rel = 1.0e-8,
static reg : on, ϵ1 = 1.0e-8, ϵ2 = 4.9e-32
dynamic reg: on, ϵ = 1.0e-13, δ = 2.0e-7
iter refine: on, reltol = 1.0e-13, abstol = 1.0e-12,
max iter = 10, stop ratio = 5.0
equilibrate: on, min_scale = 1.0e-4, max_scale = 1.0e4
max iter = 10
iter pcost dcost gap pres dres k/t μ step
---------------------------------------------------------------------------------------------
0 +7.1341e-05 -2.3333e+00 2.33e+00 2.40e-01 7.10e-01 1.00e+00 1.88e+00 ------
1 +7.1330e-05 -1.5432e-01 1.54e-01 1.64e-02 1.12e-01 3.28e-02 1.80e-01 9.61e-01
2 +6.9440e-05 -1.7266e-03 1.80e-03 2.10e-04 1.89e-03 5.65e-04 2.43e-03 9.87e-01
3 +2.2053e-05 -5.0002e-05 7.21e-05 5.56e-06 4.83e-05 1.98e-05 7.80e-05 9.87e-01
4 +1.1119e-06 -2.5543e-05 2.67e-05 1.80e-07 1.41e-06 3.21e-06 1.04e-05 9.90e-01
5 +7.9250e-09 -2.0484e-06 2.06e-06 2.06e-09 1.59e-08 2.11e-07 6.43e-07 9.90e-01
6 +2.2543e-09 -2.8390e-08 3.06e-08 2.07e-11 1.59e-10 3.12e-09 9.42e-09 9.90e-01
7 +1.8480e-09 -8.5733e-10 2.71e-09 1.79e-12 2.36e-10 3.19e-10 9.40e-10 9.15e-01
---------------------------------------------------------------------------------------------
Terminated with status = Solved
solve time = 113.289µs
-------------------------------------------------------------------------------
Summary
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Problem status: optimal
(CVXPY) Sep 05 09:28:07 AM: Optimal value: 1.848e-09
(CVXPY) Sep 05 09:28:07 AM: Compilation took 8.657e-03 seconds
(CVXPY) Sep 05 09:28:07 AM: Solver (including time spent in interface) took 8.504e-04 seconds
[9]:
1.8479854472132682e-09
[10]:
coeff_approx_clarabel = coeff_var.value
print(coeff_approx_clarabel)
print(np.linalg.norm(coeff_approx_clarabel, ord=1))
[-0.21284352 -0.00792683 1.22077035]
1.4415406931636157
Indeed, the outcome has changed. But now we have one major contribution from the largest \(k_j\) value which is also not ideal.
An exhaustive search of the various available solvers is beyond the scope of this guide. Instead we encourage you to try this yourself.
Changing the convergence thresholds¶
Rather than changing the solver altogether, we can also tweak the convergence parameters. For the OSQP solver the relevant ones are called eps_abs and eps_rel and they are both set to 1e-5 by default.
[11]:
model.solve(verbose=True, eps_abs=1e-8, eps_rel=1e-8)
===============================================================================
CVXPY
v1.5.3
===============================================================================
(CVXPY) Sep 05 09:28:07 AM: Your problem has 3 variables, 2 constraints, and 0 parameters.
(CVXPY) Sep 05 09:28:07 AM: It is compliant with the following grammars: DCP, DQCP
(CVXPY) Sep 05 09:28:07 AM: (If you need to solve this problem multiple times, but with different data, consider using parameters.)
(CVXPY) Sep 05 09:28:07 AM: CVXPY will first compile your problem; then, it will invoke a numerical solver to obtain a solution.
(CVXPY) Sep 05 09:28:07 AM: Your problem is compiled with the CPP canonicalization backend.
-------------------------------------------------------------------------------
Compilation
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Compiling problem (target solver=OSQP).
(CVXPY) Sep 05 09:28:07 AM: Reduction chain: CvxAttr2Constr -> Qp2SymbolicQp -> QpMatrixStuffing -> OSQP
(CVXPY) Sep 05 09:28:07 AM: Applying reduction CvxAttr2Constr
(CVXPY) Sep 05 09:28:07 AM: Applying reduction Qp2SymbolicQp
(CVXPY) Sep 05 09:28:07 AM: Applying reduction QpMatrixStuffing
(CVXPY) Sep 05 09:28:07 AM: Applying reduction OSQP
(CVXPY) Sep 05 09:28:07 AM: Finished problem compilation (took 8.835e-03 seconds).
-------------------------------------------------------------------------------
Numerical solver
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Invoking solver OSQP to obtain a solution.
-----------------------------------------------------------------
OSQP v0.6.3 - Operator Splitting QP Solver
(c) Bartolomeo Stellato, Goran Banjac
University of Oxford - Stanford University 2021
-----------------------------------------------------------------
problem: variables n = 9, constraints m = 11
nnz(P) + nnz(A) = 33
settings: linear system solver = qdldl,
eps_abs = 1.0e-08, eps_rel = 1.0e-08,
eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
rho = 1.00e-01 (adaptive),
sigma = 1.00e-06, alpha = 1.60, max_iter = 10000
check_termination: on (interval 25),
scaling: on, scaled_termination: off
warm start: on, polish: on, time_limit: off
iter objective pri res dua res rho time
1 0.0000e+00 1.00e+00 2.00e+02 1.00e-01 5.02e-05s
200 4.0661e-09 7.52e-06 1.44e-08 1.22e-05 2.10e-04s
400 3.8209e-09 4.95e-06 2.88e-09 2.09e-06 3.72e-04s
600 3.3607e-09 6.09e-06 1.96e-09 2.09e-06 4.33e-04s
800 2.8975e-09 6.24e-06 1.84e-09 2.09e-06 4.87e-04s
1000 2.4704e-09 6.08e-06 1.70e-09 2.09e-06 5.40e-04s
1200 2.0921e-09 5.78e-06 1.56e-09 2.09e-06 5.92e-04s
1400 1.7645e-09 5.40e-06 1.44e-09 2.09e-06 6.44e-04s
1600 1.4845e-09 5.01e-06 1.32e-09 2.09e-06 7.01e-04s
1800 1.2471e-09 4.63e-06 1.21e-09 2.09e-06 7.53e-04s
2000 1.0467e-09 4.26e-06 1.11e-09 2.09e-06 8.05e-04s
2200 8.7801e-10 3.91e-06 1.01e-09 2.09e-06 8.57e-04s
2400 7.3890e-10 3.25e-06 7.74e-10 2.09e-06 9.10e-04s
2600 6.4003e-10 2.47e-06 7.11e-10 2.09e-06 9.62e-04s
2800 5.6692e-10 2.03e-06 6.68e-10 2.09e-06 1.01e-03s
3000 5.0872e-10 1.76e-06 6.33e-10 2.09e-06 1.07e-03s
3200 4.5991e-10 1.58e-06 6.02e-10 2.09e-06 1.12e-03s
3400 4.1754e-10 1.46e-06 5.74e-10 2.09e-06 1.17e-03s
3600 3.7999e-10 1.37e-06 5.47e-10 2.09e-06 1.22e-03s
3800 3.4628e-10 1.29e-06 5.22e-10 2.09e-06 1.32e-03s
4000 3.1580e-10 1.22e-06 4.99e-10 2.09e-06 1.37e-03s
4200 2.8813e-10 1.16e-06 4.76e-10 2.09e-06 1.43e-03s
4400 2.6295e-10 1.11e-06 4.55e-10 2.09e-06 1.48e-03s
4600 2.4000e-10 1.06e-06 4.35e-10 2.09e-06 1.53e-03s
4800 2.1907e-10 1.01e-06 4.15e-10 2.09e-06 1.58e-03s
5000 1.9997e-10 9.67e-07 3.97e-10 2.09e-06 1.63e-03s
5200 1.8255e-10 9.23e-07 3.79e-10 2.09e-06 1.69e-03s
5400 1.6664e-10 8.82e-07 3.62e-10 2.09e-06 1.74e-03s
5600 1.5219e-10 1.01e-06 2.60e-10 2.09e-06 1.79e-03s
5800 1.4049e-10 6.59e-07 2.47e-10 2.09e-06 1.84e-03s
6000 1.3096e-10 5.75e-07 2.39e-10 2.09e-06 1.90e-03s
6200 1.2274e-10 5.24e-07 2.31e-10 2.09e-06 1.95e-03s
6400 1.1538e-10 4.90e-07 2.24e-10 2.09e-06 2.00e-03s
6600 1.0863e-10 4.67e-07 2.17e-10 2.09e-06 2.05e-03s
6800 1.0237e-10 4.48e-07 2.11e-10 2.09e-06 2.10e-03s
7000 9.6519e-11 4.32e-07 2.05e-10 2.09e-06 2.16e-03s
7200 9.1024e-11 4.19e-07 1.99e-10 2.09e-06 2.21e-03s
7400 8.5854e-11 4.06e-07 1.93e-10 2.09e-06 2.29e-03s
7600 8.0984e-11 3.94e-07 1.88e-10 2.09e-06 2.36e-03s
7800 7.6393e-11 3.82e-07 1.82e-10 2.09e-06 2.41e-03s
8000 7.2065e-11 3.71e-07 1.77e-10 2.09e-06 2.46e-03s
8200 6.7982e-11 3.60e-07 1.72e-10 2.09e-06 2.51e-03s
8400 6.4131e-11 3.50e-07 1.67e-10 2.09e-06 2.56e-03s
8600 6.0499e-11 3.40e-07 1.62e-10 2.09e-06 2.62e-03s
8800 5.7072e-11 3.30e-07 1.57e-10 2.09e-06 2.67e-03s
9000 5.3859e-11 8.77e-06 8.55e-12 2.09e-06 2.72e-03s
9200 5.1785e-11 5.78e-07 2.18e-13 2.09e-06 2.79e-03s
9400 5.0736e-11 7.72e-08 8.32e-14 2.09e-06 2.87e-03s
9550 5.0303e-11 4.69e-08 8.73e-14 2.09e-06 2.96e-03s
plsh 4.9640e-11 1.94e-14 1.50e-12 -------- 3.00e-03s
status: solved
solution polish: successful
number of iterations: 9550
optimal objective: 0.0000
run time: 3.00e-03s
optimal rho estimate: 2.96e-06
-------------------------------------------------------------------------------
Summary
-------------------------------------------------------------------------------
(CVXPY) Sep 05 09:28:07 AM: Problem status: optimal
(CVXPY) Sep 05 09:28:07 AM: Optimal value: 4.964e-11
(CVXPY) Sep 05 09:28:07 AM: Compilation took 8.835e-03 seconds
(CVXPY) Sep 05 09:28:07 AM: Solver (including time spent in interface) took 3.988e-03 seconds
[11]:
4.9640211694779584e-11
[12]:
coeff_approx_tight = coeff_var.value
print(coeff_approx_tight)
print(np.linalg.norm(coeff_approx_tight, ord=1))
[ 0.10925065 -1. 1.89074935]
3.0
As you can see, we have now obtained a comparatively higher L1-norm (which is still constrained by our choice of max_l1_norm=3). While the coefficient of the smallest \(k_j\) value is still fairly small in magnitude, at least the other two provide a reasonable amount of mixing.
Comparing our results¶
We will now briefly repeat the computations from the Getting started tutorial in order to assess the quality of the various coefficients.
[13]:
from qiskit.synthesis import SuzukiTrotter
from qiskit_addon_utils.problem_generators import generate_time_evolution_circuit
time = 8.0
circuits = []
for k in trotter_steps:
circ = generate_time_evolution_circuit(
hamiltonian,
synthesis=SuzukiTrotter(order=2, reps=k),
time=time,
)
circuits.append(circ)
[14]:
from qiskit.quantum_info import SparsePauliOp
L = 10
observable = SparsePauliOp.from_sparse_list([("Z", [i], 1 / L / 2) for i in range(L)], num_qubits=L)
print(observable)
SparsePauliOp(['IIIIIIIIIZ', 'IIIIIIIIZI', 'IIIIIIIZII', 'IIIIIIZIII', 'IIIIIZIIII', 'IIIIZIIIII', 'IIIZIIIIII', 'IIZIIIIIII', 'IZIIIIIIII', 'ZIIIIIIIII'],
coeffs=[0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j,
0.05+0.j, 0.05+0.j, 0.05+0.j])
[15]:
from qiskit.primitives import StatevectorEstimator
estimator = StatevectorEstimator()
job = estimator.run([(circ, observable) for circ in circuits])
result = job.result()
[16]:
evs = [res.data.evs for res in result]
print(evs)
[array(0.23799162), array(0.35754312), array(0.38649906)]
[17]:
from scipy.linalg import expm
exp_H = expm(-1j * time * hamiltonian.to_matrix())
initial_state = np.zeros(exp_H.shape[0])
initial_state[0] = 1.0
time_evolved_state = exp_H @ initial_state
exact_obs = time_evolved_state.conj() @ observable.to_matrix() @ time_evolved_state
[18]:
print("Analytical coeff solution:", evs @ coeff_analytical)
print("Plain approx. solution:", evs @ coeff_approx_plain)
print("CLARABEL approx. solution:", evs @ coeff_approx_clarabel)
print("Tight approx. solution:", evs @ coeff_approx_tight)
print("Exact solution:", exact_obs.real)
Analytical coeff solution: 0.39548478559794153
Plain approx. solution: 0.42928990901998165
CLARABEL approx. solution: 0.41833743748333746
Tight approx. solution: 0.3992304700751579
Exact solution: 0.40060242487900255
We can clearly see now, that the various approximate solutions can yield quite different results. Luckily, testing out different expansion coefficients for a fixed choice of \(k_j\) values does not require us to rerun our quantum experiments. Therefore, we encourage you to play around and tweak your expansion coefficients when using the approximate solutions.