Skip to content

Commit

Permalink
Final cleanup Radau.
Browse files Browse the repository at this point in the history
  • Loading branch information
JonasBreuling committed Aug 23, 2024
1 parent bcb4ad9 commit 0716efc
Showing 1 changed file with 3 additions and 6 deletions.
9 changes: 3 additions & 6 deletions scipy_dae/integrate/_dae/radau.py
Original file line number Diff line number Diff line change
Expand Up @@ -220,7 +220,7 @@ def solve_collocation_system(fun, t, y, h, Z0, scale, tol,
return converged, k + 1, Y, Yp, Z, rate


def solve_collocation_system2(fun, t, y, h, Y0, Yp0, scale, tol,
def solve_collocation_system2(fun, t, y, h, Yp0, scale, tol,
LU_real, LU_complex, solve_lu,
C, T, TI, A, newton_max_iter):
s, n = Yp0.shape
Expand All @@ -229,8 +229,6 @@ def solve_collocation_system2(fun, t, y, h, Y0, Yp0, scale, tol,

Yp = Yp0
Y = y + h * A.dot(Yp)
# TODO: Why is this bad?
# Y = Y0
V_dot = TI.dot(Yp)

F = np.empty((s, n))
Expand Down Expand Up @@ -623,10 +621,9 @@ def _step_impl(self):
# note: this only works when we exrapolate the derivative
# of the collocation polynomial but do not use the sth order
# collocation polynomial for the derivatives as well.
# Yp0 = self.sol(t + h * C)[1].T
Yp0 = self.sol(t + h * C)[1].T
# # Z0 = self.sol(t + h * C)[0].T - y
# # Yp0 = (1 / h) * A_inv @ Z0
Y0, Yp0 = map(np.transpose, self.sol(t + h * C))
else:
Z0 = self.sol(t + h * C)[0].T - y
scale = atol + np.abs(y) * rtol
Expand All @@ -643,7 +640,7 @@ def _step_impl(self):

if UNKNOWN_VELOCITIES:
converged, n_iter, Y, Yp, Z, rate = solve_collocation_system2(
self.fun, t, y, h, Y0, Yp0, scale, newton_tol,
self.fun, t, y, h, Yp0, scale, newton_tol,
LU_real, LU_complex, self.solve_lu,
C, T, TI, A, newton_max_iter)
else:
Expand Down

0 comments on commit 0716efc

Please sign in to comment.