Modeling and Validation of Diamagnetic Rotor Levitated by Permanent Magnetics

As an innovative, low-power consuming, and low-stiffness suspension approach, the diamagnetic levitation technique has attracted considerable interest because of its potential applicability in miniaturized mechanical systems. The foundation of a diamagnetic levitation system is mathematical modeling, which is essential for operating performance optimization and stability prediction. However, few studies on systematic mathematical modeling have been reported. In this study, a systematic mathematical model for a disc-shaped diamagnetically levitated rotor on a permanent magnet array is proposed. Based on the proposed model, the magnetic field distribution characteristics, diamagnetic levitation force characteristics (i


Introduction
Magnetic levitation technology has been increasingly applied in industrial applications such as bearings, sensors, and magnetic transportation systems because of its ability to avoid mechanical wear and lubrication [1][2][3][4][5][6][7].Generally, friction levitation in most magnetic levitation systems is achieved by employing active-feedback-controlled electromagnetics.However, the use of feedbackcontrolled electromagnetics is power-consuming and requires complicated controller designs.Moreover, given that leakage flux dramatically increases with scale downsizing, it poses an obstacle to realizing the miniaturization of mechanical systems with the utilization of feedback-controlled electromagnetics [8].In this regard, the recently developed passive levitated suspension, represented by diamagnetic levitation, has emerged as an alternative suspension approach for mechanical systems, particularly those that require low power consumption and are small.Diamagnetic levitation is the free levitation of diamagnetic objects in a magnetic field.Typical diamagnetic materials include water, bismuth, and graphite [9].In the absence of a control system, the suspension structure of a diamagnetic levitation system is simple [10].Statically and dynamically stable diamagnetic levitation can easily be achieved at room temperature without the need for an energy supply.These features make diamagnetic levitation appealing for modern micromechanical systems, such as microscale separation or manipulation [11] and micro rotor bearing systems [12][13][14][15][16][17][18][19].Additionally, the inherent characterization of low stiffness in diamagnetic levitation systems permits their potential application in high-sensitivity sensors [20][21][22][23] and ultralow-frequency energy harvesters [24,25].
For instance, Liu et al. [14,15] invented a diamagnetically levitated variable-capacitance micromotor and successfully realized its continuous rotational operation.Bleuler et al. [16] demonstrated the rotation of a diamagnetically suspended rotor by an electrostatic glass motor.Su et al. [17] reported a micromachined rotor in which highly oriented pyrolytic graphite was diamagnetically levitated by magnetics and rotated by gas flow.Our previous work [18,19] realized a rotating diamagnetically levitating rotor system driven by an electrostatic field with a maximum rotation speed of 300 r/min.
Besides the effort in designing specific micromechanical systems based on diamagnetic levitation mechanisms, attempts have been made to construct theoretical and empirical models of diamagnetic levitation systems for design and performance optimization.For instance, Chow et al. [26] theoretically analyzed the effect of the shape of magnetic sources formed by standard coil and ring magnet elements on vertical diamagnetic levitation.Pasquale et al. [27] presented a diamagnetic suspension model based on a finite element software and validated it experimentally.Cansiz et al. [12] developed a unique numerical calculation method to analyze the static and dynamic characterizations of a diamagnetic bearing.Chen et al. [13] developed a numerical model to analyze the performance and stability of diamagnetic bearings for microscale applications.
The proposed models and analyses consolidate the theoretical foundation of diamagnetic levitation.However, to date, systematic mathematical modeling of diamagnetic levitation has been rarely reported.In this study, a systematic mathematical model of a disc-shaped diamagnetically levitated rotor on a permanent magnet array was proposed and developed.Based on the proposed mathematical model, the magnetic field distribution characteristics, diamagnetic levitation force characteristic (i.e., levitation height and stiffness) expressions, and optimized theoretical conditions for realizing stable levitation were achieved.Experiments were also conducted to determine the relationship between the levitation height and force.The operating stability of the rotor was evaluated under various conditions.The experimental results agree well with the theoretical predictions, confirming the feasibility of the proposed mathematical model.
The remainder of this paper is organized as follows.Section 2 describes the diamagnetic levitation test rig.Section 3 presents the modeling and analysis of the magnetic field generated by the magnetic array.Section 4 describes the modeling and analysis of the diamagnetic force.The experimental validation and discussion are presented in Section 5, and the conclusions are presented in Section 6.

Description of Diamagnetic Levitation Test Rig
In our previous study, we developed a diamagnetic contactless suspension rotor driven by an electrostatic glass motor [18,19].In this study, a more compact test rig consisting of a rotor, electrostatic glass motor, and array of permanent magnets was designed.The disc-shaped rotor (φ10.38 mm) was comprised of pyrolytic graphite and borosilicate glass (Figure 1(a)).The thickness of pyrolytic graphite and borosilicate glass were 500 μm and 120 μm, respectively.
For the array of permanent magnets, a concentric ring-shaped magnet (φ18.68×φ7.96×12mm) encircling a cylinder magnet (φ7.94×12 mm) with opposite axial magnetization pattern was adopted, as shown in Figure 1(b).The diamagnetic rotor was passively and stably levitated at room temperature by the diamagnetic force between the magnets and rotor.The levitation gap between the rotor and magnets was 0.65 mm.An electrostatic glass-motor-driven system was designed to ensure rotation.The stator of the electrostatic glass motor was fabricated using a flexible printed circuit board (PCB) film (thickness 100 μm) with 12 poles (36 electrodes), which can be placed at the gap between the magnets and rotor.
The proposed model and corresponding analysis are directed toward static conditions.Therefore, the borosilicate glass of the rotor can be eliminated because no electrostatic force is required for the rotation of the rotor.Thus, in the proposed model and experimental validation, the rotor was composed only of pyrolytic graphite, and the measured levitation height between the rotor and magnet was 0.9 mm.

Magnetic Feld Generated by Cylinder Permanent Magnet
First, a current model was adopted to analyze the permanent magnets.According to this model, the magnet was reduced to a distribution of equivalent current.This was then input into the magnetostatic field equations as a source term, and the field was obtained using standard methods for steady currents.Mathematically, a magnet with magnetization M is equivalent to the following current distribution: Volume current density

Surface current density
where n is the outward unit normal to the magnet surface.
For a better analysis, cylindrical coordinates were adopted considering the cylindrical geometry of the magnets.The unit vectors of the cylindrical coordinates r , φ , and z were r , φ , and ẑ , respectively.As shown in Fig- ure 2, the polarization direction of the cylindrical magnet was along its axis.Therefore, the volume current density J m = ∇ × M = ∇ × M s was zero, where M s represents the magnetization.
Regarding the surface current density, given that the normal of the top and bottom surface was either parallel or antiparallel to the magnetization, j m = 0 at the top and bottom surfaces.Consequently, we considered the outer surface as follows: where R is the outer radius of the cylindrical magnet surface and z 1 and z 2 denote the axial positions of the bottom and top of the cylindrical magnet surface, respectively. (1) (2) The unit normal for the outer surface is Therefore, the surface current density is According to the current model, the vector potential A generated by the magnets is [28] where µ 0 is the magnetic permeability of vacuum, x is the observation point, and x ′ is the source point.x − x ′ represents the distance between these two points and is given in cylindrical coordinates: Because the volume current density J m = 0 , Eq. ( 6) can be written as The unit vectors x and ŷ for Cartesian coordinates x and y have the following relationship with the cylindrical coordinate unit vectors: Substituting Eq. (9a) into Eq.( 8) yields ( 4) n = r.

Magnetic Field Generated by Array of Permanent Magnets
The array of permanent magnets for the diamagnetic levitation system employed was a concentric ring-shaped magnet encircling a cylindrical magnet with an opposite axial magnetization pattern.Figure 3 shows a brief illustration of the array structure of the permanent magnets.
For the ring-shaped magnet, the following two remaining surfaces were considered: and The unit normal for these surfaces is Therefore, the surface current densities are According to Eq. ( 22), a ring-shaped magnet is expected to be equivalent to two cylindrical magnets assembled with opposite axial magnetization patterns, as illustrated in Figure 4. Therefore, the arrayed permanent magnets (18b) (

Figure 3
Array of permanent magnets structure are equivalent to three cylindrical magnets assembled with different axial magnetization patterns (see Figure 5).
Therefore, the magnetic flux density B provided by the array of permanent magnets can be written as where B R2 , B R11 , and B R12 are the magnetic flux densities generated by each cylindrical magnet.Similarly, radius B r and axial field component B z generated by the array of permanent magnets can be obtained.
Table 1 lists the parameters of the magnetic array.Based on the above deduction, the magnetic flux density generated by an array of permanent magnets can be calculated.Figure 6 shows the obtained magnetic flux density distribution in the plane parallel the top surface of the permanent magnets wherein the distance between the plane and top surface is fixed to 0.9 mm.From Figure 6, it can be observed that radial component of the magnetic flux density B r first exhibits an up-trend with increasing radial distance from the magnet center and reaches the maximum at the point that is close to the interface between the ring-shaped and cylinder magnets.Afterwards, the magnitude of B r tended to gradually (23) decline to 0, followed by another upward trend, and then a downward trend within a radical distance of 10 mm.In the case of the axial field component of the magnetic flux density B z , the changing-trend appeared to be the opposite of that of B r within a radical distance of 10 mm.Specifically, the magnitude of B z first exhibited a downward trend with an increasing radial distance and then declined to 0 at the point close to the interface between the ring-shaped and cylindrical magnets.The points corresponding to the maximum B r and minimum B z were not located at the interface between the ring-shaped and cylindrical magnets, which can be explained by the different magnetizations of the two magnets.
Figure 7 shows the obtained magnetic flux density B z and B r distribution (in the plane parallel to the top surface of the permanent magnets) versus height from the top surface of the permanent magnets.Both the magnitude of B z and B r tended to decrease with increasing height.

Magnetic Field Generated by Array of Permanent Magnets
An external magnetic field can cause a diamagnetic subject to induce an opposite magnetic field.The unit volume potential energy for a diamagnetic material with induced magnetization M d under an external magnetic flux density B is Therefore, the elementary diamagnetic force on the unit volume is (24)    where χ m is the magnetic susceptibility.For pyrolytic graphite, because the χ m is very small, given Eqs.(26,27), the induced magnetization M d can be expressed and simplified as a function of B, This simplification implies that the influence of the diamagnetic material on the external magnetic field could almost be neglected.Therefore, the expression of the diamagnetic force per unit volume as a function of the magnetic field density becomes The total magnetic force acting on the diamagnetic body can be obtained by integrating the unit force on the entire volume: (25) When the magnetic susceptibility χ m is negative, the reaction force according to Eq. (30) provided by the magnetic field will push the diamagnetic body out of the field.However, because the magnetic susceptibility χ m for most diamagnetic materials at room temperature is very small, the diamagnetic force is relatively weak.Thus, diamagnetic levitation is a promising method for micromechanical applications.
From Eq. (30), the diamagnetic force in each orthogonal direction is According to the Ostrogradsky theorem, the above equations can be expressed as For the disc-shaped rotor levitation system, because the theoretical diamagnetic susceptibility of pyrolytic graphite in the radial direction is lower than that in the axis direction and the thickness of the disc-shaped rotor is small, the resultant F y is small and can be almost neglected.Consequently, the vertical diamagnetic force is mainly responsible for overcoming the (31a) gravitational force and supporting levitation.Therefore, our primary concern is the vertical diamagnetic force property.For the disc-shaped rotor adopted in this study, Eq. (32b) is simplified as where r d denotes the radius of the disc-shaped rotor.The theoretical diamagnetic susceptibility χ m of the pyrolytic graphite material is χ m = −450 × 10 −6 .However, because of unavoidable defects, it is difficult to ensure the good quality of commercial pyrolytic graphite and thus the real χ m of pyrolytic graphite is generally lower than that of the theoretical value.Additionally, the different processing histories and rare materials of pyrolytic graphite may cause variations in the real χ m .In this work, we deduce that χ m = −320 × 10 −6 is consistent with the experimental data.Thus, it is adopted for modeling.Using Eq. ( 33), we calculated the diamagnetic levitation force applied to the rotor at different levitation heights, and the results are presented in Figure 8.As shown, the diamagnetic levitation force and levitation height are nonlinearly correlated.As the levitation height increases, the diamagnetic levitation force gradually decreases.

Stiffness
Stiffness is an important parameter for calculating force generation and resonant frequency as well as predicting system behavior.The stiffness of the diamagnetic levitation force of a diamagnetic rotor at different vertical air gaps (levitation height) can be obtained from the derivation of the magnetic force given by Eq. ( 34): (33) With a fixed air gap, the stiffness of the diamagnetic levitation force for the diamagnetic rotor at different radial distances can be obtained from the derivation of Eq. ( 33), which is Figure 9(a) shows the theoretical stiffness of the diamagnetic levitation force versus air gap (levitation height).A nonlinear relationship between the theoretical stiffness and levitation height can be clearly observed.Specifically, the theoretical stiffness of the diamagnetic levitation force first sharply decreased and then slowly decreased to zero. Figure 9(b) shows the theoretical stiffness of the diamagnetic levitation force versus the radial position of the rotor at a fixed levitation height of 0.9 mm.With a fixed levitation height, the theoretical stiffness of the diamagnetic levitation force first increased consistently with increasing radial distance and reached its peak.A further increase in the radial distance decreased the theoretical stiffness.When the radial distance increased to the point corresponding to the interface between the ring-shaped and cylindrical magnets, the theoretical stiffness decreased to zero, followed by a further upward trend.

Stability
Stability is another critical concern for magnetic levitation systems.In a diamagnetic levitation system, steady levitation cannot be ensured if the diamagnetic levitation force can only overcome gravity.As mentioned above, F y can be neglected; thus, the stability of the rotor is determined by the torque.Here, the rotor mass center was selected as the centroid for the torque calculation.For better analysis, the cylindrical coordinates expressed in Eq. (33) were transformed into Cartesian coordinates.Assuming that the rotor had an eccentricity e along the x direction, the torque induced by the diamagnetic levitation force was expressed as Eq.(36): where Figure 10 shows the lateral displacement of the rotor versus torque obtained from Eq. (36) with a fixed rotor radius of 5.19 mm and levitation height of 0.9 mm. (34)

Figure 8 Theoretical diamagnetic levitation force versus levitation height
Figure 10 shows that when the rotor shifts towards the positive x direction and is located in the stable region, it is subjected to an anticlockwise torque in the y direction and then returns to the center of the permanent magnet.However, when the rotor shifts beyond the stable region, it is subjected to a clockwise torque in the y direction and then drops off the permanent magnet.A similar phenomenon was observed when the rotor shifted towards the negative x direction.Based on Eq. (36), we analyzed the lateral displacement of the rotor versus torque under different levitation heights, and the results are shown in Figure 11.From Figure 11, it can be concluded that the levitation height can strongly affect the levitation stability.At a levitation height of 0 mm, the rotor could not stably levitate.When the levitation height increased to 0.3 mm, the rotor had three equilibrium points at which it could levitate stably.With a further increase in the levitation height, the number of equilibrium points for the rotor became only one.Moreover, from Figure 11, it can be deduced that increasing the levitation height resulted in an enlarged stable region.
To better illustrate the effect of the levitation height on the rotor stability, the torque of the rotor versus the levitation height and lateral displacement is shown in Figure 12.For clarity, the regions wherein the torque is beyond 0 are marked with dots in Figure 12.From Figure 12, it can be concluded that the number of equilibrium points for the rotor varies at different levitation heights, and that increasing the levitation height is advantageous for broadening the stable region, which corroborates the analysis derived from Figure 11.
The radius of the rotor is another critical factor that determines the levitation stability.To elucidate the effect of the rotor size on the levitation stability, we analyzed the torque of the rotor versus the rotor radius and rotor lateral displacement based on Eq. ( 36) with a fixed levitation height.For distinct graphing, the numerical values on the z-axis were non-dimensionally processed, as described below: The resultant graph is shown in Figure 13, which shows the torque of the rotor versus the rotor radius and rotor lateral displacement.Apparently, with increasing rotor radius, the equilibrium line that supported stable (38) z = M e /mgr.levitation changed from the interface between the ringshaped and cylindrical magnets to the center of the magnet.However, with a further increase in the rotor radius, the stable regions of the rotor gradually diminished and even vanished.Thus, it can be concluded that with fixed levitation, a moderate rotor radius can permit stable levitation.
The results shown in Figure 13 were achieved at a fixed levitation height of 0.9 mm.A change in the rotor radius caused a change in the weight and area.Consequently, the levitation height may vary.Additionally, when the rotor radius is small, unstable levitation may occur at the interface between the ring-shaped and cylindrical magnets.This may be explained by the fact that the corresponding B z reaches a minimum that is not sufficiently large to overcome gravity.

Experimental Setup Description
Figure 14(a) shows a photograph of the magnetic field distribution tester employed for the magnetic flux density measurement.The magnet array was placed on the XY linear motor stage and the teslameter was fixed onto the precision Z-stage.The air gap between the teslameter and magnets was fixed at 0.9 mm, which is in accordance with the levitation height of the rotor.Considering that the teslameter can only measure the magnetic flux density in the vertical direction, we performed experiments to characterize B z .
Figure 14(b) shows a schematic drawing of the setup for the levitation force measurement.The precision scale was fixed to the base to measure the force.The permanent magnet array was fixed to the stepper motor to modulate the air gap between the magnets and diamagnetic rotor.To eliminate the effect of the magnets on the precision scale, a plastic rod was used to support the diamagnetic rotor and provide a sufficient distance between the magnets and precision scale.The mass of the plastic rotor was 10.0033 g and that of the diamagnetic rotor was 0.0903 g.The stepper motor was set to descend with a distance of 0.1 mm for each step and the corresponding data were recorded until the data shown in the scale were stable.

Experimental Setup Description
Figure 15(a) presents a comparison of the magnetic flux density B z obtained from the model and experiment with a fixed levitation height of 0.9 mm.Clearly, the theoretical values were in good agreement with the experimental results, despite some small errors.These errors may be induced by several factors, such as insufficient magnetization of the magnets, material issues of the magnets, or measurement errors.
Figure 15(b) shows a comparison of the diamagnetic force obtained from the model and experiment.The theoretical diamagnetic force fitted well with the experimental data, indicating the feasibility and applicability of the proposed model.It should also be mentioned that with a levitation height of 0.9 mm, the diamagnetic force obtained in the experiment was 0.881 mN, which was comparable to the gravity of the rotor (0.885 mN, 0.0903 g).This confirmed the feasibility of the model.
The levitation stability of the rotors with various radii (5.19, 6.25, and 8.50 mm) was also evaluated.Pyrolytic graphite with a thickness of 0.5 mm acted as a rotor in the experiments.The results revealed that rotors with radii of 5.19 and 6.25 mm could realize stable levitation, as shown in Figure 16.This corroborated the aforementioned stability prediction that a moderate radius can ensure stable levitation (Figure 13).For the rotor with a radius of 8.5 mm, unstable levitation occurred without the application of an external force.This phenomenon verified the stability prediction that a rotor with a large radius undergoes unstable levitation.
The levitation experiment for the rotor with a radius of 1.5 mm was also attempted.It was found that the rotor was unable to levitate fully, which appeared to contradict the theoretical analysis.However, this could be rationalized by the fact that the B z of a rotor with a small radius was too weak to overcome gravity.Thus, even when the rotor was located at the equilibrium line (the interface between the ring-shaped and cylindrical magnets) and the torque was zero, unstable levitation still occurred.

Conclusions
(1) A systematic mathematical model of a disc-shaped diamagnetically levitated rotor on a permanent magnet array was proposed.The radial and axial magnetic-field distribution characteristics were determined using the proposed model.(2) The expression of the diamagnetic levitation force was derived, and the nonlinear relationship between the diamagnetic levitation force and levitation height was analyzed.As the levitation height increased, the diamagnetic levitation force gradually decreased.Expressions for the axial and radial stiffness were obtained, and the effects of the rotor levitation height and radial position on the stiffness were analyzed.(3) The rotor stability of the diamagnetic levitation system was analyzed.Increasing the levitation height increased the stability of the region.Moreover, with a further increase in the rotor radius, the stable regions of the rotor would gradually diminish and even vanish.Thus, when the levitation height is fixed, a moderate rotor radius permits stable levitation.(4) The relationship between the levitation height and diamagnetic levitation force was verified through experiments, and the operational stability of the rotor under different conditions was evaluated.The experimental results were in good agreement with the theoretical predictions, confirming the feasibility of the mathematical model.

Figure 1
Figure 1 Diamagnetic levitation test rig: (a) Diamagnetic rotor levitated by permanent magnets array, (b) Array of permanent magnets

Figure 6
Figure 6Theoretical magnetic flux density distribution in the plane parallel to the top surface of permanent magnets wherein the distance between the plane and top surface is 0.9 mm

Figure 7
Figure 7 Height affection for the magnetic flux density generated by the array of permanent magnets: (a) B z and (b) B r

Figure 9 Figure 10 Figure 11
Figure 9 Stiffness property of the diamagnetic levitation force (a) at different heights, and (b) at different rotor radial positions with a height of 0.9 mm

Figure 12 Figure 13
Figure 12 Torque of rotor versus levitation height and lateral displacement (rotor radius is 5.19 mm)

Figure 16 Figure 17
Figure 16 Levitation stability test with radii of (a) 5.19 mm and (b) 6.25 mm

Table 1
Parameters for the magnets array