Numerical Integration of PBH Spin-Up Equation Over Extended Time Span

Hello Everyone,

I am trying to solve the spin-up equation for a primordial black hole (PBH) system numerically over a longer integration time. The spin-up equation is given by:
fdot = C_GWf^{11/3} + C_DFf^{3/2}
At t=0, f=f_isco​. Using time reversal, my equation becomes:
fdot = -C_GWf^{11/3} - C_DFf^{3/2}
with the initial condition t=0, f=f_isco​.

I am solving this equation using solve_ivp with a step size based on the Nyquist theorem, i.e., step_size = 1/2*f_isco. With this step size, my time vector is t= (0,1e+3,step_size).This setup allows my system to evolve up to 200 Hz.

However, I want to integrate over a longer time span, such as t = (0,1e+8,step_size), to observe the effects of the second term in the spin-up equation, which are usually noticeable below 50 Hz. The problem is that this requires a huge number of steps (around 1e+14), leading to memory errors.

Can anyone suggest a method to solve this equation for a longer time span without running into memory issues?

Thanks!

Thank you for the question.

I’m not really an expert on this, but I can offer some suggestions. Are you only trying to work out the relationship between f and fdot, or is it something more complicated? For example, are you keeping track of the phase information?

One thing you could try is using a different step size for near merger, and for times far away from merger. You could evolve back to (say) 1e3, and then use those conditions as a starting point to evolve farther back with a different step size.

Good luck!

Hi Jonah, thanks for the suggestion. I used a similar method, breaking down my system into smaller parts and using short integration times and small steps. I’m mainly focused on how the frequency changes over time.

Thanks again.