irksome package

Submodules

irksome.ButcherTableaux module

class irksome.ButcherTableaux.Alexander

Bases: ButcherTableau

Third-order, diagonally implicit, 3-stage, L-stable scheme from Diagonally Implicit Runge-Kutta Methods for Stiff O.D.E.’s, R. Alexander, SINUM 14(6): 1006-1021, 1977.

class irksome.ButcherTableaux.BackwardEuler

Bases: RadauIIA

The rock-solid first-order implicit method.

class irksome.ButcherTableaux.ButcherTableau(A, b, btilde, c, order, embedded_order, gamma0)

Bases: object

Top-level class representing a Butcher tableau encoding

a Runge-Kutta method. It has members

Parameters:

• A – a 2d array containing the Butcher matrix

• b – a 1d array giving weights assigned to each stage when computing the solution at time n+1.

• btilde – If present, a 1d array giving weights for an embedded lower-order method (used in estimating temporal truncation error.)

• c – a 1d array containing weights at which time-dependent terms are evaluated.

• order – the (integer) formal order of accuracy of the method

• embedded_order – If present, the (integer) formal order of accuracy of the embedded method

• gamma0 – If present, the weight on the explicit term in the embedded lower-order method

property is_diagonally_implicit

property is_explicit

property is_fully_implicit

property is_implicit

property is_stiffly_accurate

Determines whether the method is stiffly accurate.

property num_stages

Return the number of stages the method has.

class irksome.ButcherTableaux.CollocationButcherTableau(L, order)

Bases: ButcherTableau

When an RK method is based on collocation with point sets present in FIAT, we have a general formula for producing the Butcher tableau.

Parameters:

• L – a one-dimensional class FIAT.FiniteElement of Lagrange type – the degrees of freedom must all be point evaluation.

• order – the order of the resulting RK method.

class irksome.ButcherTableaux.GaussLegendre(num_stages)

Bases: CollocationButcherTableau

Collocation method based on the Gauss-Legendre points. The order of accuracy is 2 * num_stages. GL methods are A-stable, B-stable, and symplectic.

Parameters:

num_stages – The number of stages (1 or greater)

class irksome.ButcherTableaux.LobattoIIIA(num_stages)

Bases: CollocationButcherTableau

Collocation method based on the Gauss-Lobatto points. The order of accuracy is 2 * num_stages - 2. LobattoIIIA methods are A-stable but not B- or L-stable.

Parameters:

num_stages – The number of stages (2 or greater)

class irksome.ButcherTableaux.LobattoIIIC(num_stages)

Bases: ButcherTableau

Discontinuous collocation method based on the Lobatto points. The order of accuracy is 2 * num_stages - 2. LobattoIIIC methods are A-, L-, algebraically, and B- stable.

Parameters:

num_stages – The number of stages (2 or greater)

class irksome.ButcherTableaux.PareschiRusso(x)

Bases: ButcherTableau

Second order, diagonally implicit, 2-stage. A-stable if x >= 1/4 and L-stable iff x = 1 plus/minus 1/sqrt(2).

class irksome.ButcherTableaux.QinZhang

Bases: PareschiRusso

Symplectic Pareschi-Russo DIRK

class irksome.ButcherTableaux.RadauIIA(num_stages, variant='embed_Radau5')

Bases: CollocationButcherTableau

Collocation method based on the Gauss-Radau points. The order of accuracy is 2 * num_stages - 1. RadauIIA methods are algebraically (hence B-) stable.

Parameters:

• num_stages – The number of stages (2 or greater)

• variant – Indicate whether to use the Radau5 style of embedded scheme (with variant = “embed_Radau5”) or the simple collocation type (with variant = “embed_colloc”)

irksome.bcs module

irksome.bcs.BCStageData(V, gcur, u0, u0_mult, i, t, dt)

class irksome.bcs.BoundsConstrainedBC(V, g, sub_domain, bounds, solver_parameters=None)

Bases: DirichletBC

A DirichletBC with bounds-constrained data.

property function_arg

The value of this boundary condition.

irksome.bcs.EmbeddedBCData(bc, t, dt, num_fields, butcher_tableau, ws, u0)

irksome.bcs.bc2space(bc, V)

irksome.bcs.get_sub(u, indices)

irksome.bcs.stage2spaces4bc(bc, V, Vbig, i)

used to figure out how to apply Dirichlet BC to each stage

irksome.deriv module

irksome.deriv.Dt(f)

Short-hand function to produce a TimeDerivative of the input.

class irksome.deriv.TimeDerivative(f)

Bases: Derivative

UFL node representing a time derivative of some quantity/field. Note: Currently form compilers do not understand how to process these nodes. Instead, Irksome pre-processes forms containing TimeDerivative nodes.

Initalise.

property ufl_free_indices

property ufl_index_dimensions

property ufl_shape

class irksome.deriv.TimeDerivativeRuleDispatcher(t, timedep_coeffs)

Bases: MultiFunction

Initialise.

coefficient_derivative(o)

coordinate_derivative(o)

derivative(o)

div(o)

expr(o, *ops)

Reuse object if operands are the same objects.

Use in your own subclass by setting e.g.

expr = MultiFunction.reuse_if_untouched

as a default rule.

grad(o)

reference_grad(o)

terminal(o)

time_derivative(o)

class irksome.deriv.TimeDerivativeRuleset(t, timedep_coeffs)

Bases: GenericDerivativeRuleset

Apply AD rules to time derivative expressions. WIP

Initialise.

coefficient(o)

irksome.deriv.apply_time_derivatives(expression, t, timedep_coeffs=[])

irksome.dirk_stepper module

class irksome.dirk_stepper.DIRKTimeStepper(F, butcher_tableau, t, dt, u0, bcs=None, solver_parameters=None, appctx=None, nullspace=None)

Bases: object

Front-end class for advancing a time-dependent PDE via a diagonally-implicit Runge-Kutta method formulated in terms of stage derivatives.

advance()

solver_stats()

update_bc_constants(i, c)

irksome.dirk_stepper.getFormDIRK(F, ks, butch, t, dt, u0, bcs=None)

irksome.explicit_stepper module

class irksome.explicit_stepper.ExplicitTimeStepper(F, butcher_tableau, t, dt, u0, bcs=None, solver_parameters=None, appctx=None)

Bases: DIRKTimeStepper

irksome.getForm module

irksome.getForm.getForm(F, butch, t, dt, u0, bcs=None, bc_type=None, splitting=<function AI>, nullspace=None)

Given a time-dependent variational form and a ButcherTableau, produce UFL for the s-stage RK method.

Parameters:

• F – UFL form for the semidiscrete ODE/DAE

• butch – the ButcherTableau for the RK method being used to advance in time.

• t – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time.

• dt – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time step. The user may adjust this value between time steps.

• splitting – a callable that maps the (floating point) Butcher matrix a to a pair of matrices A1, A2 such that butch.A = A1 A2. This is used to vary between the classical RK formulation and Butcher’s reformulation that leads to a denser mass matrix with block-diagonal stiffness. Some choices of function will assume that butch.A is invertible.

• u0 – a Function referring to the state of the PDE system at time t

• bcs – optionally, a DirichletBC object (or iterable thereof) containing (possible time-dependent) boundary conditions imposed on the system.

• bc_type – How to manipulate the strongly-enforced boundary conditions to derive the stage boundary conditions. Should be a string, either “DAE”, which implements BCs as constraints in the style of a differential-algebraic equation, or “ODE”, which takes the time derivative of the boundary data and evaluates this for the stage values

• nullspace – A list of tuples of the form (index, VSB) where index is an index into the function space associated with u and VSB is a :class: firedrake.VectorSpaceBasis instance to be passed to a firedrake.MixedVectorSpaceBasis over the larger space associated with the Runge-Kutta method

On output, we return a tuple consisting of four parts:

• Fnew, the Form

• k, the firedrake.Function holding all the stages. It lives in a firedrake.FunctionSpace corresponding to the s-way tensor product of the space on which the semidiscrete form lives.

• bcnew, a list of firedrake.DirichletBC objects to be posed on the stages,

• ‘nspnew’, the firedrake.MixedVectorSpaceBasis object that represents the nullspace of the coupled system

irksome.imex module

class irksome.imex.RadauIIAIMEXMethod(F, Fexp, butcher_tableau, t, dt, u0, bcs=None, it_solver_parameters=None, prop_solver_parameters=None, splitting=<function AI>, appctx=None, nullspace=None, num_its_initial=0, num_its_per_step=0)

Bases: object

Class for advancing a time-dependent PDE via a polynomial IMEX/RadauIIA method. This requires one to split the PDE into an implicit and explicit part. The class sets up two methods – advance and iterate. The former is used to move the solution forward in time, while the latter is used both to start the method (filling up the initial stage values) and can be used at each time step to increase the accuracy/stability. In the limit as the iterator is applied many times per time step, one expects convergence to the solution that would have been obtained from fully-implicit RadauIIA method.

Parameters:

• F – A ufl.Form instance describing the implicit part of the semi-discrete problem F(t, u; v) == 0, where u is the unknown firedrake.Function and `v is the :class:firedrake.TestFunction`.

• Fexp – A ufl.Form instance describing the part of the PDE that is explicitly split off.

• butcher_tableau – A ButcherTableau instance giving the Runge-Kutta method to be used for time marching. Only RadauIIA is allowed here (but it can be any number of stages).

• t – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time.

• dt – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time step. The user may adjust this value between time steps.

• u0 – A firedrake.Function containing the current state of the problem to be solved.

• bcs – An iterable of firedrake.DirichletBC containing the strongly-enforced boundary conditions. Irksome will manipulate these to obtain boundary conditions for each stage of the RK method.

• it_solver_parameters – A dict of solver parameters that will be used in solving the algebraic problem associated with the iterator.

• prop_solver_parameters – A dict of solver parameters that will be used in solving the algebraic problem associated with the propagator.

• splitting – A callable used to factor the Butcher matrix, currently, only AI is supported.

• appctx – An optional dict containing application context.

• nullspace – An optional null space object.

advance()

iterate()

Called 1 or more times to set up the initial state of the system before time-stepping. Can also be called after each call to advance

propagate()

Moves the solution forward in time, to be followed by 0 or more calls to iterate.

solver_stats()

irksome.imex.getFormExplicit(Fexp, butch, u0, UU, t, dt, splitting=None)

Processes the explicitly split-off part for a RadauIIA-IMEX method. Returns the forms for both the iterator and propagator, which really just differ by which constants are in them.

irksome.imex.riia_explicit_coeffs(k)

Computes the coefficients needed for the explicit part of a RadauIIA-IMEX method.

irksome.manipulation module

Manipulation of expressions containing TimeDerivative terms.

These can be used to do some basic checking of the suitability of a Form for use in Irksome (via check_integrals), and splitting out terms in the Form that contain a time derivative from those that don’t (via extract_terms).

class irksome.manipulation.SplitTimeForm(time: Form, remainder: Form)

Bases: NamedTuple

A container for a form split into time terms and a remainder.

Create new instance of SplitTimeForm(time, remainder)

remainder: Form

Alias for field number 1

time: Form

Alias for field number 0

irksome.manipulation.check_integrals(integrals: List[Integral], expect_time_derivative: bool = True) → List[Integral]

Check a list of integrals for linearity in the time derivative.

Parameters:

• integrals – list of integrals.

• expect_time_derivative – Are we expecting to see a time derivative?

Raises:

ValueError – if we are expecting a time derivative and don’t see one, or time derivatives are applied nonlinearly, to more than one coefficient, or more than first order.

irksome.manipulation.extract_terms(form: Form) → SplitTimeForm

Extract terms from a Form.

This splits a form (a sum of integrals) into those integrals which do contain a TimeDerivative and those that don’t.

Parameters:

form – The form to split.

Returns:

a SplitTimeForm tuple.

Raises:

ValueError – if the form does not apply anything other than first-order time derivatives to a single coefficient.

irksome.pc module

class irksome.pc.IRKAuxiliaryOperatorPC

Bases: AuxiliaryOperatorPC

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

• initialize

• update

• apply

• applyTranspose

form(pc, test, trial)

Implements the interface for AuxiliaryOperatorPC.

abstract getNewForm(pc, u0, test)

class irksome.pc.RanaBase

Bases: AuxiliaryOperatorPC

Base class for methods out of Rana, Howle, Long, Meek, & Milestone. It inherits from Firedrake’s AuxiliaryOperatorPC class and provides the preconditioning bilinear form associated with an approximation to the Butcher matrix (which is provided by subclasses).

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

• initialize

• update

• apply

• applyTranspose

form(pc, test, trial)

Implements the interface for AuxiliaryOperatorPC.

abstract getAtilde(A)

Derived classes produce a typically structured approximation to A.

class irksome.pc.RanaDU

Bases: RanaBase

Implements Rana-type preconditioner using Atilde = DU where A=LDU.

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

• initialize

• update

• apply

• applyTranspose

getAtilde(A)

Derived classes produce a typically structured approximation to A.

class irksome.pc.RanaLD

Bases: RanaBase

Implements Rana-type preconditioner using Atilde = LD where A=LDU.

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

• initialize

• update

• apply

• applyTranspose

getAtilde(A)

Derived classes produce a typically structured approximation to A.

irksome.pc.ldu(A)

irksome.pep_explicit_rk module

class irksome.pep_explicit_rk.PEPRK(ns, order, peporder)

Bases: ButcherTableau

irksome.stage module

class irksome.stage.StageValueTimeStepper(F, butcher_tableau, t, dt, u0, bcs=None, solver_parameters=None, update_solver_parameters=None, bc_constraints=None, splitting=<function AI>, basis_type=None, nullspace=None, appctx=None)

Bases: object

advance(bounds=None)

solver_stats()

irksome.stage.getBits(num_stages, num_fields, u0, UU, v, VV)

irksome.stage.getFormStage(F, butch, u0, t, dt, bcs=None, splitting=None, vandermonde=None, bc_constraints=None, nullspace=None)

Given a time-dependent variational form and a ButcherTableau, produce UFL for the s-stage RK method.

Parameters:

• F – UFL form for the semidiscrete ODE/DAE

• butch – the ButcherTableau for the RK method being used to advance in time.

• u0 – a Function referring to the state of the PDE system at time t

• t – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time.

• dt – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time step. The user may adjust this value between time steps.

• splitting – a callable that maps the (floating point) Butcher matrix a to a pair of matrices A1, A2 such that butch.A = A1 A2. This is used to vary between the classical RK formulation and Butcher’s reformulation that leads to a denser mass matrix with block-diagonal stiffness. Only AI and IA are currently supported.

• vandermonde – a numpy array encoding a change of basis to the Lagrange polynomials associated with the collocation nodes from some other (e.g. Bernstein or Chebyshev) basis. This allows us to solve for the coefficients in some basis rather than the values at particular stages, which can be useful for satisfying bounds constraints. If none is provided, we assume it is the identity, working in the Lagrange basis.

• bcs – optionally, a DirichletBC object (or iterable thereof) containing (possible time-dependent) boundary conditions imposed on the system.

• bc_constraints – optionally, a dictionary mapping (some of) the boundary conditions in bcs to triples of the form (params, lower, upper) indicating solver parameters to use and lower and upper bounds to provide in doing a bounds-constrained projection of the boundary data. Note: if these bounds change over time, the user is responsible for maintaining a handle on them and updating them between time steps.

• nullspace – A list of tuples of the form (index, VSB) where index is an index into the function space associated with u and VSB is a :class: firedrake.VectorSpaceBasis instance to be passed to a firedrake.MixedVectorSpaceBasis over the larger space associated with the Runge-Kutta method

On output, we return a tuple consisting of several parts:

• Fnew, the Form

• possibly a 4-tuple containing information needed to solve a mass matrix to update the solution (this is empty for RadauIIA methods for which there is a trivial update function.

• UU, the firedrake.Function holding all the stage time values. It lives in a firedrake.FunctionSpace corresponding to the s-way tensor product of the space on which the semidiscrete form lives.

• bcnew, a list of firedrake.DirichletBC objects to be posed on the stages,

• ‘nspnew’, the firedrake.MixedVectorSpaceBasis object that represents the nullspace of the coupled system

irksome.stage.isiterable(x)

irksome.stage.split_field(num_fields, u)

irksome.stage.split_stage_field(num_stages, num_fields, UU)

irksome.stepper module

class irksome.stepper.AdaptiveTimeStepper(F, butcher_tableau, t, dt, u0, bcs=None, appctx=None, solver_parameters=None, bc_type='DAE', splitting=<function AI>, nullspace=None, tol=0.001, dtmin=1e-15, dtmax=1.0, KI=0.06666666666666667, KP=0.13, max_reject=10, onscale_factor=1.2, safety_factor=0.9, gamma0_params=None)

Bases: StageDerivativeTimeStepper

Front-end class for advancing a time-dependent PDE via an adaptive Runge-Kutta method.

Parameters:

• F – A ufl.Form instance describing the semi-discrete problem F(t, u; v) == 0, where u is the unknown firedrake.Function and `v is the :class:firedrake.TestFunction`.

• butcher_tableau – A ButcherTableau instance giving the Runge-Kutta method to be used for time marching.

• t – A firedrake.Constant instance that always contains the time value at the beginning of a time step

• dt – A firedrake.Constant containing the size of the current time step. The user may adjust this value between time steps; however, note that the adaptive time step controls may adjust this before the step is taken.

• u0 – A firedrake.Function containing the current state of the problem to be solved.

• tol – The temporal truncation error tolerance

• dtmin – Minimal acceptable time step. An exception is raised if the step size drops below this threshhold.

• dtmax – Maximal acceptable time step, imposed as a hard cap; this can be adjusted externally once the time-stepper is instantiated, by modifying stepper.dt_max

• KI – Integration gain for step-size controller. Should be less than 1/p, where p is the expected order of the scheme. Larger values lead to faster (attempted) increases in time-step size when steps are accepted. See Gustafsson, Lundh, and Soderlind, BIT 1988.

• KP – Proportional gain for step-size controller. Controls dependence on ratio of (error estimate)/(step size) in determining new time-step size when steps are accepted. See Gustafsson, Lundh, and Soderlind, BIT 1988.

• max_reject – Maximum number of rejected timesteps in a row that does not lead to a failure

• onscale_factor – Allowable tolerance in determining initial timestep to be “on scale”

• safety_factor – Safety factor used when shrinking timestep if a proposed step is rejected

• gamma0_params – Solver parameters for mass matrix solve when using an embedded scheme with explicit first stage

• bcs – An iterable of firedrake.DirichletBC containing the strongly-enforced boundary conditions. Irksome will manipulate these to obtain boundary conditions for each stage of the RK method.

• solver_parameters – A dict of solver parameters that will be used in solving the algebraic problem associated with each time step.

• nullspace – A list of tuples of the form (index, VSB) where index is an index into the function space associated with u and VSB is a :class: firedrake.VectorSpaceBasis instance to be passed to a firedrake.MixedVectorSpaceBasis over the larger space associated with the Runge-Kutta method

advance()

Attempts to advances the system from time t to time t + dt. If the error threshhold is exceeded, will adaptively decrease the time step until the step is accepted. Also predicts new time step once the step is accepted. Note: overwrites the value u0.

class irksome.stepper.StageDerivativeTimeStepper(F, butcher_tableau, t, dt, u0, bcs=None, solver_parameters=None, splitting=<function AI>, appctx=None, nullspace=None, bc_type='DAE')

Bases: object

Front-end class for advancing a time-dependent PDE via a Runge-Kutta method formulated in terms of stage derivatives.

Parameters:

• F – A ufl.Form instance describing the semi-discrete problem F(t, u; v) == 0, where u is the unknown firedrake.Function and `v is the :class:firedrake.TestFunction`.

• butcher_tableau – A ButcherTableau instance giving the Runge-Kutta method to be used for time marching.

• t – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time.

• dt – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time step. The user may adjust this value between time steps.

• u0 – A firedrake.Function containing the current state of the problem to be solved.

• bcs – An iterable of firedrake.DirichletBC containing the strongly-enforced boundary conditions. Irksome will manipulate these to obtain boundary conditions for each stage of the RK method.

• bc_type – How to manipulate the strongly-enforced boundary conditions to derive the stage boundary conditions. Should be a string, either “DAE”, which implements BCs as constraints in the style of a differential-algebraic equation, or “ODE”, which takes the time derivative of the boundary data and evaluates this for the stage values

• solver_parameters – A dict of solver parameters that will be used in solving the algebraic problem associated with each time step.

• splitting – An callable used to factor the Butcher matrix

• appctx – An optional dict containing application context. This gets included with particular things that Irksome will pass into the nonlinear solver so that, say, user-defined preconditioners have access to it.

• nullspace – A list of tuples of the form (index, VSB) where index is an index into the function space associated with u and VSB is a :class: firedrake.VectorSpaceBasis instance to be passed to a firedrake.MixedVectorSpaceBasis over the larger space associated with the Runge-Kutta method

advance()

Advances the system from time t to time t + dt. Note: overwrites the value u0.

solver_stats()

irksome.stepper.TimeStepper(F, butcher_tableau, t, dt, u0, **kwargs)

Helper function to dispatch between various back-end classes

for doing time stepping. Returns an instance of the appropriate class.

Parameters:

• F – A ufl.Form instance describing the semi-discrete problem F(t, u; v) == 0, where u is the unknown firedrake.Function and `v iss the :class:firedrake.TestFunction`.

• butcher_tableau – A ButcherTableau instance giving the Runge-Kutta method to be used for time marching.

• t – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time.

• dt – a Function on the Real space over the same mesh as u0. This serves as a variable referring to the current time step. The user may adjust this value between time steps.

• u0 – A firedrake.Function containing the current state of the problem to be solved.

• bcs – An iterable of firedrake.DirichletBC containing the strongly-enforced boundary conditions. Irksome will manipulate these to obtain boundary conditions for each stage of the RK method.

• nullspace – A list of tuples of the form (index, VSB) where index is an index into the function space associated with u and VSB is a :class: firedrake.VectorSpaceBasis instance to be passed to a firedrake.MixedVectorSpaceBasis over the larger space associated with the Runge-Kutta method

• stage_type – Whether to formulate in terms of a stage derivatives or stage values.

• splitting – An callable used to factor the Butcher matrix

• bc_type – For stage derivative formulation, how to manipulate the strongly-enforced boundary conditions.

• solver_parameters – A dict of solver parameters that will be used in solving the algebraic problem associated with each time step.

• update_solver_parameters – A dict of parameters for inverting the mass matrix at each step (only used if stage_type is “value”)

• adaptive_parameters – A dict of parameters for use with adaptive time stepping (only used if stage_type is “deriv”)

irksome.tools module

irksome.tools.AI(A)

irksome.tools.ConstantOrZero(x, MC)

irksome.tools.IA(A)

class irksome.tools.MeshConstant(msh)

Bases: object

Constant(val=0.0)

class irksome.tools.MyReplacer(mapping)

Bases: MultiFunction

Initialise.

expr(o)

irksome.tools.getNullspace(V, Vbig, butch, nullspace)

Computes the nullspace for a multi-stage method.

Parameters:

• V – The FunctionSpace on which the original time-dependent PDE is posed.

• Vbig – The multi-stage FunctionSpace for the stage problem

• butch – The ButcherTableau defining the RK method

• nullspace – The nullspace for the original problem.

On output, we produce a MixedVectorSpaceBasis defining the nullspace for the multistage problem.

irksome.tools.is_ode(f, u)

Given a form defined over a function u, checks if (each bit of) u appears under a time derivative.

irksome.tools.replace(e, mapping)

Replace subexpressions in expression.

@param e:

An Expr or Form.

@param mapping:

A dict with from:to replacements to perform.

irksome.wso_dirk_tableaux module

class irksome.wso_dirk_tableaux.WSODIRK(ns, order, ws_order)

Bases: ButcherTableau

Module contents