Waves

Waves: vibrating strings, chains of springs, and quantum wave packets, free and charged.

ChargedWaveFunctionElement(label, visual=None, *, u=0.0, w=0.0, D=0.5, q=1.0, A_x=0.0, A_t=0.0, E=0.0, duration=1.0)

Bases: Element

Wave function psi = u + i w of a charge q in a uniform electromagnetic potential.

The wave function follows (d_t - i q A_t) psi = i D (d_x - i q A_x)^2 psi. The vector potential A_x(t) = A_x - E min(t, duration) holds a uniform electric field E from time 0 to duration, and the scalar potential A_t is constant. Expanded, the equation is psi_t = i D psi_xx + 2 D q A_x psi_x + i (q A_t - D q^2 A_x^2) psi.

Parameters:

Name	Type	Description	Default
label	str	The label of the site or the field.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
u	Annotated[float, StateSpec()]	The initial real part.	0.0
w	Annotated[float, StateSpec()]	The initial imaginary part.	0.0
D	Annotated[float, ParameterSpec(domain=open(0, oo))]	The coefficient of the Laplacian, ħ / (2 m).	0.5
q	Annotated[float, ParameterSpec()]	The electric charge.	1.0
A_x	Annotated[float, ParameterSpec()]	The vector potential at time 0.	0.0
A_t	Annotated[float, ParameterSpec()]	The scalar potential.	0.0
E	Annotated[float, ParameterSpec()]	The electric field.	0.0
duration	Annotated[float, ParameterSpec(domain=open(0, oo))]	The time the electric field is on.	1.0

Examples:

>>> from skcomplex.models.physics.waves import ChargedWaveFunctionElement
>>> ChargedWaveFunctionElement('site').q
1.0

Source code in src/skcomplex/models/physics/waves/_charged.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    u: Annotated[float, StateSpec()] = 0.0,
    w: Annotated[float, StateSpec()] = 0.0,
    D: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.5,
    q: Annotated[float, ParameterSpec()] = 1.0,
    A_x: Annotated[float, ParameterSpec()] = 0.0,
    A_t: Annotated[float, ParameterSpec()] = 0.0,
    E: Annotated[float, ParameterSpec()] = 0.0,
    duration: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.u = u
    self.w = w
    self.D = D
    self.q = q
    self.A_x = A_x
    self.A_t = A_t
    self.E = E
    self.duration = duration

DrivenAnchorInteraction(label, *, target, k=1.0, center=0.0, amplitude=0.0, omega=1.0)

Bases: HarmonicPotentialInteraction

Spring that ties a particle to an anchor point, which oscillates transversally as amplitude sin(omega t).

The energy is k |x - center - (0, amplitude sin(omega t))|² / 2, and a zero amplitude gives a fixed anchor.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
target	Annotated[str, RoleSpec()]	The path of the particle.	required
k	Annotated[float, ParameterSpec(domain=open(0, oo))]	The stiffness of the spring.	1.0
center	Annotated[ArrayLike, ParameterSpec(vector=True)]	The rest position of the anchor.	0.0
amplitude	Annotated[float, ParameterSpec()]	The amplitude of the transverse motion of the anchor.	0.0
omega	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The angular frequency of the anchor.	1.0

Examples:

>>> from skcomplex.models.physics.waves import DrivenAnchorInteraction
>>> DrivenAnchorInteraction('end', target='mass_0', center=[0.0, 0.0], amplitude=0.1).amplitude
0.1

Source code in src/skcomplex/models/physics/waves/_driven_anchor.py

def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    k: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    center: Annotated[npt.ArrayLike, ParameterSpec(vector=True)] = 0.0,
    amplitude: Annotated[float, ParameterSpec()] = 0.0,
    omega: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, target=target, k=k, center=center)
    self.amplitude = amplitude
    self.omega = omega

LatticeHamiltonianInteraction(label, *, target, D=0.5)

Bases: Interaction

Hamiltonian -D laplacian + V of a particle on the sites of a lattice.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
target	Annotated[str, RoleSpec()]	The path of the particle, whose Hilbert space holds the operators laplacian and potential.	required
D	Annotated[float, ParameterSpec(domain=open(0, oo))]	The coefficient of the Laplacian, ħ / (2 m).	0.5

Examples:

>>> from skcomplex.models.physics.waves import LatticeHamiltonianInteraction
>>> LatticeHamiltonianInteraction('kinetic', target='packet').D
0.5

Source code in src/skcomplex/models/physics/waves/_lattice_hamiltonian.py

def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    D: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.5,
) -> None:
    super().__init__(label)
    self.target = target
    self.D = D

PositionDephasingInteraction(label, *, target, gamma=0.1, n_sites=1)

Bases: Interaction

Dephasing of a particle in the position basis, one dissipator per site projector at the rate gamma.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
target	Annotated[str, RoleSpec()]	The path of the particle, whose Hilbert space holds the site projectors site_<k>.	required
gamma	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The dephasing rate.	0.1
n_sites	int	The number of sites.	1

Examples:

>>> from skcomplex.models.physics.waves import PositionDephasingInteraction
>>> PositionDephasingInteraction('dephasing', target='packet', gamma=0.5).gamma
0.5

Source code in src/skcomplex/models/physics/waves/_position_dephasing.py

def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    n_sites: int = 1,
) -> None:
    super().__init__(label)
    self.target = target
    self.gamma = gamma
    self.n_sites = n_sites

WaveFunctionElement(label, visual=None, *, u=0.0, w=0.0, D=0.5, V=0.0)

Bases: Element

Site of the wave function psi = u + i w of a particle with the potential V, by its real and imaginary parts.

The parts follow the Schrödinger equation i psi' = -D laplacian(psi) + V psi, so u' = -D laplacian(w) + V w and w' = D laplacian(u) - V u.

Parameters:

Name	Type	Description	Default
label	str	The label of the site.	required
visual	VisualSpec | None	How the site is drawn, or None for none.	None
u	Annotated[float, StateSpec()]	The initial real part.	0.0
w	Annotated[float, StateSpec()]	The initial imaginary part.	0.0
D	Annotated[float, ParameterSpec(domain=open(0, oo))]	The coefficient of the kinetic energy, hbar / (2 m).	0.5
V	Annotated[float, ParameterSpec()]	The potential at the site.	0.0

Examples:

>>> from skcomplex.models.physics.waves import WaveFunctionElement
>>> WaveFunctionElement('site', V=1.0).V
1.0

Source code in src/skcomplex/models/physics/waves/_wave_function.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    u: Annotated[float, StateSpec()] = 0.0,
    w: Annotated[float, StateSpec()] = 0.0,
    D: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.5,
    V: Annotated[float, ParameterSpec()] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.u = u
    self.w = w
    self.D = D
    self.V = V

WaveStringElement(label, visual=None, *, u=0.0, u_d1=0.0, c=1.0)

Bases: Element

Site of a string, with the transverse displacement u of order 2 and the wave speed c.

Parameters:

Name	Type	Description	Default
label	str	The label of the site.	required
visual	VisualSpec | None	How the site is drawn, or None for none.	None
u	Annotated[float, StateSpec(order=2)]	The initial displacement.	0.0
u_d1	float	The initial transverse velocity.	0.0
c	Annotated[float, ParameterSpec(domain=open(0, oo))]	The wave speed.	1.0

Examples:

>>> from skcomplex.models.physics.waves import WaveStringElement
>>> WaveStringElement('site', c=2.0).c
2.0

Source code in src/skcomplex/models/physics/waves/_wave_string.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    u: Annotated[float, StateSpec(order=2)] = 0.0,
    u_d1: float = 0.0,
    c: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.u = u
    self.u_d1 = u_d1
    self.c = c

charged_schrodinger(n_sites=256, *, length=80.0, D=0.5, q=1.0, momentum=0.0, width=2.0, A_x=0.0, A_t=0.0, E=1.0, duration=2.0)

Return a Gaussian wave packet of a charge in a prescribed electromagnetic potential, after Schwichtenberg 4.3.

The wave function follows the Schrödinger equation with the covariant derivatives d - i q A, (d_t - i q A_t) psi = i D (d_x - i q A_x)^2 psi, on a periodic line with spectral derivatives. A mode of wave number k turns at the frequency D (k - q A_x)^2 - q A_t, so the packet moves at the group velocity 2 D (k - q A_x). A constant A_x is a pure gauge that only shifts the canonical momentum, while the electric field E switched on for duration changes A_x and accelerates the packet at the rate 2 D q E.

Parameters:

Name	Type	Description	Default
n_sites	int	The number of sites of the periodic lattice.	256
length	float	The length of the periodic line, centered on 0.	80.0
D	float	The coefficient of the Laplacian, ħ / (2 m).	0.5
q	float	The electric charge.	1.0
momentum	float	The mean wave number k of the packet.	0.0
width	float	The width of the packet.	2.0
A_x	float	The vector potential at time 0.	0.0
A_t	float	The scalar potential.	0.0
E	float	The electric field.	1.0
duration	float	The time the electric field is on.	2.0

Returns:

Name	Type	Description
system	System	The model, with the real part packet.u and the imaginary part packet.w of the wave function.

Examples:

>>> from skcomplex.models.physics.waves import charged_schrodinger
>>> charged_schrodinger(64, E=0.5).elements[0].values['E']
0.5

Source code in src/skcomplex/models/physics/waves/_charged.py

@register_model(
    field='physics/waves',
    tags=('quantum', 'gauge field', 'electromagnetism', 'wave packet', 'spectral', 'pde'),
    references=('Schwichtenberg 2020',),
    defaults={'until': 8.0, 'dt': 0.1},
)
def charged_schrodinger(
    n_sites: int = 256,
    *,
    length: float = 80.0,
    D: float = 0.5,
    q: float = 1.0,
    momentum: float = 0.0,
    width: float = 2.0,
    A_x: float = 0.0,
    A_t: float = 0.0,
    E: float = 1.0,
    duration: float = 2.0,
) -> System:
    """Return a Gaussian wave packet of a charge in a prescribed electromagnetic potential, after Schwichtenberg 4.3.

    The wave function follows the Schrödinger equation with the covariant derivatives `d - i q A`,
    `(d_t - i q A_t) psi = i D (d_x - i q A_x)^2 psi`, on a periodic line with spectral derivatives. A mode of wave
    number `k` turns at the frequency `D (k - q A_x)^2 - q A_t`, so the packet moves at the group velocity
    `2 D (k - q A_x)`. A constant `A_x` is a pure gauge that only shifts the canonical momentum, while the electric
    field `E` switched on for `duration` changes `A_x` and accelerates the packet at the rate `2 D q E`.

    Args:
        n_sites:
            The number of sites of the periodic lattice.

        length:
            The length of the periodic line, centered on 0.

        D:
            The coefficient of the Laplacian, `ħ / (2 m)`.

        q:
            The electric charge.

        momentum:
            The mean wave number `k` of the packet.

        width:
            The width of the packet.

        A_x:
            The vector potential at time 0.

        A_t:
            The scalar potential.

        E:
            The electric field.

        duration:
            The time the electric field is on.

    Returns:
        system:
            The model, with the real part `packet.u` and the imaginary part `packet.w` of the wave function.

    Examples:
        >>> from skcomplex.models.physics.waves import charged_schrodinger
        >>> charged_schrodinger(64, E=0.5).elements[0].values['E']
        0.5
    """
    spacing = length / n_sites
    positions = spacing * np.arange(n_sites) - length / 2
    envelope = np.exp(-(positions**2) / (4 * width**2)) / (2 * np.pi * width**2) ** 0.25
    psi = envelope * np.exp(1j * momentum * positions)
    values = {
        'u': psi.real,
        'w': psi.imag,
        'D': D,
        'q': q,
        'A_x': A_x,
        'A_t': A_t,
        'E': E,
        'duration': duration,
    }
    space = LatticeSpace(n_sites, spacing=spacing)
    return System(
        elements=[
            FieldElement(
                'packet',
                element=ChargedWaveFunctionElement,
                space=space,
                values=values,
                discretization='spectral',
            ),
        ],
    )

schrodinger(n_sites=256, *, length=40.0, D=0.5, momentum=2.0, width=1.0, omega=0.0, formulation='fields', gamma=0.0)

Return a Gaussian wave packet that follows the Schrödinger equation, as its real and imaginary parts.

The wave function psi = u + i w follows psi_t = i D laplacian(psi) - i V psi, written as u_t = -D laplacian(w) + V w and w_t = D laplacian(u) - V u, on a periodic lattice with spectral derivatives. In units with ħ = 1, D = 1 / (2 m). The potential is the harmonic well V = omega² x² / (4 D), or zero for a free packet. A free packet of momentum k moves at the group velocity 2 D k and spreads, and its norm sum(u² + w²) is conserved.

The operator formulation builds the same particle as one quantum subsystem whose basis states are the sites, with the Hamiltonian -D laplacian + V of the same spectral Laplacian, normalised so that |psi|² is the probability of every site. It runs with the quantum solvers, and it accepts position dephasing at the rate gamma, which turns the state into a density matrix.

Parameters:

Name	Type	Description	Default
n_sites	int	The number of sites of the periodic lattice.	256
length	float	The length of the periodic box, centered on 0.	40.0
D	float	The coefficient of the Laplacian, ħ / (2 m).	0.5
momentum	float	The mean wave number of the packet.	2.0
width	float	The width of the packet.	1.0
omega	float	The angular frequency of the harmonic well.	0.0
formulation	str	'fields' for two real fields, or 'operators' for a particle on the sites of the lattice.	'fields'
gamma	float	The rate of position dephasing of the operator formulation.	0.0

Returns:

Name	Type	Description
system	System	The model, with the real part packet.u and the imaginary part packet.w of the wave function, or the quantum state schrodinger.state of the operator formulation.

Raises:

Type	Description
ValidationError	The formulation is unknown, or the field formulation has a dephasing rate.

Examples:

>>> from skcomplex.models.physics.waves import schrodinger
>>> schrodinger(64).elements[0].check().prepare().size_
64

Source code in src/skcomplex/models/physics/waves/_quantum.py

@register_model(
    field='physics/waves',
    tags=('waves', 'field', 'pde', 'quantum', 'spectral'),
    references=('Schrödinger 1926', 'Schwichtenberg 2020'),
    defaults={'until': 4.0, 'dt': 0.1, 'solver': ETDRK4Solver(dt=0.01)},
)
def schrodinger(
    n_sites: int = 256,
    *,
    length: float = 40.0,
    D: float = 0.5,
    momentum: float = 2.0,
    width: float = 1.0,
    omega: float = 0.0,
    formulation: str = 'fields',
    gamma: float = 0.0,
) -> System:
    """Return a Gaussian wave packet that follows the Schrödinger equation, as its real and imaginary parts.

    The wave function `psi = u + i w` follows `psi_t = i D laplacian(psi) - i V psi`, written as `u_t = -D laplacian(w)
    + V w` and `w_t = D laplacian(u) - V u`, on a periodic lattice with spectral derivatives. In units with `ħ = 1`,
    `D = 1 / (2 m)`. The potential is the harmonic well `V = omega² x² / (4 D)`, or zero for a free packet. A free
    packet of momentum `k` moves at the group velocity `2 D k` and spreads, and its norm `sum(u² + w²)` is conserved.

    The operator formulation builds the same particle as one quantum subsystem whose basis states are the sites, with
    the Hamiltonian `-D laplacian + V` of the same spectral Laplacian, normalised so that `|psi|²` is the probability of
    every site. It runs with the quantum solvers, and it accepts position dephasing at the rate `gamma`, which turns
    the state into a density matrix.

    Args:
        n_sites:
            The number of sites of the periodic lattice.

        length:
            The length of the periodic box, centered on 0.

        D:
            The coefficient of the Laplacian, `ħ / (2 m)`.

        momentum:
            The mean wave number of the packet.

        width:
            The width of the packet.

        omega:
            The angular frequency of the harmonic well.

        formulation:
            `'fields'` for two real fields, or `'operators'` for a particle on the sites of the lattice.

        gamma:
            The rate of position dephasing of the operator formulation.

    Returns:
        system:
            The model, with the real part `packet.u` and the imaginary part `packet.w` of the wave function, or the
            quantum state `schrodinger.state` of the operator formulation.

    Raises:
        ValidationError:
            The formulation is unknown, or the field formulation has a dephasing rate.

    Examples:
        >>> from skcomplex.models.physics.waves import schrodinger
        >>> schrodinger(64).elements[0].check().prepare().size_
        64
    """
    spacing = length / n_sites
    positions = spacing * np.arange(n_sites) - length / 2
    envelope = np.exp(-(positions**2) / (4 * width**2)) / (2 * np.pi * width**2) ** 0.25
    psi = envelope * np.exp(1j * momentum * positions)
    potential = omega**2 * positions**2 / (4 * D)
    space = LatticeSpace(n_sites, spacing=spacing)
    if formulation not in FORMULATIONS or (formulation == 'fields' and gamma):
        error_msg = (
            f'Unknown formulation {formulation!r}, or dephasing in the field formulation. Use one of '
            f'{", ".join(FORMULATIONS)}, and gamma > 0 only with the operator formulation.'
        )
        raise ValidationError(error_msg)
    if formulation == 'operators':
        amplitudes = psi / np.linalg.norm(psi)
        interactions: list[BaseInteraction] = [
            LatticeHamiltonianInteraction('kinetic', target='packet', D=D),
        ]
        if gamma > 0:
            interactions.append(
                PositionDephasingInteraction('dephasing', target='packet', gamma=gamma, n_sites=n_sites),
            )
        return System(
            'schrodinger',
            elements=[SiteElement('packet', amplitudes, hilbert=site_space(space, potential, projectors=gamma > 0))],
            interactions=interactions,
        )
    packet = FieldElement(
        'packet',
        element=WaveFunctionElement,
        space=space,
        values={'u': psi.real, 'w': psi.imag, 'D': D, 'V': potential},
        discretization='spectral',
    )
    return System(elements=[packet])

spring_chain(n_masses=15, *, length=1.0, k=1.0, mass=0.01, pluck=0.05, pluck_at=0.5, amplitude=0.0, omega=2.0)

Return masses joined in a row by springs of zero rest length between two anchors, a discrete string.

The masses start evenly spaced between the anchors at 0 and length, stretched, so the chain carries tension. With pluck, the masses start on a triangle of that height peaking at the fraction pluck_at of the length, the pluck of a guitar string. With amplitude, the left anchor moves transversally as amplitude sin(omega t) and sends a transverse wave along the chain. Springs of zero rest length are linear, so the transverse modes are independent, with frequencies 2 sqrt(k / m) sin(n pi / (2 (N + 1))).

Parameters:

Name	Type	Description	Default
n_masses	int	The number of masses between the anchors.	15
length	float	The distance between the anchors.	1.0
k	float	The stiffness of every spring.	1.0
mass	float	The mass of every mass.	0.01
pluck	float	The initial transverse displacement of the plucked point.	0.05
pluck_at	float	The plucked point, as a fraction of the length.	0.5
amplitude	float	The amplitude of the transverse motion of the left anchor.	0.0
omega	float	The angular frequency of the left anchor.	2.0

Returns:

Name	Type	Description
system	System	The model, with the masses mass_<i> and their positions mass_<i>.x in two dimensions.

Examples:

>>> from skcomplex.models.physics.waves import spring_chain
>>> len(spring_chain(4).interactions)
5

Source code in src/skcomplex/models/physics/waves/_strings.py

@register_model(
    field='physics/waves',
    tags=('waves', 'ode', 'springs', 'normal modes', 'mechanics'),
    references=('Argun et al. 2021',),
    defaults={'until': 20.0, 'dt': 0.05, 'solver': VelocityVerletSolver(dt=0.01)},
)
def spring_chain(
    n_masses: int = 15,
    *,
    length: float = 1.0,
    k: float = 1.0,
    mass: float = 0.01,
    pluck: float = 0.05,
    pluck_at: float = 0.5,
    amplitude: float = 0.0,
    omega: float = 2.0,
) -> System:
    """Return masses joined in a row by springs of zero rest length between two anchors, a discrete string.

    The masses start evenly spaced between the anchors at `0` and `length`, stretched, so the chain carries tension.
    With `pluck`, the masses start on a triangle of that height peaking at the fraction `pluck_at` of the length, the
    pluck of a guitar string. With `amplitude`, the left anchor moves transversally as `amplitude sin(omega t)` and
    sends a transverse wave along the chain. Springs of zero rest length are linear, so the transverse modes are
    independent, with frequencies `2 sqrt(k / m) sin(n pi / (2 (N + 1)))`.

    Args:
        n_masses:
            The number of masses between the anchors.

        length:
            The distance between the anchors.

        k:
            The stiffness of every spring.

        mass:
            The mass of every mass.

        pluck:
            The initial transverse displacement of the plucked point.

        pluck_at:
            The plucked point, as a fraction of the length.

        amplitude:
            The amplitude of the transverse motion of the left anchor.

        omega:
            The angular frequency of the left anchor.

    Returns:
        system:
            The model, with the masses `mass_<i>` and their positions `mass_<i>.x` in two dimensions.

    Examples:
        >>> from skcomplex.models.physics.waves import spring_chain
        >>> len(spring_chain(4).interactions)
        5
    """
    spacing = length / (n_masses + 1)
    positions = spacing * np.arange(1, n_masses + 1)
    heights = _pluck(positions, length, pluck, pluck_at) if pluck else np.zeros(n_masses)
    masses = [
        MassPointElement(f'mass_{index}', x=[float(positions[index]), float(heights[index])], x_d1=[0.0, 0.0], m=mass)
        for index in range(n_masses)
    ]
    springs = [
        SpringInteraction(f'spring_{index}', first=f'mass_{index}', second=f'mass_{index + 1}', k=k)
        for index in range(n_masses - 1)
    ]
    anchors = [
        DrivenAnchorInteraction('left_end', target='mass_0', k=k, center=[0.0, 0.0], amplitude=amplitude, omega=omega),
        DrivenAnchorInteraction('right_end', target=f'mass_{n_masses - 1}', k=k, center=[length, 0.0]),
    ]
    return System(elements=masses, interactions=[*anchors, *springs], space=ContinuousSpace(n_dim=2))

standing_waves(n_sites=99, *, length=2.0, c=1.0, mode=1, amplitude=1.0)

Return a string with fixed ends on [-length / 2, length / 2] that vibrates as a standing wave.

The displacement follows the wave equation u'' = c² u_xx, written as u' = v and v' = c² laplacian(u) on a lattice whose edges hold the string at zero. The string starts at rest in the mode amplitude sin(mode pi (x + length / 2) / length), which oscillates in place with the period 2 length / (mode c), up to the discretization.

Parameters:

Name	Type	Description	Default
n_sites	int	The number of sites inside the string.	99
length	float	The length of the string.	2.0
c	float	The speed of the waves.	1.0
mode	int	The number of half wavelengths along the string.	1
amplitude	float	The initial amplitude.	1.0

Returns:

Name	Type	Description
system	System	The model, with the displacement string.u and the velocity string.v.

Examples:

>>> from skcomplex.models.physics.waves import standing_waves
>>> standing_waves(9).elements[0].check().prepare().size_
9

Source code in src/skcomplex/models/physics/waves/_strings.py

@register_model(
    field='physics/waves',
    tags=('waves', 'field', 'pde', 'normal modes'),
    references=('Starkovich 2021',),
    defaults={'until': 4.0, 'dt': 0.05},
)
def standing_waves(
    n_sites: int = 99,
    *,
    length: float = 2.0,
    c: float = 1.0,
    mode: int = 1,
    amplitude: float = 1.0,
) -> System:
    """Return a string with fixed ends on `[-length / 2, length / 2]` that vibrates as a standing wave.

    The displacement follows the wave equation `u'' = c² u_xx`, written as `u' = v` and `v' = c² laplacian(u)` on a
    lattice whose edges hold the string at zero. The string starts at rest in the mode
    `amplitude sin(mode pi (x + length / 2) / length)`, which oscillates in place with the period
    `2 length / (mode c)`, up to the discretization.

    Args:
        n_sites:
            The number of sites inside the string.

        length:
            The length of the string.

        c:
            The speed of the waves.

        mode:
            The number of half wavelengths along the string.

        amplitude:
            The initial amplitude.

    Returns:
        system:
            The model, with the displacement `string.u` and the velocity `string.v`.

    Examples:
        >>> from skcomplex.models.physics.waves import standing_waves
        >>> standing_waves(9).elements[0].check().prepare().size_
        9
    """
    spacing = length / (n_sites + 1)
    positions = spacing * np.arange(1, n_sites + 1)
    shape = amplitude * np.sin(mode * math.pi * positions / length)
    space = LatticeSpace(n_sites, spacing=spacing)
    edges = ('fixed', 0.0)
    string = FieldElement(
        'string',
        element=WaveStringElement,
        space=space,
        boundary={'u': edges},
        values={'u': shape, 'c': c},
    )
    tension = RateInteraction(
        'tension',
        target='string',
        state='u',
        rate=lambda s: s.target.c**2 * laplacian(s.target.u),
    )
    return System(elements=[string], interactions=[tension])