PDF
Magnetic Brake
Introduction
A magnet brake in its simplest form consists of a disc of conductive material and a permanent magnet. The magnet generates a constant magnetic field, in which the disc is rotating. When a conductor moves in a magnetic field it induces currents, and the Lorentz forces from the currents slow the disc.
Model Definition
This 3D problem is solved using a stationary formulation for the electromagnetic part coupled to an ordinary differential equation for the rotational rigid body dynamics. It also illustrates the proper use of a Lorentz type induced current density term. An ungauged A-V formulation is used to solve the electromagnetic part in a fast and memory efficient way.
For a disc rotating with angular velocity ω about the z-axis, the velocity v at a point (x, y) is given by
Maxwell-Ampère’s law, expressed using a magnetic vector potential A, a scalar electric potential V, and an induced current density term of Lorentz type is the fundamental field equation used in this model. It is complemented by a current balance but no explicit gauge is provided.
Such a formulation is inherently singular because the divergence of the magnetic vector potential is not uniquely determined and there is an infinite number of possible solutions to the problem. All such solutions, however, yield the same magnetic flux and current densities. It is further known that the above formulation, when solved using an iterative solver, converges quickly and robustly provided that any explicitly supplied source current density (in this case, zero) is divergence free and no direct coarse grid solver is involved in the solution scheme.
The induced Lorentz current density term is a common source of confusion in electromagnetic modeling. In situations when the moving domain is of bounded extent in the direction of the motion or varies in this direction or contains magnetic sources that also move, Lorentz terms cannot be used. This is because the part of the magnetic flux generated by moving sources must not be included in the Lorentz term. In this situation, the induced current distribution is stationary (it stays where the magnet is and does not move with the spinning disc). Thus, the moving domain does not contain any magnetic sources that move along with it and it is unbounded and invariant in the direction of the motion.
The magnetic and electric boundary conditions on external boundaries are
Now consider how the system evolves over time. The induced torque slows the disc down, described by an ordinary differential equation (ODE) for the angular velocity ω.
Where the torque Tz is obtained as the z component of the vector
and the moment of inertia I for a disc with radius r and unit thickness equals
Here m is the disc mass and ρ is the density.
Results and Discussion
The COMSOL Multiphysics model is set up for a 1 cm thick copper disc with a radius of 10 cm and an initial angular speed of 1000 rpm. The magnet consists of a 1 T, hard (μr = 1) permanent magnet connected to an iron yoke with a 1.5 cm air gap in which the copper disc spins. The figures below show the induced eddy current density and the time evolution of the angular velocity, braking torque and dissipated power respectively.
Figure 1: The eddy current magnitude and direction at t = 0 s.
.
Figure 2: The time evolution of the angular velocity.
Figure 3: The time evolution of the braking torque.
Figure 4: The time evolution of the dissipated power.
Application Library path: ACDC_Module/Motors_and_Actuators/magnetic_brake
Modeling Instructions
From the File menu, choose New.
New
In the New window, click Model Wizard.
Model Wizard
1
In the Model Wizard window, click 3D.
2
In the Select Physics tree, select AC/DC>Electromagnetic Fields>Vector Formulations>Magnetic and Electric Fields (mef).
3
Click Add.
4
In the Select Physics tree, select Mathematics>ODE and DAE Interfaces>Global ODEs and DAEs (ge).
5
Click Add.
6
Click Study.
7
In the Select Study tree, select General Studies>Stationary.
8
Click Done.
Global Definitions
Parameters 1
1
In the Model Builder window, under Global Definitions click Parameters 1.
2
In the Settings window for Parameters, locate the Parameters section.
3
In the table, enter the following settings:
Name
Expression
Value
Description
dt
1[cm]
0.01 m
Disc thickness
dr
10[cm]
0.1 m
Disc radius
mh
12[cm]
0.12 m
Magnet height
mw
2[cm]
0.02 m
Magnet width
ml
8[cm]
0.08 m
Magnet length
mt
2[cm]
0.02 m
Magnet thickness
mg
1.5[cm]
0.015 m
Magnet air gap
mB
1[T]
1 T
Magnet flux
ymur
4000
4000
Yoke relative permeability
dV0
1000[rpm]
16.667 1/s
Disc initial angular frequency
Geometry 1
Sphere 1 (sph1)
1
In the Geometry toolbar, click Sphere.
2
In the Settings window for Sphere, locate the Size section.
3
In the Radius text field, type dr*3.
Cylinder 1 (cyl1)
1
In the Geometry toolbar, click Cylinder.
2
Right-click Cylinder 1 (cyl1) and choose Rename.
3
In the Rename Cylinder dialog box, type Disc in the New label text field.
4
Click OK.
5
In the Settings window for Cylinder, locate the Size and Shape section.
6
In the Radius text field, type dr.
7
In the Height text field, type dt.
8
Locate the Position section. In the z text field, type -dt/2.
9
Locate the Selections of Resulting Entities section. Select the Resulting objects selection check box.
A selection named Disc is then automatically created for the disc.
Now, design the magnet and the yoke on a plane, and then extrude them.
Work Plane 1 (wp1)
1
In the Geometry toolbar, click Work Plane.
2
In the Settings window for Work Plane, locate the Plane Definition section.
3
From the Plane list, choose yz-plane.
4
In the x-coordinate text field, type -mt/2.
5
Click Show Work Plane.
Work Plane 1 (wp1)>Rectangle 1 (r1)
1
In the Work Plane toolbar, click Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type mw.
4
In the Height text field, type mh-2*mw.
5
Locate the Position section. From the Base list, choose Center.
6
In the xw text field, type dr+ml/2-mw/2.
Work Plane 1 (wp1)>Rectangle 2 (r2)
1
In the Work Plane toolbar, click Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type ml-mw.
4
In the Height text field, type mw.
5
Locate the Position section. In the xw text field, type dr-ml/2+mw.
6
In the yw text field, type mh/2-mw.
Work Plane 1 (wp1)>Rectangle 3 (r3)
1
In the Work Plane toolbar, click Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type mw.
4
In the Height text field, type mh/2-mg/2.
5
Locate the Position section. In the xw text field, type dr-ml/2.
6
In the yw text field, type mg/2.
Work Plane 1 (wp1)>Mirror 1 (mir1)
1
In the Work Plane toolbar, click Transforms and choose Mirror.
2
Click the Zoom Extents button in the Graphics toolbar.
3
Select the objects r2 and r3 only.
4
In the Settings window for Mirror, locate the Input section.
5
Select the Keep input objects check box.
6
Locate the Normal Vector to Line of Reflection section. In the xw text field, type 0.
7
In the yw text field, type 1.
8
Right-click Mirror 1 (mir1) and choose Build Selected.
The geometry on the work plane should now look like the figure below.
9
In the Model Builder window, click Geometry 1.
Extrude 1 (ext1)
1
In the Geometry toolbar, click Extrude.
2
In the Settings window for Extrude, locate the Distances section.
3
In the table, enter the following settings:
Distances (m)
mt
4
Click Build All Objects.
Hide the exterior boundaries to show the yoke and the disc only.
5
In the Graphics toolbar, click the Click and Hide button.
6
In the Graphics toolbar, click the Select Boundaries button.
7
Select the exterior boundaries of the sphere (boundaries 1–8 on the object sph1).
8
In the Graphics toolbar, click the Click and Hide button again to deactivate it.
9
Click the Go to Default View button in the Graphics toolbar.
The model geometry is now complete.
Define an integration operator on the disc to compute the total dissipated power.
Definitions
Integration 1 (intop1)
1
In the Definitions toolbar, click Nonlocal Couplings and choose Integration.
2
In the Settings window for Integration, type intdisc in the Operator name text field.
3
Select Domain 2 only.
4
Right-click Integration 1 (intop1) and choose Rename.
5
In the Rename Integration dialog box, type Integration over Disc in the New label text field.
6
Click OK.
Add Material
1
In the Home toolbar, click Add Material to open the Add Material window.
2
Go to the Add Material window.
3
In the tree, select Built-in>Air.
4
Click Add to Component in the window toolbar.
5
In the tree, select Built-in>Copper.
6
Click Add to Component in the window toolbar.
7
In the Home toolbar, click Add Material to close the Add Material window.
Materials
Copper (mat2)
1
In the Settings window for Material, locate the Geometric Entity Selection section.
2
From the Selection list, choose Disc.
Material 3 (mat3)
1
In the Model Builder window, right-click Materials and choose Blank Material.
2
Select Domains 3–6 only.
3
In the Settings window for Material, locate the Material Contents section.
4
In the table, enter the following settings:
Property
Variable
Value
Unit
Property group
Relative permeability
mur_iso ; murii = mur_iso, murij = 0
ymur
1
Basic
Electrical conductivity
sigma_iso ; sigmaii = sigma_iso, sigmaij = 0
0
S/m
Basic
Relative permittivity
epsilonr_iso ; epsilonrii = epsilonr_iso, epsilonrij = 0
1
1
Basic
5
In the Label text field, type Yoke.
6
In the Home toolbar, click Windows and choose Add Material from Library.
Add Material
1
Go to the Add Material window.
2
In the tree, select AC/DC>Hard Magnetic Materials>Sintered NdFeB Grades (Chinese Standard)>N50 (Sintered NdFeB).
3
Click Add to Component 1 (comp1).
4
In the Home toolbar, click Add Material to close the Add Material window.
Materials
N50 (Sintered NdFeB) (mat4)
1
Select Domain 7 only.
2
In the Settings window for Material, locate the Material Contents section.
3
In the table, enter the following settings:
Property
Variable
Value
Unit
Property group
Electrical conductivity
sigma_iso ; sigmaii = sigma_iso, sigmaij = 0
1
S/m
Basic
Recoil permeability
murec_iso ; murecii = murec_iso, murecij = 0
1
1
Remanent flux density
Remanent flux density norm
normBr
mB
T
Remanent flux density
4
In the Label text field, type Generic magnet.
Add a Mass Properties node to compute the moment of inertia of the disc.
5
Click the Show More Options button in the Model Builder toolbar.
6
In the Show More Options dialog box, in the tree, select the check box for the node General>Variable Utilities.
7
Click OK.
Definitions
Mass Properties 1 (mass1)
1
In the Model Builder window, under Component 1 (comp1) right-click Definitions and choose Variable Utilities>Mass Properties.
2
In the Settings window for Mass Properties, locate the Source Selection section.
3
Click Clear Selection.
4
From the Selection list, choose Disc.
5
Locate the Variables section. From the Frame list, choose Spatial  (x, y, z).
Use the density of the material (in this case, copper) to compute the moment of inertia.
6
Locate the Density section. In the Density expression text field, type mat2.def.rho.
7
Clear all check boxes except Create moment of inertia variables.
The Mass Properties node will define the moment of inertia variables mass1.Ixx, mass1.Ixy, etc.
Magnetic and Electric Fields (mef)
1
In the Model Builder window, under Component 1 (comp1) click Magnetic and Electric Fields (mef).
2
In the Settings window for Magnetic and Electric Fields, click to expand the Equation section.
3
From the Equation form list, choose Stationary to solve a stationary problem for the magnetic field.
Ampère’s Law and Current Conservation 1
1
In the Model Builder window, under Component 1 (comp1)>Magnetic and Electric Fields (mef) click Ampère’s Law and Current Conservation 1.
2
In the Settings window for Ampère’s Law and Current Conservation, locate the Constitutive Relation Jc-E section.
3
From the σ list, choose User defined. In the associated text field, type 1.
Ampère’s Law and Current Conservation 2
1
In the Physics toolbar, click Domains and choose Ampère’s Law and Current Conservation.
2
Select Domain 7 only.
3
In the Settings window for Ampère’s Law and Current Conservation, locate the Constitutive Relation B-H section.
4
From the Magnetization model list, choose Remanent flux density.
5
Specify the e vector as
0
x
0
y
1
z
6
In the Label text field, type Permanent Magnet.
Ampère’s Law and Current Conservation 3
1
In the Physics toolbar, click Domains and choose Ampère’s Law and Current Conservation.
2
In the Settings window for Ampère’s Law and Current Conservation, locate the Domain Selection section.
3
From the Selection list, choose Disc.
4
Right-click Ampère’s Law and Current Conservation 3 and choose Rename.
5
In the Rename Ampère’s Law and Current Conservation dialog box, type Disc in the New label text field.
6
Click OK.
Velocity (Lorentz Term) 1
1
In the Physics toolbar, click Domains and choose Velocity (Lorentz Term).
2
In the Settings window for Velocity (Lorentz Term), locate the Domain Selection section.
3
From the Selection list, choose Disc.
4
Locate the Velocity (Lorentz Term) section. Specify the v vector as
-y*W
x
x*W
y
0
z
Specifying the variable W as angular velocity (that will be used later as a state variable in the ODE) couples the Magnetic and Electric Fields interface with the ODE.
Global ODEs and DAEs (ge)
Global Equations 1
1
In the Model Builder window, under Component 1 (comp1)>Global ODEs and DAEs (ge) click Global Equations 1.
2
In the Settings window for Global Equations, locate the Global Equations section.
3
Specify the differential equation that describes the dynamics of the disc in the table. Torque is computed by integrating the Lorentz force contribution.
Name
f(u,ut,utt,t) (1)
Initial value (u_0) (1)
Description
W
Wt-intdisc(x*mef.FLtzy-y*mef.FLtzx)/mass1.Izz
2*pi*dV0
Angular Velocity
4
Locate the Units section. Click Select Dependent Variable Quantity.
5
In the Physical Quantity dialog box, type angularfrequency in the text field.
6
Click Filter.
7
In the tree, select General>Angular frequency (rad/s).
8
Click OK.
9
In the Settings window for Global Equations, locate the Units section.
10
Click Select Source Term Quantity.
11
In the Physical Quantity dialog box, type angularacceleration in the text field.
12
Click Filter.
13
In the tree, select General>Angular acceleration (rad/s^2).
14
Click OK.
This is the differential equation that describes the dynamics of the disc.
Definitions
Variables 1
1
In the Model Builder window, under Component 1 (comp1) right-click Definitions and choose Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
axialTorque
-intdisc(x*mef.FLtzy-y*mef.FLtzx)
N·m
Axial torque
totalHeating
intdisc(mef.Qh)
W
Total heating
Global Variable Probe 1 (var1)
1
In the Definitions toolbar, click Probes and choose Global Variable Probe.
2
In the Settings window for Global Variable Probe, click to expand the Table and Window Settings section.
3
From the Plot window list, choose New window.
Global Variable Probe 2 (var2)
1
In the Definitions toolbar, click Probes and choose Global Variable Probe.
2
In the Settings window for Global Variable Probe, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1)>Definitions>Variables>axialTorque - Axial torque - N·m.
3
Locate the Table and Window Settings section. From the Plot window list, choose New window.
Global Variable Probe 3 (var3)
1
In the Definitions toolbar, click Probes and choose Global Variable Probe.
2
In the Settings window for Global Variable Probe, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1)>Definitions>Variables>totalHeating - Total heating - W.
3
Locate the Table and Window Settings section. From the Plot window list, choose New window.
Mesh 1
Size 1
In the Model Builder window, under Component 1 (comp1) right-click Mesh 1 and choose Size.
Size
1
In the Settings window for Size, locate the Element Size section.
2
From the Predefined list, choose Coarser.
Size 1
1
In the Model Builder window, click Size 1.
2
In the Settings window for Size, locate the Element Size section.
3
From the Predefined list, choose Extremely fine.
4
Locate the Geometric Entity Selection section. From the Geometric entity level list, choose Domain.
5
Select Domains 2–4 only.
Free Tetrahedral 1
1
In the Model Builder window, right-click Mesh 1 and choose Free Tetrahedral.
2
In the Settings window for Free Tetrahedral, click to expand the Scale Geometry section.
3
In the z-direction scale text field, type 1.1.
4
Click Build All.
Study 1
Step 1: Stationary
The stationary solver will compute initial values for the time dependent solver.
1
In the Model Builder window, under Study 1 click Step 1: Stationary.
2
In the Settings window for Stationary, locate the Physics and Variables Selection section.
3
In the table, clear the Solve for check box for the Global ODEs and DAEs interface.
Time Dependent
1
In the Study toolbar, click Study Steps and choose Time Dependent>Time Dependent.
2
In the Settings window for Time Dependent, locate the Study Settings section.
3
In the Times text field, type range(0,0.5,25).
4
In the Model Builder window, click Study 1.
5
In the Settings window for Study, locate the Study Settings section.
6
Clear the Generate default plots check box.
Solution 1 (sol1)
1
In the Study toolbar, click Show Default Solver.
For a time-dependent problem, the default error estimate will be too conservative and may result in stagnation of the iterative solver. To overcome this, solver tollerances are tuned.
2
In the Model Builder window, expand the Solution 1 (sol1) node, then click Time-Dependent Solver 1.
3
In the Settings window for Time-Dependent Solver, click to expand the Absolute Tolerance section.
4
In the Variables list, select Angular Velocity (comp1.ODE1).
5
From the Method list, choose Unscaled.
6
From the Tolerance method list, choose Manual.
7
In the Tolerance text field, type 0.1.
8
Click to expand the Time Stepping section. From the Steps taken by solver list, choose Intermediate.
9
Find the Algebraic variable settings subsection. From the Error estimation list, choose Exclude algebraic.
10
In the Model Builder window, expand the Study 1>Solver Configurations>Solution 1 (sol1)>Time-Dependent Solver 1 node.
11
Right-click Time-Dependent Solver 1 and choose Fully Coupled.
12
In the Settings window for Iterative, click to expand the Error section.
13
In the Factor in error estimate text field, type 1.
14
In the Study toolbar, click Compute.
Results
During the computation, the angular velocity as a function of time has been plotted together with torque and dissipated power. Give a proper label and compare with those reported in documentation.
Probe Plot Group 1
1
In the Model Builder window, under Results click Probe Plot Group 1.
2
In the Settings window for 1D Plot Group, type Angular velocity VS time in the Label text field.
Probe Plot Group 2
1
In the Model Builder window, under Results click Probe Plot Group 2.
2
In the Settings window for 1D Plot Group, type Axial torque VS time in the Label text field.
Probe Plot Group 3
1
In the Model Builder window, under Results click Probe Plot Group 3.
2
In the Settings window for 1D Plot Group, type Total heating VS time in the Label text field.
3D Plot Group 4
1
In the Home toolbar, click Add Plot Group and choose 3D Plot Group.
2
In the Settings window for 3D Plot Group, locate the Plot Settings section.
3
Clear the Plot dataset edges check box.
Volume 1
1
Right-click 3D Plot Group 4 and choose Volume.
2
In the Settings window for Volume, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1>Magnetic and Electric Fields>Currents and charge>mef.normJ - Current density norm - A/m².
Selection 1
1
Right-click Volume 1 and choose Selection.
2
Select Domain 2 only.
Volume 2
1
In the Model Builder window, right-click 3D Plot Group 4 and choose Volume.
2
In the Settings window for Volume, locate the Expression section.
3
In the Expression text field, type 1.
4
Locate the Coloring and Style section. From the Coloring list, choose Uniform.
5
From the Color list, choose Gray.
Selection 1
1
Right-click Volume 2 and choose Selection.
2
Select Domains 3–7 only.
Arrow Surface 1
1
In the Model Builder window, right-click 3D Plot Group 4 and choose Arrow Surface.
2
In the Settings window for Arrow Surface, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1>Magnetic and Electric Fields>Currents and charge>mef.Jx,mef.Jy,mef.Jz - Current density.
3
Locate the Arrow Positioning section. In the Number of arrows text field, type 1000.
4
Locate the Coloring and Style section. From the Color list, choose White.
3D Plot Group 4
1
In the Model Builder window, click 3D Plot Group 4.
2
In the Settings window for 3D Plot Group, locate the Data section.
3
From the Time (s) list, choose 0.
4
In the 3D Plot Group 4 toolbar, click Plot.
5
Click the Zoom Extents button in the Graphics toolbar.
The plot shows the current density on the surface of the disc at t = 0 (when the disc is still spinning).
6
From the Time (s) list, choose 25.
7
In the 3D Plot Group 4 toolbar, click Plot.
Now, the plot shows the current density when the disc has almost stopped.