qmat.solvers.generic.diffops ============================ .. py:module:: qmat.solvers.generic.diffops .. autoapi-nested-parse:: Contains various specialized implementation of :class:`DiffOp` classes. Attributes ---------- .. autoapisummary:: qmat.solvers.generic.diffops.T qmat.solvers.generic.diffops.DIFFOPS Classes ------- .. autoapisummary:: qmat.solvers.generic.diffops.Dahlquist qmat.solvers.generic.diffops.Lorenz qmat.solvers.generic.diffops.ProtheroRobinson Functions --------- .. autoapisummary:: qmat.solvers.generic.diffops.registerDiffOp Module Contents --------------- .. py:data:: T .. py:data:: DIFFOPS :type: dict[str, type[qmat.solvers.generic.DiffOp]] Dictionary containing all specialized :class:`DiffOp` classes .. py:function:: registerDiffOp(cls: type[T]) -> type[T] Class decorator to register a specialized :class:`DiffOp` class in `qmat` .. py:class:: Dahlquist(lam=1j) .. autoapi-inheritance-diagram:: qmat.solvers.generic.diffops.Dahlquist :parts: 1 Implements a Dahlquist differential operator .. math:: f(u,t) = \lambda u .. note:: This class is implemented for illustration and testing purposes. To solve with many :math:`\lambda` values, consider using the :class:`qmat.solvers.dahlquist.Dahlquist` class instead. :Parameters: **lam** (*complex, optional*) -- The :math:`\lambda` value. The default is 1j. .. py:attribute:: lam :value: 1j .. py:method:: evalF(u, t, out) Evaluate :math:`f(u,t)` and store the result into `out`. :Parameters: * **u** (*np.ndarray*) -- Input solution for the evaluation. * **t** (*float*) -- Time for the evaluation. * **out** (*np.ndarray*) -- Output array in which is stored the evaluation. .. py:class:: Lorenz(sigma=10, rho=28, beta=8 / 3, nativeFSolve=False) .. autoapi-inheritance-diagram:: qmat.solvers.generic.diffops.Lorenz :parts: 1 RHS of the Lorenz system, which can be written : .. math:: \frac{dx}{dt} = \sigma (y-x), \; \frac{dy}{dt} = x (\rho - z) - y, \; \frac{dz}{dt} = xy - \beta z, with starting initial solution :math:`u_0=(x_0,y_0,z_0)=(5, -5, 20)`. Considering the three dimensional vector :math:`u=(x,y,z)`, the formal expression of :math:`f` is then .. math:: f(u,t) = [ \sigma (y-x), x (\rho - z) - y, xy - \beta z ] :Parameters: * **sigma** (*float, optional*) -- The :math:`\sigma` parameter (default=10). * **rho** (*float, optional*) -- The :math:`\rho` parameter (default=28). * **beta** (*float, optional*) -- The :math:`\beta` parameter (default=8/3). * **nativeFSolve** (*bool, optional*) -- Wether or not using the native fSolve method (default is False). .. py:attribute:: params List containing :math:`\sigma`, :math:`\rho` and :math:`\beta` .. py:attribute:: newton Parameters for the Newton iteration used in native fSolve .. py:attribute:: gemv Level-2 blas gemv function used in the native solver (just for flex, very small speedup) .. py:method:: test() :classmethod: Class method to test the `DiffOp` implementation. :Parameters: * **t0** (*float, optional*) -- Evaluation time to test the instance. The default is 0. * **dt** (*float, optional*) -- Time-step to test the `fSolve` method. The default is 1e-1. * **eps** (*float, optional*) -- Perturbation added in the expected solution to test the `fSolve` method. The default is 1e-3. * **instance** (`DiffOp`, optional) -- Instance to be tested. If not provided (`None`), an instance is created using the default constructor. .. py:method:: evalF(u, t, out) Evaluate :math:`f(u,t)` and store the result into `out`. :Parameters: * **u** (*np.ndarray*) -- Input solution for the evaluation. * **t** (*float*) -- Time for the evaluation. * **out** (*np.ndarray*) -- Output array in which is stored the evaluation. .. py:method:: fSolve_NATIVE(a, rhs, t, out) Solve :math:`u-\alpha f(u,t)=rhs` for given :math:`u,t,rhs`, using a Newton iteration with exact Jacobian of :math:`f(u,t)`. :Parameters: * **a** (*float*) -- The :math:`\alpha` coefficient. * **rhs** (*np.ndarray*) -- The right hand side. * **t** (*float*) -- Time for the evaluation. * **out** (*np.ndarray*) -- Input-output array used as initial guess, in which is stored the solution. .. py:class:: ProtheroRobinson(epsilon=0.001, nonLinear=False, nativeFSolve=True) .. autoapi-inheritance-diagram:: qmat.solvers.generic.diffops.ProtheroRobinson :parts: 1 Implement the Prothero-Robinson problem: .. math:: \frac{du}{dt} = -\frac{u-g(t)}{\epsilon} + \frac{dg}{dt}, \quad u(0) = g(0), with :math:`\epsilon` a stiffness parameter, that makes the problem more stiff the smaller it is (usual taken value is :math:`\epsilon=1e^{-3}`). Exact solution is given by :math:`u(t)=g(t)`, and this implementation uses :math:`g(t)=\cos(t)`. Implement also the non-linear form of this problem: .. math:: \frac{du}{dt} = -\frac{u^3-g(t)^3}{\epsilon} + \frac{dg}{dt}, \quad u(0) = g(0). To use an other exact solution, one just have to derivate this class and overload the `g` and `dg` methods. For instance, to use :math:`g(t)=e^{-0.2*t}`, define and use the following class: >>> class MyProtheroRobinson(ProtheroRobinson): >>> >>> def g(self, t): >>> return np.exp(-0.2 * t) >>> >>> def dg(self, t): >>> return (-0.2) * np.exp(-0.2 * t) Reference --------- A. Prothero and A. Robinson, *On the stability and accuracy of one-step methods for solving stiff systems of ordinary differential equations*, Mathematics of Computation, **28** (1974), pp. 145–162. :Parameters: * **epsilon** (*float, optional*) -- Stiffness parameter. The default is 1e-3. * **nonLinear** (*bool, optional*) -- Wether or not to use the non-linear form of the problem. The default is False. * **nativeFSolve** (*bool, optional*) -- Wether or not use the native fSolver using exact Jacobian. The default is True. .. py:attribute:: epsilon :value: 0.001 Value used for :math:`\epsilon`. .. py:attribute:: newton Parameters used for the Newton iteration in `fSolve`. .. py:attribute:: evalF Evaluate :math:`f(u,t)` and store the result into `out`. :Parameters: * **u** (*np.ndarray*) -- Input solution for the evaluation. * **t** (*float*) -- Time for the evaluation. * **out** (*np.ndarray*) -- Output array in which is stored the evaluation. .. py:attribute:: jac .. py:method:: test() :classmethod: Test both linear and non-linear version of this differential operator. .. py:property:: nonLinear Wether the current operator is non-linear .. py:method:: g(t) .. py:method:: dg(t) .. py:method:: evalF_LIN(u, t, out) .. py:method:: evalF_NONLIN(u, t, out) .. py:method:: jac_LIN(u, t) .. py:method:: jac_NONLIN(u, t) .. py:method:: fSolve_NATIVE(a, rhs, t, out) Solve :math:`u-\alpha f(u,t)=rhs` for given :math:`u,t,rhs`, using a Newton iteration with exact Jacobian (derivative) of :math:`f(u,t)`. :Parameters: * **a** (*float*) -- The :math:`\alpha` coefficient. * **rhs** (*np.ndarray*) -- The right hand side. * **t** (*float*) -- Time for the evaluation. * **out** (*np.ndarray*) -- Input-output array used as initial guess, in which is stored the solution.