# Any hints on implementing this optimization problem solver?

Started by May 9, 2012
Dear all,

It's been a while I've last posted to this group. Nevertheless, I've always appreciated the quality of answers due to due richness of experience of its members.

I have the following problem:

Given a digital signal from a sensor, I need to extract some properties of this signal. I cannot give further details on the sensor, nor on the data or on the properties to extract.

I could extract the wanted parameters by using quite simple maths. The accuracy of extracted parameter was acceptable. However, now the requirements on the accuracy have changed and they are much more severe than before. I have to use another method: solve a nonlinear optimization problem.

There is no simpler way to extract the parameters with the required accuracy. My maths program shows me that the accuracy using the Levenberg-Marquardt algorithm is good, as well as the convergence. With a simple gradient descent or a line search I have slow convergence and a poor dependency on the initial guess of the parameters.

So far so good.

The main problem is that I have to implement this solver on a *fixed-point* DSP and I am afraid of the following things:
(A) computing the objective function, its Jacobian and its Hessian matrix
(B) solving a linear equation system

With the gradient descent approach, there is no (B), and (A) is simplified. However, there are the disadvantages stated above.

Does anyone have practical experience with this? In particular on a fixed-point DSP? Which linear equation solver to choose? Or does someone have even some bunch of code to inspire myself?
Any practical hints or considerations?

Thanks a lot

Andreas

_____________________________________
Hi Andreas,

If you are constrained to use a fixed point processor, you would need to
simulate floating point by software subroutines. It might be quite
questionable that a nonlinear optimization can be implemented in fixed point.
My guess is that even a Q64.32 format would not have enough range and, for
small values, enough significant digits to ensure stability in general
case.

It might happen though that your problem is well conditioned and can be scaled
so that both the function and derivatives would have enough significant
digits in a fixed point format. But again, this means that you will have
to recode all the subroutines to work with fixed point numbers instead
of floating point, at the price of loosing generality of floating point
calculations.

Regarding an implementation, there are quite a few solvers for money, like
NAG or IMSL. They do not provide with the sources, and they most probably do
not have a C6x DSP port.

There is an open source GNU Scientific library that have an implementation
of the LM algorithm in Nonlinear Least-Squares Fitting.

Another option is to look at MINPACK on netlib.org, but it is written in
Fortran.

For the B you might want to use the same gsl, or obtain Linpack package at
http://www.sundancedsp.com. There is clapack on netlib.org, which is free.

Rgds,

Andrew

P.S. A steepest descend method is used just to illustrate how things work, for
real problems at least a conjugated gradient need to be used to ensure
practical convergence rate.

> Subject: Any hints on implementing this optimization problem solver?
> Posted by: "a...@gmail.com" a...@gmail.com
> Date: Tue May 8, 2012 8:43 pm ((PDT))
>
> Dear all,
>
> It's been a while I've last posted to this group. Nevertheless, I've always
> appreciated the quality of answers due to due richness of experience of its
> members.
>
> I have the following problem:
>
> Given a digital signal from a sensor, I need to extract some properties of
> this signal. I cannot give further details on the sensor, nor on the data or
> on the properties to extract.
>
> I could extract the wanted parameters by using quite simple maths. The
> accuracy of extracted parameter was acceptable. However, now the requirements
> on the accuracy have changed and they are much more severe than before. I
> have to use another method: solve a nonlinear optimization problem.
>
> There is no simpler way to extract the parameters with the required
> accuracy. My maths program shows me that the accuracy using the
> Levenberg-Marquardt algorithm is good, as well as the convergence. With a
> simple gradient descent or a line search I have slow convergence and a poor
> dependency on the initial guess of the parameters.
>
> So far so good.
>
> The main problem is that I have to implement this solver on a *fixed-point*
> DSP and I am afraid of the following things:
> (A) computing the objective function, its Jacobian and its Hessian matrix
> (B) solving a linear equation system
>
> With the gradient descent approach, there is no (B), and (A) is simplified.
> However, there are the disadvantages stated above.
>
> Does anyone have practical experience with this? In particular on a
> fixed-point DSP? Which linear equation solver to choose? Or does someone have
> even some bunch of code to inspire myself?
> Any practical hints or considerations?
>
> Thanks a lot
>
> Andreas

_____________________________________
*079CCE2A* MVK .S2 0x388C ,SP

*0780006A* MVKH .S2 0x0000 ,SP

*07BF07A2 *AND .S2 -8,SP,SP

*071ACE2A* MVK .S2 0x359C,DP

*0700006A* MVKH .S2 0x0000 ,DP

*020000FA* ZERO .L2 B4

*091003A2 *MVC .S2 B4 , FADCR

*0A1003A2 *MVC .S2 B4,FMCR

*020C902A *MVK .S2 0x1920,B4

*0200006A* MVKH .S2 0x0000,B4

*00100362 *B.S2 B4
Andreas-

Agree with Andrew -- get it working with floating-point software emulation first. It sounds like in your situation
results are key: if you can verify your results are correct, then you know "how much faster" you have to be, and can
start making an optimization plan. For example, some code sections may not need floating-point, some sections might
be able to use longer-precision intrinsics, or maybe can otherwise take advantage of SIMD/pipeline features of the
device.

Speaking of that: which DSP are you using? What speed?

-Jeff

> Hi Andreas,
>
> If you are constrained to use a fixed point processor, you would need to
> simulate floating point by software subroutines. It might be quite
> questionable that a nonlinear optimization can be implemented in fixed point.
> My guess is that even a Q64.32 format would not have enough range and, for
> small values, enough significant digits to ensure stability in general
> case.
>
> It might happen though that your problem is well conditioned and can be scaled
> so that both the function and derivatives would have enough significant
> digits in a fixed point format. But again, this means that you will have
> to recode all the subroutines to work with fixed point numbers instead
> of floating point, at the price of loosing generality of floating point
> calculations.
>
> Regarding an implementation, there are quite a few solvers for money, like
> NAG or IMSL. They do not provide with the sources, and they most probably do
> not have a C6x DSP port.
>
> There is an open source GNU Scientific library that have an implementation
> of the LM algorithm in Nonlinear Least-Squares Fitting.
>
> Another option is to look at MINPACK on netlib.org, but it is written in
> Fortran.
>
> For the B you might want to use the same gsl, or obtain Linpack package at
> http://www.sundancedsp.com. There is clapack on netlib.org, which is free.
>
> Rgds,
>
> Andrew
>
> P.S. A steepest descend method is used just to illustrate how things work, for
> real problems at least a conjugated gradient need to be used to ensure
> practical convergence rate.
>
>> Subject: Any hints on implementing this optimization problem solver?
>> Posted by: "a...@gmail.com" a...@gmail.com
>> Date: Tue May 8, 2012 8:43 pm ((PDT))
>>
>> Dear all,
>>
>> It's been a while I've last posted to this group. Nevertheless, I've always
>> appreciated the quality of answers due to due richness of experience of its
>> members.
>>
>> I have the following problem:
>>
>> Given a digital signal from a sensor, I need to extract some properties of
>> this signal. I cannot give further details on the sensor, nor on the data or
>> on the properties to extract.
>>
>> I could extract the wanted parameters by using quite simple maths. The
>> accuracy of extracted parameter was acceptable. However, now the requirements
>> on the accuracy have changed and they are much more severe than before. I
>> have to use another method: solve a nonlinear optimization problem.
>>
>> There is no simpler way to extract the parameters with the required
>> accuracy. My maths program shows me that the accuracy using the
>> Levenberg-Marquardt algorithm is good, as well as the convergence. With a
>> simple gradient descent or a line search I have slow convergence and a poor
>> dependency on the initial guess of the parameters.
>>
>> So far so good.
>>
>> The main problem is that I have to implement this solver on a *fixed-point*
>> DSP and I am afraid of the following things:
>> (A) computing the objective function, its Jacobian and its Hessian matrix
>> (B) solving a linear equation system
>>
>> With the gradient descent approach, there is no (B), and (A) is simplified.
>> However, there are the disadvantages stated above.
>>
>> Does anyone have practical experience with this? In particular on a
>> fixed-point DSP? Which linear equation solver to choose? Or does someone have
>> even some bunch of code to inspire myself?
>> Any practical hints or considerations?
>>
>> Thanks a lot
>>
>> Andreas

_____________________________________
Dear Andrew, dear Jeff,

Thanks for the valuable hints. When I first thought about the implementation, the only thing I had in mind was variable scaling for each single computation. This is quite a nightmare. I didn't think of floating-point simulation.

Thanks to you, I think this is the way to go - if in the end I will really need to go the way of nonlinear optimization. In the mean time I found a way to improve results a lot by keeping the computational overhead low. However, those results are not as good as for nonlinear optimization but they seem to be sufficient.

The other thing is that the computational resources are low and the requirements on execution time are high:
- the processor is of TI's 16-bit low-power family. therefore I think going beyond 32-bits is not an option. anyway, it seems that the problem is still well conditioned for 32-bit float.
- the number of input samples for the optimization problem is rather low, the number of parameters as well (it is not a large-scale optimization problem)
- however, I can only afford to have 50us computation time until convergence
- with my customized Levenberg-Marquardt implementation convergence is reached at <8 iterations, which makes approx 6us per iteration - or 1000 cycles on the DSP I use.
. I think this is pretty tight to compute the model, its derivatives, Hessian matrix, residual, the parameter update etc.

But what's good is that I now know a lot of possibilities to solve the problem and the trade-offs they imply. I'm sure I'll find a good way.

Thanks a lot.
Have a nice week end,

Andreas
>
>It's been a while I've last posted to this group. Nevertheless, I've always appreciated the quality of answers due to due richness of experience of its members.
>
>I have the following problem:
>
>Given a digital signal from a sensor, I need to extract some properties of this signal. I cannot give further details on the sensor, nor on the data or on the properties to extract.
>
>I could extract the wanted parameters by using quite simple maths. The accuracy of extracted parameter was acceptable. However, now the requirements on the accuracy have changed and they are much more severe than before. I have to use another method: solve a nonlinear optimization problem.
>
>There is no simpler way to extract the parameters with the required accuracy. My maths program shows me that the accuracy using the Levenberg-Marquardt algorithm is good, as well as the convergence. With a simple gradient descent or a line search I have slow convergence and a poor dependency on the initial guess of the parameters.
>
>So far so good.
>
>The main problem is that I have to implement this solver on a *fixed-point* DSP and I am afraid of the following things:
>(A) computing the objective function, its Jacobian and its Hessian matrix
>(B) solving a linear equation system
>
>With the gradient descent approach, there is no (B), and (A) is simplified. However, there are the disadvantages stated above.
>
>Does anyone have practical experience with this? In particular on a fixed-point DSP? Which linear equation solver to choose? Or does someone have even some bunch of code to inspire myself?
>Any practical hints or considerations?
>
>Thanks a lot
>
>Andreas
>
>_____________________________________

_____________________________________