001// Copyright (c) FIRST and other WPILib contributors.
002// Open Source Software; you can modify and/or share it under the terms of
003// the WPILib BSD license file in the root directory of this project.
004
005package org.wpilib.math.system;
006
007import java.util.Arrays;
008import org.ejml.simple.SimpleMatrix;
009import org.wpilib.math.linalg.Matrix;
010import org.wpilib.math.numbers.N1;
011import org.wpilib.math.system.proto.LinearSystemProto;
012import org.wpilib.math.system.struct.LinearSystemStruct;
013import org.wpilib.math.util.Nat;
014import org.wpilib.math.util.Num;
015import org.wpilib.util.protobuf.Protobuf;
016import org.wpilib.util.protobuf.ProtobufSerializable;
017import org.wpilib.util.struct.Struct;
018import org.wpilib.util.struct.StructSerializable;
019
020/**
021 * A plant defined using state-space notation.
022 *
023 * <p>A plant is a mathematical model of a system's dynamics.
024 *
025 * <p>For more on the underlying math, read
026 * https://file.tavsys.net/control/controls-engineering-in-frc.pdf.
027 *
028 * @param <States> Number of states.
029 * @param <Inputs> Number of inputs.
030 * @param <Outputs> Number of outputs.
031 */
032public class LinearSystem<States extends Num, Inputs extends Num, Outputs extends Num>
033    implements ProtobufSerializable, StructSerializable {
034  /** Continuous system matrix. */
035  private final Matrix<States, States> m_A;
036
037  /** Continuous input matrix. */
038  private final Matrix<States, Inputs> m_B;
039
040  /** Output matrix. */
041  private final Matrix<Outputs, States> m_C;
042
043  /** Feedthrough matrix. */
044  private final Matrix<Outputs, Inputs> m_D;
045
046  /**
047   * Construct a new LinearSystem from the four system matrices.
048   *
049   * @param A The system matrix A.
050   * @param B The input matrix B.
051   * @param C The output matrix C.
052   * @param D The feedthrough matrix D.
053   * @throws IllegalArgumentException if any matrix element isn't finite.
054   */
055  public LinearSystem(
056      Matrix<States, States> A,
057      Matrix<States, Inputs> B,
058      Matrix<Outputs, States> C,
059      Matrix<Outputs, Inputs> D) {
060    for (int row = 0; row < A.getNumRows(); ++row) {
061      for (int col = 0; col < A.getNumCols(); ++col) {
062        if (!Double.isFinite(A.get(row, col))) {
063          throw new IllegalArgumentException(
064              "Elements of A aren't finite. This is usually due to model implementation errors.");
065        }
066      }
067    }
068    for (int row = 0; row < B.getNumRows(); ++row) {
069      for (int col = 0; col < B.getNumCols(); ++col) {
070        if (!Double.isFinite(B.get(row, col))) {
071          throw new IllegalArgumentException(
072              "Elements of B aren't finite. This is usually due to model implementation errors.");
073        }
074      }
075    }
076    for (int row = 0; row < C.getNumRows(); ++row) {
077      for (int col = 0; col < C.getNumCols(); ++col) {
078        if (!Double.isFinite(C.get(row, col))) {
079          throw new IllegalArgumentException(
080              "Elements of C aren't finite. This is usually due to model implementation errors.");
081        }
082      }
083    }
084    for (int row = 0; row < D.getNumRows(); ++row) {
085      for (int col = 0; col < D.getNumCols(); ++col) {
086        if (!Double.isFinite(D.get(row, col))) {
087          throw new IllegalArgumentException(
088              "Elements of D aren't finite. This is usually due to model implementation errors.");
089        }
090      }
091    }
092
093    this.m_A = A;
094    this.m_B = B;
095    this.m_C = C;
096    this.m_D = D;
097  }
098
099  /**
100   * Returns the system matrix A.
101   *
102   * @return the system matrix A.
103   */
104  public Matrix<States, States> getA() {
105    return m_A;
106  }
107
108  /**
109   * Returns an element of the system matrix A.
110   *
111   * @param row Row of A.
112   * @param col Column of A.
113   * @return the system matrix A at (i, j).
114   */
115  public double getA(int row, int col) {
116    return m_A.get(row, col);
117  }
118
119  /**
120   * Returns the input matrix B.
121   *
122   * @return the input matrix B.
123   */
124  public Matrix<States, Inputs> getB() {
125    return m_B;
126  }
127
128  /**
129   * Returns an element of the input matrix B.
130   *
131   * @param row Row of B.
132   * @param col Column of B.
133   * @return The value of the input matrix B at (i, j).
134   */
135  public double getB(int row, int col) {
136    return m_B.get(row, col);
137  }
138
139  /**
140   * Returns the output matrix C.
141   *
142   * @return Output matrix C.
143   */
144  public Matrix<Outputs, States> getC() {
145    return m_C;
146  }
147
148  /**
149   * Returns an element of the output matrix C.
150   *
151   * @param row Row of C.
152   * @param col Column of C.
153   * @return the double value of C at the given position.
154   */
155  public double getC(int row, int col) {
156    return m_C.get(row, col);
157  }
158
159  /**
160   * Returns the feedthrough matrix D.
161   *
162   * @return the feedthrough matrix D.
163   */
164  public Matrix<Outputs, Inputs> getD() {
165    return m_D;
166  }
167
168  /**
169   * Returns an element of the feedthrough matrix D.
170   *
171   * @param row Row of D.
172   * @param col Column of D.
173   * @return The feedthrough matrix D at (i, j).
174   */
175  public double getD(int row, int col) {
176    return m_D.get(row, col);
177  }
178
179  /**
180   * Computes the new x given the old x and the control input.
181   *
182   * <p>This is used by state observers directly to run updates based on state estimate.
183   *
184   * @param x The current state.
185   * @param clampedU The control input.
186   * @param dt Timestep for model update in seconds.
187   * @return the updated x.
188   */
189  public Matrix<States, N1> calculateX(
190      Matrix<States, N1> x, Matrix<Inputs, N1> clampedU, double dt) {
191    var discABpair = Discretization.discretizeAB(m_A, m_B, dt);
192
193    return discABpair.getFirst().times(x).plus(discABpair.getSecond().times(clampedU));
194  }
195
196  /**
197   * Computes the new y given the control input.
198   *
199   * <p>This is used by state observers directly to run updates based on state estimate.
200   *
201   * @param x The current state.
202   * @param clampedU The control input.
203   * @return the updated output matrix Y.
204   */
205  public Matrix<Outputs, N1> calculateY(Matrix<States, N1> x, Matrix<Inputs, N1> clampedU) {
206    return m_C.times(x).plus(m_D.times(clampedU));
207  }
208
209  /**
210   * Returns the LinearSystem with the outputs listed in outputIndices.
211   *
212   * <p>This is used by state observers such as the Kalman Filter.
213   *
214   * @param outputIndices the list of output indices to include in the sliced system.
215   * @return the sliced LinearSystem with outputs set to row vectors of LinearSystem.
216   * @throws IllegalArgumentException if any outputIndices are outside the range of system outputs.
217   * @throws IllegalArgumentException if number of outputIndices exceeds the number of system
218   *     outputs.
219   * @throws IllegalArgumentException if duplication exists in outputIndices.
220   */
221  public LinearSystem<States, Inputs, ? extends Num> slice(int... outputIndices) {
222    for (int index : outputIndices) {
223      if (index < 0 || index >= m_C.getNumRows()) {
224        throw new IllegalArgumentException(
225            "Output indices out of range. This is usually due to model implementation errors.");
226      }
227    }
228
229    if (outputIndices.length >= m_C.getNumRows()) {
230      throw new IllegalArgumentException(
231          "More outputs requested than available. This is usually due to model implementation "
232              + "errors.");
233    }
234
235    int[] outputIndicesList = Arrays.stream(outputIndices).distinct().sorted().toArray();
236
237    if (outputIndices.length != outputIndicesList.length) {
238      throw new IllegalArgumentException(
239          "Duplicate indices exist.  This is usually due to model implementation " + "errors.");
240    }
241
242    SimpleMatrix new_C_Storage = new SimpleMatrix(outputIndices.length, m_C.getNumCols());
243    SimpleMatrix new_D_Storage = new SimpleMatrix(outputIndices.length, m_D.getNumCols());
244    int row = 0;
245    for (var index : outputIndicesList) {
246      new_C_Storage.setRow(row, 0, m_C.extractRowVector(index).getData());
247      new_D_Storage.setRow(row, 0, m_D.extractRowVector(index).getData());
248      row++;
249    }
250
251    return new LinearSystem<>(m_A, m_B, new Matrix<>(new_C_Storage), new Matrix<>(new_D_Storage));
252  }
253
254  @Override
255  public String toString() {
256    return String.format(
257        "Linear System: A\n%s\n\nB:\n%s\n\nC:\n%s\n\nD:\n%s\n",
258        m_A.toString(), m_B.toString(), m_C.toString(), m_D.toString());
259  }
260
261  /**
262   * Creates an implementation of the {@link Protobuf} interface for linear systems.
263   *
264   * @param <States> The number of states of the linear systems this serializer processes.
265   * @param <Inputs> The number of inputs of the linear systems this serializer processes.
266   * @param <Outputs> The number of outputs of the linear systems this serializer processes.
267   * @param states The number of states of the linear systems this serializer processes.
268   * @param inputs The number of inputs of the linear systems this serializer processes.
269   * @param outputs The number of outputs of the linear systems this serializer processes.
270   * @return The protobuf implementation.
271   */
272  public static <States extends Num, Inputs extends Num, Outputs extends Num>
273      LinearSystemProto<States, Inputs, Outputs> getProto(
274          Nat<States> states, Nat<Inputs> inputs, Nat<Outputs> outputs) {
275    return new LinearSystemProto<>(states, inputs, outputs);
276  }
277
278  /**
279   * Creates an implementation of the {@link Struct} interface for linear systems.
280   *
281   * @param <States> The number of states of the linear systems this serializer processes.
282   * @param <Inputs> The number of inputs of the linear systems this serializer processes.
283   * @param <Outputs> The number of outputs of the linear systems this serializer processes.
284   * @param states The number of states of the linear systems this serializer processes.
285   * @param inputs The number of inputs of the linear systems this serializer processes.
286   * @param outputs The number of outputs of the linear systems this serializer processes.
287   * @return The struct implementation.
288   */
289  public static <States extends Num, Inputs extends Num, Outputs extends Num>
290      LinearSystemStruct<States, Inputs, Outputs> getStruct(
291          Nat<States> states, Nat<Inputs> inputs, Nat<Outputs> outputs) {
292    return new LinearSystemStruct<>(states, inputs, outputs);
293  }
294}