From 76efa6aaf9a6e022d831f86867ed5825085eca31 Mon Sep 17 00:00:00 2001 From: lisawim Date: Tue, 31 Oct 2023 08:44:52 +0100 Subject: [PATCH] Moved solve_system to generic ProblemDAE class --- pySDC/projects/DAE/misc/ProblemDAE.py | 61 +++++++++-- pySDC/projects/DAE/problems/simple_DAE.py | 103 +----------------- .../DAE/problems/synchronous_machine.py | 35 +----- .../DAE/problems/transistor_amplifier.py | 69 +----------- 4 files changed, 56 insertions(+), 212 deletions(-) diff --git a/pySDC/projects/DAE/misc/ProblemDAE.py b/pySDC/projects/DAE/misc/ProblemDAE.py index 97c64b9d87..91cc6bc148 100644 --- a/pySDC/projects/DAE/misc/ProblemDAE.py +++ b/pySDC/projects/DAE/misc/ProblemDAE.py @@ -1,25 +1,64 @@ import numpy as np +from scipy.optimize import root -from pySDC.core.Problem import ptype +from pySDC.core.Problem import ptype, WorkCounter from pySDC.implementations.datatype_classes.mesh import mesh class ptype_dae(ptype): - """ - Interface class for DAE problems. Ensures that all parameters are passed that are needed by DAE sweepers + r""" + This class implements a generic DAE class and illustrates the interface class for DAE problems. + It ensures that all parameters are passed that are needed by DAE sweepers. + + Parameters + ---------- + nvars : int + Number of unknowns of the problem class. + newton_tol : float + Tolerance for the nonlinear solver. + + Attributes + ---------- + work_counters : WorkCounter + Counts the work, here the number of function calls during the nonlinear solve is logged. """ dtype_u = mesh dtype_f = mesh def __init__(self, nvars, newton_tol): - """ - Initialization routine - - Args: - problem_params (dict): custom parameters for the example - dtype_u: mesh data type (will be passed parent class) - dtype_f: mesh data type (will be passed parent class) - """ + """Initialization routine""" super().__init__((nvars, None, np.dtype('float64'))) self._makeAttributeAndRegister('nvars', 'newton_tol', localVars=locals(), readOnly=True) + + self.work_counters['newton'] = WorkCounter() + + def solve_system(self, impl_sys, u0, t): + r""" + Solver for nonlinear implicit system (defined in sweeper). + + Parameters + ---------- + impl_sys : callable + Implicit system to be solved. + u0 : dtype_u + Initial guess for solver. + t : float + Current time :math:`t`. + + Returns + ------- + me : dtype_u + Numerical solution. + """ + + me = self.dtype_u(self.init) + opt = root( + impl_sys, + u0, + method='hybr', + tol=self.newton_tol, + ) + me[:] = opt.x + self.work_counters['newton'].niter += opt.nfev + return me diff --git a/pySDC/projects/DAE/problems/simple_DAE.py b/pySDC/projects/DAE/problems/simple_DAE.py index 859dbf7c5f..21ff17be26 100644 --- a/pySDC/projects/DAE/problems/simple_DAE.py +++ b/pySDC/projects/DAE/problems/simple_DAE.py @@ -1,7 +1,6 @@ import warnings import numpy as np from scipy.interpolate import interp1d -from scipy.optimize import root from pySDC.projects.DAE.misc.ProblemDAE import ptype_dae from pySDC.core.Problem import WorkCounter @@ -38,8 +37,7 @@ class pendulum_2d(ptype_dae): ---------- work_counters : WorkCounter Counts the work, i.e., number of function calls of right-hand side is called and stored in - ``work_counters['rhs']``, the number of function calls of the Newton-like solver is stored in - ``work_counters['newton']``. + ``work_counters['rhs']``. References ---------- @@ -57,7 +55,6 @@ def __init__(self, nvars, newton_tol): # solution = data[:, 1:] # self.u_ref = interp1d(t, solution, kind='cubic', axis=0, fill_value='extrapolate') self.t_end = 0.0 - self.work_counters['newton'] = WorkCounter() self.work_counters['rhs'] = WorkCounter() def eval_f(self, u, du, t): @@ -86,36 +83,6 @@ def eval_f(self, u, du, t): self.work_counters['rhs']() return f - def solve_system(self, impl_sys, u0, t): - r""" - Solver for nonlinear implicit system (defined in sweeper). - - Parameters - ---------- - impl_sys : callable - Implicit system to be solved. - u0 : dtype_u - Initial guess for solver. - t : float - Current time :math:`t`. - - Returns - ------- - me : dtype_u - Numerical solution. - """ - - me = self.dtype_u(self.init) - opt = root( - impl_sys, - u0, - method='hybr', - tol=self.newton_tol, - ) - me[:] = opt.x - self.work_counters['newton'].niter += opt.nfev - return me - def u_exact(self, t): """ Approximation of the exact solution generated by spline interpolation of an extremely accurate numerical reference solution. @@ -177,8 +144,7 @@ class simple_dae_1(ptype_dae): ---------- work_counters : WorkCounter Counts the work, i.e., number of function calls of right-hand side is called and stored in - ``work_counters['rhs']``, the number of function calls of the Newton-like solver is stored in - ``work_counters['newton']``. + ``work_counters['rhs']``. References ---------- @@ -190,7 +156,6 @@ def __init__(self, nvars, newton_tol): """Initialization routine""" super().__init__(nvars, newton_tol) - self.work_counters['newton'] = WorkCounter() self.work_counters['rhs'] = WorkCounter() def eval_f(self, u, du, t): @@ -222,36 +187,6 @@ def eval_f(self, u, du, t): self.work_counters['rhs']() return f - def solve_system(self, impl_sys, u0, t): - r""" - Solver for nonlinear implicit system (defined in sweeper). - - Parameters - ---------- - impl_sys : callable - Implicit system to be solved. - u0 : dtype_u - Initial guess for solver. - t : float - Current time :math:`t`. - - Returns - ------- - me : dtype_u - Numerical solution. - """ - - me = self.dtype_u(self.init) - opt = root( - impl_sys, - u0, - method='hybr', - tol=self.newton_tol, - ) - me[:] = opt.x - self.work_counters['newton'].niter += opt.nfev - return me - def u_exact(self, t): """ Routine for the exact solution. @@ -297,8 +232,7 @@ class problematic_f(ptype_dae): Specific parameter of the problem. work_counters : WorkCounter Counts the work, i.e., number of function calls of right-hand side is called and stored in - ``work_counters['rhs']``, the number of function calls of the Newton-like solver is stored in - ``work_counters['newton']``. + ``work_counters['rhs']`` References ---------- @@ -312,7 +246,6 @@ def __init__(self, nvars, newton_tol, eta=1): self._makeAttributeAndRegister('eta', localVars=locals()) self.work_counters['rhs'] = WorkCounter() - self.work_counters['newton'] = WorkCounter() def eval_f(self, u, du, t): r""" @@ -340,36 +273,6 @@ def eval_f(self, u, du, t): self.work_counters['rhs']() return f - def solve_system(self, impl_sys, u0, t): - r""" - Solver for nonlinear implicit system (defined in sweeper). - - Parameters - ---------- - impl_sys : callable - Implicit system to be solved. - u0 : dtype_u - Initial guess for solver. - t : float - Current time :math:`t`. - - Returns - ------- - me : dtype_u - Numerical solution. - """ - - me = self.dtype_u(self.init) - opt = root( - impl_sys, - u0, - method='hybr', - tol=self.newton_tol, - ) - me[:] = opt.x - self.work_counters['newton'].niter += opt.nfev - return me - def u_exact(self, t): """ Routine for the exact solution. diff --git a/pySDC/projects/DAE/problems/synchronous_machine.py b/pySDC/projects/DAE/problems/synchronous_machine.py index 1f99f5644a..a3309ab369 100644 --- a/pySDC/projects/DAE/problems/synchronous_machine.py +++ b/pySDC/projects/DAE/problems/synchronous_machine.py @@ -1,7 +1,6 @@ import numpy as np import warnings from scipy.interpolate import interp1d -from scipy.optimize import root from pySDC.projects.DAE.misc.ProblemDAE import ptype_dae from pySDC.core.Problem import WorkCounter @@ -147,8 +146,7 @@ class synchronous_machine_infinite_bus(ptype_dae): Defines the mechanical torque applied to the rotor shaft. work_counters : WorkCounter Counts the work, i.e., number of function calls of right-hand side is called and stored in - ``work_counters['rhs']``, the number of function calls of the Newton-like solver is stored in - ``work_counters['newton']``. + ``work_counters['rhs']``. References ---------- @@ -191,7 +189,6 @@ def __init__(self, nvars, newton_tol): self.T_m = 0.854 self.work_counters['rhs'] = WorkCounter() - self.work_counters['newton'] = WorkCounter() def eval_f(self, u, du, t): r""" @@ -269,36 +266,6 @@ def eval_f(self, u, du, t): self.work_counters['rhs']() return f - def solve_system(self, impl_sys, u0, t): - r""" - Solver for nonlinear implicit system (defined in sweeper). - - Parameters - ---------- - impl_sys : callable - Implicit system to be solved. - u0 : dtype_u - Initial guess for solver. - t : float - Current time :math:`t`. - - Returns - ------- - me : dtype_u - Numerical solution. - """ - - me = self.dtype_u(self.init) - opt = root( - impl_sys, - u0, - method='hybr', - tol=self.newton_tol, - ) - me[:] = opt.x - self.work_counters['newton'].niter += opt.nfev - return me - def u_exact(self, t): """ Approximation of the exact solution generated by spline interpolation of an extremely accurate numerical reference solution. diff --git a/pySDC/projects/DAE/problems/transistor_amplifier.py b/pySDC/projects/DAE/problems/transistor_amplifier.py index b160be8ef1..039678ebbe 100644 --- a/pySDC/projects/DAE/problems/transistor_amplifier.py +++ b/pySDC/projects/DAE/problems/transistor_amplifier.py @@ -1,6 +1,5 @@ import warnings import numpy as np -from scipy.optimize import root from scipy.interpolate import interp1d from pySDC.projects.DAE.misc.ProblemDAE import ptype_dae @@ -58,8 +57,7 @@ class one_transistor_amplifier(ptype_dae): The end time at which the reference solution is determined. work_counters : WorkCounter Counts the work, i.e., number of function calls of right-hand side is called and stored in - ``work_counters['rhs']``, the number of function calls of the Newton-like solver is stored in - ``work_counters['newton']``. + ``work_counters['rhs']``. References ---------- @@ -79,7 +77,6 @@ def __init__(self, nvars, newton_tol): self.t_end = 0.0 self.work_counters['rhs'] = WorkCounter() - self.work_counters['newton'] = WorkCounter() def eval_f(self, u, du, t): r""" @@ -116,36 +113,6 @@ def eval_f(self, u, du, t): self.work_counters['rhs']() return f - def solve_system(self, impl_sys, u0, t): - r""" - Solver for nonlinear implicit system (defined in sweeper). - - Parameters - ---------- - impl_sys : callable - Implicit system to be solved. - u0 : dtype_u - Initial guess for solver. - t : float - Current time :math:`t`. - - Returns - ------- - me : dtype_u - Numerical solution. - """ - - me = self.dtype_u(self.init) - opt = root( - impl_sys, - u0, - method='hybr', - tol=self.newton_tol, - ) - me[:] = opt.x - self.work_counters['newton'].niter += opt.nfev - return me - def u_exact(self, t): """ Approximation of the exact solution generated by spline interpolation of an extremely accurate numerical @@ -229,8 +196,7 @@ class two_transistor_amplifier(ptype_dae): The end time at which the reference solution is determined. work_counters : WorkCounter Counts the work, i.e., number of function calls of right-hand side is called and stored in - ``work_counters['rhs']``, the number of function calls of the Newton-like solver is stored in - ``work_counters['newton']``. + ``work_counters['rhs']``. References ---------- @@ -250,7 +216,6 @@ def __init__(self, nvars, newton_tol): self.t_end = 0.0 self.work_counters['rhs'] = WorkCounter() - self.work_counters['newton'] = WorkCounter() def eval_f(self, u, du, t): r""" @@ -290,36 +255,6 @@ def eval_f(self, u, du, t): self.work_counters['rhs']() return f - def solve_system(self, impl_sys, u0, t): - r""" - Solver for nonlinear implicit system (defined in sweeper). - - Parameters - ---------- - impl_sys : callable - Implicit system to be solved. - u0 : dtype_u - Initial guess for solver. - t : float - Current time :math:`t`. - - Returns - ------- - me : dtype_u - Numerical solution. - """ - - me = self.dtype_u(self.init) - opt = root( - impl_sys, - u0, - method='hybr', - tol=self.newton_tol, - ) - me[:] = opt.x - self.work_counters['newton'].niter += opt.nfev - return me - def u_exact(self, t): """ Dummy exact solution that should only be used to get initial conditions for the problem. This makes