shenfun.utilities package

Submodules

shenfun.utilities.integrators module

Module for some integrators.

  • RK4: Runge-Kutta fourth order

  • ETD: Exponential time differencing Euler method

  • ETDRK4: Exponential time differencing Runge-Kutta fourth order

See, e.g., H. Montanelli and N. Bootland “Solving periodic semilinear PDEs in 1D, 2D and 3D with exponential integrators”, https://arxiv.org/pdf/1604.08900.pdf

Integrators are set up to solve equations like

\[\frac{\partial u}{\partial t} = L u + N(u)\]

where \(u\) is the solution, \(L\) is a linear operator and \(N(u)\) is the nonlinear part of the right hand side.

class shenfun.utilities.integrators.RK4(T, L=<function RK4.<lambda>>, N=<function RK4.<lambda>>, update=<function RK4.<lambda>>, **params)[source]

Bases: shenfun.utilities.integrators.IntegratorBase

Regular 4’th order Runge-Kutta integrator

Parameters
  • T (TensorProductSpace)

  • L (function) – To compute linear part of right hand side

  • N (function) – To compute nonlinear part of right hand side

  • update (function) – To be called at the end of a timestep

  • params (dictionary) – Any relevant keyword arguments

setup(dt)[source]

Set up RK4 ODE solver

solve(u, u_hat, dt, trange)[source]

Integrate forward in end_time

Parameters
  • u (array) – The solution array in physical space

  • u_hat (array) – The solution array in spectral space

  • dt (float) – Timestep

  • trange (two-tuple) – Time and end time

class shenfun.utilities.integrators.ETDRK4(T, L=<function ETDRK4.<lambda>>, N=<function ETDRK4.<lambda>>, update=<function ETDRK4.<lambda>>, **params)[source]

Bases: shenfun.utilities.integrators.IntegratorBase

Exponential time differencing Runge-Kutta 4’th order method

H. Montanelli and N. Bootland “Solving periodic semilinear PDEs in 1D, 2D and 3D with exponential integrators”, https://arxiv.org/pdf/1604.08900.pdf

Parameters
  • T (TensorProductSpace)

  • L (function) – To compute linear part of right hand side

  • N (function) – To compute nonlinear part of right hand side

  • update (function) – To be called at the end of a timestep

  • params (dictionary) – Any relevant keyword arguments

setup(dt)[source]

Set up ETDRK4 ODE solver

solve(u, u_hat, dt, trange)[source]

Integrate forward in time

Parameters
  • u (array) – The solution array in physical space

  • u_hat (array) – The solution array in spectral space

  • dt (float) – Timestep

  • trange (two-tuple) – Time and end time

class shenfun.utilities.integrators.ETD(T, L=<function ETD.<lambda>>, N=<function ETD.<lambda>>, update=<function ETD.<lambda>>, **params)[source]

Bases: shenfun.utilities.integrators.IntegratorBase

Exponential time differencing Euler method

H. Montanelli and N. Bootland “Solving periodic semilinear PDEs in 1D, 2D and 3D with exponential integrators”, https://arxiv.org/pdf/1604.08900.pdf

Parameters
  • T (TensorProductSpace)

  • L (function) – To compute linear part of right hand side

  • N (function) – To compute nonlinear part of right hand side

  • update (function) – To be called at the end of a timestep

  • params (dictionary) – Any relevant keyword arguments

setup(dt)[source]

Set up ETD ODE solver

solve(u, u_hat, dt, trange)[source]

Integrate forward in time

Parameters
  • u (array) – The solution array in physical space

  • u_hat (array) – The solution array in spectral space

  • dt (float) – Timestep

  • trange (two-tuple) – Time and end time

shenfun.utilities.lagrangian_particles module

class shenfun.utilities.lagrangian_particles.LagrangianParticles(points, dt, u_hat)[source]

Bases: object

Class for tracking Lagrangian particles

Parameters
  • points (array) – Initial location of particles. (D, N) array, with N particles in D dimensions

  • dt (float) – Time step

  • u_hat (Function) – Spectral Galerkin Function for the Eulerian velocity

rhs()[source]
step()[source]

Module contents

Module for implementing helper functions.

shenfun.utilities.inheritdocstrings(cls)[source]

Method used for inheriting docstrings from parent class

Use as decorator:

@inheritdocstrings
class Child(Parent):

and Child will use the same docstrings as parent even if a method is overloaded. The Child class may overload the docstring as well and a new docstring defined for a method in Child will overload the Parent.

shenfun.utilities.dx(u)[source]

Compute integral of u over domain

\[\int_{\Omega} u dx\]
Parameters

u (Array) – The Array to integrate

shenfun.utilities.clenshaw_curtis1D(u, quad='GC')[source]

Clenshaw-Curtis integration in 1D

class shenfun.utilities.CachedArrayDict[source]

Bases: collections.abc.MutableMapping

Dictionary for caching Numpy arrays (work arrays)

Example

>>> import numpy as np
>>> from shenfun.utilities import CachedArrayDict
>>> work = CachedArrayDict()
>>> a = np.ones((3, 4), dtype=int)
>>> w = work[(a, 0, True)] # create work array with shape as a
>>> print(w.shape)
(3, 4)
>>> print(w)
[[0 0 0 0]
 [0 0 0 0]
 [0 0 0 0]]
>>> w2 = work[(a, 1, True)] # Get different(note 1!) array of same shape/dtype
values() → an object providing a view on D's values[source]
shenfun.utilities.outer(a, b, c)[source]

Return outer product $c_{i,j} = a_i b_j$

Parameters
  • a (Array of shape (N, …))

  • b (Array of shape (N, …))

  • c (Array of shape (N*N, …))

  • The outer product is taken over the first index of a and b,

  • for all remaining indices.

shenfun.utilities.apply_mask(u_hat, mask)[source]