A Halbach array is a specialized arrangement of permanent magnets designed to generate a strong, uniform magnetic field in the designated region. This configuration has gained increasing attention in recent years due to its potential applications in high-efficiency electric machines, magnetic levitation systems, and advanced electromagnetic devices. Halbach arrays can be classified into two main types: linear Halbach arrays and circular Halbach arrays. A linear Halbach array consists of a sequence of permanent magnets with a rotating magnetization pattern that enhances the magnetic field on one side while canceling it on the other. A circular Halbach array arranges the magnets in a cylindrical configuration, not only concentrating the magnetic field within the inner region while minimizing the external field but also generating a desired magnetic field in the target region. The circular Halbach array has a wide range of potential applications across various fields, including Magnetic Resonance Imaging (MRI) and Nuclear Magnetic Resonance (NMR) spectroscopy [23], particle trapping and manipulation, wireless power transfer (WPT) and energy applications, high-precision sensors and metrology, as well as magnetic levitation and transport. The structure of the electromagnetic system is usually so complex and sophisticated that its design process has been highly dependent on the engineer's experience and intuition [27], the permanent magnet arrangement with different magnetization directions also presents a challenge.
To address these challenges, researchers have increasingly turned to computational methods for the systematic design of these complex systems. Topology optimization has emerged as a powerful tool for designing structures with optimized material distribution to achieve desired performance objectives [5,6], enabling the creation of multi-material, high-performance, and multifunctional electromagnetic structures. With the growing emphasis on electromagnetic development, topology optimization has been increasingly applied to a wide range of electromagnetic applications. The design of high-performance and powerdense electric machines necessitates a thorough exploration of the design space to identify optimal configurations [12]. Topology optimization has been widely utilized in electric machine design. Credo et al. [9] employed topology optimization to design a high-speed synchronous reluctance motor for electric vehicles. Tian et al. [32] applied a level set-based topology optimization method to design a generator, aiming to maximize magnetic energy while minimizing deviations from the target magnetic field. Additionally, Tian et al.
[31] utilized level set-based topology optimization to design a synchronous reluctance motor with the objective of maximizing average torque while minimizing torque ripple. Jung et al. [19] employed a density-based multi-material topology optimization approach for the design of magnetic systems, with the material interpolation scheme to distinguish different materials-permanent magnets (PM), ferromagnetic materials, and air-based on their respective material properties.
The design of ferromagnetic materials can also be utilized in the development of flexible magnetism-sensitive robots, where the magnetization direction in soft active materials is tailored to achieve the desired performance [21,30,38]. Additionally, ferromagnetic material design plays a crucial role in the advancement of compliant electronics [10,41] and bionic medical devices [8,42].
The rational design of Halbach arrays and the magnetic flux control present significant challenges due to the complex interplay between material properties and magnetization. Choi et al. [7] employed a parameter optimization method based on finite element analysis, defining the magnetization direction of each element as a design variable. Lee et al. [23] proposed an isoparametric projection method to optimize both the permanent magnet strength and magnetization directions in a Halbach array simultaneously. Their approach incorporates penalization and extrusion schemes to achieve different desired angle sets. Furthermore, Lee et al. [22] extended the three-dimensional design of segmented permanent magnet arrays by applying topology optimization with a density-based method.
Structural optimization-based methods leverage advanced computational techniques to achieve optimal designs. These methods include the ground structure method [11], the homogenization method [26], the solid isotropic material with penalization (SIMP) method [29], and level set-based approaches [4,36]. Recent advancements in this domain have introduced innovative techniques such as the flexible building block method [20], the unit cell approach [33], the radial basis function (RBF) levelset method [37], and the cardinal basis function (CBF) level set method [18].
Compared to the Solid Isotropic Material with Penalization (SIMP) method, the level-set approach enables the direct evolution of the design boundary, ensuring a well-defined and sharp boundary in the final design while reducing the need for postprocessing. The ability of the level-set method to produce distinct boundaries without transitional grey regions between different materials makes it particularly advantageous for multimaterial topology optimization problems [34]. However, a major issue of conventional level-set approaches is the necessity of a reinitialization scheme to maintain numerical stability [17,18]. During boundary evolution, the level-set function can become excessively flat or steep, leading to numerical instability. To solve this issue, the level set function must be periodically regularized into a signed distance function (SDF) to enhance numerical stability and improve the accuracy of material property mapping. This reinitialization process disrupts the optimization procedure at regular intervals, increasing computational costs and requiring additional iterations before convergence to the final design [37]. To address this issue, we propose a Cardinal Basis Function (CBF)-based level-set method for the design of a circular Halbach array. By parameterizing the level-set function using a CBF representation, the original Hamilton-Jacobi partial differ-ential equation (PDE) is reformulated into a system of ordinary differential equations (ODEs), reducing computational costs and accelerating convergence. In this framework, the conventional reinitialization process is replaced by minimizing a distance regularization energy function, thereby eliminating the need for periodic reinitialization steps. The optimization is carried out using the Method of Moving Asymptotes (MMA), a gradient-based optimization algorithm that is good for topology optimization. The proposed CBF-based level set approach significantly reduces the number of convergence steps and improves the computational efficiency in the topology optimization process.
In the rest of this paper, section 2 focuses on the modeling of the circular Halbach array, while Section 3 is dedicated to the topology optimization of the circular Halbach array. This includes an introduction to the conventional level-set method, the Cardinal Basis Function (CBF) parametric level-set method, problem formulation, design sensitivity analysis, and the solution of the adjoint equations. Section 4 presents numerical results, followed by discussions and conclusions in the final section.
A circular Halbach array is an arrangement of permanent magnets in a cylindrical configuration with a central void, as illustrated in Figure 1. In this figure, the red regions represent the permanent magnets, while the light blue regions represent the air. The permanent magnet section is typically composed of multiple PM segments, each with different magnetization directions. Increasing the number of magnet segments provides additional control over the magnetic flux distribution within the designated region.
Figure 2 presents two classical configurations of circular Halbach arrays: the dipole and quadrupole arrangements. In the dipole configuration, the magnetic flux within the central region follows a dipole pattern, primarily aligned in the vertical direction. In contrast, the quadrupole configuration generates a quadrupole field, producing a strong magnetic field with alternating north and south poles around the circumference. These two configurations demonstrate how the arrangement of permanent magnets in a circular Halbach array influences the magnetic flux distribution within the designated central region.
In this design example, the permanent magnet used is the N54 sintered NdFeB (neodymium-iron-boron) magnet, with a remanent flux density norm of 1.47 T and a recoil permeability of 1.05. The geometric parameters of the design, as depicted in Figure 1, include a central hole with a diameter of 0.1 m and a surrounding permanent magnet ring with a diameter of 0.2 m. The outer air region is modeled as a square with a side length of 0.4 m, subject to a Dirichlet boundary condition, where the magnetic vector potential is set to A = 0. This model is formulated as a 2.5-D representation, with a thickness of 0.05 m. The central hole, outlined by a dashed line in Figure 1, represents the measurement region used to evaluate the total magnetic flux. Additionally, the shaded area in the upper-right portion of the model in Figure 1 defines the design domain, where topology optimization is applied. The simulation and finite element analysis (FEA) in this study are conducted using COMSOL Multiphysics. The CBF-based level-set topology optimization is implemented using MATLAB.
Topology optimization, a shape optimization method, employs algorithmic models to optimize material distribution within a predefined design domain, considering specified objective functions, constraints, and boundary conditions. In recent years, topology optimization has garnered increasing popularity and attention within engineering design circles. Its scope has expanded significantly to address a wide array of challenges involving multiphysics coupling, spanning electromagnetics, thermodynamics, acoustics, solid mechanics, and fluid mechanics, among others. The level set method, devised by S. Osher and J.A. Sethian, employs a high-dimensional function to implicitly represent the 2-D contour. Pioneered by Sethian and Wiegmann [1] and further refined by Wang [35, 37] and Allaire [4], this method has emerged as a promising approach for shape and topology optimization. It ensures a clear boundary between phases without a grey region, significantly enhancing precision and optimization accuracy. In this framework, the structural boundary is implicitly represented as the 2-D contour of a level set function with one higher dimension. Implicitly embedded within the level set function Φ(x,t). Depending on the sign of the level set function, the design domain can be partitioned into three distinct regions, representing the material, the interface, and the void, respectively, as follows:
where D denotes the design domain. The evolution of the boundary dynamics is governed by the Hamilton-Jacobi equation:
The normal velocity field V n can be determined through shape sensitivity analysis. Solving the Hamilton-Jacobi equa-tion outlined above enables the updating of the normal velocity field, which subsequently governs the evolution of the structural boundary.
Set Method This study presents a parametric level set method formulated using Cardinal Basis Functions (CBFs). The CBF is constructed via a radial basis function (RBF) partition of unity (RBF-POU) collocation approach, ensuring that each basis function assumes a value of one at its center node and zero at all other nodes within the design domain [28]. An RBF is radially symmetric and can form a basis set of kernel functions for describing the level set function. Given N different nodes x 1 , . . . , x N ∈ R d , a level set function Φ(x,t) can be interpolated as selected CBF kernel function in the Equation 3, the details of formulate the parametric level set function can be found in citation [18].
In Equation ( 3), where µ j (t) is actually the value of the level set function on the j th node. Ψ j (x) is the constructed CBF, which is equal to 1 at the center node and 0 at other nodes. The CBFs Ψ j (x), j = 1, . . . , N have the Kronecker delta property as described as follows:
A distance-regularized level set function is preferable throughout the optimization process to maintain numerical stability and avoid the need for frequent reinitialization. For a given level set function Φ, the distance regularization energy functional R is defined as:
where P(|∇Φ|) is the regularization energy potential density, ∇Φ is the gradient of the level set function, and |∇Φ| is the norm of ∇Φ.
To achieve a signed distance property, the norm of the gradient |∇Φ| of the level set function should reach 1. To maintain a flat surface, |∇Φ| needs to be equal to 0. Based on this idea, the "double-well" regularization energy potential can be formulated as follows [24]:
After parameterizing the level-set function using a CBF representation, the original Hamilton-Jacobi partial differential equation (PDE) is reformulated into a system of ordinary differential equations (ODEs) [18]:
Where the normal velocity V n can be obtained as: In this topology optimization problem, the objective is to maximize the area-averaged vertical directional magnetic flux within the central hole, which serves as the measurement region (depicted by the dashed line in Figure 1). By maximizing the area-averaged vertical directional magnetic flux in this region, the design ensures the generation of a uniform magnetic field in the vertical direction. In this design example, the magnetization direction is not considered a design variable. Instead, two fixed sets of magnetization angles are predefined, and the optimization focuses solely on shaping the permanent magnet distribution within the designated design domain. The objective function J is formulated as follows:
The region Ω 3 is shown in Figure 4. The regions Ω 1 and Ω 2 represent the permanent magnets within the design domain. The B represents the magnetic flux density, and the A represents the vector potential. The U represents the space of admissible vector potentials, where the arbitrary vector potential Ā belongs. The weak form governing equation of the magnetostatic system is given by: a(A, Ā) = l( Ā). The energy bilinear form a(A, Ā) and source linear form l( Ā), which will be discussed in Section 3.4.
The Γ is the Dirichlet essential boundary, and H 1 (Ω) is the Sobolev space of order one [3]. The magnetic field intensity is written with the magnetic flux density and the permanent magnetization by the constitutive relation;
Where the H is the magnetic field intensity, the µ is magnetic permeability. The M 0 is the permanent magnetization of the permanent magnet. The magnetization model of the permanent magnet used in this example is the remanent flux density. The equation of magnetic flux density B and the remanent flux density B r is shown below:
The remanent flux density B r is a vector defined as the product of the remanent flux density norm and the remanent flux direction. The unit vector e is the remanent flux direction. The parameter µ 0 represents the vacuum magnetic permeability constant. Similarly, µ rec is a material constant that characterizes the magnetic properties of permanent magnets.
This section presents the methodology for conducting shape sensitivity analysis. The material time derivative is introduced for computing shape sensitivity.
A general form of the objective function is defined as:
Where the Ω 1 and Ω 2 represent the PM region, and the Ω 3 represents the Air region. The k is the localizing factor used to define the integral domain for the objective function J. The schematic of the magnetostatic system interface is shown in Figure 5.
The weak form of the governing equation of the magnetostatic system is shown below:
Where the J represents the current density and M 0 represents the permanent magnetization in each material.
To address this problem, the Lagrange multiplier method is employed to reformulate the PDE-constrained problem as an unconstrained optimization problem. This is achieved by defining the Lagrangian functional L as follows:
The material time derivative is employed to derive the shape sensitivity formulation:
The derivative of the Lagrangian function is presented directly as follows. The interface of two permanent magnets with different magnetization when the permanent magnet 1 is Ω 1 and the permanent magnet 2 is Ω 2 .
Where the B r2 and B r1 are the remanent flux density of the two permanent magnet, the γ represents the interface between the permanent magnet 1 and permanent magnet 2, the Ā2 represents the adjoint variable can b solving the adjoint equation which will be discussed in section 3.5.
By taking the equation 12 into the equation 18, we can get:
Since the remanent flux density B r are the same for different permanent magnet, the equation 19 can be rewritten as:
By combining with the equation 8, the equation 21 can rewrite as:
(22) With the chain rule, we have:
The sensitivity of the objective function F can be formulated by combining the equation 22 and 23 as follows:
) Since e is a unit vector representing the remanent flux direction in each permanent magnet, it is also predefined in this design example. The general form of e 1 and e 2 can be written as:
In this design example, since the magnetization angle of the two permanent magnets is redefined, θ 1 is defined as π 3 and θ 2 is defined as 5π 3 . By plugging Equation 25 into Equation 24, we have:
It is noted that the boundary integration in the above equation can be converted into a domain integration by using the Dirac delta function δ as:
In addition, the sensitivity of the distance regularization energy functional R can be derived as follows:
where d p is defined as [24]:
where p ′ (s) denotes the first derivative of the regularization energy potential density, as defined in Equation 6. Consequently, the shape sensitivity of the total objective function J, which couples the distance regularization energy function R with the original objective function F, can be expressed as a domain integral in the following form:
The weak form of the adjoint equation calculated is shown as follows:
where A m is the area of the measurement region, which is a scalar.
So the equation 32 can be rewritten as the final adjoint equation:
Copyright © 2025 by ASME This final adjoint equation can solve the adjoint variable Ā2 .
Considering symmetry, only a quarter of the circular Halbach array is used as the design domain. In this design example, the objective is to maximize the area-averaged vertical magnetic flux within the central hole, aiming to generate a uniform vertical magnetic field in that region. The permanent magnet employed is an N54 sintered NdFeB (neodymium-iron-boron) magnet with a remanent flux density norm of 1.47 T and a recoil permeability of 1.05. FIGURE 6. Conformal mapping from the irregular design domain to a standard rectangular redesign domain [14,16]. This transformation simplifies the level set function representation, facilitating efficient topology optimization.
The entire design domain is discretized into a 76×101 grid. To facilitate level set representation, conformal mapping theory is applied to parameterize the irregular 2D triangular-meshed design domain onto a structured quadrilateral-meshed rectangular domain [2,13,15,25,30,39,40], as illustrated in Figure 6. This transformation simplifies the definition of the level set function. The optimized design is subsequently mapped back to the original irregular domain for finite element analysis and verification.
Two distinct permanent magnets are defined within the design domain, differentiated by a Heaviside function, with each possessing a different magnetization direction. The topology optimization process is capped at 100 iterations; however, convergence is achieved within just 10 iterations using the CBF-based level set method. Figure 7 presents the topology optimization evolution, where the blue region represents permanent magnet 1, the orange region represents permanent magnet 2, and the black and white arrows represent their corresponding magnetization direction. Figure 8 illustrates the evolution of the level set function. Figure 9 provides the objective function curve, the regularization energy plot, and the volume ratio of the two permanent magnets. The optimization follows the objective function given in Equation 9, resulting in a final volume ratio of 48.241% for permanent magnet 1 and 51.759% for permanent magnet 2.
The reference design is depicted in Figure 10, where each permanent magnet occupies an equal volume and standard shape. In contrast, Figure 11 illustrates the topology-optimized design of the circular Halbach array. In this optimized configuration, the individual permanent magnets do not maintain equal volumes; however, the total volume remains consistent with the reference design. The optimized design consists of four distinct types of permanent magnets, each possessing the same remanent flux density norm but differing in magnetization direction.
After the topology optimization, numerical verification is performed to evaluate the effectiveness of the optimized design. vertical direction) within the central hole compared to the reference design.
By evaluating the performance, as formulated in Equation 9, of the optimized circular Halbach array, the results are summarized in Table 1. The total magnetic flux density in the vertical direction within the central hole is computed as 1.4199 Tesla for the reference design, whereas for the optimized design, it increases to 1.5118 Tesla. This corresponds to a 6.47% improvement over the reference design, demonstrating that the optimization successfully meets the design objective.
The final extruded 2.5D prototype of the optimized circular Halbach array is shown in Figure 13 with a thickness of 0.05m. 5118 T within the central hole, representing a 6.47% improvement compared to 1.4199 T in the reference design, as summarized in Table 1.
In this paper, we propose a cardinal basis function (CBF)based level set method for designing a circular Halbach array with the objective of generating a uniform vertical magnetic field within the central hole. The cardinal basis function (CBF)-based topology optimization significantly reduces computational cost and accelerates convergence. The reinitialization scheme of the conventional level set method is replaced by incorporating the regularization energy function, which enhances optimization stability.
The core contribution of this work lies in integrating level set-based topology optimization with Halbach array design. Traditionally, engineers design Halbach arrays primarily through parameter optimization, which heavily relies on experience and intuition. In contrast, topology optimization provides a rigorous and mathematically grounded approach, enabling the systematic design of electromagnetic systems with high accuracy and precision.
Current state-of-the-art methods for Halbach array optimization primarily utilize density-based approaches and parameter optimization. However, a review of the existing literature revealed no prior work employing level set-based topology optimization for Halbach array design. The level set method offers distinct advantages, including the ability to define clear and welldefined boundaries, effectively handle topological changes, and facilitate multi-material optimization. In future work, we plan to incorporate additional materials into the design process, allowing for greater control over the magnetic flux distribution.
In this study, the design variable is restricted to the shape of the permanent magnets, without consideration of their magnetization direction. The optimized circular Halbach array maintains the same total volume of permanent magnet material as the reference design, yet achieves a 6.47% improvement in performance. Future research will extend the approach to include magnetization direction as an additional design variable, enabling simultaneous optimization of both shape and magnetization. This enhancement will provide greater design flexibility and further improve performance.
Copyright © 2025 by ASME electromagnetic systems. Springer, 2019. [28] Ali Safdari-Vaighani, Alfa Heryudono, and Elisabeth Larsson. A radial basis function partition of unity collocation method for convection-diffusion equations arising in financial applications. Journal of Scientific Computing, 64(2):341-367, 2015. [29] Ole Sigmund. On the design of compliant mechanisms using topology optimization. Journal of Structural Mechanics, 25(4):493-524, 1997. [30] Jiawei Tian, Manqi Li, Zhonghao Han, Yong Chen, Xianfeng David Gu, QJ Ge, and Shikui Chen. Conformal topology optimization of multi-material ferromagnetic soft active structures using an extended level set method. Computer Methods in Applied Mechanics and Engineering, 389:114394, 2022. [31] Jiawei Tian, David Torrey, Fang Luo, Jon Longtin, and Shikui Chen. Robust topology optimization of synchronous reluctance motors using cardinal basis function based level set method. In International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, volume 87318, page V03BT03A041. American Society of Mechanical Engineers, 2023. [32] Jiawei Tian, Ran Zhuang, Juan Cilia, Arvind Rangarajan, Fang Luo, Jon Longtin, and Shikui Chen. Topology optimization of permanent magnets for generators using level set methods. In International design engineering technical conferences and computers and information in engineering conference, volume 86236, page V03BT03A037. American Society of Mechanical Engineers, 2022. [33] Hongqing Vincent Wang. A unit cell approach for lightweight structure and compliant mechanism. Georgia Institute of Technology, 2005. [34] Michael Yu Wang and Xiaoming Wang. "color" level sets: a multi-phase method for structural topology optimization with multiple materials. Computer methods in applied mechanics and engineering, 193(6-8):469-496, 2004. [35] Michael Yu Wang and Xiaoming Wang. "color" level sets: a multi-phase method for structural topology optimization with multiple materials. Computer Methods in Applied Mechanics and Engineering, 193(6):469-496, 2004. [36] Michael Yu Wang, Xiaoming Wang, and Dongming Guo. A level set method for structural topology optimization. Computer methods in applied mechanics and engineering, 192(1-2):227-246, 2003. [37] Shengyin Wang and Michael Yu Wang. Radial basis functions and level set method for structural topology optimization. International journal for numerical methods in engineering, 65(12):2060-2090, 2006. [38] Michael Wehner, Ryan L Truby, Daniel J Fitzgerald, Bobak Mosadegh, George M Whitesides, Jennifer A Lewis, and Robert J Wood. An integrated design and fabrication strategy for entirely soft, autonomous robots. nature,