|
| 1 | +"""Riemannian Adam Optimizer for Variational Quantum Circuits. |
| 2 | +
|
| 3 | +This module implements Adam [1]_ with a manifold-aware retraction step, |
| 4 | +for optimizing variational quantum circuit parameters on Riemannian |
| 5 | +manifolds (Bloch sphere). |
| 6 | +
|
| 7 | +Standard Adam, as shipped in ``qiskit_algorithms.optimizers.ADAM``, treats |
| 8 | +parameters as living in unconstrained Euclidean space: its ``bounds`` |
| 9 | +support level is ``ignored``, so nothing stops iterates from drifting |
| 10 | +outside a bounded interval such as ``[0, pi]``. |
| 11 | +
|
| 12 | +Gate-rotation angles instead live on a flat manifold (a product of circles |
| 13 | +for full-period parameters, a bounded interval for others). Because that |
| 14 | +manifold has zero curvature, the moment estimates need no parallel |
| 15 | +transport between tangent spaces — only a retraction step after each |
| 16 | +update: wrap angles with a ``2*pi`` period back into range, and clip |
| 17 | +angles with a non-periodic bound (e.g. QIOCE's ``[0, pi]``) instead of |
| 18 | +letting them escape the valid domain. |
| 19 | +
|
| 20 | +Gradients are approximated via central finite differences, so no |
| 21 | +analytical gradient is required. |
| 22 | +
|
| 23 | +References |
| 24 | +---------- |
| 25 | +.. [1] Kingma, D. P., & Ba, J. (2015). Adam: A method for stochastic |
| 26 | + optimization. arXiv:1412.6980. |
| 27 | +.. [2] Becigneul, G., & Ganea, O. E. (2019). Riemannian adaptive |
| 28 | + optimization methods. ICLR. |
| 29 | +""" |
| 30 | + |
| 31 | +import numpy as np |
| 32 | +from qiskit_algorithms.optimizers import ( |
| 33 | + Optimizer, |
| 34 | + OptimizerResult, |
| 35 | + OptimizerSupportLevel, |
| 36 | +) |
| 37 | + |
| 38 | + |
| 39 | +class RiemannianAdamOptimizer(Optimizer): |
| 40 | + """Adam optimizer with manifold-aware retraction for VQC parameters, |
| 41 | + inspired by [1]_. |
| 42 | +
|
| 43 | + Parameters |
| 44 | + ---------- |
| 45 | + maxiter : int, default=100 |
| 46 | + Maximum number of iterations. |
| 47 | + lr : float, default=0.1 |
| 48 | + Learning rate. |
| 49 | + beta1 : float, default=0.9 |
| 50 | + Exponential decay rate for the first moment estimate. |
| 51 | + beta2 : float, default=0.999 |
| 52 | + Exponential decay rate for the second moment estimate. |
| 53 | + eps : float, default=1e-8 |
| 54 | + Term added to the denominator for numerical stability. |
| 55 | + fd_epsilon : float, default=1e-5 |
| 56 | + Finite-difference step size for gradient approximation. |
| 57 | + Should be small (~1e-5 to 1e-7) for accurate numerical gradients. |
| 58 | + tol : float, default=1e-6 |
| 59 | + Convergence tolerance on the gradient norm. |
| 60 | +
|
| 61 | + Attributes |
| 62 | + ---------- |
| 63 | + trajectory_ : list of ndarray |
| 64 | + Optimization trajectory (parameter history). |
| 65 | + loss_history_ : list of float |
| 66 | + Loss function values at each iteration. |
| 67 | +
|
| 68 | + Examples |
| 69 | + -------- |
| 70 | + >>> from pyriemann_qiskit.optimization.riemannian_adam import ( |
| 71 | + ... RiemannianAdamOptimizer |
| 72 | + ... ) |
| 73 | + >>> optimizer = RiemannianAdamOptimizer(maxiter=100, lr=0.1) |
| 74 | + >>> # Use with QuanticNCH or other quantum classifiers |
| 75 | + >>> # qaoa_optimizer=optimizer |
| 76 | +
|
| 77 | + Notes |
| 78 | + ----- |
| 79 | + .. versionadded:: 0.7.0 |
| 80 | +
|
| 81 | + References |
| 82 | + ---------- |
| 83 | + .. [1] Becigneul, G., & Ganea, O. E. (2019). Riemannian adaptive |
| 84 | + optimization methods. ICLR. |
| 85 | + """ |
| 86 | + |
| 87 | + def __init__( |
| 88 | + self, |
| 89 | + maxiter=100, |
| 90 | + lr=0.1, |
| 91 | + beta1=0.9, |
| 92 | + beta2=0.999, |
| 93 | + eps=1e-8, |
| 94 | + fd_epsilon=1e-5, |
| 95 | + tol=1e-6, |
| 96 | + ): |
| 97 | + super().__init__() |
| 98 | + self._maxiter = maxiter |
| 99 | + self._lr = lr |
| 100 | + self._beta1 = beta1 |
| 101 | + self._beta2 = beta2 |
| 102 | + self._eps = eps |
| 103 | + self._fd_epsilon = fd_epsilon |
| 104 | + self._tol = tol |
| 105 | + self.trajectory_ = [] |
| 106 | + self.loss_history_ = [] |
| 107 | + |
| 108 | + def minimize(self, fun, x0, jac=None, bounds=None): |
| 109 | + """Minimize the objective function using Riemannian Adam. |
| 110 | +
|
| 111 | + Parameters |
| 112 | + ---------- |
| 113 | + fun : callable |
| 114 | + Objective function to minimize. Should accept a 1D array |
| 115 | + and return a scalar. |
| 116 | + x0 : ndarray |
| 117 | + Initial parameter vector. |
| 118 | + jac : callable, optional |
| 119 | + Gradient function (not used, included for compatibility). |
| 120 | + bounds : list of tuples, optional |
| 121 | + Parameter bounds. A bound spanning exactly ``2*pi`` is |
| 122 | + treated as periodic (wrapped); any other bound is clipped. |
| 123 | +
|
| 124 | + Returns |
| 125 | + ------- |
| 126 | + OptimizerResult |
| 127 | + Result object containing: |
| 128 | +
|
| 129 | + - x: optimal parameters |
| 130 | + - fun: objective value at optimal parameters |
| 131 | + - nfev: number of function evaluations |
| 132 | + - nit: number of iterations |
| 133 | + """ |
| 134 | + x = np.array(x0, dtype=float) |
| 135 | + n = len(x) |
| 136 | + nfev = 0 |
| 137 | + |
| 138 | + f_current = fun(x) |
| 139 | + nfev += 1 |
| 140 | + |
| 141 | + self.trajectory_ = [x.copy()] |
| 142 | + self.loss_history_ = [f_current] |
| 143 | + |
| 144 | + # Precompute vectorised bounds arrays once (avoids repeated per-element |
| 145 | + # work). Only treat as periodic if the span is a full 2*pi rotation |
| 146 | + # period. Bounds like [0, pi] (QIOCE) must be clipped, not wrapped — |
| 147 | + # wrapping with period pi teleports parameters to the wrong end of |
| 148 | + # the range. |
| 149 | + if bounds is not None: |
| 150 | + lower_b = np.array([b[0] if b[0] is not None else np.nan for b in bounds]) |
| 151 | + upper_b = np.array([b[1] if b[1] is not None else np.nan for b in bounds]) |
| 152 | + span = upper_b - lower_b |
| 153 | + periodic = ~(np.isnan(lower_b) | np.isnan(upper_b)) & np.isclose( |
| 154 | + span, 2 * np.pi |
| 155 | + ) |
| 156 | + period_b = np.where(periodic, span, 1.0) # dummy for non-periodic |
| 157 | + has_lower = ~np.isnan(lower_b) & ~periodic |
| 158 | + has_upper = ~np.isnan(upper_b) & ~periodic |
| 159 | + |
| 160 | + m = np.zeros(n) |
| 161 | + v = np.zeros(n) |
| 162 | + |
| 163 | + iteration = 0 |
| 164 | + for iteration in range(self._maxiter): |
| 165 | + # Estimate gradient using central differences. |
| 166 | + # Modify x[i] in-place and restore — avoids 2n array allocations |
| 167 | + # per iter. Clamp perturbations to stay within bounds so loss is |
| 168 | + # never evaluated at an illegal parameter value. |
| 169 | + grad = np.zeros(n) |
| 170 | + for i in range(n): |
| 171 | + orig = x[i] |
| 172 | + step = self._fd_epsilon |
| 173 | + |
| 174 | + x_plus = orig + step |
| 175 | + x_minus = orig - step |
| 176 | + if bounds is not None: |
| 177 | + lo, hi = bounds[i] |
| 178 | + if lo is not None: |
| 179 | + x_minus = max(x_minus, lo) |
| 180 | + if hi is not None: |
| 181 | + x_plus = min(x_plus, hi) |
| 182 | + actual_step = (x_plus - x_minus) / 2.0 |
| 183 | + |
| 184 | + x[i] = x_plus |
| 185 | + f_plus = fun(x) |
| 186 | + nfev += 1 |
| 187 | + x[i] = x_minus |
| 188 | + f_minus = fun(x) |
| 189 | + nfev += 1 |
| 190 | + x[i] = orig |
| 191 | + if actual_step != 0.0: |
| 192 | + grad[i] = (f_plus - f_minus) / (2.0 * actual_step) |
| 193 | + |
| 194 | + grad_norm = np.linalg.norm(grad) |
| 195 | + if grad_norm < self._tol: |
| 196 | + break |
| 197 | + |
| 198 | + # Adam moment updates (Euclidean — the manifold here is flat, |
| 199 | + # so no parallel transport is needed between iterates). |
| 200 | + t = iteration + 1 |
| 201 | + m = self._beta1 * m + (1 - self._beta1) * grad |
| 202 | + v = self._beta2 * v + (1 - self._beta2) * grad**2 |
| 203 | + m_hat = m / (1 - self._beta1**t) |
| 204 | + v_hat = v / (1 - self._beta2**t) |
| 205 | + |
| 206 | + x_new = x - self._lr * m_hat / (np.sqrt(v_hat) + self._eps) |
| 207 | + |
| 208 | + # Riemannian retraction back onto the parameter manifold — |
| 209 | + # vectorised, same convention as AndersonAccelerationOptimizer. |
| 210 | + if bounds is not None: |
| 211 | + x_new = np.where( |
| 212 | + periodic, |
| 213 | + lower_b + np.mod(x_new - lower_b, period_b), |
| 214 | + x_new, |
| 215 | + ) |
| 216 | + x_new = np.where(has_lower, np.maximum(x_new, lower_b), x_new) |
| 217 | + x_new = np.where(has_upper, np.minimum(x_new, upper_b), x_new) |
| 218 | + |
| 219 | + x = x_new |
| 220 | + f_current = fun(x) |
| 221 | + nfev += 1 |
| 222 | + |
| 223 | + self.trajectory_.append(x.copy()) |
| 224 | + self.loss_history_.append(f_current) |
| 225 | + |
| 226 | + result = OptimizerResult() |
| 227 | + result.x = x |
| 228 | + result.fun = f_current |
| 229 | + result.nfev = nfev |
| 230 | + result.nit = iteration + 1 if iteration >= 0 else 0 |
| 231 | + return result |
| 232 | + |
| 233 | + @property |
| 234 | + def settings(self): |
| 235 | + """Return optimizer settings.""" |
| 236 | + return { |
| 237 | + "maxiter": self._maxiter, |
| 238 | + "lr": self._lr, |
| 239 | + "beta1": self._beta1, |
| 240 | + "beta2": self._beta2, |
| 241 | + "eps": self._eps, |
| 242 | + "fd_epsilon": self._fd_epsilon, |
| 243 | + "tol": self._tol, |
| 244 | + } |
| 245 | + |
| 246 | + def get_support_level(self): |
| 247 | + """Return support level dictionary.""" |
| 248 | + return { |
| 249 | + "gradient": OptimizerSupportLevel.ignored, |
| 250 | + "bounds": OptimizerSupportLevel.supported, |
| 251 | + "initial_point": OptimizerSupportLevel.required, |
| 252 | + } |
![pre-commit-ci[bot]](https://avatars.githubusercontent.com/in/68672?v=4&size=40){kind=avatar}
0 commit comments