Right now we are spending a lot of time copying from the current Jacobian to a smaller matrix to use BLAS to multiply the kepler+drift Jacobian (jac_ij).
The BLAS routine seems to be faster than this!
So, it's not clear why the copying is so slow. It might be faster to simply insert the kepler+drift array into a full 7n x 7n Jacobian array and carry out the multiplication of the full matrix with BLAS. Not sure if sparse matrices can be used - the bodies which are not being touched will have the identity matrix for that sub-step.
Anyways, we need to see whether we can experiment with this and figure out how to speed it up.
To do:
[x] Edit kepler_driftij to put the jac_ij into a 7n x 7n matrix.
I tried this out for the outer solar system problem, and this added 67% to the run time. So, the BLAS multiplications are taking significant time, and it does pay to copy back and forth.
Right now we are spending a lot of time copying from the current Jacobian to a smaller matrix to use BLAS to multiply the kepler+drift Jacobian (jac_ij).
The BLAS routine seems to be faster than this!
So, it's not clear why the copying is so slow. It might be faster to simply insert the kepler+drift array into a full 7n x 7n Jacobian array and carry out the multiplication of the full matrix with BLAS. Not sure if sparse matrices can be used - the bodies which are not being touched will have the identity matrix for that sub-step.
Anyways, we need to see whether we can experiment with this and figure out how to speed it up.
To do:
kepler_driftij
to put thejac_ij
into a 7n x 7n matrix.jac_ij
back to the identity matrix.