commons-dev mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From Evan Ward <evan.w...@nrl.navy.mil>
Subject Re: [math] refactoring least squares
Date Mon, 24 Feb 2014 19:41:34 GMT

On 02/24/2014 01:23 PM, Gilles wrote:
> On Mon, 24 Feb 2014 11:49:26 -0500, Evan Ward wrote:
>> One way to improve performance would be to provide pre-allocated space
>> for the Jacobian and reuse it for each evaluation.
>
> Do you have actual data to back this statement?

I did some tests with the CircleVectorial problem from the test cases.
The Jacobian is 1000x2, and I ran it 1000 times until hotspot stopped
making it run faster. The first column is the current state of the code.
The second column is with one less matrix allocation each problem
evaluation.

       trunk    -1 alloc  %change
lu     0.90 s   0.74      -17%
chol   0.90     0.74      -17%
qr     0.87     0.70      -20%


I also see similar reductions in runtime using 1e6 observations and 3
observations. We could save 3-4 allocations per evaluation, which
extrapolates to 60%-80% in run time.

>> The
>> LeastSquaresProblem interface would then be:
>>
>> void evaluate(RealVector point, RealVector resultResiduals, RealVector
>> resultJacobian);
>>
>> I'm interested in hearing your ideas on other approaches to solve this
>> issue. Or even if this is an issue worth solving.
>
> Not before we can be sure that in-place modification (rather than
> reallocation) always provides a performance benefit.

I would like to hear other ideas for improving the performance.

>
>
> Best Regards,
> Gilles
>
>> Evan
>

On 02/24/2014 12:09 PM, Luc Maisonobe wrote:
> Hi Evan,
>
> Le 24/02/2014 17:49, Evan Ward a écrit :
>> I've looked into improving performance further, but it seems any further
>> improvements will need big API changes for memory management.
>>
>> Currently using Gauss-Newton with Cholesky (or LU) requires 4 matrix
>> allocations _each_ evaluation. The objective function initially
>> allocates the Jacobian matrix. Then the weights are applied through
>> matrix multiplication, allocating a new matrix. Computing the normal
>> equations allocates a new matrix to hold the result, and finally the
>> decomposition allocates it's own matrix as a copy. With QR there are 3
>> matrix allocations each model function evaluation, since there is no
>> need to compute the normal equations, but the third allocation+copy is
>> larger. Some empirical sampling data I've collected with the jvisualvm
>> tool indicates that matrix allocation and copying takes 30% to 80% of
>> the execution time, depending on the dimension of the Jacobian.
>>
>> One way to improve performance would be to provide pre-allocated space
>> for the Jacobian and reuse it for each evaluation. The
>> LeastSquaresProblem interface would then be:
>>
>> void evaluate(RealVector point, RealVector resultResiduals, RealVector
>> resultJacobian);
>>
>> I'm interested in hearing your ideas on other approaches to solve this
>> issue. Or even if this is an issue worth solving.
> Yes, I think this issue is worth solving, especially since we are going
> to ship 3.3 and need to fix as much as possible before the release, thus
> avoiding future problems. Everything spotted now is worth fixing now.
>
> Your approach seems reasonable, as long as the work arrays are really
> allocated at the start of the optimization and shared only through a few
> documented methods like the one you propose. This would mean we can say
> in the javadoc that these area should be used only to fulfill the API
> requirements and not copied elsewhere, as they *will* be modified as the
> algorithm run, and are explicitly devoted to avoid reallocation. I guess
> this kind of problems is more important when lots of observations are
> performed, which correspond to very frequent use case (at least in the
> fields I know about).
>
> For the record, what you propose seems similar to what is done in the
> ODE package, as the state vector and its first derivatives are also kept
> in preallocated arrays which are reused throughout the integration and
> are used to exchange data between the Apache Commons Math algorithm and
> the user problem to be solved. So it is somehting we already do elsewhere.

OK. We could keep the Evaluation interface, which would just reference
the pre-allocated residuals and matrix. If the result parameters are
null the LSP could allocate a matrix of the correct size automatically.
So then the interface would look like:

Evaluation evaluate(RealVector point, RealVector resultResiduals,
RealVector resultJacobian);

>
> best regards,
> Luc
>
>> Best Regards,
>> Evan
>>
>>
>> ---------------------------------------------------------------------
>> To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
>> For additional commands, e-mail: dev-help@commons.apache.org
>>
>>
>>
> ---------------------------------------------------------------------
> To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
> For additional commands, e-mail: dev-help@commons.apache.org
>


On 02/24/2014 01:16 PM, Gilles wrote:
> On Mon, 24 Feb 2014 18:09:26 +0100, Luc Maisonobe wrote:
>> Hi Evan,
>>
>> Le 24/02/2014 17:49, Evan Ward a écrit :
>>> I've looked into improving performance further, but it seems any
>>> further
>>> improvements will need big API changes for memory management.
>>>
>>> Currently using Gauss-Newton with Cholesky (or LU) requires 4 matrix
>>> allocations _each_ evaluation. The objective function initially
>>> allocates the Jacobian matrix. Then the weights are applied through
>>> matrix multiplication, allocating a new matrix. Computing the normal
>>> equations allocates a new matrix to hold the result, and finally the
>>> decomposition allocates it's own matrix as a copy. With QR there are 3
>>> matrix allocations each model function evaluation, since there is no
>>> need to compute the normal equations, but the third allocation+copy is
>>> larger. Some empirical sampling data I've collected with the jvisualvm
>>> tool indicates that matrix allocation and copying takes 30% to 80% of
>>> the execution time, depending on the dimension of the Jacobian.
>>>
>>> One way to improve performance would be to provide pre-allocated space
>>> for the Jacobian and reuse it for each evaluation. The
>>> LeastSquaresProblem interface would then be:
>>>
>>> void evaluate(RealVector point, RealVector resultResiduals, RealVector
>>> resultJacobian);
>>>
>>> I'm interested in hearing your ideas on other approaches to solve this
>>> issue. Or even if this is an issue worth solving.
>>
>> Yes, I think this issue is worth solving, especially since we are going
>> to ship 3.3 and need to fix as much as possible before the release, thus
>> avoiding future problems. Everything spotted now is worth fixing now.
>>
>> Your approach seems reasonable, as long as the work arrays are really
>> allocated at the start of the optimization and shared only through a few
>> documented methods like the one you propose. This would mean we can say
>> in the javadoc that these area should be used only to fulfill the API
>> requirements and not copied elsewhere, as they *will* be modified as the
>> algorithm run, and are explicitly devoted to avoid reallocation. I guess
>> this kind of problems is more important when lots of observations are
>> performed, which correspond to very frequent use case (at least in the
>> fields I know about).
>>
>> For the record, what you propose seems similar to what is done in the
>> ODE package, as the state vector and its first derivatives are also kept
>> in preallocated arrays which are reused throughout the integration and
>> are used to exchange data between the Apache Commons Math algorithm and
>> the user problem to be solved. So it is somehting we already do
>> elsewhere.
>
> If I understand correctly what is being discussed, I do not agree with
> this approach.
>
> The optimization/fitting algorithms must use matrix abstractions.
> If performance improvements can achieved, they must happen at the level
> of the appropriate matrix implementations.
>

The matrix abstractions will still be used in the interface. As far as I
can tell none of the optimizers or linear algebra classes use the matrix
abstractions internally. For example LU, QR, and Cholesky all copy the
matrix data to an internal double[][]. I tried computing the normal
equation in GaussNewton as j.transpose().multiply(j), but the
performance was bad because j.transpose() creates a copy of the matrix.
That's why we have the current ugly for loop implementation with
getEntry() and setEntry(). Maybe matrix "views" could help solve the issue.

>
> Best regards,
> Gilles
>
>
> ---------------------------------------------------------------------
> To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
> For additional commands, e-mail: dev-help@commons.apache.org
>

Regards,
Evan



---------------------------------------------------------------------
To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
For additional commands, e-mail: dev-help@commons.apache.org


Mime
View raw message