# Solving for an operator by minimization

Please note that I am looking for numerical algorithms that will tell me what the operator is that minimizes a problem.

I have a 2x2 complex hermitian operator that is a function of two variables, so $\tilde{O} = \tilde{O}(x,y)$. I will not have this operator in a closed form, but rather as the output of a numberical simulation who's inputs are just $x$ and $y$. I want to find out which values of $x$ and $y$ will allow me to map a given vector $\psi_{in}$ to another given vector $\phi_{1}$. I feel that some kind of gradient descent algorithm is in order. For instance, I can define

$d = | \tilde{O}(x,y)\psi_{in} - \phi_{1} |$

Now I can minimize this value over the two input values $x$ and $y$. So I would have an equation like this:

$\frac{\partial^2 d}{\partial_x \partial_y} = \frac{\partial^2 | \tilde{O}(x,y)\psi_{in} - \phi_{1} |}{\partial_x \partial_y}$

We might do some minimization from here but I am not sure how to procede and I feel there is a simpler numerical recipe that someone might know of.

Any help will be most appreciated.

-
You gave the answer yourself: Numerical Recipes is the name of the book you want to look for your answers in. Just to give you some keywords, you are trying to find roots of a multivariate nonlinear equation, and presumably you don't want to have to calculate gradients. – Yoav Kallus May 5 '12 at 0:04

Why not use the standard numerical methods for solving a system of equations, available in Maple, Matlab, Mathematica etc?

-
I think the problem is that $\tilde O(x,y)$ is not available in a form where one could analytically calculate derivatives of it. So these will have to be approximated by finite differences, regardless whether the solution is direct (e.g. by Newton's method) or by minimization. However, packages that can do an optimization problem in two variables, even if analytical derivatives of the function to be minimized are unavailable, are not hard to find. – Michael Renardy May 4 '12 at 18:25
There are many pitfalls in approximating derivatives by finite differences, and there is a growing literature on alternatives; please see my answer below. – Pait May 5 '12 at 19:04

You haven't said so, but I'm assuming that $\psi$ and $\phi$ are vectors. These could more generally be functions in some function space, and you would typically discretize those functions to work with vectors in your numerical computations.

I'm also assuming that you have one or more $\psi$, $\phi$ examples containing measured values or values that have been computed by some numerical model. Call these examples $\psi_{1}$, $\phi_{1}$, $\psi_{2}$, $\phi_{2}$, $\ldots$, $\psi_{n}$,$\phi_{n}$.

There may well be no exact solution to $O(x,y)\psi=\phi$ that works for all of these examples, so a typical approach would be to turn this into a nonlinear least squares problem of the form

$\min_{x,y} \sum_{k=1}^{n} \| O(x,y)\psi_{k}-\phi_{k} \|_{2}^{2}$

Your data and parameters are complex valued, but these are most easily dealt with by splitting the real and imaginary parts to produce a nonlinear least squares problems involving only real parameters and measured values.

Once you've reduced the problem to a conventional real nonlinear least squares problem, you can use standard methods to solve it. The Levenberg-Marquadt method is most commonly used in practice. In the LM method, you can use finite difference approximations to get the required derivatives- most software for doing this has the ability to do this finite differencing for you.

-
Hi, Yes $\phi$ and $\psi$ are 2 element complex vectors. Elements of a 2-d Hilbert space. – Ben Sprott May 4 '12 at 19:08
Thanks very much for the suggestion! – Ben Sprott May 4 '12 at 19:10
One important issue with the LM method is that it a local search method. In some circumstances, if your least squares problem has local minima it may return a solution that is locally but not globally optimal. It's possible in that case that you'd miss a solution that really did satisfy $O(x,y)\psi=\phi$. A simple approach that can help with this problem is to start LM running from many different initial guesses and then take the best solution found. – Brian Borchers May 4 '12 at 19:13
I should mention that I do not need to minimize for more than one pair at a time. I will want to take a single pair and find x,y that form a solution, and then a new pair will come along and new x,y values will do. Also, $\psi$ is fixed. – Ben Sprott May 4 '12 at 20:04
Will gradient descent be fine in the case I have mentioned? – Ben Sprott May 4 '12 at 20:05

You are interested in what is called direct, or derivative-free optimization. The algorithm can only use a "zero-order oracle" which, when asked, supplies the value of the function at a point. It used to be a not-very-popular subject with a few ad hoc algorithms. The title of a 1996 paper by Margaret Wright, "Direct Search Methods: Once Scorned, Now Respectable" says a lot. It is available at

http://cm.bell-labs.com/cm/cs/doc/96/4-02.ps.gz

A recent book which gives a readable overview of the field is

http://www.amazon.com/Introduction-Derivative-Free-Optimization-Mps-Siam-Series/dp/0898716683

Another book that mentions direct optimization is

http://www.amazon.com/Introductory-Lectures-Convex-Optimization-Applied/dp/1402075537/ref=sr_1_1?s=books&ie=UTF8&qid=1336182876&sr=1-1

-
What happens when I am not evaluating a "function at a point" but an operator at a point and in particular, an operator that can be written as a noncommuting product of operators each defined on different variables or components of the "point". The "point" is (x,y) and the operator is $\tilde{O_1}(x) \tilde{O_1}(y)$. If you have any thoughts, they would be helpful. – Ben Sprott May 17 '12 at 14:54
If each operator can be parametrized by a finite vector (essentially, a list of real-valued parameters), then the direct optimization algorithms are worth exploring. Otherwise, I do not know. – Pait Aug 7 '12 at 0:14