Maximum $d$ you may get is a total number $N$ of $y_n$'s you have - 1. Larger $d$'s yield curvatures that can't be represented by N points.

Solve a matrix equation $A\vec{t}=\vec{y}$ where $\vec{y}$ is a measurements vector (size $N \times 1$) and $A$ is a normalized version of matrix $X$ (size $N \times N$), where each column $A_j$ is calculated as:

$A_{ij} = \frac {X_{ij}} {\sum \limits_j (X_{ij})^2}$

$X = \left [ \begin{matrix} x_0^0 & \cdots &x_0^{N-1} \\ \vdots & \ddots & \vdots \\ x_{N-1}^0 & \cdots & x_{N-1}^{N-1} \end{matrix} \right ]$

d you are interested in is a position of largest coefficient in a vector $\vec{t}$ (size $N \times 1$ starting from zero). Your approximation would be $y = t_d \cdot x^d$.

This is a sort of exhaustive search, yes.

Maximum $d$ you may get is a total number $N$ of $y_n$'s you have - 1. Larger $d$'s yield curvatures that can't be represented by N points.

Solve a matrix equation $A\vec{t}=\vec{y}$ where $\vec{y}$ is a measurements vector (size $N \times 1$) and $A$ is a normalized version of matrix $X$ (size $N \times N$), where each column $A_j$ is calculated as:

$A_{ij} = \frac {X_{ij}} {\sum \limits_j (X_{ij})^2}$

$X = \left [ \begin{matrix} x_0^0 & \cdots &x_0^{N-1} \\ \vdots & \ddots & \vdots \\ x_{N-1}^0 & \cdots & x_{N-1}^{N-1} \end{matrix} \right ]$

d you are interested in is a position of largest coefficient in a vector $\vec{t}$ (size $N \times 1$ starting from zero). Your approximation would be $y = t_d \cdot x^d$.

You can achieve better results if your solver finds solution optimal in $L^1$ sense.

This is a sort of exhaustive search, yes.