Automated Controller Design for the PMSM using Dynamic Mode Decomposition

This work presents a method to automatically generate a high performance controller for the permanent magnet synchronous motor (PMSM). The method consists of two components, a nominal system identification and a harmonic component identification. Both identification methods are based on dynamic mode decomposition (DMD). The nominal system identification is used to assign the feedback gaines by matching the desired closed-loop eigenvalues and eigenvectors. The harmonic system identification is used to generate vectors that are multiplied by a delay embedding of the current to predict harmonic components at the next time-step. The method is applied to two experimental test setups, one interior permanent magnet (IPM) and one surface mount permanent magnet (SPM) motor. It is shown that the automatically generated feedback controller is able to achieve a more precise transient response than the traditional rule-based PI controller. It is also shown the harmonic compensation method is able to reduce total demand distortion (TDD) on phase currents better than a traditional adaptive filter approach without the need for gain tuning. This work shows a novel approach to using DMD for the complete system identification of the PMSM, and lays the foundation for using DMD with delay embeddings to analyze and manipulate harmonic signals in a predictive way.


I. INTRODUCTION
T HIS work provides a generalizable methodology to automate the design of a current controller for the PMSM. The controller consists of two components, one nominal feedback controller, containing proportional and integral feedback, and the other, predictive current harmonic compensation. The system is first identified using techniques from DMD with some necessary modifications. Many theoretical papers have been written on the DMD algorithm for identifying linear representations from high dimensional and perhaps non-linear systems. However, this work addresses some of the practical concerns that must be dealt with when applying DMD to motor controls. This includes, removing the bias due to noisy measurements, handling the timevarying frequency of the system and actuation commands. The algorithmic implementation is provided on Github for the interested reader [1].
Traditional motor drives with parameter identification capabilities typically identify individual motor parameters such as stator resistance, d and q axis inductances and magnetic flux [2]. These parameters are used to determine proportional integral (PI) gains based on rules obtained from transfer function analysis. The problem with this approach is that there are many non-idealities that are not included in the individual parameter identification. This may include sensor measurement errors, delay caused by the inverter, dead-band time, difference between voltage command and actual phase voltage. With these non-idealities the drive manufacturers assume that the user will further tune the PI gains from the recommended values to get the desired performance.
The purpose of this work is to eliminate the need for finetuning controller gains and get predictable, high-level performance. The method is designed for field oriented control (FOC) of PMSM currents and could potentially be extended to speed control with some modifications. The methodology presented here identifies the system between the voltage commands provided by the controller and the current measurements fed back into the processor. In this way, the non-idealities are better accounted for in the lumped system identification than individual parameter identification. The nominal feedback controller design first identifies a discrete state-space model that is free from bias due to measurement noise [3]. The changing speed of the motor creates a timevarying system that cannot be identified by the standard DMD algorithm. To handle the changing speed during system identification, augmented states are created so that the timevarying component is extracted from the identified statetransition matrix. The identified system model is then used to set the closed-loop feedback gains by assigning the desired closed-loop eigenvalues and eigenvectors. In addition to the nominal feedback controller a predictive current harmonic compensation voltage is added to the final voltage command to reduce the overall TDD on the phase currents.
In PMSM machines, current harmonics lead to additional losses, non-sinusoidal back-emf and torque ripple [4] [5]. Current harmonics can be caused by non-uniform flux linkage distribution, current and position measurement errors, and insufficient inverter frequency for a given motor fundamental frequency. It is often desirable to eliminate the harmonics during the design of the motor, however this comes with increased manufacturing complexity and cost [6]. Also, eliminating harmonics through machine design does not account for harmonics caused by controller artifacts or measurement errors. This work develops a new methodology for harmonic identification and compensation based on DMD [7]. The proposed method has the advantages of not requiring any hand-tuning and it compensates harmonics in a predictive way. The proposed method is compared with a traditional harmonic compensation algorithm, adaptive selective harmonic elimination (ASHE), on two experimental motors [8].
The fundamental frequency of an electric motor is determined by the number of pole-pairs and the rated speed of the motor. If the fundamental frequency is too high for a given switching frequency, the inverter will not be able to produce a smooth sinusoid current. It is often desirable to have a high number of pole-pairs in the design of a PMSM machine, which increases the fundamental frequency. Increasing the number of pole-pairs increases torque density, reduces torque ripple and requires less magnetic material [9] [10]. In this work, it is shown that adding the proposed harmonic compensation method to a nominal feedback controller reduces the TDD at a given switching frequency.
DMD with delay embeddings is perfectly suited to identify PMSM current harmonics because it is able to decompose a signal into its most dominant modes, and determine their amplitudes and the frequencies at which they occur. The DMD algorithm is constrained in that it only works on systems with constant frequencies. This is fine for harmonic identification but the harmonic compensation must work at all operating speeds of the motor. In order to use DMD for harmonic compensation on the PMSM, an algorithm was developed to automatically generate a look-up table of the DMD compensation vectors for each operating frequency of the motor. The algorithm is based on the fact that the DMD eigenvectors produce a Fourier basis on sinusoidal signals. This paper is organized as follows. The first section covers works related to this one, including parameter identification methods for the PMSM and other applications of DMD to electrical systems. The DMD algorithm is described in section III. The nominal system identification is described in section IV and includes the necessary modifications to use FBDMD on a system with actuation commands and augmented states. Section V describes how to use the model identified in section IV to automatically generate the feedback gains. The experimental results for the nominal feedback controller are then shown in section V-B. Section VI describes the relationship between DMD and the discrete Fourier transform (DFT). The harmonic identification is described in section VII-A with experimental results comparing DMD to the DFT. The harmonic compensation is described in section VII-B. This section describes the algorithm that was developed to use the DMD for harmonic compensation at all operating speeds of the motor and is based on the similarities between DMD and the DFT as described in section VI.

II. RELEVANT WORKS
In [2] a comprehensive review of parameter identification techniques for the PMSM is given. Some involve online parameter identification to deal with changing parameter values for control and condition monitoring. Other methods attempt to identify the parameters without spinning the motor. The proposed method differs in that the complete system is identified as opposed to individual motor parameters. This allows for measurement errors and other non-idealities to be captured in the model of the system that is ultimately used for controller design.
The works most similar to this one are the recursive leastsquares (RLS) based finite-control-set model predictive control (FCS-MPC) of [11] and [12]. In [12], the authors extend the FCS-MPC method of [11] by including interlocking time of the inverter, and propose a DMD method for harmonic compensation. These approaches use data to identify seven models for each of the switching states of a two-stage inverter at a set operating point. The identified model is accurate at only that operating point so a learning rate is used to adapt the model to changing speed and other parameter variations. These approaches are limited in their practicality for the following reasons. First, the initialization of the models must be done with data from the motor operating at a constant speed. This is not possible in applications where the motor speed cannot be controlled by an external source. Another limitation of the approach is that the FCS scheme used does not allow for space-vector pulse-width modulation (SV-PWM). Only one of the seven voltages of the inverter can be applied at a given time-step making the current regulation performance quite poor. The data collection is also quite complex, requiring data to be collected for each of the seven switching states in order to create the seven models for prediction. The method proposed in this article requires only one model to be identified and allows for data collected at varying motor speeds. It also allows for traditional SV-PWM techniques to be used, greatly improving the regulation performance.
The proposed harmonic compensation technique in [12] uses a lifted dictionary of cosine and sine terms to identify and compensate harmonics in a predictive way. The lifted feature space of cosine and sine terms is intended to compensate the 6 th , 12 th and 18 th order harmonics of the electrical speed. Of the three models compared in the paper, the model with harmonic compensation performed the worst. Lifted features of cosine and sine terms are not robust to disturbances, which is why adaptive filters have traditionally been used for harmonic compensation [13]. Another limitation of this approach is that the harmonic amplitude estimate must be adapted as speed and current change. This is because the amplitude of current harmonics in PMSMs are determined by a combination of speed and current [14]. The delay embedding approach proposed here has several advantages over the method used in [12]. The harmonic amplitude information is contained within the delay embeddings of current measurements; thus the amplitude prediction does not need to adapt to changing speed or current. Also, with a sufficient embedding length, our approach can handle signals with a low signal-to-noise ratio.
For harmonic compensation in the PMSM, there are several works which use a Fourier decomposition to determine compensating current components to reduce torque ripple [15]- [17]. In [16] Espinosa et al. adaptively tune Fourier coefficients for periodic torque ripple minimization. In [18] a method to reduce torque ripple by updating the current reference is presented. This has the drawback in that the current controller is unable to compensate harmonics at high motor speed due to the bandwidth of the controller. In [19] the authors focus on eliminating the zero-axis current of a six-leg three-phase PMSM. In this way the harmonic suppression of the current is able to reduce torque ripple and peak currents.
Several works have used DMD to identify harmonics in power grids. In [20] the authors use Total DMD (TDMD) to identify harmonics on a power system in a way that is debiased towards noisy data. They show the method is able to more accurately identify harmonics than DFT based methods on simulated signals with noise and real-time measurements. In [21] the authors use DMD to estimate frequency and amplitudes of harmonics and compare with the DFT on simulated signals. The method presented here improves on previous harmonic compensation methods in that it is completely automated, requiring no human in the loop for gain tuning, and it compensates harmonics in a predictive way without needing to adapt to changing speed or current. To our knowledge, this is the first work to extend the use of DMD with delay embeddings beyond harmonic identification to harmonic compensation.

III. DYNAMIC MODE DECOMPOSITION
DMD was developed by the fluid dynamics community as a way to identify spacio-temporal structures from high dimensional data [7]. While motor drive systems are not high dimensional, in order to represent the harmonics of the PMSM in a linear way and reduce the effect of noise, the dimensionality is increased using delay embeddings of the current, as will be described later. DMD is an algorithm that identifies the best-fit linear dynamical system that advances measurements forward in time. DMD is able to identify the characteristic frequencies of dynamic systems as well as their growth and decay rates. Additionally, DMD is able to provide linear representations of some non-linear systems.
Many variants of DMD have been developed to handle different use cases. For instance, in [3] Dawson et al. show analytically that DMD is biased to noisy data and present methods for eliminating the bias, and in [22] DMD was extended to handle nonlinear systems by allowing for augmented states. Both of these techniques are applied in this work to get a reliable model of the PMSM system including harmonic components. For more information on different use cases of DMD with practical code implementations see the references by Kutz and Brunton [23] [24].

A. THE DMD ALGORITHM
The FBDMD method used by Dawson et al. in [3] was used in this work to eliminate the effect of noisy data. The FBDMD approach solves for the state-transition matrix in both the forward and backward direction and then averages the two matrices to get the unbiased estimate.
In order to identify sinusoidal components in a signal using DMD, delay embeddings must be used [25]. The lower the signal to noise ratio, the more delay embeddings are needed to extract the signal. To form the delay embeddings, the data is stacked in the so called Hankel matrix with the most recent time-step being at the top.
The data matrix Z is passed into Algorithm 1 as X. VOLUME 4, 2021 Algorithm 1 FBDMD by Dawson et al in [3] 1) Take the singular value decomposition (SVD) of data matrix X, letting X = SΣV * 2) Truncate SVD by only considering the first r modes.
The matrices containing the first r columns of S and V * , and the first r rows and columns of Σ are denoted as S r , V r and Σ r , respectively. 3) Project data X onto first r modes, i.e. X r = S r X 4) Shift projected data X r to get one time-step advanced projected data G 1 = X r (t 0 : t m−1 ) G 2 = X r (t 1 : t m ) 5) Using the FBDMD approach to remove bias solve for forward and backward state-transition matrix in projected space Given the initial conditions in a column vector z 0 , the future state at time-step k is determined using the DMD eigenvectors and eigenvalues as: where the columns of Φ contain the eigenvectors, individually denoted as φ φ φ j , for each mode j. The diagonal matrix Λ contains the eigenvalues, individually denoted as λ j . The inverse of the eigenvectors is denoted Γ, and the rows of Γ are denoted γ γ γ j . The signal can be represented as as sum of the individual modes as: where r denote number of modes used in the SVD truncation in step 2. of Algorithm 1. It is this ability to separate modes in the data that allows us to identify and compensate harmonics of a given frequency. As an example, given a multi-modal signal with a DC offset as: x(t) = 0.1 sin(300t) + sin(50t) + 1.5 The DMD algorithm is able to use four delay embeddings plus the current time-step to identify and separate the five modes of the system as shown in Fig. (1). In this work, only harmonics of a specific harmonic order and magnitude are extracted. The modal magnitudes are calculated as: only the first row of M is of interest since this row corresponds to stepping forward one time-step and the other rows are required only for the delayed embeddings. The actual magnitude of the sinusoid is two times the individual value in M due to the fact each magnitude is repeated because of the complex conjugate eigenvalues. Once the harmonic modes with significant amplitude are identified, a compensation vector is created for each operating speed which is then used in real-time to extract and predict harmonic components from the delay embedding of the signal.

IV. NOMINAL PMSM SYSTEM IDENTIFICATION
The proposed controller design requires that a complete statespace model first be identified from data. This model is then used to set the closed-loop feedback gains by assigning the closed-loop eigenvalues and eigenvectors. One of the difficulties with this approach is that DMD does not work on systems with changing frequencies since the system eigenvalues must remain constant. To deal with this, augmented states are created, for the system identification that include the currents multiplied by the motor speed. This eliminates the changing frequency in the model identified. The other difficulty is that noisy measurements will cause a bias on the identified model that cannot be eliminated by collecting more data [3]. To deal with this, the FBDMD algorithm, which eliminates the bias due to noisy measurements, is modified to handle a system with actuation commands and augmented states.

A. PMSM MODEL
The continuous time state-variable model of the PMSM in the dq coordinate frame can be represented in state-space form as: where R s is the stator resistance, ω e is the electrical speed in rad/s, λ f is the permanent magnet flux, and L q and L d are the q and d axis inductances, respectively. The system identified from data is actually the discretetime representation of (6). The discrete-time state-space matrices will be denoted as A d , B d and D d . The back-emf term caused by the flux is represented as a disturbance D d . The flux term D d is identified during the DMD system identification and compensated directly. That is to say, it does not play a role in determining the closed-loop feedback gains.

B. ARRANGING ACTUATION INPUTS FOR FBDMD
As discussed earlier, DMD has a bias towards noisy data that shifts the eigenvalues of the identified system. In order to get an unbiased estimate of the system, the model can be solved in both the forward and backwards direction. The forward direction is the standard solution, solving for the matrix that transitions from X to Y . The backward solution, solves for the matrix that transitions from Y to X. The inverse of the backwards estimate is then averaged with the forward estimate to get the unbiased estimate. Without the need for dimensionality reduction, DMD is equivalent to a least-squares solution. The forward estimate is calculated as: The backward estimate is calculated as: Combining the forward and backward estimations the unbiased estimate is: Typically, Y is equivalent to X but shifted one time-step into the future. However, the PMSM system identification includes augmented states, actuation commands and measured speed, which do not transition from one time-step to the next linearly. This causes a misidentification of the parameters if the data is stacked in the traditional way. To handle the non-linear transition of actuation commands and augmented states, only the d and q axis currents, are advanced forward in time in the Y data matrix. The augmented states, actuation commands and measured speed are kept at the same time step as in the X data matrix.
The DMD matrix should then solve as: The parameters of interest are only in the first to rows of A.
The augmented states, actuation commands and measured speed, solve to the identity matrix. The relevant parameters are extracted from A to form the discrete-time state-space model.

C. INITIALIZATION DATA
It is important that the system be excited in such a way that the resistance, inductance and flux parameters can be seen in the measured currents. To do this, the motor v q and v d voltage commands are actuated in such a way that the motor ramps up to close to the rated speed. This ensures that the back-emf produced by the magnetic flux is sufficiently high. A sinusoid signal with an amplitude of 10% the rated voltage is added to the ramping v q voltage to excite the q axis inductance. Similarly, v d is excited with a sinusoid with an amplitude of 20% the rated voltage to excite the d axis inductance. Both sinusoid signals were set to a frequency of 500 Hz. This frequency is large enough to generate some impedance in the inductance but not too large that the current measurements become too low. Once the desired speed is reached, v q is held constant and v d is held at negative 10% the rated voltage. The same procedure is repeated with v q going in the opposite direction. The initialization voltage wave-forms for the IPM are shown in Fig. 2a. The measured response of i d , i q and ω e parameters for the IPM can be seen in Fig. 2b. The above initialization routine takes approximately one minute and uses about 80 MB of memory with a samplingtime of 50 µs. The memory requirement for DMD and FBDMD is the same; however, the computation time for FBDMD is double that of the DMD algorithm. Using an Intel i7-7700 with clock at 2.80 GHZ the DMD solution computed offline takes 89 ms while the FBDMD solution takes 194 ms. The FBDMD approach could potentially be solved in an online fashion to adapt to changing motor parameters; however, a weighting parameter would need to be tuned. More research is needed to understand how to automate the selection of this weighting parameter for the online application of FBDMD to motor identification.

V. FEEDBACK CONTROLLER DESIGN
The state-space model identified in Sec. IV-B is used to design the feedback gains by assigning the closed-loop  2: Data collection for IPM motor system identification eigenvalues and eigenvectors. In this work, the eigenvectors and eigenvalues are assigned in a least-squares since. The interested reader may refer to the following references for additional eigenstructure assignment techniques [26] [27] [28].

A. METHODOLOGY
Once the A d , B d and D d system matrices of (12) have been identified, the system is converted to a form that includes integral error.
where ∆T is the inverter sampling time, u(t k ) contains the q and d axis voltages, i r contains the q and d axis current references and x i indicates the integral error terms. The integral error terms are calculated as: The modified system matrices with integral error are defined as: The closed-closed loop system is then calculated as: where K d is the feedback gain matrix. A desired transient response is determined in simulation and the simulated closed-loop system matrix is stored as A cl_s . The stability of the controller is guaranteed if the closed-loop eigenvalues of the nominal feedback controller are sufficiently less than one. Some margin should be given in case there are minor errors in the identified model.
The goal is to assign the feedback gain matrix, K d , so the actual and desired closed-loop systems match.
The feedback gain K d is then set such that the actual and desired closed-loop systems are equivalent in the leastsquares since. Plugging in (16) for A cl and solving for K d , the feedback gain is set as: This approach is similar to pole-placement, except the eigenvectors are set in addition to the poles.

B. EXPERIMENTAL RESULTS FOR TRANSIENT RESPONSE
The method used to determine controller gains outlined in section V-A was performed on both an SPM and an IPM in hardware. The motor parameters can be seen in Table 1. The experimental setup used a dSPACE MicroLabBox for control and data collection along with an APS two-stage inverter. The hardware used for testing can be seen in Fig. (3).  For testing purposes, three different transient responses with slow, medium and fast rise times were determined in simulation. The closed-loop system matrices for the three transient responses were stored in A cl_s . The proposed method was compared to the traditional rule based PI controller without further gain tuning. With the motor spinning freely a step in i q was given with i d commanded to remain at zero. The response to the step commands were recorded and can be seen in Fig. (4) for the IPM and Fig. (5) for the SPM.
As seen in the top rows of Figs. (4) and (5), the transient responses of the proposed approach match almost exactly to the desired transient response without any hand-tuning. The rule-based PI controller on the other hand, has significant overshoot. This discrepancy between desired transient response as determined by the closed-loop transfer function and the actual transient response in hardware is likely due to a combination of the following not accounted for in the closedloop transfer function: system delays, interlocking time and differences between the commanded voltage and actual voltage caused by space-vector modulation and inverter voltage drop. While there are methods available to identify these nonidealities such as described in [29] and [30], the complexity of the identification and controller design is greatly increased compared to the proposed method. These results show clearly the proposed method provides better results in automatic controller generation than the traditional rule-based PI controller.

VI. RELATIONSHIP OF DMD TO DISCRETE FOURIER TRANSFORM
The key differences between the Fourier transform and DMD can be summarized as follows. While the DMD algorithm identifies a linear model of a system from data, the Fourier transform identifies the frequencies contained within a signal. The DMD algorithm, since it operates on system data with multiple inputs, is able to distinguish the contribution of each signal onto other signals. The DMD algorithm can also identify the growth or decay rates of particular modes while the Fourier transform only considers constant or sinusoidal signals. Also, the modes identified by the DMD algorithm can be used to filter signals and the eigenvalues are used to predict the signal at future time-steps. With these differences in mind, it is important to understand the similarities when applied to sinusoidal signals.
The discrete Fourier transform (DFT) is defined as: The k th frequency is a relative frequency based on the number of time-steps N and the sampling period. The k th frequency is converted to rad/s as: where ω k is the k th frequency in rad/s, N are the number of time-steps in the signal and ∆T is the sampling period. To If the IDFT is written in matrix format the connection between DMD on a sinusoidal signal can be clearly seen. The IDFT can be calculated using a matrix dot product as: where F is the vector of frequency components. The k th column of the M F matrix is written as: In the right vector of (23) the relative frequency k is replaced with the actual frequency ω k in rad/s by solving (20) for k.
Similarly, the vector produced by DMD for a purely sinusoid signal with N − 1 delay embeddings is: where ω h is the frequency of the sinusoid in rad/s. The similarity between the IDFT and DMD with delay embeddings is made clear by comparing (23) with (25). The frequency of the DMD eigenvector ω h is determined precisely by the SVD truncation, while the frequency of the IDFT ω k is discretized based on the sampling time and the number of samples. The DMD algorithm is able to determine an exact

VII. IDENTIFYING AND COMPENSATING MOTOR HARMONICS
The harmonic analysis is divided into two parts, identification and compensation vector generation. The first part, requires the motor be run at constant speed so that the most dominant harmonic components can be identified. The second part, uses the previously identified harmonics to automatically generate vectors that are used in a look-up table to compensate harmonics at all operating speeds of the motor. An algorithm was developed to generate these compensation vectors at all operating speeds of the motor, and is based on the concept of a Fourier basis. These vectors are what the DMD algorithm would have produced on sinusoidal signals generated at each operating frequency of the motor, but without the need for more data. The harmonic analysis was performed on both the SPM and IPM test beds. DMD is compared with the DFT for harmonic identification on both motors. For harmonic compensation the DMD based approach is compared with a traditional adaptive filtering approach called ASHE. While fairly robust, the ASHE method takes time to converge after transients and requires parameter tuning.

A. MOTOR HARMONIC IDENTIFICATION
In order to understand the effect of the number of delay embeddings in estimating current harmonics, DMD, with 49 and 299 delay embeddings, was compared in predicting harmonic components at future time-steps. Fig. (6) shows the comparison of DMD with 49 delay embeddings and 299 delay embeddings on the IPM. The top plot in Fig.  (6) shows the uncompensated current prediction using DMD with 49 delay embeddings. The bottom plot of Fig. (6) shows the prediction using 299 delay embeddings. Both prediction methods used eight modes in the SVD truncation. Clearly, the 299 delay embedding is able to predict more accurately than the 49 delay embedding. The MAE for predicting harmonics using 49 and 299 delay embeddings was 6.6 mA and 2.7 mA, respectively. With 299 delays used in the embedding the current harmonic prediction takes 18µs. This is well within the half switching period limit of 25µs. The memory requirement for the harmonic identification for the q and daxis currents is about 640 kB for two seconds of data.
While the dSPACE MicroLabBox used for testing has a more powerful CPU than a typical digital signal processor (DSP), the DMD with delay embedding approach is computed most efficiently with parallel processing devices such as field programmable gate arrays (FPGAs), graphical processing units (GPUs) or other chips designed for matrix multiplication. FPGAs have been used in motor control for many years, and GPUs continue to decrease in price for a given throughput. Thus, the computation time for DMD with delay embeddings is not seen as prohibitive for real-world applications. It should also be noted that the computation time is determined by the number of delay embeddings and not the number of harmonics compensated, unlike traditional approaches. Given a sufficient number of harmonics, the delay-embedding approach would be faster than traditional methods, even on a serial processor. The harmonic identification was done by running the motors at close to the rated speed and recording the i q and i d currents. The currents were stacked into a Hankel matrix as in (1) with 299 delay embeddings and passed through the DMD algorithm. The most dominant modes identified by DMD will not necessarily be harmonic multiples of the electrical speed. Only frequencies that are harmonics of the electrical speed, and whose magnitudes are greater than a threshold value, are considered for compensation. To check if an identified harmonic is a harmonic multiple of the electrical speed the DMD eigenvalues must be converted to continuous time.
where ∆T is the sampling period of the inverter. The imaginary component of the continuous time eigenvalue, λ i_c , indicates the harmonic frequency in rad/s. The modulus of λ i_c is taken with the electrical speed to identify those modes for which frequencies are within 5% of a speed harmonic. The magnitudes of these modes are calculated using (5) and only those over a set threshold are considered for compensation. A comparison was done between DMD with delay embeddings and the FFT for harmonic identification on both the SPM and IPM motors. The harmonic identification for the SPM is shown in Fig. (7). The FFT frequency spectrum is shown in blue and the most dominant DMD modes are shown For the SPM, current harmonics with amplitudes greater than 5e−3 were identified for compensation. With this threshold amplitude, the 6 th order harmonic was identified to be compensated on both the d and q-axis currents. Fig. (8b) shows the comparison of the FFT with DMD with delay embeddings for current harmonic identification on the IPM. For the q-axis current, no harmonics were identified for compensation, as the lowest harmonic frequency was of the 70 th order, with very small magnitudes as shown in Fig (8a). The IPM current harmonics on the d-axis were found to be the 6 th , 12 th , 18 th and 48 th order as shown in Fig. (8b). The 48 th order harmonic was not considered for compensation as it had a magnitude below 1e−3.
The experimental results show DMD with delay embeddings produces a more distinctive representation of the most dominant harmonics than the FFT. Harmonic identification with the FFT, suffers from the picket fence effect, which causes the energy of the signal to be spread out across

B. HARMONIC COMPENSATION
The current harmonic compensation is done by predicting the harmonics at the next time-step using the DMD eigenvectors and eigenvalues and applying a corresponding compensation voltage to cancel the predicted harmonic. Only harmonic modes identified to be harmonic multiples of the electrical speed and of sufficient magnitude are included in the compensation vector. To predict the harmonic component at the next time step, the dot product is taken between the delay embedding of the current measurements and the compensation vector. i where I d contains the delay embeddings of current measurements. In this work the most recent current measurement is stored in the first element of the delay embedding vector and the n th previous measurement is stored in the last element. The compensation vector is taken from the first row of A h as this predicts the harmonic at the next time-step.
A h is formed by the dot product of Φ c , Λ c and Γ c which are the DMD eigenvectors, eigenvalues and inverse eigenvectors that have been identified for compensation, respectively. This compensation vector is only accurate at the speed at which the harmonics were originally identified. The next section shows the algorithm that was developed to automatically generate these vectors at all operating speeds of the motor without the need to collect more data.

1) Automatically generating compensation vectors for all operating speeds
The only inputs needed to generate the compensation vector are the number of delay embeddings and the harmonic frequencies to be compensated. A compensation vector is generated for each speed in one rad/s increments and stored in a lookup table. The first step in creating the compensation vector is to determine the continuous-time eigenvalue for each order harmonic to be compensated. Each continuoustime eigenvalue is determined by where o is the harmonic order and ω e is the motor electrical speed. The harmonic orders to be compensated were identified in Sec. VII. Each eigenvalue has a respective complex conjugate. The eigenvalues are stacked into a diagonal matrix.
where j is the number of harmonics to compensate. They are then converted to discrete-time.
The DMD eigenvectors are constructed as in (25), where for each operating electrical speed the DMD vector is created with ω h = λ c . Each eigenvector has a respective complex conjugate vector. The complete set of eigenvectors is then where there is a complex-conjugate pair of eigenvectors for each of the j harmonics identified for compensation. The current measurements have a DC component that must be extracted from the harmonic prediction at the next timestep. To do this, an eigenvector representing the DC mode is concatenated to Φ before obtaining the inverse of the eigenvectors. The DC eigenvector is a vector of constants VOLUME 4, 2021 with magnitude equal to one and length N + 1. The inverse of the eigenvectors is then calculated as: Finally the compensation matrix can be calculated as Only the top row of A h is important as it corresponds to predicting the harmonics at the next time-step. A h is calculated for all operating speeds of the motor and the top row of each A h is stored in a lookup table as the compensation vector for a given speed. The potential uses of the proposed method extend well beyond motor controls. The DMD eigenvectors are able extract all relevant information about sinusoidal signals. Some potential use cases are listed here: • The vectors can be used as filters to remove the DC or harmonic components of a signal in real-time. • Harmonic components can be predicted multiple timesteps in the future. • The amplitude of a sinusoid can be determined from only a few measurements of a signal. • The phase angle of a signal can be determined in a predictive way. For example, phase-locked loops generate an output signal whose phase is related to the phase of an input signal. They are used in radio and telecommunications as well as power systems and micro-grids. The DMD vectors with delay embeddings of the signal can output the phase information in a predictive way, potentially having an advantage over phaselocked-loops which take time to converge. The approach is also highly robust to noise when using DMD vectors with a large number of delay embeddings. A simulation was done using DMD vectors with 299 delay embeddings to predict the phase angle of a sinusoid. Gaussian noise was added to a sinusoid signal with amplitude of one. The Gaussian noise had an average magnitude three times that of the sinusoid. The predicted phase angle is shown in purple in Fig. 9 and the actual phase angle is shown in yellow. The original sinusoid is shown in red and the noisy signal used for predicting the phase is shown in blue. As can be seen in the figure, even with a high level of noise, DMD is able to accurately predict the phase angle. This makes a strong case for considering DMD with delay embeddings as a new method to replace PLLs for sinusoid signals with low signal to noise ratios.

2) Compensating harmonics
The harmonic prediction at the next time-step is calculated as . . .
where i h (t k+1 ) is the predicted harmonic current for the next time-step, i(t k ) is the current measurement at the current  It is important to add the predicted harmonic component to the current measurements for the delay embedding because once the harmonic is compensated, the harmonic information would otherwise be lost in the embedding. In order to avoid instabilities during large disturbances, no compensation is done if the actual current is not within two percent of the reference current and if the predicted harmonic magnitude is greater than four times the originally identified harmonic magnitude.
The proposed method uses a feed-forward gain to go from harmonic current prediction to compensation voltage. This differs from other harmonic compensation methods which subtract the estimated harmonic from the current reference and allow the feedback controller to eliminate the harmonic. The feedback controller is however, limited by its bandwidth and cannot compensate higher order harmonics.
The feed-forward gain is calculated using the system matrices determined from data in (12).The feed-forward gain is calculated as: Only the diagonal components of H are used to compute the compensation voltage. The harmonic compensation voltages for the d and q-axis are calculated The compensation voltages v hq and v h d are added to the voltage commands generated by the nominal feedback controller. The block diagram of complete controller is shown in Fig.  (10).

3) Traditional Harmonic Compensation
The ASHE approach used for benchmark comparison is based on traditional signal processing has been used in a vari- The drawbacks of the method are that it requires parameter tuning and takes time to converge after transients. The DMD approach developed in this work, on the other hand, does not require parameter tuning and does not need time to adapt after transients. Fig. (11) shows ASHE on a simulated signal. The signal being tracked is shown in red. The weighted cosine and sine signals that are summed to give the prediction are shown in blue and green, respectively. Their sum, which is the prediction, is shown in cyan. It is seen that by summing the weighted sinusoidals we get a new sinusoid of a different phase and amplitude but the same frequency. Over-time this sinusoid converges to the sinusoid in the input signal.

4) Harmonic Compensation Experimental Results
The harmonic compensation uses a vector dot product which can be computed efficiently on an FPGA or GPU device. In these experiments, a dSPACE MicroLabBox was used with only the CPU activated. The run-time of the control loop without harmonic compensation was 10.8 µs. With a delay embedding of 299 the compensation had a run-time of 18 µs. This is within acceptable limits for a 20 kHz switching frequency. The approach presented here was successfully able to compensate the desired harmonics on both the SPM and the IPM. For the SPM, it can be seen in Fig. (12) that the 6 th order harmonic has been eliminated. For the IPM, only the d axis   Fig. (13b), much of the harmonic energy of the signal has been eliminated, including the 6 th order harmonic. The TDD reduction with a nominal current of one amp can be seen in Table 2. Ideally, the delay embeddings would contain multiple cycles of the sinusoid signal; however, if the fundamental frequency of the motor is too low this may not be the case. In this case the harmonic prediction will not be as robust to noise. However, for low fundamental frequencies, the feedback controller is likely to have enough bandwidth to eliminate harmonics without the need for additional compensation.
In the experiments it can be seen that the uncompensated harmonic amplitude of the IPM was about ten times less than that of the SPM. This is because the fundamental frequency for the IPM was 200 rad/s compared to 1200 rad/s for the SPM. The PI controller at 20 kHz has a high enough bandwidth to eliminate most of the harmonics at such a low fundamental frequency.
The predictive approach presented here was compared to ASHE during transient conditions. While the ASHE algorithm performs well in the steady state, it takes time to adapt to transients in current and speed. Fig. (14) shows the comparison of ASHE with our DMD based approach after a step in the i q current. It can be seen in the plot of the i q current error that our predictive approach compensates harmonics just after the transient (red), while the ASHE algorithm takes time to adapt (yellow). A comparison was also done with the proposed method and the ASHE method with a varying speed. The current was held constant at 1A while the speed was allowed to ramp up freely. The DMD method was able to compensate better during the speed change as shown in Fig. (15). From the figure it can be seen that the i q current error from the DMD compensation (green) has a lower harmonic amplitude than the ASHE method. Both the ASHE and DMD compensation methods reduced the mean absolute error (MAE) between the current reference and the actual current compared to the uncompensated PI controller. The proposed DMD method had the greatest reduction in MAE. The MAE for the uncompensated, ASHE compensated and DMD compensated were 0.0337, 0.028 and 0.0172 respectively. This work presents a novel method to automate the controller design of the PMSM, providing high-level performance without the need for human intervention. The controller design, based on the DMD algorithm, automatically generates nominal feedback gains and provides predictive harmonic compensation.
A PMSM system model is first identified using DMD and then used to set the feedback gains by matching a predetermined closed-loop set of eigenvalues and eigenvectors. The method is able to identify a complete model for the PMSM from control inputs to measured current values. This is an improvement to traditional approaches in that many of the non-idealities of the system such as, inverter delay, inverter voltage drop, wiring and measurement errors are included in the model during identification. The approach  is able to identify the PMSM model even with changing speed. In addition, the bias from noisy data is removed by extending the FBDMD algorithm to a system with actuation commands and augmented states. The feedback gains were set by assigning the closed-loop eigenstructure of the system by matching a desired set of eigenvalues and eigenvectors in a least-square manor. On low dimensional systems such as this one, the results are similar to pole-placement. However, the authors feel it is important to make the connection between DMD and eigenstrucure assignment since DMD identifies the most dominant eigenvalues and eigenvectors of a system. The approach was shown to give a more predictable transient response than a rule-based PI controller on two experimental motors.
In addition to the nominal feedback controller, DMD was also used to identify and compensate PMSM current har- monics in a linear predictive way. Having smooth sinusoid currents is critical for efficient energy conversion, smooth mechanical operation and long-life of the equipment. It was shown through experiments that the proposed method is able to identify dominant harmonic frequencies and amplitudes with more clarity than FFT based methods. The harmonic compensation method developed in this work has been shown to reduce the phase TDD for a given switching frequency. To our knowledge, this is the first extension of DMD with delay embeddings beyond harmonic identification to harmonic compensation. A method based on the concept of a Fourier basis was developed to generate the DMD compensation vectors at all operating frequencies of the motor. The experiments show that the proposed method performs better than the traditional ASHE harmonic compensation method under transient conditions.