Improving Excel's solver speed for differential equations

355 Views Asked by At

I have an Excel spreadsheet that contains a linear differential 5 DOF vehicle dynamics model. I use the built-in solver in a loop to numerically solve the equations. I used it qualitatively for a while, but then in the search for increased accuracy I needed to dramatically reduce the time step. This highlighted the need to increase the macro speed and review what was a hastily written piece of code.

By using the Timer function I could establish the solver setup and solve times so I then set about trying to reduce the time taken for each loop. The improvements I implemented were:

  • Switching screen updating off
  • Setting sheet calculation to manual and only updating necessary rows
  • Switching the status bar update off
  • Switching off enable events
  • Using the Simplex LP engine
  • Reducing the required solving precision
  • Using the result of the previous time step as a start point for the next row
  • Moving all derived calculations to another sheet so only the necessary calculation are updated
  • Reconfiguring the spreadsheet so that the solving row is constant (the results are then copied to rows below), meaning that the solver only needs to be setup once.

The improvements were very substantial, dropping the solving time per row from about 1.3 to 0.35 seconds. The current time of 0.35 seconds is almost entirely solving time - the remaining few milliseconds is copying time.

However, at a minimum of 1000 rows per run, the total run time is still over 16 minutes. For a full run, it would be 85 minutes.

Although I'm pleased with the improvements, the model still seems to be slow for only 5 DOF. So, my question is does anyone know if I can improve the solving time or can comment on whether I've hit the bottom limit for 5 such equations? Has anyone done a back to back for a similar kind of model between Excel and Matlab (or Octave, perhaps) to compare solving times?

The equations are:

  1. Lateral forces:

$$(m \cdot \dot{v})+(m \cdot rad(r) \cdot u)+(m_s \cdot h' \cdot rad(\ddot{\varphi}))+(deg(C_{\alpha f}) \cdot rad(\alpha_f))+(deg(C_{\alpha r}) \cdot rad(\alpha_r)) = 0$$

  1. Yaw moments:

$$(I_z \cdot rad(\dot{r}))+(I_{z,s} \cdot rad(\Theta_r) - I_{xz,s}) \cdot rad(\ddot{\varphi}))+((a-t_{pf}) \cdot deg(C_{\alpha f}) \cdot rad(\alpha_f))-((b-t_{pr}) \cdot deg(C_{\alpha r)} \cdot rad(\alpha_r)) = 0$$

  1. Roll moments:

$$((I_{x,s}+(m_s \cdot h'^2)) \cdot rad(\ddot{\varphi}))+(m_s \cdot h' \cdot \dot{v})+(m_s \cdot h' \cdot rad(r) \cdot u)+(((I_z \cdot rad(\Theta_r))-I_{xz}) \cdot rad(\dot{r}))+((deg(C_\varphi)-(m_s \cdot g \cdot h')) \cdot rad(\varphi))+(deg(K_{\varphi}) \cdot rad(\dot\varphi)) = 0$$

  1. Front slip angular velocity balance:

$$rad(\dot\alpha_f)-(v/\sigma_f)-((a \cdot rad(r))/\sigma_f)+((u \cdot \varepsilon_f \cdot rad(\varphi))/\sigma_f)+(((c_f \cdot C_{\alpha f})+1) \cdot ((u \cdot rad(\alpha_f))/\sigma_f))+((u \cdot rad(\delta))/(\sigma_f \cdot i_s)) = 0$$

  1. Rear slip angular velocity balance:

$$rad(\dot\alpha_r)-(v/\sigma_r)+((b \cdot rad(r))/\sigma_r)+((u \cdot \varepsilon_r \cdot rad(\varphi))/\sigma_r)+(((c_r \cdot C_{\alpha r})+1) \cdot ((u \cdot rad(\alpha_r))/\sigma_r)) = 0$$

Where (symbol: meaning [unit]; value_typical ; value_simplified):

$m$: mass [kg]; 2692; 950

$m_s$: sprung mass [kg]; 2412; 0

$u$: forward speed [m/s]; 27.8; 35.76311111

$a$: CG to front axle [m]; 1.608; 0.88

$b$: CG to rear axle [m]; 1.335; 1.82

$C_{\alpha f}$: front cornering stiffness [N/deg]; 3225; 1066

$C_{\alpha r}$: rear cornering stiffness [N/deg]; 3400; 900

$I_z$: yaw moment of inertia [kgm^2]; 4130; 1200

$i_s$: steering ratio [-]; 14.3; 1

$\varepsilon_f$: front roll steer coefficient [-]; 0.052; 0

$\varepsilon_r$: rear roll steer coefficient [-]; -0.18; 0

$\Theta_r$: roll axis inclination [deg]; 0.914941467; 0

$\sigma_f$: front lateral tyre relaxation length [m]; 0.6408; 0.001

$\sigma_r$: rear lateral tyre relaxation length [m]; 0.3952; 0.001

$C_\varphi$: total roll stiffness [Nm/deg]; 2650; 0

$c_f$: front compliance steer [deg/N]; 0.000064; 0

$c_r$: rear compliance steer [deg/N]; -0.000041; 0

$h'$: sprung mass CG height wrt roll axis [m]; 0.4435; 0

$I_{x,s}$: sprung mass roll moment of inertia [kgm^2]; 786 0

$I_{xz,s}$: sprung mass roll-yaw product of inertia [kgm^2]; 107; 0

$I_{z,s}$: sprung mass yaw moment of inertia [kgm^2]; 3320; 0

$K_{\varphi}$: roll damping [Nms/deg]; 200; 0

$t_{pf}$: front pneumatic trail [m]; 0.0364; 0

$t_{pr}$: rear pneumatic trail [m]; 0.0236; 0

$I_{xz}$: roll-yaw product of inertia [kgm^2]; 150; 0

and finally, the input time history:

$\delta$: steer angle [deg]; 22.92; 22.92

In the above I have included a set of values that are typical for a modern saloon and a set of values that represent a highly simplified model that I have used to correlate against a different model with further reduced DOF. The two sets are included to see if there is a significant difference in solving speed.

The input time history can be a simple step input, with a value of 0 deg from 0-0.1 seconds, then 22.92 deg from there onwards for a few seconds.

The equations are a little clunky because I wanted the angular parameters in degrees rather than radians.

Many thanks,

Simon.