Source code for pathsim.solvers.rkdp87

########################################################################################
##
##                EXPLICIT ADAPTIVE TIMESTEPPING RUNGE-KUTTA INTEGRATORS
##                                 (solvers/rkf78.py)
##
##                                 Milan Rother 2024
##
########################################################################################

# IMPORTS ==============================================================================

from ._rungekutta import ExplicitRungeKutta


# SOLVERS ==============================================================================

[docs] class RKDP87(ExplicitRungeKutta): """13-stage 8-th order embedded Runge-Kutta method from Dormand and Prince with embedded 7-th order method for 8-th order truncation error estimate that can be used to adaptively control the timestep. This solver is a great choice if extremely high accuracy is required. """ def __init__(self, *solver_args, **solver_kwargs): super().__init__(*solver_args, **solver_kwargs) #number of stages in RK scheme self.s = 13 #order of scheme and embedded method self.n = 8 self.m = 7 #flag adaptive timestep solver self.is_adaptive = True #intermediate evaluation times self.eval_stages = [0.0, 1/18, 1/12, 1/8, 5/16, 3/8, 59/400, 93/200, 5490023248/9719169821, 13/20, 1201146811/1299019798, 1.0, 1.0] #extended butcher table self.BT = {0:[1/18], 1:[1/48, 1/16], 2:[1/32, 0, 3/32], 3:[5/16, 0, -75/64, 75/64], 4:[3/80, 0, 0, 3/16, 3/20], 5:[29443841/614563906, 0, 0, 77736538/692538347, -28693883/1125000000, 23124283/1800000000], 6:[16016141/946692911, 0, 0, 61564180/158732637, 22789713/633445777, 545815736/2771057229, -180193667/1043307555], 7:[39632708/573591083, 0, 0, -433636366/683701615, -421739975/2616292301, 100302831/723423059, 790204164/839813087, 800635310/3783071287], 8:[246121993/1340847787, 0, 0, -37695042795/15268766246, -309121744/1061227803, -12992083/490766935, 6005943493/2108947869, 393006217/1396673457, 123872331/1001029789], 9:[-1028468189/846180014, 0, 0, 8478235783/508512852, 1311729495/1432422823, -10304129995/1701304382, -48777925059/3047939560, 15336726248/1032824649, -45442868181/3398467696, 3065993473/597172653], 10:[185892177/718116043, 0, 0, -3185094517/667107341, -477755414/1098053517, -703635378/230739211, 5731566787/1027545527, 5232866602/850066563, -4093664535/808688257, 3962137247/1805957418, 65686358/487910083], 11:[403863854/491063109, 0, 0, -5068492393/434740067, -411421997/543043805, 652783627/914296604, 11173962825/925320556, -13158990841/6184727034, 3936647629/1978049680, -160528059/685178525, 248638103/1413531060, 0], 12:[14005451/335480064, 0, 0, 0, 0, -59238493/1068277825, 181606767/758867731, 561292985/797845732, -1041891430/1371343529, 760417239/1151165299, 118820643/751138087, -528747749/2220607170, 1/4]} #coefficients for lower order solution evaluation bh = [13451932/455176623, 0, 0, 0, 0, -808719846/976000145, 1757004468/5645159321, 656045339/265891186, -3867574721/1518517206, 465885868/322736535, 53011238/667516719, 2/45, 0] #coefficients for truncation error self.TR = [a-b for a, b in zip(self.BT[12], bh)]