Skip to content

Controllers

Controllers in OpenSMC integrate sliding surfaces and reaching laws to compute the actual control input \(u\) applied to the plant. OpenSMC includes 17 different controller types, from classical SMC to adaptive, fuzzy, and higher-order variants.

Usage Example

from opensmc.controllers import AdaptiveSMC
from opensmc.surfaces import LinearSurface
from opensmc.reaching import ConstantRate

# Compose an Adaptive SMC controller
surface = LinearSurface(c=5.0)
reaching = ConstantRate(k=2.0)
controller = AdaptiveSMC(surface, reaching, adaptation_gain=0.1)

# Compute control input for a given state
u = controller.compute(state=[1.0, 0.0], reference=[0.0, 0.0])

Available Controllers

controllers

OpenSMC Controllers — 15 SMC + 2 baselines (PID, LQR).

Classes

AdaptiveSMC

Bases: Controller

Adaptive SMC with online gain adaptation.

Parameters

surface : SlidingSurface or None Pluggable surface. If None, uses s = c*e + edot. reaching : ReachingLaw or None Pluggable reaching law (overrides internal logic if provided). c : float Surface slope. gamma : float Adaptation rate (eta_dot = gamma * |s|). eta0 : float Initial switching gain. eta_max : float Upper bound on eta (projection). lam : float Proportional sliding gain. dt : float Integration timestep for eta update (forward Euler). threshold : float Adaptation dead-zone: eta only adapts when |s| > threshold.

Source code in opensmc/controllers/adaptive.py
class AdaptiveSMC(Controller):
    """Adaptive SMC with online gain adaptation.

    Parameters
    ----------
    surface : SlidingSurface or None
        Pluggable surface. If None, uses s = c*e + edot.
    reaching : ReachingLaw or None
        Pluggable reaching law (overrides internal logic if provided).
    c : float
        Surface slope.
    gamma : float
        Adaptation rate (eta_dot = gamma * |s|).
    eta0 : float
        Initial switching gain.
    eta_max : float
        Upper bound on eta (projection).
    lam : float
        Proportional sliding gain.
    dt : float
        Integration timestep for eta update (forward Euler).
    threshold : float
        Adaptation dead-zone: eta only adapts when |s| > threshold.
    """

    def __init__(self, surface=None, reaching=None,
                 c=10.0, gamma=100.0, eta0=0.1, eta_max=100.0,
                 lam=5.0, dt=1e-4, threshold=0.0):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.gamma = gamma
        self.eta0 = eta0
        self.eta_max = eta_max
        self.lam = lam
        self.dt = dt
        self.threshold = threshold
        self._eta = eta0

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute adaptive SMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2]
        xref : array [xref1, xref2]
        plant : Plant or None

        Returns
        -------
        u : float
        info : dict with keys 's', 'eta', 'u_eq', 'u_rob'
        """
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Sliding variable
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = self.c * e + edot

        # Adaptation law: eta_dot = gamma * |s| (with projection)
        if np.abs(s) > self.threshold:
            self._eta += self.dt * self.gamma * np.abs(s)
        self._eta = min(self._eta, self.eta_max)

        # Equivalent control: from sdot = c*edot + eddot, with eddot = -u
        # for double integrator (xddot = u), sdot = c*edot - u + u_eq
        # Setting sdot = 0 => u_eq = c*edot = -c*x[1] when xref = 0
        # But we want u = u_eq - sdot_desired, and for eddot = -u:
        #   sdot = c*edot - u, so u = c*edot - sdot
        # If we want sdot = -eta*sign(s) - lam*s, then:
        #   u = c*edot + eta*sign(s) + lam*s
        u_eq = -self.c * x[1]

        # Robust control (positive sign drives sdot negative when s > 0)
        if self.reaching is not None:
            u_rob = -self.reaching.compute(s, t=t, **kwargs)
        else:
            u_rob = self._eta * np.sign(s) + self.lam * s

        u = u_eq + u_rob

        info = {
            's': s,
            'eta': self._eta,
            'u_eq': u_eq,
            'u_rob': u_rob,
            'e': e,
            'edot': edot,
        }
        return u, info

    def reset(self):
        """Reset adaptive gain to initial value."""
        super().reset()
        self._eta = self.eta0
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute adaptive SMC control.

Parameters

t : float x : array [x1, x2] xref : array [xref1, xref2] plant : Plant or None

Returns

u : float info : dict with keys 's', 'eta', 'u_eq', 'u_rob'

Source code in opensmc/controllers/adaptive.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute adaptive SMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2]
    xref : array [xref1, xref2]
    plant : Plant or None

    Returns
    -------
    u : float
    info : dict with keys 's', 'eta', 'u_eq', 'u_rob'
    """
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Sliding variable
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = self.c * e + edot

    # Adaptation law: eta_dot = gamma * |s| (with projection)
    if np.abs(s) > self.threshold:
        self._eta += self.dt * self.gamma * np.abs(s)
    self._eta = min(self._eta, self.eta_max)

    # Equivalent control: from sdot = c*edot + eddot, with eddot = -u
    # for double integrator (xddot = u), sdot = c*edot - u + u_eq
    # Setting sdot = 0 => u_eq = c*edot = -c*x[1] when xref = 0
    # But we want u = u_eq - sdot_desired, and for eddot = -u:
    #   sdot = c*edot - u, so u = c*edot - sdot
    # If we want sdot = -eta*sign(s) - lam*s, then:
    #   u = c*edot + eta*sign(s) + lam*s
    u_eq = -self.c * x[1]

    # Robust control (positive sign drives sdot negative when s > 0)
    if self.reaching is not None:
        u_rob = -self.reaching.compute(s, t=t, **kwargs)
    else:
        u_rob = self._eta * np.sign(s) + self.lam * s

    u = u_eq + u_rob

    info = {
        's': s,
        'eta': self._eta,
        'u_eq': u_eq,
        'u_rob': u_rob,
        'e': e,
        'edot': edot,
    }
    return u, info
reset()

Reset adaptive gain to initial value.

Source code in opensmc/controllers/adaptive.py
def reset(self):
    """Reset adaptive gain to initial value."""
    super().reset()
    self._eta = self.eta0

AggregatedHSMC

Bases: Controller

Aggregated Hierarchical SMC (Qian & Yi 2015, Sec 4.2).

Surfaces

s1 = c1e1 + e2 (trolley subsystem) s2 = c2e3 + e4 (swing subsystem) S = alpha*s1 + s2 (hierarchical combination)

Control

ueq1 = -(c1xdot + f1) / b1 ueq2 = -(c2thetadot + f2) / b2 den = alphab1 + b2 u = (alphab1ueq1 + b2ueq2 + kappaS + etasign(S)) / den

Parameters

c1 : float — trolley surface slope c2 : float — swing surface slope alpha : float — hierarchical weighting kappa : float — proportional reaching gain eta : float — switching gain u_max : float — control saturation

Source code in opensmc/controllers/hsmc.py
class AggregatedHSMC(Controller):
    """Aggregated Hierarchical SMC (Qian & Yi 2015, Sec 4.2).

    Surfaces:
      s1 = c1*e1 + e2     (trolley subsystem)
      s2 = c2*e3 + e4     (swing subsystem)
      S  = alpha*s1 + s2   (hierarchical combination)

    Control:
      ueq1 = -(c1*xdot + f1) / b1
      ueq2 = -(c2*thetadot + f2) / b2
      den  = alpha*b1 + b2
      u = (alpha*b1*ueq1 + b2*ueq2 + kappa*S + eta*sign(S)) / den

    Parameters
    ----------
    c1 : float — trolley surface slope
    c2 : float — swing surface slope
    alpha : float — hierarchical weighting
    kappa : float — proportional reaching gain
    eta : float — switching gain
    u_max : float — control saturation
    """

    def __init__(self, surface=None, reaching=None,
                 c1=0.7, c2=8.2, alpha=-2.3, kappa=3.0, eta=0.1,
                 u_max=50.0):
        super().__init__(surface=surface, reaching=reaching)
        self.c1 = c1
        self.c2 = c2
        self.alpha = alpha
        self.kappa = kappa
        self.eta = eta
        self.u_max = u_max

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute aggregated HSMC control.

        Parameters
        ----------
        t : float
        x : array [x_trolley, xdot, theta, thetadot]
        xref : array [xref, xdot_ref, theta_ref, thetadot_ref]
        plant : Plant with get_fb(x) method or None

        Returns
        -------
        u : float
        info : dict with 's', 's1', 's2'
        """
        e = xref - x
        s1 = self.c1 * e[0] + e[1]
        s2 = self.c2 * e[2] + e[3]
        S = self.alpha * s1 + s2

        f1, f2, b1, b2 = _get_fb(x, plant)

        ueq1 = -(self.c1 * x[1] + f1) / b1
        ueq2 = -(self.c2 * x[3] + f2) / b2
        den = self.alpha * b1 + b2

        u = (self.alpha * b1 * ueq1 + b2 * ueq2
             + self.kappa * S + self.eta * np.sign(S)) / den
        u = np.clip(u, -self.u_max, self.u_max)

        info = {
            's': S,
            's1': s1,
            's2': s2,
            'e': e,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute aggregated HSMC control.

Parameters

t : float x : array [x_trolley, xdot, theta, thetadot] xref : array [xref, xdot_ref, theta_ref, thetadot_ref] plant : Plant with get_fb(x) method or None

Returns

u : float info : dict with 's', 's1', 's2'

Source code in opensmc/controllers/hsmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute aggregated HSMC control.

    Parameters
    ----------
    t : float
    x : array [x_trolley, xdot, theta, thetadot]
    xref : array [xref, xdot_ref, theta_ref, thetadot_ref]
    plant : Plant with get_fb(x) method or None

    Returns
    -------
    u : float
    info : dict with 's', 's1', 's2'
    """
    e = xref - x
    s1 = self.c1 * e[0] + e[1]
    s2 = self.c2 * e[2] + e[3]
    S = self.alpha * s1 + s2

    f1, f2, b1, b2 = _get_fb(x, plant)

    ueq1 = -(self.c1 * x[1] + f1) / b1
    ueq2 = -(self.c2 * x[3] + f2) / b2
    den = self.alpha * b1 + b2

    u = (self.alpha * b1 * ueq1 + b2 * ueq2
         + self.kappa * S + self.eta * np.sign(S)) / den
    u = np.clip(u, -self.u_max, self.u_max)

    info = {
        's': S,
        's1': s1,
        's2': s2,
        'e': e,
    }
    return u, info

ClassicalSMC

Bases: Controller

Classical SMC — generic surface + reaching law composition.

Parameters

surface : SlidingSurface or None Pluggable sliding surface. If None, uses s = ce + edot. reaching : ReachingLaw or None Pluggable reaching law. If None, uses -ks - etasign(s). c : float Surface slope (used when surface is None). k : float Exponential reaching gain (used when reaching is None). eta : float Switching gain (used when reaching is None). J : float Inertia / plant gain (u_eq = J(-c*edot + xref_ddot)). Set J=1 for a unit double integrator.

Source code in opensmc/controllers/classical.py
class ClassicalSMC(Controller):
    """Classical SMC — generic surface + reaching law composition.

    Parameters
    ----------
    surface : SlidingSurface or None
        Pluggable sliding surface. If None, uses s = c*e + edot.
    reaching : ReachingLaw or None
        Pluggable reaching law. If None, uses -k*s - eta*sign(s).
    c : float
        Surface slope (used when surface is None).
    k : float
        Exponential reaching gain (used when reaching is None).
    eta : float
        Switching gain (used when reaching is None).
    J : float
        Inertia / plant gain (u_eq = J*(-c*edot + xref_ddot)).
        Set J=1 for a unit double integrator.
    """

    def __init__(self, surface=None, reaching=None,
                 c=10.0, k=10.0, eta=1.1, J=1.0):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.k = k
        self.eta = eta
        self.J = J

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute classical SMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2] — current state
        xref : array [xref1, xref2] or [xref1, xref2, xref_ddot]
        plant : Plant or None
        **kwargs :
            xref_ddot : float — reference acceleration (default 0)

        Returns
        -------
        u : float
        info : dict with keys 's', 'u_eq', 'u_rob', 'e', 'edot'
        """
        # Error: e = xref - x (desired minus actual)
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Reference acceleration
        xref_ddot = 0.0
        if len(xref) > 2:
            xref_ddot = xref[2]
        xref_ddot = kwargs.get('xref_ddot', xref_ddot)

        # Sliding variable
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = self.c * e + edot

        # Equivalent control: from sdot = eddot + c*edot = (xref_ddot - u/J) + c*edot
        # Setting sdot = 0: u_eq = J*(xref_ddot + c*edot)
        if self.surface is not None:
            # When using pluggable surface, u_eq = 0 (user handles model-based terms)
            u_eq = 0.0
        else:
            u_eq = self.J * (xref_ddot + self.c * edot)

        # Reaching / robust control
        # For sdot = -reaching(s), we need u_rob = +reaching output
        # since reaching laws return NEGATIVE values (e.g., -k*sign(s))
        # and we need sdot = reaching(s) < 0 when s > 0
        # From sdot = (xref_ddot + c*edot - u/J), sdot = reaching(s)
        # => u = J*(xref_ddot + c*edot - reaching(s)) = u_eq - J*reaching(s)
        if self.reaching is not None:
            u_rob = -self.J * self.reaching.compute(s, t=t, **kwargs)
        else:
            u_rob = self.J * (self.k * s + self.eta * np.sign(s))

        u = u_eq + u_rob

        info = {
            's': s,
            'u_eq': u_eq,
            'u_rob': u_rob,
            'e': e,
            'edot': edot,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute classical SMC control.

Parameters

t : float x : array [x1, x2] — current state xref : array [xref1, xref2] or [xref1, xref2, xref_ddot] plant : Plant or None **kwargs : xref_ddot : float — reference acceleration (default 0)

Returns

u : float info : dict with keys 's', 'u_eq', 'u_rob', 'e', 'edot'

Source code in opensmc/controllers/classical.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute classical SMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2] — current state
    xref : array [xref1, xref2] or [xref1, xref2, xref_ddot]
    plant : Plant or None
    **kwargs :
        xref_ddot : float — reference acceleration (default 0)

    Returns
    -------
    u : float
    info : dict with keys 's', 'u_eq', 'u_rob', 'e', 'edot'
    """
    # Error: e = xref - x (desired minus actual)
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Reference acceleration
    xref_ddot = 0.0
    if len(xref) > 2:
        xref_ddot = xref[2]
    xref_ddot = kwargs.get('xref_ddot', xref_ddot)

    # Sliding variable
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = self.c * e + edot

    # Equivalent control: from sdot = eddot + c*edot = (xref_ddot - u/J) + c*edot
    # Setting sdot = 0: u_eq = J*(xref_ddot + c*edot)
    if self.surface is not None:
        # When using pluggable surface, u_eq = 0 (user handles model-based terms)
        u_eq = 0.0
    else:
        u_eq = self.J * (xref_ddot + self.c * edot)

    # Reaching / robust control
    # For sdot = -reaching(s), we need u_rob = +reaching output
    # since reaching laws return NEGATIVE values (e.g., -k*sign(s))
    # and we need sdot = reaching(s) < 0 when s > 0
    # From sdot = (xref_ddot + c*edot - u/J), sdot = reaching(s)
    # => u = J*(xref_ddot + c*edot - reaching(s)) = u_eq - J*reaching(s)
    if self.reaching is not None:
        u_rob = -self.J * self.reaching.compute(s, t=t, **kwargs)
    else:
        u_rob = self.J * (self.k * s + self.eta * np.sign(s))

    u = u_eq + u_rob

    info = {
        's': s,
        'u_eq': u_eq,
        'u_rob': u_rob,
        'e': e,
        'edot': edot,
    }
    return u, info

CombiningHSMC

Bases: Controller

Combining Hierarchical SMC (Qian & Yi 2015, Sec 4.4).

Intermediate variable with sign-switching (Eq 4.58): c_eff = |c| if e1*e3 >= 0, else -|c| z = e1 + c_eff * e3 zdot = e2 + c_eff * e4 s = alpha * z + zdot

Control

ueq = (alpha(e2 + c_effe4) - (f1 + c_efff2)) / (b1 + c_effb2) usw = (kappas + etasign(s)) / (b1 + c_eff*b2) u = ueq + usw

Tuning note: c must satisfy c << l (cable length) to avoid denominator discontinuity. For l=0.305, c=0.01 is stable.

Parameters

c : float — combining coefficient (|c| << cable length) alpha : float — surface parameter kappa : float — reaching gain eta : float — switching gain u_max : float — control saturation

Source code in opensmc/controllers/hsmc.py
class CombiningHSMC(Controller):
    """Combining Hierarchical SMC (Qian & Yi 2015, Sec 4.4).

    Intermediate variable with sign-switching (Eq 4.58):
      c_eff = |c| if e1*e3 >= 0, else -|c|
      z    = e1 + c_eff * e3
      zdot = e2 + c_eff * e4
      s    = alpha * z + zdot

    Control:
      ueq = (alpha*(e2 + c_eff*e4) - (f1 + c_eff*f2)) / (b1 + c_eff*b2)
      usw = (kappa*s + eta*sign(s)) / (b1 + c_eff*b2)
      u = ueq + usw

    Tuning note: c must satisfy c << l (cable length) to avoid
    denominator discontinuity. For l=0.305, c=0.01 is stable.

    Parameters
    ----------
    c : float — combining coefficient (|c| << cable length)
    alpha : float — surface parameter
    kappa : float — reaching gain
    eta : float — switching gain
    u_max : float — control saturation
    """

    def __init__(self, surface=None, reaching=None,
                 c=0.010, alpha=0.280, kappa=1.0, eta=0.025,
                 u_max=50.0):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.alpha = alpha
        self.kappa = kappa
        self.eta = eta
        self.u_max = u_max

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute combining HSMC control.

        Parameters
        ----------
        t : float
        x : array [x_trolley, xdot, theta, thetadot]
        xref : array [xref, xdot_ref, theta_ref, thetadot_ref]
        plant : Plant with get_fb(x) or None

        Returns
        -------
        u : float
        info : dict with 's', 'z', 'c_eff'
        """
        e = xref - x
        e1, e2, e3, e4 = e[0], e[1], e[2], e[3]
        c_abs = np.abs(self.c)

        # Sign-switching (Eq 4.58)
        c_eff = c_abs if e1 * e3 >= 0 else -c_abs

        # Intermediate variable
        z = e1 + c_eff * e3
        zdot = e2 + c_eff * e4
        s = self.alpha * z + zdot

        # Plant dynamics decomposition
        f1, f2, b1, b2 = _get_fb(x, plant)
        den = b1 + c_eff * b2

        # Avoid division by zero
        if np.abs(den) < 1e-10:
            den = np.sign(den) * 1e-10 if den != 0 else 1e-10

        # Equivalent control
        ueq = (self.alpha * (e2 + c_eff * e4) - (f1 + c_eff * f2)) / den

        # Switching control
        usw = (self.kappa * s + self.eta * np.sign(s)) / den

        u = np.clip(ueq + usw, -self.u_max, self.u_max)

        info = {
            's': s,
            'z': z,
            'c_eff': c_eff,
            'e': e,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute combining HSMC control.

Parameters

t : float x : array [x_trolley, xdot, theta, thetadot] xref : array [xref, xdot_ref, theta_ref, thetadot_ref] plant : Plant with get_fb(x) or None

Returns

u : float info : dict with 's', 'z', 'c_eff'

Source code in opensmc/controllers/hsmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute combining HSMC control.

    Parameters
    ----------
    t : float
    x : array [x_trolley, xdot, theta, thetadot]
    xref : array [xref, xdot_ref, theta_ref, thetadot_ref]
    plant : Plant with get_fb(x) or None

    Returns
    -------
    u : float
    info : dict with 's', 'z', 'c_eff'
    """
    e = xref - x
    e1, e2, e3, e4 = e[0], e[1], e[2], e[3]
    c_abs = np.abs(self.c)

    # Sign-switching (Eq 4.58)
    c_eff = c_abs if e1 * e3 >= 0 else -c_abs

    # Intermediate variable
    z = e1 + c_eff * e3
    zdot = e2 + c_eff * e4
    s = self.alpha * z + zdot

    # Plant dynamics decomposition
    f1, f2, b1, b2 = _get_fb(x, plant)
    den = b1 + c_eff * b2

    # Avoid division by zero
    if np.abs(den) < 1e-10:
        den = np.sign(den) * 1e-10 if den != 0 else 1e-10

    # Equivalent control
    ueq = (self.alpha * (e2 + c_eff * e4) - (f1 + c_eff * f2)) / den

    # Switching control
    usw = (self.kappa * s + self.eta * np.sign(s)) / den

    u = np.clip(ueq + usw, -self.u_max, self.u_max)

    info = {
        's': s,
        'z': z,
        'c_eff': c_eff,
        'e': e,
    }
    return u, info

DiscreteSMC

Bases: Controller

Discrete-time SMC with Gao reaching law.

Parameters

surface : SlidingSurface or None Pluggable surface. If None, uses s = ce + edot. reaching : ReachingLaw or None If provided, overrides internal Gao reaching law. c : float Surface slope. q_rate : float Exponential reaching rate (must satisfy qTs < 1). epsilon : float Constant reaching gain. Ts : float Sampling period. u_max : float Control saturation.

Source code in opensmc/controllers/discrete.py
class DiscreteSMC(Controller):
    """Discrete-time SMC with Gao reaching law.

    Parameters
    ----------
    surface : SlidingSurface or None
        Pluggable surface. If None, uses s = c*e + edot.
    reaching : ReachingLaw or None
        If provided, overrides internal Gao reaching law.
    c : float
        Surface slope.
    q_rate : float
        Exponential reaching rate (must satisfy q*Ts < 1).
    epsilon : float
        Constant reaching gain.
    Ts : float
        Sampling period.
    u_max : float
        Control saturation.
    """

    def __init__(self, surface=None, reaching=None,
                 c=10.0, q_rate=10.0, epsilon=5.0, Ts=0.01,
                 u_max=500.0):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.q_rate = q_rate
        self.epsilon = epsilon
        self.Ts = Ts
        self.u_max = u_max
        # Input-to-surface gain for double integrator with ZOH:
        #   g = Ts + c * Ts^2 / 2
        self._g = Ts + c * Ts**2 / 2.0

    @property
    def quasi_sliding_band(self):
        """Theoretical quasi-sliding band width."""
        return self.epsilon * self.Ts / (2.0 - self.q_rate * self.Ts)

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute discrete SMC control.

        Parameters
        ----------
        t : float — current discrete time
        x : array [x1, x2] — current state
        xref : array [xref1, xref2] — reference state
        plant : Plant or None

        Returns
        -------
        u : float
        info : dict with 's', 'band'
        """
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Sliding variable
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = self.c * e + edot

        # Discrete reaching law control:
        #   s(k+1) = (1-q*Ts)*s(k) - eps*Ts*sign(s(k))
        # Requires: u = -(q*Ts*s + eps*Ts*sign(s) + c*Ts*x2) / g
        if self.reaching is not None:
            u_rob = -self.reaching.compute(s, t=t, **kwargs)
            u = u_rob
        else:
            # For double integrator xddot = u, the discrete reaching law gives:
            #   s(k+1) - s(k) = c*Ts*edot + Ts*eddot ≈ c*Ts*edot - g*u
            # Setting s(k+1) = (1-q*Ts)*s - eps*Ts*sign(s):
            #   u = (q*Ts*s + eps*Ts*sign(s) + c*Ts*edot) / g
            u = (self.q_rate * self.Ts * s
                 + self.epsilon * self.Ts * np.sign(s)
                 + self.c * self.Ts * edot) / self._g

        u = np.clip(u, -self.u_max, self.u_max)

        info = {
            's': s,
            'band': self.quasi_sliding_band,
            'e': e,
            'edot': edot,
        }
        return u, info
Attributes
quasi_sliding_band property

Theoretical quasi-sliding band width.

Functions
compute(t, x, xref, plant=None, **kwargs)

Compute discrete SMC control.

Parameters

t : float — current discrete time x : array [x1, x2] — current state xref : array [xref1, xref2] — reference state plant : Plant or None

Returns

u : float info : dict with 's', 'band'

Source code in opensmc/controllers/discrete.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute discrete SMC control.

    Parameters
    ----------
    t : float — current discrete time
    x : array [x1, x2] — current state
    xref : array [xref1, xref2] — reference state
    plant : Plant or None

    Returns
    -------
    u : float
    info : dict with 's', 'band'
    """
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Sliding variable
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = self.c * e + edot

    # Discrete reaching law control:
    #   s(k+1) = (1-q*Ts)*s(k) - eps*Ts*sign(s(k))
    # Requires: u = -(q*Ts*s + eps*Ts*sign(s) + c*Ts*x2) / g
    if self.reaching is not None:
        u_rob = -self.reaching.compute(s, t=t, **kwargs)
        u = u_rob
    else:
        # For double integrator xddot = u, the discrete reaching law gives:
        #   s(k+1) - s(k) = c*Ts*edot + Ts*eddot ≈ c*Ts*edot - g*u
        # Setting s(k+1) = (1-q*Ts)*s - eps*Ts*sign(s):
        #   u = (q*Ts*s + eps*Ts*sign(s) + c*Ts*edot) / g
        u = (self.q_rate * self.Ts * s
             + self.epsilon * self.Ts * np.sign(s)
             + self.c * self.Ts * edot) / self._g

    u = np.clip(u, -self.u_max, self.u_max)

    info = {
        's': s,
        'band': self.quasi_sliding_band,
        'e': e,
        'edot': edot,
    }
    return u, info

DynamicSMC

Bases: Controller

Dynamic SMC — continuous control via integrated udot.

Parameters

surface : SlidingSurface or None Pluggable surface. If None, uses s = ce + edot. reaching : ReachingLaw or None If provided, used for the sigma-level reaching. c : float Surface slope. K : float Switching gain on sigma (lives in udot, not u). lam : float Dynamic surface parameter (sigma = sdot + lams). q : float Proportional gain on sigma. u_max : float Control saturation limit. dt : float Integration timestep for u update.

Source code in opensmc/controllers/dynamic.py
class DynamicSMC(Controller):
    """Dynamic SMC — continuous control via integrated udot.

    Parameters
    ----------
    surface : SlidingSurface or None
        Pluggable surface. If None, uses s = c*e + edot.
    reaching : ReachingLaw or None
        If provided, used for the sigma-level reaching.
    c : float
        Surface slope.
    K : float
        Switching gain on sigma (lives in udot, not u).
    lam : float
        Dynamic surface parameter (sigma = sdot + lam*s).
    q : float
        Proportional gain on sigma.
    u_max : float
        Control saturation limit.
    dt : float
        Integration timestep for u update.
    """

    def __init__(self, surface=None, reaching=None,
                 c=10.0, K=20.0, lam=10.0, q=5.0,
                 u_max=500.0, dt=1e-4):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.K = K
        self.lam = lam
        self.q = q
        self.u_max = u_max
        self.dt = dt
        self._u = 0.0

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute dynamic SMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2]
        xref : array [xref1, xref2]
        plant : Plant or None

        Returns
        -------
        u : float
        info : dict with keys 's', 'sigma', 'udot'
        """
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Sliding variable
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = self.c * e + edot

        # Dynamic sliding variable:
        #   sigma = u + (c + lam)*x[1] + lam*c*x[0]
        # This is derived from sigma = sdot + lam*s where
        #   sdot = c*edot + eddot = c*(-x[1]) + (u + d - xref_ddot)
        # For regulation (xref=0): sigma = u + (c+lam)*x[1] + lam*c*x[0]
        sigma = self._u + (self.c + self.lam) * x[1] + self.lam * self.c * x[0]

        # udot: discontinuity lives here, not in u
        #   udot = -K*sign(sigma) - q*sigma - (c+lam)*u - lam*c*x[1]
        udot = (-self.K * np.sign(sigma)
                - self.q * sigma
                - (self.c + self.lam) * self._u
                - self.lam * self.c * x[1])

        # Integrate u (forward Euler)
        self._u = np.clip(self._u + self.dt * udot,
                          -self.u_max, self.u_max)

        info = {
            's': s,
            'sigma': sigma,
            'udot': udot,
            'e': e,
            'edot': edot,
        }
        return self._u, info

    def reset(self):
        """Reset internal control state."""
        super().reset()
        self._u = 0.0
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute dynamic SMC control.

Parameters

t : float x : array [x1, x2] xref : array [xref1, xref2] plant : Plant or None

Returns

u : float info : dict with keys 's', 'sigma', 'udot'

Source code in opensmc/controllers/dynamic.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute dynamic SMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2]
    xref : array [xref1, xref2]
    plant : Plant or None

    Returns
    -------
    u : float
    info : dict with keys 's', 'sigma', 'udot'
    """
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Sliding variable
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = self.c * e + edot

    # Dynamic sliding variable:
    #   sigma = u + (c + lam)*x[1] + lam*c*x[0]
    # This is derived from sigma = sdot + lam*s where
    #   sdot = c*edot + eddot = c*(-x[1]) + (u + d - xref_ddot)
    # For regulation (xref=0): sigma = u + (c+lam)*x[1] + lam*c*x[0]
    sigma = self._u + (self.c + self.lam) * x[1] + self.lam * self.c * x[0]

    # udot: discontinuity lives here, not in u
    #   udot = -K*sign(sigma) - q*sigma - (c+lam)*u - lam*c*x[1]
    udot = (-self.K * np.sign(sigma)
            - self.q * sigma
            - (self.c + self.lam) * self._u
            - self.lam * self.c * x[1])

    # Integrate u (forward Euler)
    self._u = np.clip(self._u + self.dt * udot,
                      -self.u_max, self.u_max)

    info = {
        's': s,
        'sigma': sigma,
        'udot': udot,
        'e': e,
        'edot': edot,
    }
    return self._u, info
reset()

Reset internal control state.

Source code in opensmc/controllers/dynamic.py
def reset(self):
    """Reset internal control state."""
    super().reset()
    self._u = 0.0

FixedTimeSMC

Bases: Controller

Fixed-time SMC with bi-power reaching law.

Parameters

surface : SlidingSurface or None Pluggable surface. If None, uses s = c*e + edot. reaching : ReachingLaw or None If provided, overrides bi-power reaching. c : float Surface slope. alpha : float Sub-linear reaching gain (handles small |s|). beta : float Super-linear reaching gain (handles large |s|). p : float Sub-linear power, 0 < p < 1. q : float Super-linear power, q > 1.

Source code in opensmc/controllers/fixed_time.py
class FixedTimeSMC(Controller):
    """Fixed-time SMC with bi-power reaching law.

    Parameters
    ----------
    surface : SlidingSurface or None
        Pluggable surface. If None, uses s = c*e + edot.
    reaching : ReachingLaw or None
        If provided, overrides bi-power reaching.
    c : float
        Surface slope.
    alpha : float
        Sub-linear reaching gain (handles small |s|).
    beta : float
        Super-linear reaching gain (handles large |s|).
    p : float
        Sub-linear power, 0 < p < 1.
    q : float
        Super-linear power, q > 1.
    """

    def __init__(self, surface=None, reaching=None,
                 c=10.0, alpha=5.0, beta=5.0, p=0.5, q=1.5):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.alpha = alpha
        self.beta = beta
        self.p = p
        self.q = q

    @property
    def T_max(self):
        """Upper bound on settling time (independent of IC)."""
        return 1.0 / (self.alpha * (1.0 - self.p)) + \
               1.0 / (self.beta * (self.q - 1.0))

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute fixed-time SMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2]
        xref : array [xref1, xref2]
        plant : Plant or None

        Returns
        -------
        u : float
        info : dict with 's', 'u_eq', 'u_ft', 'T_max'
        """
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Sliding variable
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = self.c * e + edot

        # Equivalent control: from sdot = c*edot + eddot, with eddot = -u
        # for double integrator xddot = u => u_eq = c*edot
        u_eq = self.c * edot

        # Bi-power reaching: sdot = -alpha*|s|^p*sign(s) - beta*|s|^q*sign(s)
        # u = c*edot - sdot = c*edot + alpha*|s|^p*sign(s) + beta*|s|^q*sign(s)
        # Clip |s| to avoid overflow in |s|^q when s is large
        s_clipped = np.clip(np.abs(s), 0.0, 1e6)
        if self.reaching is not None:
            u_ft = -self.reaching.compute(s, t=t, **kwargs)
        else:
            u_ft = (self.alpha * s_clipped ** self.p * np.sign(s)
                    + self.beta * s_clipped ** self.q * np.sign(s))

        u = u_eq + u_ft

        info = {
            's': s,
            'u_eq': u_eq,
            'u_ft': u_ft,
            'T_max': self.T_max,
            'e': e,
            'edot': edot,
        }
        return u, info
Attributes
T_max property

Upper bound on settling time (independent of IC).

Functions
compute(t, x, xref, plant=None, **kwargs)

Compute fixed-time SMC control.

Parameters

t : float x : array [x1, x2] xref : array [xref1, xref2] plant : Plant or None

Returns

u : float info : dict with 's', 'u_eq', 'u_ft', 'T_max'

Source code in opensmc/controllers/fixed_time.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute fixed-time SMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2]
    xref : array [xref1, xref2]
    plant : Plant or None

    Returns
    -------
    u : float
    info : dict with 's', 'u_eq', 'u_ft', 'T_max'
    """
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Sliding variable
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = self.c * e + edot

    # Equivalent control: from sdot = c*edot + eddot, with eddot = -u
    # for double integrator xddot = u => u_eq = c*edot
    u_eq = self.c * edot

    # Bi-power reaching: sdot = -alpha*|s|^p*sign(s) - beta*|s|^q*sign(s)
    # u = c*edot - sdot = c*edot + alpha*|s|^p*sign(s) + beta*|s|^q*sign(s)
    # Clip |s| to avoid overflow in |s|^q when s is large
    s_clipped = np.clip(np.abs(s), 0.0, 1e6)
    if self.reaching is not None:
        u_ft = -self.reaching.compute(s, t=t, **kwargs)
    else:
        u_ft = (self.alpha * s_clipped ** self.p * np.sign(s)
                + self.beta * s_clipped ** self.q * np.sign(s))

    u = u_eq + u_ft

    info = {
        's': s,
        'u_eq': u_eq,
        'u_ft': u_ft,
        'T_max': self.T_max,
        'e': e,
        'edot': edot,
    }
    return u, info

FuzzySMC

Bases: Controller

Fuzzy SMC — smooth gain modulation via Mamdani inference.

Parameters

surface : SlidingSurface or None Pluggable surface. If None, uses s = c*e + edot. reaching : ReachingLaw or None If provided, overrides fuzzy reaching. c : float Surface slope. k : float Maximum switching gain. s_range : float Normalization range for s in fuzzy inference. sdot_range : float Normalization range for sdot in fuzzy inference. dt : float Timestep for numerical sdot estimation.

Source code in opensmc/controllers/fuzzy.py
class FuzzySMC(Controller):
    """Fuzzy SMC — smooth gain modulation via Mamdani inference.

    Parameters
    ----------
    surface : SlidingSurface or None
        Pluggable surface. If None, uses s = c*e + edot.
    reaching : ReachingLaw or None
        If provided, overrides fuzzy reaching.
    c : float
        Surface slope.
    k : float
        Maximum switching gain.
    s_range : float
        Normalization range for s in fuzzy inference.
    sdot_range : float
        Normalization range for sdot in fuzzy inference.
    dt : float
        Timestep for numerical sdot estimation.
    """

    def __init__(self, surface=None, reaching=None,
                 c=10.0, k=15.0, s_range=1.0, sdot_range=5.0,
                 dt=1e-4):
        super().__init__(surface=surface, reaching=reaching)
        self.c = c
        self.k = k
        self.s_range = s_range
        self.sdot_range = sdot_range
        self.dt = dt
        self._s_prev = 0.0
        self._first_call = True

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute fuzzy SMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2]
        xref : array [xref1, xref2]
        plant : Plant or None

        Returns
        -------
        u : float
        info : dict with 's', 'k_fuzzy', 'sdot'
        """
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Sliding variable
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = self.c * e + edot

        # Estimate sdot numerically
        if self._first_call:
            sdot = 0.0
            self._first_call = False
        else:
            sdot = (s - self._s_prev) / self.dt
        self._s_prev = s

        # Fuzzy inference: k_fuzzy in [-1, 1]
        k_fuzzy = _fuzzy_inference(s, sdot,
                                   s_range=self.s_range,
                                   sdot_range=self.sdot_range)

        # Equivalent control: for xddot=u with s = c*e + edot,
        # sdot = c*edot + eddot = c*edot - u => u_eq = c*edot
        u_eq = self.c * edot

        # Fuzzy reaching control: u_rob replaces sign(s) with k_fuzzy
        if self.reaching is not None:
            u_rob = -self.reaching.compute(s, t=t, **kwargs)
        else:
            u_rob = self.k * k_fuzzy

        u = u_eq + u_rob

        info = {
            's': s,
            'k_fuzzy': k_fuzzy,
            'sdot': sdot,
            'e': e,
            'edot': edot,
        }
        return u, info

    def reset(self):
        """Reset internal state."""
        super().reset()
        self._s_prev = 0.0
        self._first_call = True
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute fuzzy SMC control.

Parameters

t : float x : array [x1, x2] xref : array [xref1, xref2] plant : Plant or None

Returns

u : float info : dict with 's', 'k_fuzzy', 'sdot'

Source code in opensmc/controllers/fuzzy.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute fuzzy SMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2]
    xref : array [xref1, xref2]
    plant : Plant or None

    Returns
    -------
    u : float
    info : dict with 's', 'k_fuzzy', 'sdot'
    """
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Sliding variable
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = self.c * e + edot

    # Estimate sdot numerically
    if self._first_call:
        sdot = 0.0
        self._first_call = False
    else:
        sdot = (s - self._s_prev) / self.dt
    self._s_prev = s

    # Fuzzy inference: k_fuzzy in [-1, 1]
    k_fuzzy = _fuzzy_inference(s, sdot,
                               s_range=self.s_range,
                               sdot_range=self.sdot_range)

    # Equivalent control: for xddot=u with s = c*e + edot,
    # sdot = c*edot + eddot = c*edot - u => u_eq = c*edot
    u_eq = self.c * edot

    # Fuzzy reaching control: u_rob replaces sign(s) with k_fuzzy
    if self.reaching is not None:
        u_rob = -self.reaching.compute(s, t=t, **kwargs)
    else:
        u_rob = self.k * k_fuzzy

    u = u_eq + u_rob

    info = {
        's': s,
        'k_fuzzy': k_fuzzy,
        'sdot': sdot,
        'e': e,
        'edot': edot,
    }
    return u, info
reset()

Reset internal state.

Source code in opensmc/controllers/fuzzy.py
def reset(self):
    """Reset internal state."""
    super().reset()
    self._s_prev = 0.0
    self._first_call = True

ITSMC

Bases: Controller

Integral Terminal SMC (Al Ghanimi 2026).

Parameters

surface : SlidingSurface or None If provided, overrides the built-in ITSMC surface. reaching : ReachingLaw or None If provided, overrides the built-in reaching law. c1 : float Linear error gain. c2 : float Integral terminal gain. p : int Numerator of fractional power (p < q). q : int Denominator of fractional power. K : float Saturation-based switching gain. lambda_s : float Proportional sliding gain. phi : float Boundary layer thickness for saturation. dt : float Integration timestep for integral state E.

Source code in opensmc/controllers/itsmc.py
class ITSMC(Controller):
    """Integral Terminal SMC (Al Ghanimi 2026).

    Parameters
    ----------
    surface : SlidingSurface or None
        If provided, overrides the built-in ITSMC surface.
    reaching : ReachingLaw or None
        If provided, overrides the built-in reaching law.
    c1 : float
        Linear error gain.
    c2 : float
        Integral terminal gain.
    p : int
        Numerator of fractional power (p < q).
    q : int
        Denominator of fractional power.
    K : float
        Saturation-based switching gain.
    lambda_s : float
        Proportional sliding gain.
    phi : float
        Boundary layer thickness for saturation.
    dt : float
        Integration timestep for integral state E.
    """

    def __init__(self, surface=None, reaching=None,
                 c1=5.0, c2=2.0, p=5, q=7,
                 K=5.0, lambda_s=3.0, phi=0.05, dt=1e-4):
        super().__init__(surface=surface, reaching=reaching)
        self.c1 = c1
        self.c2 = c2
        self.p = p
        self.q = q
        self.K = K
        self.lambda_s = lambda_s
        self.phi = phi
        self.dt = dt
        self._E = 0.0  # integral state: int(|e|^(p/q)*sign(e)) dt

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute ITSMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2]
        xref : array [xref1, xref2] or [xref1, xref2, xref_ddot]
        plant : Plant or None
        **kwargs :
            xref_ddot : float — reference acceleration (default 0)

        Returns
        -------
        u : float
        info : dict with 's', 'E', 'u_eq', 'u_rob', 'e', 'edot'
        """
        e = xref[0] - x[0]
        edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

        # Reference acceleration
        xref_ddot = 0.0
        if len(xref) > 2:
            xref_ddot = xref[2]
        xref_ddot = kwargs.get('xref_ddot', xref_ddot)

        # ITSMC surface: s = edot + c1*e + c2*E
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t,
                                     E=self._E, **kwargs)
        else:
            s = edot + self.c1 * e + self.c2 * self._E

        # Equivalent control:
        #   From sdot = (xref_ddot - u - d) + c1*edot + c2*|e|^(p/q)*sign(e) = 0
        #   u_eq = xref_ddot + c1*edot + c2*|e|^(p/q)*sign(e)
        frac_e = _fractional_sign(e, self.p, self.q)
        u_eq = xref_ddot + self.c1 * edot + self.c2 * frac_e

        # Robust control: u_rob = K*sat(s/phi) + lambda_s*s
        if self.reaching is not None:
            u_rob = self.reaching.compute(s, t=t, **kwargs)
        else:
            u_rob = self.K * _sat(s / self.phi) + self.lambda_s * s

        u = u_eq + u_rob

        # Update integral state (forward Euler)
        self._E += frac_e * self.dt

        info = {
            's': s,
            'E': self._E,
            'u_eq': u_eq,
            'u_rob': u_rob,
            'e': e,
            'edot': edot,
        }
        return u, info

    def reset(self):
        """Reset integral state."""
        super().reset()
        self._E = 0.0
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute ITSMC control.

Parameters

t : float x : array [x1, x2] xref : array [xref1, xref2] or [xref1, xref2, xref_ddot] plant : Plant or None **kwargs : xref_ddot : float — reference acceleration (default 0)

Returns

u : float info : dict with 's', 'E', 'u_eq', 'u_rob', 'e', 'edot'

Source code in opensmc/controllers/itsmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute ITSMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2]
    xref : array [xref1, xref2] or [xref1, xref2, xref_ddot]
    plant : Plant or None
    **kwargs :
        xref_ddot : float — reference acceleration (default 0)

    Returns
    -------
    u : float
    info : dict with 's', 'E', 'u_eq', 'u_rob', 'e', 'edot'
    """
    e = xref[0] - x[0]
    edot = xref[1] - x[1] if len(xref) > 1 else -x[1]

    # Reference acceleration
    xref_ddot = 0.0
    if len(xref) > 2:
        xref_ddot = xref[2]
    xref_ddot = kwargs.get('xref_ddot', xref_ddot)

    # ITSMC surface: s = edot + c1*e + c2*E
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t,
                                 E=self._E, **kwargs)
    else:
        s = edot + self.c1 * e + self.c2 * self._E

    # Equivalent control:
    #   From sdot = (xref_ddot - u - d) + c1*edot + c2*|e|^(p/q)*sign(e) = 0
    #   u_eq = xref_ddot + c1*edot + c2*|e|^(p/q)*sign(e)
    frac_e = _fractional_sign(e, self.p, self.q)
    u_eq = xref_ddot + self.c1 * edot + self.c2 * frac_e

    # Robust control: u_rob = K*sat(s/phi) + lambda_s*s
    if self.reaching is not None:
        u_rob = self.reaching.compute(s, t=t, **kwargs)
    else:
        u_rob = self.K * _sat(s / self.phi) + self.lambda_s * s

    u = u_eq + u_rob

    # Update integral state (forward Euler)
    self._E += frac_e * self.dt

    info = {
        's': s,
        'E': self._E,
        'u_eq': u_eq,
        'u_rob': u_rob,
        'e': e,
        'edot': edot,
    }
    return u, info
reset()

Reset integral state.

Source code in opensmc/controllers/itsmc.py
def reset(self):
    """Reset integral state."""
    super().reset()
    self._E = 0.0

IncrementalHSMC

Bases: Controller

Incremental Hierarchical SMC (Qian & Yi 2015, Sec 4.3).

Three-layer incremental hierarchical sliding surface

Layer 1: s1 = c1e1 + e2 Layer 2: s2 = c2_effe3 + s1 Layer 3: s3 = c3*e4 + s2

Sign-switching rule (Eq 4.34): c2_eff = +|c2| if x[2]s1 <= 0 c2_eff = -|c2| if x[2]s1 > 0

Control (derived from sdot3 = 0): ueq = -(c3f2 + f1 + c2_effthetadot + c1xdot) / (c3b2 + b1) usw = +(kappas3 + etasign(s3)) / (c3*b2 + b1) u = ueq + usw

Parameters

c1 : float — layer 1 slope c2 : float — layer 2 coupling (absolute value; sign switches) c3 : float — layer 3 slope kappa : float — reaching gain eta : float — switching gain u_max : float — control saturation

Source code in opensmc/controllers/hsmc.py
class IncrementalHSMC(Controller):
    """Incremental Hierarchical SMC (Qian & Yi 2015, Sec 4.3).

    Three-layer incremental hierarchical sliding surface:
      Layer 1: s1 = c1*e1 + e2
      Layer 2: s2 = c2_eff*e3 + s1
      Layer 3: s3 = c3*e4 + s2

    Sign-switching rule (Eq 4.34):
      c2_eff = +|c2| if x[2]*s1 <= 0
      c2_eff = -|c2| if x[2]*s1 > 0

    Control (derived from sdot3 = 0):
      ueq = -(c3*f2 + f1 + c2_eff*thetadot + c1*xdot) / (c3*b2 + b1)
      usw = +(kappa*s3 + eta*sign(s3)) / (c3*b2 + b1)
      u = ueq + usw

    Parameters
    ----------
    c1 : float — layer 1 slope
    c2 : float — layer 2 coupling (absolute value; sign switches)
    c3 : float — layer 3 slope
    kappa : float — reaching gain
    eta : float — switching gain
    u_max : float — control saturation
    """

    def __init__(self, surface=None, reaching=None,
                 c1=0.7, c2=2.0, c3=0.3, kappa=5.0, eta=0.5,
                 u_max=50.0):
        super().__init__(surface=surface, reaching=reaching)
        self.c1 = c1
        self.c2 = c2
        self.c3 = c3
        self.kappa = kappa
        self.eta = eta
        self.u_max = u_max

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute incremental HSMC control.

        Parameters
        ----------
        t : float
        x : array [x_trolley, xdot, theta, thetadot]
        xref : array [xref, xdot_ref, theta_ref, thetadot_ref]
        plant : Plant with get_fb(x) or None

        Returns
        -------
        u : float
        info : dict with 's' (=s3), 's1', 's2', 'c2_eff'
        """
        e = xref - x

        # Layer 1
        s1 = self.c1 * e[0] + e[1]

        # Sign-switching (Eq 4.34)
        c2_abs = np.abs(self.c2)
        c2_eff = c2_abs if x[2] * s1 <= 0 else -c2_abs

        # Layer 2
        s2 = c2_eff * e[2] + s1

        # Layer 3
        s3 = self.c3 * e[3] + s2

        # Plant dynamics decomposition
        f1, f2, b1, b2 = _get_fb(x, plant)
        den = self.c3 * b2 + b1

        # Equivalent control (uses STATES, not errors)
        ueq = -(self.c3 * f2 + f1 + c2_eff * x[3] + self.c1 * x[1]) / den

        # Switching control (positive sign because den < 0 for crane)
        usw = (self.kappa * s3 + self.eta * np.sign(s3)) / den

        u = np.clip(ueq + usw, -self.u_max, self.u_max)

        info = {
            's': s3,
            's1': s1,
            's2': s2,
            'c2_eff': c2_eff,
            'e': e,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute incremental HSMC control.

Parameters

t : float x : array [x_trolley, xdot, theta, thetadot] xref : array [xref, xdot_ref, theta_ref, thetadot_ref] plant : Plant with get_fb(x) or None

Returns

u : float info : dict with 's' (=s3), 's1', 's2', 'c2_eff'

Source code in opensmc/controllers/hsmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute incremental HSMC control.

    Parameters
    ----------
    t : float
    x : array [x_trolley, xdot, theta, thetadot]
    xref : array [xref, xdot_ref, theta_ref, thetadot_ref]
    plant : Plant with get_fb(x) or None

    Returns
    -------
    u : float
    info : dict with 's' (=s3), 's1', 's2', 'c2_eff'
    """
    e = xref - x

    # Layer 1
    s1 = self.c1 * e[0] + e[1]

    # Sign-switching (Eq 4.34)
    c2_abs = np.abs(self.c2)
    c2_eff = c2_abs if x[2] * s1 <= 0 else -c2_abs

    # Layer 2
    s2 = c2_eff * e[2] + s1

    # Layer 3
    s3 = self.c3 * e[3] + s2

    # Plant dynamics decomposition
    f1, f2, b1, b2 = _get_fb(x, plant)
    den = self.c3 * b2 + b1

    # Equivalent control (uses STATES, not errors)
    ueq = -(self.c3 * f2 + f1 + c2_eff * x[3] + self.c1 * x[1]) / den

    # Switching control (positive sign because den < 0 for crane)
    usw = (self.kappa * s3 + self.eta * np.sign(s3)) / den

    u = np.clip(ueq + usw, -self.u_max, self.u_max)

    info = {
        's': s3,
        's1': s1,
        's2': s2,
        'c2_eff': c2_eff,
        'e': e,
    }
    return u, info

LQR

Bases: Controller

Linear Quadratic Regulator for benchmark comparison.

u = -K @ (x - xref)

where K = R^{-1} B^T P and P solves the CARE.

Parameters

A : ndarray (n, n) State matrix of the linearized system. B : ndarray (n, m) Input matrix. Q : ndarray (n, n) or None State cost matrix. If None, uses identity. R : ndarray (m, m) or None Input cost matrix. If None, uses identity. u_max : float or None Control saturation limit.

Example

For a double integrator (xddot = u): A = [[0, 1], [0, 0]] B = [[0], [1]] Q = [[10, 0], [0, 1]] (penalize position more) R = [[1]]

Source code in opensmc/controllers/lqr.py
class LQR(Controller):
    """Linear Quadratic Regulator for benchmark comparison.

    u = -K @ (x - xref)

    where K = R^{-1} B^T P and P solves the CARE.

    Parameters
    ----------
    A : ndarray (n, n)
        State matrix of the linearized system.
    B : ndarray (n, m)
        Input matrix.
    Q : ndarray (n, n) or None
        State cost matrix. If None, uses identity.
    R : ndarray (m, m) or None
        Input cost matrix. If None, uses identity.
    u_max : float or None
        Control saturation limit.

    Example
    -------
    For a double integrator (xddot = u):
        A = [[0, 1], [0, 0]]
        B = [[0], [1]]
        Q = [[10, 0], [0, 1]]  (penalize position more)
        R = [[1]]
    """

    def __init__(self, A, B, Q=None, R=None, u_max=None):
        super().__init__()
        self.A = np.atleast_2d(A).astype(float)
        self.B = np.atleast_2d(B).astype(float)
        n = self.A.shape[0]
        m = self.B.shape[1]

        if Q is None:
            self.Q = np.eye(n)
        else:
            self.Q = np.atleast_2d(Q).astype(float)

        if R is None:
            self.R = np.eye(m)
        else:
            self.R = np.atleast_2d(R).astype(float)

        self.u_max = u_max
        self._n = n
        self._m = m

        # Solve CARE and compute gain
        self.P, self.K = self._solve_care()

    def _solve_care(self):
        """Solve the continuous algebraic Riccati equation."""
        try:
            from scipy.linalg import solve_continuous_are
            P = solve_continuous_are(self.A, self.B, self.Q, self.R)
        except ImportError:
            n = self._n
            P = np.eye(n) * 10
        except Exception:
            # CARE failed (ill-conditioned, uncontrollable, etc.)
            # Fallback: rough diagonal gain
            n = self._n
            P = np.eye(n) * 10
        K = np.linalg.solve(self.R, self.B.T @ P)
        return P, K

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute LQR control: u = -K @ (x - xref).

        Parameters
        ----------
        t : float
        x : ndarray (n,)
        xref : ndarray (n,) — reference state (regulation target)

        Returns
        -------
        u : float or ndarray
        info : dict with 's' (= position error), 'K'
        """
        x_arr = np.atleast_1d(x[:self._n]).astype(float)
        xref_arr = np.atleast_1d(xref[:self._n]).astype(float)

        # Pad xref if shorter than state
        if len(xref_arr) < self._n:
            xref_arr = np.concatenate([xref_arr, np.zeros(self._n - len(xref_arr))])

        e_state = x_arr - xref_arr
        u = -self.K @ e_state

        if self.u_max is not None:
            u = np.clip(u, -self.u_max, self.u_max)

        # Return scalar for single-input
        if self._m == 1:
            u = float(u.ravel()[0])

        info = {
            's': float(xref_arr[0] - x_arr[0]),
            'e_state': e_state,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute LQR control: u = -K @ (x - xref).

Parameters

t : float x : ndarray (n,) xref : ndarray (n,) — reference state (regulation target)

Returns

u : float or ndarray info : dict with 's' (= position error), 'K'

Source code in opensmc/controllers/lqr.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute LQR control: u = -K @ (x - xref).

    Parameters
    ----------
    t : float
    x : ndarray (n,)
    xref : ndarray (n,) — reference state (regulation target)

    Returns
    -------
    u : float or ndarray
    info : dict with 's' (= position error), 'K'
    """
    x_arr = np.atleast_1d(x[:self._n]).astype(float)
    xref_arr = np.atleast_1d(xref[:self._n]).astype(float)

    # Pad xref if shorter than state
    if len(xref_arr) < self._n:
        xref_arr = np.concatenate([xref_arr, np.zeros(self._n - len(xref_arr))])

    e_state = x_arr - xref_arr
    u = -self.K @ e_state

    if self.u_max is not None:
        u = np.clip(u, -self.u_max, self.u_max)

    # Return scalar for single-input
    if self._m == 1:
        u = float(u.ravel()[0])

    info = {
        's': float(xref_arr[0] - x_arr[0]),
        'e_state': e_state,
    }
    return u, info

NFTSMC

Bases: Controller

Nonsingular Fast Terminal SMC (Al Ghanimi 2026).

Parameters

surface : SlidingSurface or None If provided, overrides the built-in NFTSMC surface. reaching : ReachingLaw or None If provided, overrides the built-in reaching law. alpha : float Linear error gain. beta : float Terminal error gain. gamma : float Terminal power (1 < gamma < 2, singularity-free). k1 : float Proportional reaching gain. k2 : float Power reaching gain. rho : float Reaching power (typically ~5/7). eta : float Continuous switching gain. delta : float Smoothing parameter for sn = s/(|s|+delta). kI : float Integral action gain (0 = disabled). g_input : float Input gain (for xddot = f + g*u, default g=1). dt : float Timestep for integral state update.

Source code in opensmc/controllers/nftsmc.py
class NFTSMC(Controller):
    """Nonsingular Fast Terminal SMC (Al Ghanimi 2026).

    Parameters
    ----------
    surface : SlidingSurface or None
        If provided, overrides the built-in NFTSMC surface.
    reaching : ReachingLaw or None
        If provided, overrides the built-in reaching law.
    alpha : float
        Linear error gain.
    beta : float
        Terminal error gain.
    gamma : float
        Terminal power (1 < gamma < 2, singularity-free).
    k1 : float
        Proportional reaching gain.
    k2 : float
        Power reaching gain.
    rho : float
        Reaching power (typically ~5/7).
    eta : float
        Continuous switching gain.
    delta : float
        Smoothing parameter for sn = s/(|s|+delta).
    kI : float
        Integral action gain (0 = disabled).
    g_input : float
        Input gain (for xddot = f + g*u, default g=1).
    dt : float
        Timestep for integral state update.
    """

    def __init__(self, surface=None, reaching=None,
                 alpha=5.0, beta=2.0, gamma=1.4,
                 k1=10.0, k2=3.0, rho=0.714,
                 eta=1.0, delta=1e-6,
                 kI=0.0, g_input=1.0, dt=1e-4):
        super().__init__(surface=surface, reaching=reaching)
        self.alpha = alpha
        self.beta = beta
        self.gamma = gamma
        self.k1 = k1
        self.k2 = k2
        self.rho = rho
        self.eta = eta
        self.delta = delta
        self.kI = kI
        self.g_input = g_input
        self.dt = dt
        self._e_int = 0.0  # integral of error for integral action

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute NFTSMC control.

        Parameters
        ----------
        t : float
        x : array [x1, x2]
        xref : array [xref1, xref2] — reference (constant or slow-varying)
        plant : Plant or None
        **kwargs :
            f_plant : float — plant nonlinearity f(x) if known (default 0)

        Returns
        -------
        u : float
        info : dict with 's', 'ds', 'sn', 'e', 'edot', 'e_int'
        """
        # Error convention: e = x - xref (for NFTSMC, per notebook)
        # The surface drives e -> 0 from above/below
        e = x[0] - xref[0]
        edot = x[1] - (xref[1] if len(xref) > 1 else 0.0)

        # Plant nonlinearity
        f_plant = kwargs.get('f_plant', 0.0)

        # NFTSMC surface: s = edot + alpha*e + beta*|e|^gamma*sign(e)
        if self.surface is not None:
            s = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            s = (edot
                 + self.alpha * e
                 + self.beta * np.abs(e) ** self.gamma * np.sign(e))

        # Surface derivative term: ds = beta*gamma*|e|^(gamma-1)*edot
        ds = (self.beta * self.gamma
              * (np.abs(e) + 1e-15) ** (self.gamma - 1.0) * edot)

        # Continuous sign approximation: sn = s / (|s| + delta)
        sn = s / (np.abs(s) + self.delta)

        # Control law:
        #   u = (-alpha*edot - ds - k1*s - k2*|s|^rho*sign(s) - eta*sn) / g
        if self.reaching is not None:
            u_rob = self.reaching.compute(s, t=t, **kwargs)
            u = (-self.alpha * edot - ds + u_rob) / self.g_input
        else:
            u = (-self.alpha * edot
                 - ds
                 - self.k1 * s
                 - self.k2 * np.abs(s) ** self.rho * np.sign(s)
                 - self.eta * sn) / self.g_input

        # Integral action (optional)
        if self.kI > 0:
            self._e_int += e * self.dt
            u -= self.kI * self._e_int

        info = {
            's': s,
            'ds': ds,
            'sn': sn,
            'e': e,
            'edot': edot,
            'e_int': self._e_int,
        }
        return u, info

    def reset(self):
        """Reset integral state."""
        super().reset()
        self._e_int = 0.0
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute NFTSMC control.

Parameters

t : float x : array [x1, x2] xref : array [xref1, xref2] — reference (constant or slow-varying) plant : Plant or None **kwargs : f_plant : float — plant nonlinearity f(x) if known (default 0)

Returns

u : float info : dict with 's', 'ds', 'sn', 'e', 'edot', 'e_int'

Source code in opensmc/controllers/nftsmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute NFTSMC control.

    Parameters
    ----------
    t : float
    x : array [x1, x2]
    xref : array [xref1, xref2] — reference (constant or slow-varying)
    plant : Plant or None
    **kwargs :
        f_plant : float — plant nonlinearity f(x) if known (default 0)

    Returns
    -------
    u : float
    info : dict with 's', 'ds', 'sn', 'e', 'edot', 'e_int'
    """
    # Error convention: e = x - xref (for NFTSMC, per notebook)
    # The surface drives e -> 0 from above/below
    e = x[0] - xref[0]
    edot = x[1] - (xref[1] if len(xref) > 1 else 0.0)

    # Plant nonlinearity
    f_plant = kwargs.get('f_plant', 0.0)

    # NFTSMC surface: s = edot + alpha*e + beta*|e|^gamma*sign(e)
    if self.surface is not None:
        s = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        s = (edot
             + self.alpha * e
             + self.beta * np.abs(e) ** self.gamma * np.sign(e))

    # Surface derivative term: ds = beta*gamma*|e|^(gamma-1)*edot
    ds = (self.beta * self.gamma
          * (np.abs(e) + 1e-15) ** (self.gamma - 1.0) * edot)

    # Continuous sign approximation: sn = s / (|s| + delta)
    sn = s / (np.abs(s) + self.delta)

    # Control law:
    #   u = (-alpha*edot - ds - k1*s - k2*|s|^rho*sign(s) - eta*sn) / g
    if self.reaching is not None:
        u_rob = self.reaching.compute(s, t=t, **kwargs)
        u = (-self.alpha * edot - ds + u_rob) / self.g_input
    else:
        u = (-self.alpha * edot
             - ds
             - self.k1 * s
             - self.k2 * np.abs(s) ** self.rho * np.sign(s)
             - self.eta * sn) / self.g_input

    # Integral action (optional)
    if self.kI > 0:
        self._e_int += e * self.dt
        u -= self.kI * self._e_int

    info = {
        's': s,
        'ds': ds,
        'sn': sn,
        'e': e,
        'edot': edot,
        'e_int': self._e_int,
    }
    return u, info
reset()

Reset integral state.

Source code in opensmc/controllers/nftsmc.py
def reset(self):
    """Reset integral state."""
    super().reset()
    self._e_int = 0.0

NestedHOSMC

Bases: Controller

Nested r-sliding mode controller (Levant 2003, Shtessel Sec 6.6.1).

Recursive definition

N_{i,r} = (Σ_{j=0}^{i-1} |σ^{(j)}|^{q/(r-j)})^{1/q} Ψ_{0,r} = sign(σ) Ψ_{i,r} = sign(σ^{(i)} + β_i · N_{i,r} · Ψ_{i-1,r}) u = -α · Ψ_{r-1,r}

where q = lcm(1, 2, ..., r) ensures integer exponents.

Special cases

r=1: u = -α · sign(σ) (relay control) r=2: u = -α · sign(σ̇ + β|σ|^½·sign(σ)) (prescribed convergence)

Parameters

order : int Sliding mode order r (= relative degree of the system). alpha : float Control magnitude. Must be large enough for convergence. betas : list of float or None Gains [β_1, ..., β_{r-1}]. Length must be r-1. If None, uses [1.0, 2.0, 3.0, ...] (heuristic defaults). c : float Surface slope for computing σ = c·e + edot when order=2 and no explicit sigma_derivs are provided. g_sign : float Sign of the input gain in σ^{(r)} = h + g·u. Default -1 for σ=e on standard plants where σ̈=-u.

Source code in opensmc/controllers/hosmc.py
class NestedHOSMC(Controller):
    """Nested r-sliding mode controller (Levant 2003, Shtessel Sec 6.6.1).

    Recursive definition:
        N_{i,r} = (Σ_{j=0}^{i-1} |σ^{(j)}|^{q/(r-j)})^{1/q}
        Ψ_{0,r} = sign(σ)
        Ψ_{i,r} = sign(σ^{(i)} + β_i · N_{i,r} · Ψ_{i-1,r})
        u = -α · Ψ_{r-1,r}

    where q = lcm(1, 2, ..., r) ensures integer exponents.

    Special cases:
        r=1: u = -α · sign(σ)               (relay control)
        r=2: u = -α · sign(σ̇ + β|σ|^½·sign(σ))  (prescribed convergence)

    Parameters
    ----------
    order : int
        Sliding mode order r (= relative degree of the system).
    alpha : float
        Control magnitude. Must be large enough for convergence.
    betas : list of float or None
        Gains [β_1, ..., β_{r-1}]. Length must be r-1.
        If None, uses [1.0, 2.0, 3.0, ...] (heuristic defaults).
    c : float
        Surface slope for computing σ = c·e + edot when order=2
        and no explicit sigma_derivs are provided.
    g_sign : float
        Sign of the input gain in σ^{(r)} = h + g·u.
        Default -1 for σ=e on standard plants where σ̈=-u.
    """

    def __init__(self, order=2, alpha=10.0, betas=None, c=1.0, g_sign=-1.0):
        super().__init__()
        if order < 1:
            raise ValueError(f"order must be >= 1, got {order}")
        self.order = order
        self.alpha = alpha
        self.c = c
        self.g_sign = np.sign(g_sign) if g_sign != 0 else -1.0

        # Default betas
        if betas is None:
            self.betas = [float(i) for i in range(1, order)]
        else:
            if len(betas) != order - 1:
                raise ValueError(
                    f"betas length {len(betas)} != order-1 = {order-1}")
            self.betas = list(betas)

        # Precompute q = lcm(1, ..., r)
        self._q = _compute_lcm(range(1, order + 1)) if order > 1 else 1

    def _compute_psi(self, sigma_derivs):
        """Compute Ψ_{r-1,r} from the derivative vector."""
        r = self.order
        q = self._q

        psi = np.sign(sigma_derivs[0])  # Ψ_{0,r}

        for i in range(1, r):
            # N_{i,r} = (Σ_{j=0}^{i-1} |σ^{(j)}|^{q/(r-j)})^{1/q}
            terms = 0.0
            for j in range(i):
                exp = q / (r - j)
                terms += np.abs(sigma_derivs[j]) ** exp
            N_i = terms ** (1.0 / q) if terms > 1e-30 else 0.0

            # Ψ_{i,r} = sign(σ^{(i)} + β_i · N_{i,r} · Ψ_{i-1,r})
            psi = np.sign(sigma_derivs[i] + self.betas[i-1] * N_i * psi)

        return psi

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute nested HOSMC control.

        Parameters
        ----------
        t : float
        x : ndarray — state vector
        xref : ndarray — reference
        sigma_derivs : list of float, optional
            [σ, σ̇, ..., σ^{(r-1)}]. If not provided, computed from
            (x, xref) assuming 2nd-order system with σ = c·e + edot.

        Returns
        -------
        u : float
        info : dict with 's', 'psi'
        """
        sigma_derivs = kwargs.get('sigma_derivs', None)

        if sigma_derivs is None:
            # Default: 2nd-order system, σ = c·e + edot, σ̇ = c·edot
            e = xref[0] - x[0]
            edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]
            sigma = self.c * e + edot
            sigma_dot = self.c * edot
            sigma_derivs = [sigma, sigma_dot]

        psi = self._compute_psi(sigma_derivs)
        v = -self.alpha * psi
        u = v / self.g_sign

        info = {
            's': sigma_derivs[0],
            'psi': psi,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute nested HOSMC control.

Parameters

t : float x : ndarray — state vector xref : ndarray — reference sigma_derivs : list of float, optional [σ, σ̇, ..., σ^{(r-1)}]. If not provided, computed from (x, xref) assuming 2nd-order system with σ = c·e + edot.

Returns

u : float info : dict with 's', 'psi'

Source code in opensmc/controllers/hosmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute nested HOSMC control.

    Parameters
    ----------
    t : float
    x : ndarray — state vector
    xref : ndarray — reference
    sigma_derivs : list of float, optional
        [σ, σ̇, ..., σ^{(r-1)}]. If not provided, computed from
        (x, xref) assuming 2nd-order system with σ = c·e + edot.

    Returns
    -------
    u : float
    info : dict with 's', 'psi'
    """
    sigma_derivs = kwargs.get('sigma_derivs', None)

    if sigma_derivs is None:
        # Default: 2nd-order system, σ = c·e + edot, σ̇ = c·edot
        e = xref[0] - x[0]
        edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]
        sigma = self.c * e + edot
        sigma_dot = self.c * edot
        sigma_derivs = [sigma, sigma_dot]

    psi = self._compute_psi(sigma_derivs)
    v = -self.alpha * psi
    u = v / self.g_sign

    info = {
        's': sigma_derivs[0],
        'psi': psi,
    }
    return u, info

PID

Bases: Controller

Standard PID controller for benchmark comparison.

u = Kp * e + Kd * edot + Ki * integral(e)

For multi-input systems (n_inputs > 1), applies independent PID per channel using diagonal gain matrices.

Parameters

Kp : float or ndarray Proportional gain. Kd : float or ndarray Derivative gain. Ki : float or ndarray Integral gain (default 0 = PD controller). dt : float Time step for integral accumulation. u_max : float or None Control saturation limit. If None, no saturation. n_channels : int Number of independent control channels (default 1).

Source code in opensmc/controllers/pid.py
class PID(Controller):
    """Standard PID controller for benchmark comparison.

    u = Kp * e + Kd * edot + Ki * integral(e)

    For multi-input systems (n_inputs > 1), applies independent PID
    per channel using diagonal gain matrices.

    Parameters
    ----------
    Kp : float or ndarray
        Proportional gain.
    Kd : float or ndarray
        Derivative gain.
    Ki : float or ndarray
        Integral gain (default 0 = PD controller).
    dt : float
        Time step for integral accumulation.
    u_max : float or None
        Control saturation limit. If None, no saturation.
    n_channels : int
        Number of independent control channels (default 1).
    """

    def __init__(self, Kp=10.0, Kd=5.0, Ki=0.0, dt=1e-4,
                 u_max=None, n_channels=1):
        super().__init__()
        self.Kp = np.atleast_1d(Kp).astype(float)
        self.Kd = np.atleast_1d(Kd).astype(float)
        self.Ki = np.atleast_1d(Ki).astype(float)
        self.dt = dt
        self.u_max = u_max
        self.n_channels = n_channels
        self._e_int = np.zeros(n_channels)

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute PID control.

        For a 2nd-order single-channel system:
            e = xref[0] - x[0], edot = xref[1] - x[1]

        For multi-channel (e.g., 2-link arm):
            e = [xref[0]-x[0], xref[2]-x[2]], edot = [xref[1]-x[1], xref[3]-x[3]]

        Returns
        -------
        u : float or ndarray
        info : dict with 's' (= e for PID), 'e', 'edot'
        """
        nc = self.n_channels

        if nc == 1:
            e = xref[0] - x[0]
            edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]

            self._e_int[0] += e * self.dt
            u = float(self.Kp[0] * e + self.Kd[0] * edot + self.Ki[0] * self._e_int[0])

            if self.u_max is not None:
                u = np.clip(u, -self.u_max, self.u_max)

            info = {'s': e, 'e': e, 'edot': edot}
            return u, info
        else:
            # Multi-channel: extract position/velocity pairs
            e = np.zeros(nc)
            edot = np.zeros(nc)
            for ch in range(nc):
                pos_idx = 2 * ch
                vel_idx = 2 * ch + 1
                if pos_idx < len(xref):
                    e[ch] = xref[pos_idx] - x[pos_idx]
                if vel_idx < len(xref):
                    edot[ch] = xref[vel_idx] - x[vel_idx]
                elif vel_idx < len(x):
                    edot[ch] = -x[vel_idx]

            self._e_int += e * self.dt

            Kp = self.Kp if len(self.Kp) == nc else np.full(nc, self.Kp[0])
            Kd = self.Kd if len(self.Kd) == nc else np.full(nc, self.Kd[0])
            Ki = self.Ki if len(self.Ki) == nc else np.full(nc, self.Ki[0])

            u = Kp * e + Kd * edot + Ki * self._e_int

            if self.u_max is not None:
                u = np.clip(u, -self.u_max, self.u_max)

            info = {'s': e[0], 'e': e, 'edot': edot}
            return u, info

    def reset(self):
        """Reset integral accumulator."""
        self._e_int = np.zeros(self.n_channels)
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute PID control.

For a 2nd-order single-channel system

e = xref[0] - x[0], edot = xref[1] - x[1]

For multi-channel (e.g., 2-link arm): e = [xref[0]-x[0], xref[2]-x[2]], edot = [xref[1]-x[1], xref[3]-x[3]]

Returns

u : float or ndarray info : dict with 's' (= e for PID), 'e', 'edot'

Source code in opensmc/controllers/pid.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute PID control.

    For a 2nd-order single-channel system:
        e = xref[0] - x[0], edot = xref[1] - x[1]

    For multi-channel (e.g., 2-link arm):
        e = [xref[0]-x[0], xref[2]-x[2]], edot = [xref[1]-x[1], xref[3]-x[3]]

    Returns
    -------
    u : float or ndarray
    info : dict with 's' (= e for PID), 'e', 'edot'
    """
    nc = self.n_channels

    if nc == 1:
        e = xref[0] - x[0]
        edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]

        self._e_int[0] += e * self.dt
        u = float(self.Kp[0] * e + self.Kd[0] * edot + self.Ki[0] * self._e_int[0])

        if self.u_max is not None:
            u = np.clip(u, -self.u_max, self.u_max)

        info = {'s': e, 'e': e, 'edot': edot}
        return u, info
    else:
        # Multi-channel: extract position/velocity pairs
        e = np.zeros(nc)
        edot = np.zeros(nc)
        for ch in range(nc):
            pos_idx = 2 * ch
            vel_idx = 2 * ch + 1
            if pos_idx < len(xref):
                e[ch] = xref[pos_idx] - x[pos_idx]
            if vel_idx < len(xref):
                edot[ch] = xref[vel_idx] - x[vel_idx]
            elif vel_idx < len(x):
                edot[ch] = -x[vel_idx]

        self._e_int += e * self.dt

        Kp = self.Kp if len(self.Kp) == nc else np.full(nc, self.Kp[0])
        Kd = self.Kd if len(self.Kd) == nc else np.full(nc, self.Kd[0])
        Ki = self.Ki if len(self.Ki) == nc else np.full(nc, self.Ki[0])

        u = Kp * e + Kd * edot + Ki * self._e_int

        if self.u_max is not None:
            u = np.clip(u, -self.u_max, self.u_max)

        info = {'s': e[0], 'e': e, 'edot': edot}
        return u, info
reset()

Reset integral accumulator.

Source code in opensmc/controllers/pid.py
def reset(self):
    """Reset integral accumulator."""
    self._e_int = np.zeros(self.n_channels)

QuasiContinuous2SMC

Bases: Controller

Quasi-continuous 2-sliding mode controller (Levant 2005, Shtessel Eq 4.30).

u = -α · (σ̇ + β|σ|^{1/2} · sign(σ)) / (|σ̇| + β|σ|^{1/2})

Parameters

alpha : float Control magnitude. Must satisfy α·Km - C > β²/2 where Km is the minimum input gain and C is the disturbance bound. beta : float Convergence parameter. Larger β gives faster convergence but requires larger α. Typical: β = 1.0 to 5.0. c : float Surface slope for computing σ = c·e + edot (when no surface is provided). surface : SlidingSurface or None Optional pluggable surface for σ computation. epsilon : float Small positive constant to avoid division by zero when both σ and σ̇ are very close to zero. Default 1e-10.

Source code in opensmc/controllers/quasi_continuous.py
class QuasiContinuous2SMC(Controller):
    """Quasi-continuous 2-sliding mode controller (Levant 2005, Shtessel Eq 4.30).

    u = -α · (σ̇ + β|σ|^{1/2} · sign(σ)) / (|σ̇| + β|σ|^{1/2})

    Parameters
    ----------
    alpha : float
        Control magnitude. Must satisfy α·Km - C > β²/2 where
        Km is the minimum input gain and C is the disturbance bound.
    beta : float
        Convergence parameter. Larger β gives faster convergence
        but requires larger α. Typical: β = 1.0 to 5.0.
    c : float
        Surface slope for computing σ = c·e + edot (when no surface
        is provided).
    surface : SlidingSurface or None
        Optional pluggable surface for σ computation.
    epsilon : float
        Small positive constant to avoid division by zero when
        both σ and σ̇ are very close to zero. Default 1e-10.
    """

    def __init__(self, alpha=10.0, beta=2.0, c=1.0,
                 surface=None, epsilon=1e-10, g_sign=-1.0):
        super().__init__(surface=surface)
        if alpha <= 0:
            raise ValueError(f"alpha must be positive, got {alpha}")
        if beta <= 0:
            raise ValueError(f"beta must be positive, got {beta}")
        self.alpha = alpha
        self.beta = beta
        self.c = c
        self.epsilon = epsilon
        # g_sign: sign of the input gain in σ̈ = h + g·u.
        # Default -1 for σ=e on DoubleIntegrator where σ̈ = -u.
        self.g_sign = np.sign(g_sign) if g_sign != 0 else -1.0

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute quasi-continuous 2-SM control.

        Parameters
        ----------
        t : float
            Current time.
        x : ndarray
            State vector. For 2nd-order: [x1, x2].
        xref : ndarray
            Reference. For 2nd-order: [xref1, xref2].
        plant : Plant or None
            Not used (model-free controller).

        Returns
        -------
        u : float
            Control signal, bounded by |u| ≤ α.
        info : dict
            Contains 's' (σ), 'sdot' (σ̇), 'e', 'edot'.
        """
        # Tracking errors
        e = xref[0] - x[0]
        edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]

        # Sliding variable σ
        if self.surface is not None:
            sigma = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            sigma = self.c * e + edot

        # σ̇ estimate (for linear surface σ = c·e + edot)
        sigma_dot = self.c * edot

        # Quasi-continuous control law (Eq 4.30)
        #   numerator = σ̇ + β|σ|^{1/2}·sign(σ)
        #   denominator = |σ̇| + β|σ|^{1/2}
        abs_sigma_half = self.beta * np.abs(sigma) ** 0.5
        numerator = sigma_dot + abs_sigma_half * np.sign(sigma)
        denominator = np.abs(sigma_dot) + abs_sigma_half + self.epsilon

        # The textbook assumes g > 0. For g < 0, flip sign.
        v = -self.alpha * numerator / denominator
        u = v / self.g_sign

        info = {
            's': sigma,
            'sdot': sigma_dot,
            'e': e,
            'edot': edot,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute quasi-continuous 2-SM control.

Parameters

t : float Current time. x : ndarray State vector. For 2nd-order: [x1, x2]. xref : ndarray Reference. For 2nd-order: [xref1, xref2]. plant : Plant or None Not used (model-free controller).

Returns

u : float Control signal, bounded by |u| ≤ α. info : dict Contains 's' (σ), 'sdot' (σ̇), 'e', 'edot'.

Source code in opensmc/controllers/quasi_continuous.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute quasi-continuous 2-SM control.

    Parameters
    ----------
    t : float
        Current time.
    x : ndarray
        State vector. For 2nd-order: [x1, x2].
    xref : ndarray
        Reference. For 2nd-order: [xref1, xref2].
    plant : Plant or None
        Not used (model-free controller).

    Returns
    -------
    u : float
        Control signal, bounded by |u| ≤ α.
    info : dict
        Contains 's' (σ), 'sdot' (σ̇), 'e', 'edot'.
    """
    # Tracking errors
    e = xref[0] - x[0]
    edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]

    # Sliding variable σ
    if self.surface is not None:
        sigma = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        sigma = self.c * e + edot

    # σ̇ estimate (for linear surface σ = c·e + edot)
    sigma_dot = self.c * edot

    # Quasi-continuous control law (Eq 4.30)
    #   numerator = σ̇ + β|σ|^{1/2}·sign(σ)
    #   denominator = |σ̇| + β|σ|^{1/2}
    abs_sigma_half = self.beta * np.abs(sigma) ** 0.5
    numerator = sigma_dot + abs_sigma_half * np.sign(sigma)
    denominator = np.abs(sigma_dot) + abs_sigma_half + self.epsilon

    # The textbook assumes g > 0. For g < 0, flip sign.
    v = -self.alpha * numerator / denominator
    u = v / self.g_sign

    info = {
        's': sigma,
        'sdot': sigma_dot,
        'e': e,
        'edot': edot,
    }
    return u, info

QuasiContinuousHOSMC

Bases: Controller

Quasi-continuous r-sliding mode controller (Levant 2005, Shtessel Sec 6.6.2).

Recursive definition

φ_{0,r} = σ, N_{0,r} = |σ|, Ψ_{0,r} = sign(σ) φ_{i,r} = σ^{(i)} + β_i · N_{i-1,r}^{(r-i)/(r-i+1)} · Ψ_{i-1,r} N_{i,r} = |σ^{(i)}| + β_i · N_{i-1,r}^{(r-i)/(r-i+1)} Ψ_{i,r} = φ_{i,r} / N_{i,r} u = -α · Ψ_{r-1,r}

The controller is continuous everywhere except at σ = σ̇ = ... = σ^{(r-1)} = 0. In practice (with noise and discrete sampling), the control is always continuous.

Special cases

r=1: u = -α · sign(σ) r=2: u = -α · (σ̇+β|σ|^½sign(σ)) / (|σ̇|+β|σ|^½) (QC 2-SM, Eq 4.30)

Properties
  • |u| ≤ α (bounded by design)
  • Significantly reduced chattering compared to nested
  • Same convergence guarantees (Theorem 6.3)
  • Preferred by textbook authors (p.228)
Parameters

order : int Sliding mode order r (= relative degree). alpha : float Control magnitude bound. betas : list of float or None Gains [β_1, ..., β_{r-1}]. Length must be r-1. If None, uses [1.0, 2.0, 3.0, ...] (heuristic defaults). c : float Surface slope for default σ computation. g_sign : float Sign of input gain. Default -1 for σ=e systems. epsilon : float Small constant to avoid division by zero.

Source code in opensmc/controllers/hosmc.py
class QuasiContinuousHOSMC(Controller):
    """Quasi-continuous r-sliding mode controller (Levant 2005, Shtessel Sec 6.6.2).

    Recursive definition:
        φ_{0,r} = σ, N_{0,r} = |σ|, Ψ_{0,r} = sign(σ)
        φ_{i,r} = σ^{(i)} + β_i · N_{i-1,r}^{(r-i)/(r-i+1)} · Ψ_{i-1,r}
        N_{i,r} = |σ^{(i)}| + β_i · N_{i-1,r}^{(r-i)/(r-i+1)}
        Ψ_{i,r} = φ_{i,r} / N_{i,r}
        u = -α · Ψ_{r-1,r}

    The controller is continuous everywhere except at σ = σ̇ = ... = σ^{(r-1)} = 0.
    In practice (with noise and discrete sampling), the control is always continuous.

    Special cases:
        r=1: u = -α · sign(σ)
        r=2: u = -α · (σ̇+β|σ|^½sign(σ)) / (|σ̇|+β|σ|^½)  (QC 2-SM, Eq 4.30)

    Properties:
        - |u| ≤ α (bounded by design)
        - Significantly reduced chattering compared to nested
        - Same convergence guarantees (Theorem 6.3)
        - Preferred by textbook authors (p.228)

    Parameters
    ----------
    order : int
        Sliding mode order r (= relative degree).
    alpha : float
        Control magnitude bound.
    betas : list of float or None
        Gains [β_1, ..., β_{r-1}]. Length must be r-1.
        If None, uses [1.0, 2.0, 3.0, ...] (heuristic defaults).
    c : float
        Surface slope for default σ computation.
    g_sign : float
        Sign of input gain. Default -1 for σ=e systems.
    epsilon : float
        Small constant to avoid division by zero.
    """

    def __init__(self, order=2, alpha=10.0, betas=None, c=1.0,
                 g_sign=-1.0, epsilon=1e-10):
        super().__init__()
        if order < 1:
            raise ValueError(f"order must be >= 1, got {order}")
        self.order = order
        self.alpha = alpha
        self.c = c
        self.g_sign = np.sign(g_sign) if g_sign != 0 else -1.0
        self.epsilon = epsilon

        if betas is None:
            self.betas = [float(i) for i in range(1, order)]
        else:
            if len(betas) != order - 1:
                raise ValueError(
                    f"betas length {len(betas)} != order-1 = {order-1}")
            self.betas = list(betas)

    def _compute_psi(self, sigma_derivs):
        """Compute Ψ_{r-1,r} from the derivative vector."""
        r = self.order
        eps = self.epsilon

        # Level 0
        phi = sigma_derivs[0]
        N = np.abs(sigma_derivs[0])
        psi = phi / (N + eps)  # ≈ sign(σ)

        for i in range(1, r):
            exp = (r - i) / (r - i + 1)
            N_pow = (N + eps) ** exp

            phi = sigma_derivs[i] + self.betas[i-1] * N_pow * psi
            N = np.abs(sigma_derivs[i]) + self.betas[i-1] * N_pow
            psi = phi / (N + eps)

        return psi

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute quasi-continuous HOSMC control.

        Parameters
        ----------
        t : float
        x : ndarray — state vector
        xref : ndarray — reference
        sigma_derivs : list of float, optional
            [σ, σ̇, ..., σ^{(r-1)}]. If not provided, defaults to
            2nd-order system with σ = c·e + edot, σ̇ = c·edot.

        Returns
        -------
        u : float, bounded by |u| ≤ α
        info : dict with 's', 'psi'
        """
        sigma_derivs = kwargs.get('sigma_derivs', None)

        if sigma_derivs is None:
            e = xref[0] - x[0]
            edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]
            sigma = self.c * e + edot
            sigma_dot = self.c * edot
            sigma_derivs = [sigma, sigma_dot]

        psi = self._compute_psi(sigma_derivs)
        v = -self.alpha * psi
        u = v / self.g_sign

        info = {
            's': sigma_derivs[0],
            'psi': psi,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute quasi-continuous HOSMC control.

Parameters

t : float x : ndarray — state vector xref : ndarray — reference sigma_derivs : list of float, optional [σ, σ̇, ..., σ^{(r-1)}]. If not provided, defaults to 2nd-order system with σ = c·e + edot, σ̇ = c·edot.

Returns

u : float, bounded by |u| ≤ α info : dict with 's', 'psi'

Source code in opensmc/controllers/hosmc.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute quasi-continuous HOSMC control.

    Parameters
    ----------
    t : float
    x : ndarray — state vector
    xref : ndarray — reference
    sigma_derivs : list of float, optional
        [σ, σ̇, ..., σ^{(r-1)}]. If not provided, defaults to
        2nd-order system with σ = c·e + edot, σ̇ = c·edot.

    Returns
    -------
    u : float, bounded by |u| ≤ α
    info : dict with 's', 'psi'
    """
    sigma_derivs = kwargs.get('sigma_derivs', None)

    if sigma_derivs is None:
        e = xref[0] - x[0]
        edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]
        sigma = self.c * e + edot
        sigma_dot = self.c * edot
        sigma_derivs = [sigma, sigma_dot]

    psi = self._compute_psi(sigma_derivs)
    v = -self.alpha * psi
    u = v / self.g_sign

    info = {
        's': sigma_derivs[0],
        'psi': psi,
    }
    return u, info

TwistingSMC

Bases: Controller

Twisting 2-sliding mode controller (Levant 1993, Shtessel Eq 4.13).

u = -(r1 · sign(σ) + r2 · sign(σ̇))

For a 2nd-order system x = [x1, x2], with σ = c·e + edot (or any pluggable surface), and σ̇ available from the state.

Parameters

r1 : float Gain on sign(σ). Must satisfy r1 > r2 > 0. r2 : float Gain on sign(σ̇). Must satisfy r2 > 0. c : float Surface slope (used when no surface is provided): σ = c·e + edot. surface : SlidingSurface or None Optional pluggable surface for σ computation.

Source code in opensmc/controllers/twisting.py
class TwistingSMC(Controller):
    """Twisting 2-sliding mode controller (Levant 1993, Shtessel Eq 4.13).

    u = -(r1 · sign(σ) + r2 · sign(σ̇))

    For a 2nd-order system x = [x1, x2], with σ = c·e + edot (or any
    pluggable surface), and σ̇ available from the state.

    Parameters
    ----------
    r1 : float
        Gain on sign(σ). Must satisfy r1 > r2 > 0.
    r2 : float
        Gain on sign(σ̇). Must satisfy r2 > 0.
    c : float
        Surface slope (used when no surface is provided): σ = c·e + edot.
    surface : SlidingSurface or None
        Optional pluggable surface for σ computation.
    """

    def __init__(self, r1=10.0, r2=5.0, c=1.0, surface=None, g_sign=-1.0):
        super().__init__(surface=surface)
        if r2 <= 0:
            raise ValueError(f"r2 must be positive, got {r2}")
        if r1 <= r2:
            raise ValueError(f"r1 must be > r2, got r1={r1}, r2={r2}")
        self.r1 = r1
        self.r2 = r2
        self.c = c
        # g_sign: sign of the input gain in σ̈ = h + g·u.
        # For σ = e = xref - x on DoubleIntegrator (xddot=u):
        #   σ̈ = -u, so g = -1 → g_sign = -1 (default).
        # For σ = x on DoubleIntegrator: σ̈ = u, g = +1 → g_sign = +1.
        self.g_sign = np.sign(g_sign) if g_sign != 0 else -1.0

    def compute(self, t, x, xref, plant=None, **kwargs):
        """Compute twisting control.

        Parameters
        ----------
        t : float
            Current time.
        x : ndarray
            State vector. For 2nd-order: [x1, x2].
        xref : ndarray
            Reference. For 2nd-order: [xref1, xref2].
        plant : Plant or None
            Not used (model-free controller).

        Returns
        -------
        u : float
            Control signal.
        info : dict
            Contains 's' (σ), 'sdot' (σ̇), 'e', 'edot'.
        """
        # Tracking errors
        e = xref[0] - x[0]
        edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]

        # Sliding variable σ
        if self.surface is not None:
            sigma = self.surface.compute(e, edot, t=t, **kwargs)
        else:
            sigma = self.c * e + edot

        # σ̇ for linear surface σ = c·e + edot:
        # σ̇ = c·edot + eddot
        # For the double integrator eddot = xref_ddot - u,
        # but the twisting controller doesn't need the exact σ̇ —
        # it uses the state-based estimate.
        # For σ = c·e + edot, σ̇ ≈ c·edot + (acceleration difference)
        # In practice, with full-state measurement on 2nd-order systems,
        # σ̇ can be computed as c·edot (neglecting acceleration term
        # which is what the controller produces).
        #
        # However, the textbook formulation (Eq 4.11) assumes σ̈ = h + g·u,
        # and the controller acts on the (σ, σ̇) plane. For a double
        # integrator with σ = e (position error), σ̇ = edot.
        #
        # Standard approach: use σ̇ = c·edot for linear surface.
        sigma_dot = self.c * edot

        # Twisting control law (Eq 4.13)
        # The textbook assumes σ̈ = h + g·u with g > 0.
        # For g < 0 (e.g., σ=e on DI where σ̈=-u), we flip the sign.
        v = -(self.r1 * np.sign(sigma) + self.r2 * np.sign(sigma_dot))
        u = v / self.g_sign  # u = v/g_sign

        info = {
            's': sigma,
            'sdot': sigma_dot,
            'e': e,
            'edot': edot,
        }
        return u, info
Functions
compute(t, x, xref, plant=None, **kwargs)

Compute twisting control.

Parameters

t : float Current time. x : ndarray State vector. For 2nd-order: [x1, x2]. xref : ndarray Reference. For 2nd-order: [xref1, xref2]. plant : Plant or None Not used (model-free controller).

Returns

u : float Control signal. info : dict Contains 's' (σ), 'sdot' (σ̇), 'e', 'edot'.

Source code in opensmc/controllers/twisting.py
def compute(self, t, x, xref, plant=None, **kwargs):
    """Compute twisting control.

    Parameters
    ----------
    t : float
        Current time.
    x : ndarray
        State vector. For 2nd-order: [x1, x2].
    xref : ndarray
        Reference. For 2nd-order: [xref1, xref2].
    plant : Plant or None
        Not used (model-free controller).

    Returns
    -------
    u : float
        Control signal.
    info : dict
        Contains 's' (σ), 'sdot' (σ̇), 'e', 'edot'.
    """
    # Tracking errors
    e = xref[0] - x[0]
    edot = (xref[1] - x[1]) if len(xref) > 1 else -x[1]

    # Sliding variable σ
    if self.surface is not None:
        sigma = self.surface.compute(e, edot, t=t, **kwargs)
    else:
        sigma = self.c * e + edot

    # σ̇ for linear surface σ = c·e + edot:
    # σ̇ = c·edot + eddot
    # For the double integrator eddot = xref_ddot - u,
    # but the twisting controller doesn't need the exact σ̇ —
    # it uses the state-based estimate.
    # For σ = c·e + edot, σ̇ ≈ c·edot + (acceleration difference)
    # In practice, with full-state measurement on 2nd-order systems,
    # σ̇ can be computed as c·edot (neglecting acceleration term
    # which is what the controller produces).
    #
    # However, the textbook formulation (Eq 4.11) assumes σ̈ = h + g·u,
    # and the controller acts on the (σ, σ̇) plane. For a double
    # integrator with σ = e (position error), σ̇ = edot.
    #
    # Standard approach: use σ̇ = c·edot for linear surface.
    sigma_dot = self.c * edot

    # Twisting control law (Eq 4.13)
    # The textbook assumes σ̈ = h + g·u with g > 0.
    # For g < 0 (e.g., σ=e on DI where σ̈=-u), we flip the sign.
    v = -(self.r1 * np.sign(sigma) + self.r2 * np.sign(sigma_dot))
    u = v / self.g_sign  # u = v/g_sign

    info = {
        's': sigma,
        'sdot': sigma_dot,
        'e': e,
        'edot': edot,
    }
    return u, info