Simulate physical magnetic collision quicker

I'm building a simulation that includes physical magnetised objects. The simulation consists of multiple rods of fixed length. On each rod, there are a number of spheres spread out evenly over the rod. Each sphere is magnetised and therefore exerts a force and torque on each other sphere. If spheres overlap, there is a force that pushes the spheres out each other to resemble the physical object. There is no force preventing the rods themselves to overlap or cross, only the spheres (for ease of modelling) I sum the forces and torques per rod to calculate the movement of each rod. I use a dynamic timestep that limits the movement in each iteration to 10% of sphere diameter. This way, the spheres cannot go through each other.

The problem: Once the spheres of two rods touch, the simulation freezes: The timestep that's taken becomes very small. This is due to the large magnetic forces once the spheres almost touch, and the large repelling force once the spheres overlap. It creates an oscillating behaviour. I can therefore not complete the simulation, as it would take too long.

Now, I dont include any code as I'm mainly wondering what is the appropriate strategy to deal with simulating this kind of physical interaction, and making the simulation go quicker.

I've tried:

modelling the outward-pushing force differently, e.g. springs of different stiffnesses, exponential functions.

Increasing sphere size (works a bit, but does not represent reality well, i.e. some spheres should be approx 8x larger than others).

I hope everything is clear and that someone can help me out!

Thank you in advance.

[Walter Roberson]

Dynamic time step? Are you using ode45 or one of the stiff solvers?

[Kaz Middelhoek]

With dynamic timestep I mean this: each iteration, I calculate the velocity of each rod due to all forces. I then choose dt such that v*dt is equal to 0.1*diameter of the smallest spheres. If velocities are high, the timestep I take is small. I don't use ode45 or any stiff solvers. I will read up on them. It does seem to be the problem I'm facing, so you've definitely pointed me in the right direction. Thanks!

[Sean Brennan]

If the forces are modeled as spring-type (or even exponential), and there are no forces related to the velocities of the spheres, then the final solution MUST oscillate. The damping effect that will give a final non-oscillatory solution will require some type of velocity penalty in the formulation and simulation.

Note that: if the forces are very high, then the simulation may predict motion larger than 1/10 the sphere size within one time step. In this case, the solver will have a very difficult time and almost certainly give poor results.

Also, time steps should be chosen by the time behavior of the system, not by the physical size. If you are seeing oscillations in your solution, this is actually helpful in finding a good time-step. A general guideline is to choose a time step at least 1/10 the size of the oscillation period, or smaller. So if you see oscillations occur every 1/100th of a second, then the time step you should use should at least be 1/1000th of a second (e.g. 0.001). If the time step is of the same size as the period, then the numeric solution will not be valid. Technically, the Nyquist theorem states that the time step can never be larger than 1/2 the period. So the 1/10 of the period (or more) is to add a much-needed factor of safety, as using the Nyquist rule will give poor results still. When doing simulations for the first time, I usually use 1/100th of the period as a starting point.