In the realm of precision motion control, the harmonic drive gear stands as a pivotal component, renowned for its exceptional torque density, minimal backlash, compact form factor, and high positional accuracy. These attributes have cemented its role in demanding applications across aerospace, robotics, and defense systems. As these systems evolve towards greater integration, modular harmonic drive gear assemblies—incorporating the motor, gearbox, sensors, and electronics into a single unit—have become increasingly prevalent. For engineers and researchers, the primary concern shifts from the isolated characteristics of the harmonic drive gear itself to the overall servo performance of the system utilizing this assembly. Consequently, developing a holistic, system-level model that accurately captures the dynamic behavior is paramount for effective controller design, performance simulation, and predictability. In this work, I present a detailed methodology for modeling such integrated precision harmonic drive gear systems. This approach synthesizes linear dynamics with dominant nonlinearities, establishes a practical simulation framework, and introduces an efficient parameter identification strategy. The ultimate goal is to furnish a model that is both physically insightful for engineering application and sufficiently accurate for performance prediction.
The core of any dynamic model for a harmonic drive gear system lies in representing the mechanical power transmission. A common and effective simplification is to treat the system as a two-mass model interconnected by a torsional spring, which represents the effective stiffness of the harmonic drive gear transmission. The first mass corresponds to the motor rotor and the input side of the gear (wave generator or flexspline, depending on configuration), while the second mass represents the load inertia reflected through the gear reduction. The governing equations of motion for this linear core are derived from torque balance. Let $J_m$ and $J_l$ denote the motor and load moments of inertia, respectively. The motor output torque is $T_m$, and the angular displacements are $\theta_m$ and $\theta_l$. The torque transmitted through the harmonic drive gear, $T_s$, is a function of the twist across its internal stiffness $K_s$ and the gear reduction ratio $K_w$. Disturbance torques on the motor and load sides are $T_{md}$ and $T_{ld}$.
The equations are:
$$ J_m \ddot{\theta}_m = T_m – T_{md} – T_s $$
$$ J_l \ddot{\theta}_l = T’_s – T_{ld} $$
$$ T_s = K_s (K_w \theta_m – \theta_l) $$
$$ T’_s = T_s / K_w $$
For the actuation, a DC servo motor model is typically employed. The electrical circuit equation and the torque generation equation are:
$$ T_m = K_m i $$
$$ u = iR + L \frac{di}{dt} + K_e \dot{\theta}_m $$
Here, $u$ is the terminal voltage, $i$ is the armature current, $R$ and $L$ are the armature resistance and inductance, $K_m$ is the torque constant, and $K_e$ is the back-EMF constant. Often, the inductance $L$ is negligible for dynamic analysis within the servo bandwidth. Combining these sets of equations and applying the Laplace transform yields the linear transfer function foundation of the harmonic drive gear system. However, this linear model is insufficient to capture the rich, nonlinear behavior that fundamentally influences precision performance.
The fidelity of a harmonic drive gear system model hinges on the incorporation of key nonlinearities. First, the drive amplifier itself imposes current saturation limits, which can be modeled with a simple piecewise function: $I(t) = I_{max}$ if $I(t) > I_{max}$, $I(t) = I_{min}$ if $I(t) < I_{min}$, otherwise $I(t) = I(t)$. Second, the feedback sensors, typically optical encoders, introduce quantization. The angular resolution is $\delta\theta = 2\pi / (N \cdot n_e)$, where $N$ is the interpolation factor and $n_e$ is the line count. Velocity derived via direct differentiation has a resolution of $\delta\omega = 2\delta\theta / t_s$, with $t_s$ as the sampling interval. Third, and most significantly, are the friction, damping, and load-dependent disturbances. We consolidate these into equivalent disturbance torques $T_{md}$ and $T_{ld}$. For the motor side, $T_{md}$ often primarily represents friction ($T_{mf}$). The load side disturbance is more complex: $T_{ld} = T_{lf} + T_{nb}$, encompassing load-side friction/damping and any mass unbalance torque.
Friction is notoriously nonlinear. A versatile model like the Stribeck curve effectively captures the transition from static to Coulomb to viscous friction. The friction torque $T_f$ as a function of velocity $\dot{\theta}$ and the net motor torque $T_m$ can be described piecewise:
$$
T_f(\dot{\theta}, T_m) =
\begin{cases}
T_m, & \text{if } \dot{\theta}=0 \text{ and } T_s^- < T_m < T_s^+ \\
T_s^+, & \text{if } \dot{\theta}=0 \text{ and } T_m \ge T_s^+ \\
T_s^-, & \text{if } \text{if } \dot{\theta}=0 \text{ and } T_m \le T_s^- \\
T^+_{Stribeck}(\dot{\theta}) = T_C^+ + (T_s^+ – T_C^+) e^{-(\dot{\theta}/\Omega^+)^\delta} + B^+\dot{\theta}, & \text{if } \dot{\theta} > 0 \\
T^-_{Stribeck}(\dot{\theta}) = T_C^- + (T_s^- – T_C^-) e^{-(\dot{\theta}/\Omega^-)^\delta} + B^-\dot{\theta}, & \text{if } \dot{\theta} < 0
\end{cases}
$$
Here, $T_s^+$, $T_s^-$ are static friction levels, $T_C^+$, $T_C^-$ are Coulomb friction levels, $\Omega$ is the Stribeck velocity, $\delta$ is an empirical exponent, and $B$ is the viscous damping coefficient.
Mass unbalance torque, relevant for rotating loads not perfectly balanced, is sinusoidal: $T_{nb} = k m g \rho \sin(\alpha_0 + \theta_l)$. The parameters are the equivalent unbalanced mass $m$, gravitational acceleration $g$, a proportionality factor $k$, eccentricity $\rho$, load angle $\theta_l$, and initial phase $\alpha_0$.
The fourth critical nonlinearity is inherent to the harmonic drive gear itself: its meshing stiffness and hysteresis loss. The torque-transmission relationship is not linear and exhibits hysteresis, leading to a stiffness that varies with torque level and direction. This can be approximated by a piecewise linear stiffness model. Let the transmitted load-side torque magnitude be $|T’_s|$. The effective stiffness $K_s$ can be defined as:
$$
K_s =
\begin{cases}
K_1, & \text{if } |T’_s| \le T_1 \\
K_2, & \text{if } T_1 < |T’_s| \le T_2 \\
K_3, & \text{if } |T’_s| > T_2
\end{cases}
$$
where $K_1, K_2, K_3$ are distinct stiffness values and $T_1, T_2$ are torque breakpoints. Furthermore, backlash or mechanical play, denoted $\Delta$, introduces a dead-zone in the position transmission. The actual load position $\theta’_l$ accounting for backlash can be modeled as: $\theta’_l(t) = \theta_l(t) – \Delta$ if $\dot{\theta}_l(t)>0$ and a reversal condition is met; $\theta’_l(t) = \theta_l(t) + \Delta$ if $\dot{\theta}_l(t)<0$ and a reversal condition is met; otherwise, it holds the previous value during the dead-zone traversal.
Synthesizing these elements, a comprehensive block diagram model of the harmonic drive gear servo system can be constructed. This model integrates the linear plant (motor electrical/mechanical equations and two-mass mechanics) with the nonlinear blocks for saturation, quantization, friction, mass unbalance, variable stiffness, and backlash. This forms the basis for implementing a high-fidelity simulation in tools like MATLAB Simulink. A representative Simulink model would contain subsystems for the motor driver (with current controller $C$ and feedback gain $K_h$), motor armature, motor rotor dynamics, load dynamics, friction function, stiffness function, backlash function, mass unbalance function, and encoder simulation.
A crucial step in making this model practical is the acquisition of accurate parameter values. Parameters fall into two categories. The first are readily available from component datasheets: motor constants ($R, L, K_m, K_e, J_m$), driver current limits ($I_{max}, I_{min}$), encoder specifications ($N, n_e$), harmonic drive gear reduction ratio $K_w$, and nominal stiffness/backlash values. The second category contains parameters difficult to measure directly or that vary with assembly and operating conditions. These include the current controller dynamics $C$, the feedback gain $K_h$, and the full set of friction, damping, mass unbalance, and load inertia parameters ($T_s^+, T_s^-, T_C^+, T_C^-, B, J_l, m g \rho$).
Traditional methods for identifying these parameters often involve elaborate static or quasi-static tests with specialized fixtures, which are time-consuming and costly. I propose an efficient alternative that leverages a single frequency sweep experiment to estimate most of these parameters simultaneously. The procedure begins with identifying $C$ and $K_h$. By injecting a known signal into the driver and measuring its current output, standard system identification techniques (e.g., basic least squares) can be applied to estimate the transfer function of the current loop. The static friction levels $T_s^+$ and $T_s^-$ can be estimated from the dead-zone width of the step or low-frequency sinusoidal response, corresponding to the maximum torque command before motion initiates.
For the remaining parameters ($T_C^+, T_C^-, B, J_l, m g \rho$), a simplified model is used to formulate a least-squares identification problem. Recognizing that the current-controlled driver, motor armature, and back-EMF can often be lumped into a single effective gain $K_a$ over the frequency range of interest, and consolidating all disturbances to the load side as $T_d$, the system reduces to a simplified form. The load velocity $\omega_l(s)$ in response to input voltage $u(s)$ and disturbance $T_d(s)$ is:
$$ \omega_l(s) = \frac{K_a K_m K_w}{J_l s + B} \left[ u(s) – \frac{T_d(s)}{K_a K_m K_w} \right] $$
Discretizing this with a zero-order hold (sample time $t_s$) yields a difference equation:
$$ \omega_l(k) = e^{-(B/J_l)t_s} \omega_l(k-1) + \frac{K_a K_m K_w}{B}(1 – e^{-(B/J_l)t_s}) \left[ u(k-1) – \frac{T_d(k-1)}{K_a K_m K_w} \right] $$
The disturbance $T_d$ encapsulates Coulomb friction and mass unbalance. Using a signum function based on a Stribeck velocity threshold $\Omega$ to separate positive and negative motion regimes, we can express $T_d$ as:
$$ T_d = P(\omega_l(k)) T_C^+ + N(\omega_l(k)) T_C^- + m g \rho \sin(\alpha_0 + \theta(k)) $$
where $P(\omega_l(k)) = 0.5 \sigma(\omega_l(k)) (1+\sigma(\omega_l(k)))$ and $N(\omega_l(k)) = -0.5 \sigma(\omega_l(k)) (1-\sigma(\omega_l(k)))$, with $\sigma(\omega_l(k))$ being 0 for $|\omega_l|<\Omega$, 1 for $\omega_l>\Omega$, and -1 for $\omega_l<-\Omega$. Substituting this into the difference equation gives a linear-in-the-parameters form:
$$ \omega_l(k) = \begin{bmatrix} \omega_l(k-1) & u(k-1) & -P(\omega_l(k-1)) & -N(\omega_l(k-1)) & -\sin(\alpha_0+\theta(k-1)) \end{bmatrix} \Theta $$
with the parameter vector $\Theta$ defined as:
$$ \Theta = \begin{bmatrix} e^{-(B/J_l)t_s} & \frac{K_a K_m K_w}{B}(1 – e^{-(B/J_l)t_s}) & \frac{T_C^+}{B}(1 – e^{-(B/J_l)t_s}) & \frac{T_C^-}{B}(1 – e^{-(B/J_l)t_s}) & \frac{m g \rho}{B}(1 – e^{-(B/J_l)t_s}) \end{bmatrix}^T $$
Collecting input-output data over $N_s$ samples into matrices $\mathbf{Y}$ and $\boldsymbol{\Phi}$, we have:
$$ \mathbf{Y} = \boldsymbol{\Phi} \Theta + \mathbf{E} $$
The least-squares estimate is $\hat{\Theta} = (\boldsymbol{\Phi}^T \boldsymbol{\Phi})^{-1} \boldsymbol{\Phi}^T \mathbf{Y}$. From $\hat{\Theta}$, the physical parameters are recovered:
$$ \hat{J}_l = \frac{K_a K_m K_w t_s (\hat{\Theta}_1 – 1)}{\hat{\Theta}_2 \ln(\hat{\Theta}_1)}, \quad \hat{B} = \frac{K_a K_m K_w (\hat{\Theta}_1 – 1)}{\hat{\Theta}_2} $$
$$ \hat{T}_C^+ = \frac{K_a K_m K_w \hat{\Theta}_3}{\hat{\Theta}_2}, \quad \hat{T}_C^- = \frac{K_a K_m K_w \hat{\Theta}_4}{\hat{\Theta}_2}, \quad \widehat{m g \rho} = \frac{K_a K_m K_w \hat{\Theta}_5}{\hat{\Theta}_2} $$
This method provides a efficient, data-driven way to populate the harmonic drive gear system model with credible parameters.
To validate the entire modeling methodology, a comparative analysis between simulation and experimental results is essential. I constructed a detailed Simulink model based on the described framework for a specific modular harmonic drive gear assembly. The model parameters were set using a combination of datasheet values and the identification procedure outlined above. For experimental validation, the actual harmonic drive gear system was tested under identical conditions.
First, open-loop time-domain responses were compared. Applying a 5V, 1Hz square wave command, both the simulated and actual system outputs were recorded. The results showed strong agreement in terms of rise time, settling behavior, and the presence of nonlinear effects like stiction and backlash. To quantify the match, a normalized degree of fit $p$ can be calculated over $n_p$ sample points:
$$ p = 1 – \frac{\sum_{i=1}^{n_p} (x_{e_i} – x_{s_i})^2}{\sum_{i=1}^{n_p} x_{e_i}^2} $$
where $x_{e_i}$ is the experimental value and $x_{s_i}$ is the simulation value. For the open-loop square wave response, the matching degree exceeded 89%.
Second, open-loop frequency-domain characteristics were analyzed via a 5V, 1-100Hz sine sweep. The magnitude and phase responses from the simulation and experiment were plotted. The magnitude plots matched closely across most of the spectrum, with a fit of about 96%. The phase plots showed good agreement in the mid-frequency range, though slight deviations were observed at very low and high frequencies, which can be attributed to simplifications in the friction model and slight inaccuracies in the assumed harmonic drive gear stiffness profile. The overall phase match was around 72%.
| Component | Parameter | Symbol | Value | Unit |
|---|---|---|---|---|
| DC Motor | Armature Resistance | $R$ | 5.6 | $\Omega$ |
| Armature Inductance | $L$ | $2.8 \times 10^{-3}$ | H | |
| Torque Constant | $K_m$ | 0.517 | N·m/A | |
| Back-EMF Constant | $K_e$ | 0.517 | V·s/rad | |
| Rotor Inertia | $J_m$ | $6.82 \times 10^{-4}$ | kg·m² | |
| Driver | Max Current | $I_{max}$ | 3.0 | A |
| Min Current | $I_{min}$ | -3.0 | A | |
| Encoder | Interpolation Factor | $N$ | 40 | – |
| Line Count | $n_e$ | 18000 | – | |
| Sampling Time | $t_s$ | 0.001 | s | |
| Harmonic Drive Gear | Reduction Ratio | $K_w$ | 80 | – |
| Stiffness (Piecewise) | $K_1, K_2, K_3$ | $5.4\times10^5, 8.8\times10^5, 9.8\times10^5$ | N·m/rad | |
| Backlash | $\Delta$ | $2.2 \times 10^{-4}$ | rad | |
| Load | Load Inertia | $J_l$ | $2.35 \times 10^{-2}$ | kg·m² |
Finally, the closed-loop performance was evaluated by implementing the same velocity controller on both the simulation model and the physical harmonic drive gear system. A frequency sweep analysis of the closed-loop system yielded the bandwidth and stability margins. The results, summarized in the table below, demonstrate that the model successfully predicts key servo performance metrics.
| Performance Metric | Experimental System | Simulation Model |
|---|---|---|
| Bandwidth (-3 dB) | 70 Hz | 75 Hz |
| Rise Time (10%-90%) | 0.012 s | 0.011 s |
| Overshoot | 25 % | 30 % |
| Settling Time (to 2%) | 0.05 s | 0.06 s |
| Steady-State Error (Position) | 6 arcseconds | 6 arcseconds |
The matching degree for the closed-loop frequency response magnitude was approximately 86%, and for the phase, it was about 70%. These results consistently show a match exceeding 70% across all tested domains, validating the correctness and practical utility of the proposed harmonic drive gear system modeling approach.
In conclusion, the modeling framework presented here addresses the need for a comprehensive, system-level representation of precision harmonic drive gear assemblies. By integrating a linear two-mass dynamic core with strategically chosen nonlinearities—including drive saturation, sensor quantization, Stribeck friction, mass unbalance, piecewise stiffness, and backlash—the model captures the essential behaviors that influence servo performance. A significant contribution is the parameter identification strategy that leverages a single frequency sweep experiment, employing a least-squares algorithm to efficiently estimate parameters like Coulomb friction, viscous damping, load inertia, and mass unbalance. This greatly enhances practicality by reducing the need for costly and time-consuming dedicated tests. The validation through extensive time-domain and frequency-domain comparisons between simulation and experiment confirms the model’s accuracy, with all key results showing a matching degree over 70%. This harmonic drive gear system model, with its parameters grounded in physical meaning and an efficient identification process, serves as a powerful tool for engineers. It enables advanced control algorithm development, system performance prediction, and design optimization for applications relying on the exceptional capabilities of the harmonic drive gear, ultimately accelerating the development cycle and improving the reliability of precision motion systems.
The journey of modeling a harmonic drive gear system is an exercise in balancing complexity with practicality. The harmonic drive gear, with its unique operating principle, introduces specific challenges that must be abstracted into a computationally tractable yet accurate form. The piecewise stiffness model, while a simplification of the true hysteresis, proves sufficient for predicting bandwidth and resonance effects. The friction model, though static, captures the dominant low-velocity phenomena that are critical for precision pointing and settling. Future refinements could incorporate dynamic friction models like the LuGre model for even greater fidelity in simulating very slow or reversing motions. Similarly, the stiffness identification could be enhanced to directly capture the hysteresis loop from dynamic data. Nevertheless, the current framework establishes a solid foundation. The iterative process of model-building, parameter identification, and validation fosters a deeper understanding of the harmonic drive gear system’s behavior. It reveals how the interaction between the motor’s electrical time constant, the mechanical resonances introduced by the harmonic drive gear’s compliance, and the various nonlinearities shapes the overall transfer function and transient response. This understanding is invaluable when tuning controllers to achieve aggressive performance specifications without exciting unstable modes. In essence, a good model transforms the harmonic drive gear from a black-box component into a transparent, predictable element within the larger system design. As harmonic drive gear technology continues to advance, with developments in materials, tooth profiles, and wave generator designs, the underlying modeling principles described here will remain relevant, adaptable to new data and characteristics. The continuous dialogue between model prediction and experimental measurement is the cornerstone of engineering precision systems, and this work provides a robust methodology for conducting that dialogue effectively for systems centered on the remarkable harmonic drive gear.