Accurate and Simple Modeling of Eddy Current Braking Torque: Analysis and Experimental Validation

A compact, reliable, and straightforward mathematical modeling of the eddy current (EC) braking torque is proposed in this article. First, the braking magnetic force, which curbs the movement of a rotating disk (RD), is evaluated using the basic laws of electromagnetism. Second, the braking torque is related to the braking force through a polynomial function model. The model parameters are evaluated by plugging in measurement samples in the equation of motion taking into account the nonlinear behavior of the aerodynamic drag forces. The proposed model neither has limitations on the braking system geometry nor the disk rotational speed as long as the measurement samples cover the whole speed range. Moreover, it is validated experimentally for constant and alternating magnetic flux profiles with best fit rates (BFRs) more than 85% and 93% for both rotational speed and braking torque, respectively. Furthermore, programmable braking application is demonstrated using the proposed approach and its applicability is verified experimentally.

E CBS are contactless technology in which magnetic forces hinder the motion of a metallic structure [1]. The contactfree feature of the ECBs widens the room of innovative applications and improves the performance of the existing mechanical systems. For instance, automotive applications use ECBs to slowdown a vehicle velocity avoiding the problems associated with mechanical braking [2], [3], [4], [5]. Besides, undesired vibrations in metallic structures can be reduced or suppressed using ECDs [6], [7]. Furthermore, some rehabilitation devices use magnetic braking to improve functional strength training of muscles [8]. In addition, ECBs are used in some haptic interfaces as actuators [9]. Thus, the ECBs are involved in various applications due to the simple construction and the efficient functionality.
The braking mechanism depends on the relative normal motion of a conductor inside a region containing magnetic flux density generated by a PM or an electromagnet. EMFs are generated due to this movement, which create ECs on the conductor body. According to Lenz's law, the magnetic flux generated by these ECs opposes the original magnetic flux, but it cannot cancel it totally. The exposure of the ECs to the magnetic flux resultant causes magnetic forces or torque, which curbs the original driving torque. Accordingly, the overall speed of the conductor is braked based on the applied magnetic field strength and the initial speed. This phenomenon has been studied intensively to provide mathematical formulations to predict its behavior under time variant or steady magnetic flux [10], [11], [12], [13], [14], [15].
The ECs induced in a disk rotating under external magnetic flux have been calculated considering the demagnetization effects using Maxwell's equations [10]. A simplified numerical analysis has used mesh distortion instead of expressing the velocity clearly in Maxwell's equations to estimate the ECs generated on a moving ring [11]. Moreover, other articles provided a parametric model of the braking torque [12]. This model considered the torque as a quadratic function of the braking current with nonlinear coefficients. These coefficients were represented by second-order polynomials of the rotational speed. On the other hand, different models were based on the analytical solutions of the ECs due to the exposure to a time-independent magnetic field [13]. The surface charge density is calculated using the Hall effect concept, while the generated electric field has been derived using Coulomb's law. However, this analysis neglected the nonlinear behavior of the viscous friction drag forces and it is limited to a specific structure of a pole shape.
A similar methodology has been presented, where the ECs generated from time-varying fields with different spatial variations are analyzed and modeled using auxiliary potentials [14]. Another analytical modeling has been introduced to investigate the interaction between a time-varying magnetic flux and a low speed RD [15]. However, both the studies have some degree of computational complexity since they depend on evaluating an infinite series which is truncated with sufficient terms for a sufficient accuracy. Besides, the high-frequency magnetic flux sources are limited to RDs with low speeds or low moment of inertia [16], [17].
In the above-mentioned models [10], [11], [12], [13], [14], [15], the EC formulation has been carried out based on exact solutions, but they are limited to specific poles' cross section. Moreover, the infinite sum-based solutions are timeconsuming, especially for high-accuracy calculations. Besides, the braking torque effect on the RD has been considered without taking into account the nonlinear characteristics of the aerodynamic drag forces. More details about rotationaltype ECB systems can be found in [18] In this article, we present a simple and reliable model for the ECB torque, which retards a diamagnetic RD regardless of the poles' geometry. Moreover, the nonlinear behavior of the aerodynamic drag forces is accommodated through the RD dynamics study. Furthermore, the proposed model is used to generate programmable braking which is required for rehabilitation devices and robotic control [8].
The article is organized as follows. In Section II, the ECB system is presented along with the explanation of the magnetic braking mechanism to deduce the induced currents' formula. This is followed by the derivation of the movement equation, which describes the disk dynamics in Section III. The experimental validation of the proposed model is discussed in Section IV. To ensure the model significance, an application is proposed in Section V. Finally, the conclusion of the article is summarized in Section VI.

II. ELECTROMAGNETIC BRAKING MECHANISM
The used ECB system consists of a non-magnetic RD inserted in an air gap, which separates a current-carrying coil and a fixed iron disk. Thus, the coil is wounded on horseshoe iron core to focus the magnetic flux on the RD as depicted in Fig. 1. The aluminum disk is rotating in the azimuth plane generating motional emf, which can be calculated as follows: where − → u , − → B , and − → dl are the vector quantities, which represent the disk linear velocity, the magnetic flux density, and the differential length, respectively. The limits of the integral are defined according to the emf direction, which can be obtained from the cross product of − → u and − → B . Since the disk rotates clockwise φ direction and the magnetic flux is aligned with the positive or negative z-axis, the emf is in the inward or outward radial direction. The disk regions going in the PPA induces EC, which creates magnetic flux opposing the original flux. On the contrary, the induced current in the regions going away from the PPA generates a magnetic flux, which supports the original flux, as indicated in Fig. 2(b) and (d). The ECs generated on an RD are simulated using ANSYS Maxwell-3-D full-wave solver for the structure shown in Fig. 1. The initial rotational velocity of the disk is 300 r/min and the coil current is 1 A. The EC density on the upper surface braked disk is shown in Fig. 2(e) after reaching steady-state rotational speed. The coil turns number and its current in this full-wave simulation does not match the fabricated coil, which will be used later for the experimental validation. This mismatch introduces a simplified visualized theory explanation without computation complexity, especially for such fabricated coil with very large turns number. Consequently, the magnetic field strength of the fabricated coil and the associated ECs generated on the RD are very large compared with the simulated results.
To propose simple forms of the induced currents on the RD, the braking force, and the braking torque, let us consider the following initial assumptions. • Both the rotational velocity and magnetic flux are constants as for the integration (1).
• The magnetic flux density is in direct proportionality with the coil current. These simplifications are reasonable and acceptable as long as the magnetic flux distribution is uniform over the PPA, and the B-H curve of the coil core is in the linear region. Therefore, induced ECs can be deduced based on (1) and expressed in the following form: where represents the angular speed in r/min and I c represents the coil current passes through the coil. The constant C eddy represents the combined effect of the material resistivity, number of turns, pole shape cross section, and core permeability assuming the gap between the RD and the core is maintained to a fixed spacing. The braking force can be calculated as follows: Applying the cross product, we can find that the magnetic force direction coincides the φ direction but opposes the original rotation. Based on the aforementioned simplifications, the magnetic flux density can come out of the integral. Therefore, we can express the braking force in terms of the coil current and the rotational velocity using another constant including C eddy and the integration path information as follows: where C f andφ represent the braking force constant and a unit vector in the counterclockwise rotational direction, respectively. F brk is proportional to the coil current square due to the double effect of the magnetic field on the RD. The associated braking torque with this rotational force can be calculated using the following integral: where γ includes both C eddy and C f in addition to other parameters related to the disk shape. To this end, relations (2), (4), and (5b) are used in a crude form without considering the magnetic flux nonuniformity on the RD and the B-H curve nonlinearity. Therefore, the simple braking torque model must be multiplied by a correction function to generalize its use through considering such actual effects in addition to taking any latent nonlinearity into account. The experimental observations assure the fact that γ is a current-dependent parameter, which decreases versus the injected current in the coil. The explanation of such dependency is that the EC generated on the RD is a volume current not a surface current, which takes different contour shapes with different coil currents. Accordingly, C eddy and C f are also dependent on the coil current, but absorbing such dependencies inside γ is sufficient to present a compact model of the braking torque. From another perspective, the increase in the coil current results in an increase in the braking torque. However, the braking torque profile does not increase linearly or quadratically versus the coil current since it also depends on the steady-state rotational speed ( ss ), which decreases with the increase in the coil current. Physically, the braking torque increase cannot stop the RD as the motion is a necessary condition for the braking torque existence. Therefore, the braking torque is a multiplication of three quantities, one of them increases quadratically, whereas the other two quantities decrease monotonically. Therefore, it is a monotone increasing function as will be shown in the following section. To consider the nonlinearities of the EC on the disk, γ is written as an nth-order polynomial in the coil current only as follows: where b k represents the polynomial coefficients, which are easily estimated through taking some measured samples.
These samples must cover the whole operating ranges of the coil current and the RD speed to express the coefficients in addition to other calculated parameters correctly. To estimate the impact of the ECs on the movement of the disk, we have to apply Newton's law related to the angular motion. This will be discussed in detail in Section III.

III. GOVERNMENTAL MOVEMENT EQUATION DERIVATION
The mechanical dynamic equation of the system is obtained by applying Newton's law of the rotational motion, which states that the total torque affecting a rotating body equals to the rate of change in the angular momentum. Therefore, we can write the movement equation as follows: where J represents the moment of inertia factor including all the moving parts not the RD only. T dev represents the developed torque by the source of motion, which may be mechanical motor, electrical motor whether dc or ac, or even the physical motion of the human organs. T brk represents the braking torque, which depends on the coil current square and the disk rotational velocity. The losses associated with the aerodynamic drag forces are represented by T drg , which is considered as a nonlinear term versus rotational speed as follows [19]: where d is a constant that includes the physical parameters of the rotation medium (air) and the geometry of all the rotating objects. The constant d is calculated from practical measurements at steady-state rotation of the RD to consider other moving parts such as the shaft and gears. Also, J is calculated by comparing the step response of the model to the measured one when the disk is rotating freely without braking. Representing aerodynamic drag torque in a nonlinear actual form is crucial to get accurate mathematical modeling of the problem as a whole for different speed ranges. On contrast, neglecting aerodynamic drag [20] or representing it as a linear function in the rotational velocity [13] limits the model accuracy to a specific speed range. We consider a cheap and a commercially available dc PM coreless motor with a simple equivalent circuit, which consists of a power supply E connected in series to a resistance R a and an inductance L a representing the motor windings. The back emf is generated on the inductance, which depends on the rotational velocity through the constant k v . Moreover, the generated torque is related to the driven current via the constant k t . Therefore, T dev can be expressed as follows [21]: where I a and I nl represent the rotor current during loading and without loading, respectively. The first term in the last line in (9) represents the peak torque that can be developed by the motor, whereas the second term expresses the energy transformation impact on such developed torque. Substituting (5), (8), and (9) in (7), the governmental movement equation can be represented in a simplified form as follows: The model in (10) approximates the effects of the ECs in the term γ I 2 c . Moreover, it considers the nonlinear aerodynamic drag effect. The constants k t , k v , and I nl depend on the motor specifications. The parameter γ depends on the core and disk shapes, their material, the number of turns of the coil, the gap between the core, and the RD and its medium. We should note that the equation of mechanical motion is derived in a general form without basic dependency on the motor type. In addition, the methodology of modeling the ECB torque does not depend on the motor type as well. In Section III-A, we propose a procedure to identify these constants using available measurements.

A. Proposed Model Parameters' Evaluation
The steady-state movement of the RD occurs when the coil carries time-invariant current or when there is no current passing through the coil. Thus, the right-hand side of (10) becomes zero in both the cases and the movement equation can be written as follows: where (11a) and (11b) correspond to I c = 0 and I c ̸ = 0, respectively, and 0 and ss represent the rotational speed of the RD before the magnetic braking and its steady-state rotational speed after that due to the steady magnetic field, respectively. 1) Aerodynamic Drag Coefficient Calculations: Referring to (11a), the viscous friction drag coefficient d is calculated experimentally according to the following relation: It is obvious that the coefficient d is calculated using the motor parameters and the steady-state rotational speed without magnetic braking. To validate the correctness of the range that estimated drag coefficient belongs to, it is calculated based on the theoretical equations found in [19] of a one side wetted disk of radius R rotating in an air medium in the normal temperature as follows:   to including the conversion constant of the rotational speed to match r/min unit. The experimental system that we use contains two disks with a radius of 3.2 cm. Therefore, the    Table II in the first row up to down, respectively.
drag coefficient is quadrupled to consider the two RDs with both sides.
2) Procedure of Calculating γ Polynomial Coefficients: Referring to (11b), we can show that The dependency of γ on the coil current I c and the steady-state rotational speed ss is clear in (14). As ss itself depends also on the coil current, we have proposed γ to be a polynomial function of the coil current only as shown in (6). By measuring some pairs of the steady-state rotational speed after braking and the coil current, the model of the correction function γ can be identified. It is noted that the dependence of γ on the aerodynamic drag torque has unremarkable effect at low speeds. However, the braking torque model is dramatically affected by the aerodynamic drag at high-speed zone. Therefore, it can be neither neglected nor considered as a constant or linear term versus rotational speed to serve wide ranges of rotational velocity.

A. System Description
The proposed model is validated using a cheap and a commercially available dc servo motor control module from LJ CREATE as illustrated in Fig. 3(a). This module has been designed to allow several experiments to control a dc motor speed and position using an analog or a digital controller. The dc motor is driven by a voltage input connected to an analog drive amplifier; and its rotational speed is measured by a tachogenerator directly coupled to the dc motor shaft. A light aluminum RD disk is fixed on the motor shaft whose motion can be braked through external load torque generated by a PM as illustrated in Fig. 3(a). The braking PM in the system is replaced by a current-carrying coil to control the magnetic flux. The coil is wounded as five layers on a horseshoe-shaped core of iron, with µ r ≈ 4000, length = 15 cm, and core diameter = 2.5 cm. The coil inductance is L c = 15 mH and its resistance is R c = 5.4 with the number of turns ≈ 2500 turns. Increasing the number of turns increases the magnetic flux density at the expense of increasing the resistance of the coil, and it may cause breakdown for the insulators of the coil wire. Definitely, using high conductive wires allows more compact coil size and more thermal stability, but it is a costly solution for just experimental validation. Therefore, we used different thicknesses of copper wires during coil wrapping to increase the number of turns and allow high current without damaging the insulators. The coil is connected to the output of a CA, which contains eight power transistors 2N3055 to serve as a heat sink and guarantee stable operation without thermal runaway. The CA transfer characteristic is nonlinear as illustrated in Fig. 4. The CA input is connected to the National Instrument USB-6251 data logger card to provide the required driving voltage.

B. Performance Evaluation of the Proposed Model
The motor parameters are listed in Table I according to [22]. The calculated aerodynamic drag coefficient d is 10.58 × 10 −6 which is very close to the theoretical value with small deviation due to other moving parts such as gears and shaft. As the J factor depends on all the rotating parts, its value is calculated by comparing the measured step response, without applying the magnetic braking, to the numerical solution of (10) as illustrated in Fig. 5.
The coil is energized by five different dc currents, and the corresponding steady-state rotational velocities ss are measured using the data logger and a built-in tachometer to evaluate the parameters of the γ function. Each steady-state velocity depends on the coil current, and consequently, γ is a function of the coil current only, which represents the independent variable (6). The measurements of ss versus I c and the corresponding γ values are listed in Table II. The measured values of γ are used to identify the coefficients of γ polynomial using MATLAB curve fitting tool. Based on the experimental observations, the suitable fit function, which covers all the changes in γ versus I c , is a fourthorder polynomial. The constant coefficients b k s in the function describing γ are calculated using the curve fitting tool based on the measured data shown in Table II; their values are listed in Table III. The comparison between the measured and calculated values of γ is shown in Fig. 6 which indicates that they are fit well. Moreover, the model responses, when the coil is injected with the different dc currents, are compared with the measured response of the same dc currents to check the model validity in the transient region as shown in Fig. 7. Besides, the proposed model is tested when the coil carries low-frequency time-variant current to investigate the dynamic behavior of the proposed model. The low-frequency ac current enables the motor to respond to such a variation, and thus, the transformer emf effect can be neglected. To assess the performance of the proposed model, the BFR (%) is calculated (15) as follows: where x meas stands for the measured quantity, either in r/min or T brk in N · m, whereas x mdl represents the quantity estimated from the proposed model and its mean value is used in the denominator in (15). The estimated profiles of the rotational speed and the braking torque are illustrated in Fig. 8. In addition, the BFRs of the rotational velocity  and the braking torque are 86.37% and 93.81%, respectively, which assures the model validity. In Section V, we propose available applications, which can use the proposed model as a design guide. Moreover, we validate the applicability of these applications using our RD system with the external coil.

V. PROGRAMMABLE BRAKING AS AN APPLICATION BASED ON THE PROPOSED MODEL
Functional strength training of the weak muscles requires a gradually increased resistance level, or it may follow variable resistance profiles according to physician recommendations. Therefore, some rehabilitation devices use intermittent braking through mechanical position offset of a PM away from an RD to provide such requirements [8]. It is possible to use programmable braking through replacing the PM with a current-carrying coil with suitable design to provide controllable resistance. The dc motor system is used to validate the proposed application without changing its basic feed voltage. Definitely, forcing a dc motor to follow a required rotational speed profile is available by changing its feed voltage, but it is used here for validation only. We assume that an arbitrary rotational speed profile consists of ten truncated cosine terms with frequencies in the range of [0 0.4 Hz]. As aforementioned, this low-frequency range simulating the actual mechanical changes permits the disk to track the variations and enables the designers to neglect the transient behavior of the braking coil. Both the cosine amplitudes and frequencies are selected randomly using MATLAB built-in functions. The overall velocity variation versus time is illustrated in Fig. 9(a), whereas the corresponding spectrum is shown in Fig. 9(b). It is required to find the coil current that meets the desired profile given the rotational velocity function and its derivative. Therefore, (10) is solved as an algebraic equation for the coil current not as a differential equation. Then, the CA transfer characteristics equation is used to get the required base current, which is injected using the data logger card NI-USB 6251. The resultant driving base current is illustrated in Fig. 10(a), whereas the measured rotational velocity is compared with the required profile and illustrated in Fig. 10(b). The similarity between both the profiles is estimated using the BFR criterion, which exceeds 85.53%.
The proposed model of the EC braking can be extended to cover a wide range of applications in control engineering. The model can be used to stabilize a moving object velocity corrupted with undesired disturbances, where the main control parameter is the current injected in the braking coil. Moreover, the velocity of a moving object can be regulated or forced to take a desired profile which represents a fundamental requirement for real-time motion tracking robots. On the contrary, the ECB can be used as an artificial latent disturbance, which can be used for comparing the performance of different control techniques in rejecting such an external disturbance [23].

VI. CONCLUSION
A direct mathematical model of the EC braking ECB torque retarding an RD has been presented based on the basic laws of electromagnetism and given related measurements to deduce the model parameters. Moreover, the nonlinear aerodynamic drag forces are considered and validated theoretically and experimentally. The ECB torque model has been tested and compared with the experimental results for both steady magnetic flux and low-frequency time-varying magnetic flux. It achieves identical curve profiles with BFRs of 86% and 95% for rotational velocity and braking torque, respectively. The programmable braking has been proposed as an application and its applicability has been discussed and verified experimentally with BFR exceeding 85%. The systematic procedure followed in this article can be generalized for different ECB systems regardless of the poles number or the pole cross section shape. Therefore, the model can be extended to serve various control applications.