How to simulate permanent magnets using the PDE Toolbox in R2025a
I am trying to model and visualize two permanent magnet spheres using the below code snippet:
function twoMagnets_femodel_demo()
clc; clear; close all;
% ——————————-
% 1) Define geometry parameters
% ——————————-
% Let’s assume "1/2-inch magnets" => diameter = 0.5 inch => radius = 0.25 in
% 1 inch ~ 0.0254 m, so 0.25 in ~ 0.00635 m
R_sphere = 0.00635; % [m]
% Bounding box of ~5 cm each side
boxDim = 0.05; % [m]
% ——————————-
% 2) Build geometry parts
% ——————————-
% (a) Air region: bounding box
boxGm = multicuboid(boxDim, boxDim, boxDim);
% (b) Sphere #1
sphere1Gm = multisphere(R_sphere);
% (c) Sphere #2
sphere2Gm = multisphere(R_sphere);
% Translate spheres so that they’re centered in the box, side-by-side
% Let the box corners go from (0,0,0) to (0.05,0.05,0.05), so center ~ (0.025,0.025,0.025)
centerBox = [0.025, 0.025, 0.025];
offset1 = centerBox + [-R_sphere, 0, 0]; % shift left in x
offset2 = centerBox + [ R_sphere, 0, 0]; % shift right in x
sphere1Gm = translate(sphere1Gm, offset1);
sphere2Gm = translate(sphere2Gm, offset2);
% ——————————-
% 3) Combine geometries
% ——————————-
% addCell(…) is the new PDE Toolbox function that
% keeps each shape as a separate "cell" in the final geometry.
gm = addCell(boxGm, sphere1Gm);
gm = addCell(gm, sphere2Gm);
% Visualize & label
figure;
pdegplot(gm, FaceAlpha=0.2, CellLabels="on");
title("Bounding Box + Two Spheres");
axis equal
% At this point we should see "Cell #1 = bounding box, #2 = first sphere, #3 = second sphere"
% ——————————-
% 4) Create femodel for magnetostatics
% ——————————-
model = femodel(AnalysisType="magnetostatic", Geometry=gm);
% Specify vacuum permeability in SI (H/m)
model.VacuumPermeability = 4*pi*1e-7; % ~ 1.2566370614e-6
% ——————————-
% 5) Set material properties
% ——————————-
% We have 3 cells:
% 1 => bounding box (air)
% 2 => sphere #1
% 3 => sphere #2
%
% For air and the spheres, assume relative permeability ~ 1
model.MaterialProperties(1) = materialProperties(RelativePermeability=1);
model.MaterialProperties(2) = materialProperties(RelativePermeability=1);
model.MaterialProperties(3) = materialProperties(RelativePermeability=1);
% ——————————-
% 6) Assign magnetization
% ——————————-
% In new PDE Toolbox, you can set magnetization on the "CellLoad" for a subdomain (cell).
% Suppose sphere #1 is magnetized in +z; sphere #2 in -z.
% For a typical neodymium magnet, you might set M ~ 1e5 or 1e6 A/m.
M1 = [0; 0; 1e5];
M2 = [0; 0;-1e5];
model.CellLoad(2) = cellLoad(Magnetization=M1); % sphere1
model.CellLoad(3) = cellLoad(Magnetization=M2); % sphere2
% ——————————-
% 7) Boundary conditions
% ——————————-
% Typically, you’d set a "zero" magnetic potential (A = 0) on the outer faces
% of the bounding box to approximate an open boundary.
% Identify the face IDs of the bounding box:
% pdegplot(gm, FaceAlpha=0.1, CellLabels="on", FaceLabels="on")
% The bounding box typically has 6 faces with IDs [1..6].
% So we do:
model.FaceBC(1:6) = faceBC(MagneticPotential=[0;0;0]);
% (If you find different face numbering, adjust accordingly.)
% ——————————-
% 8) Generate mesh
% ——————————-
% We can specify "Hcell" to refine the mesh in certain cells, e.g. the magnets,
% or keep it uniform. For example, refine in cells 2 & 3 (the spheres):
internalCells = [2 3];
model = generateMesh(model, Hcell={internalCells, R_sphere/3});
% ——————————-
% 9) Solve
% ——————————-
R = solve(model);
% R is a MagnetostaticResults object with fields:
% R.MagneticPotential, R.MagneticField, R.MagneticFluxDensity, R.Mesh, etc.
% ——————————-
% 10) Postprocessing
% ——————————-
% Extract Bx, By, Bz at the mesh nodes
Bx = R.MagneticFluxDensity.Bx;
By = R.MagneticFluxDensity.By;
Bz = R.MagneticFluxDensity.Bz;
Bmag = sqrt(Bx.^2 + By.^2 + Bz.^2);
% Simple 3D plot of |B|
figure;
pdeplot3D(R.Mesh, ColorMapData=Bmag);
title("|B| flux density magnitude");
colorbar;
% Optionally, a vector plot:
figure;
pdeplot3D(R.Mesh, FlowData=[Bx By Bz]);
title("Magnetic Field Vectors");
end
While logically complete, the above code doesn’t work. The current error I am struggling with is that addCell() only works inside the original geometry. And throws the following error here:
Error using pde.DiscreteGeometry/addCell
Added cells must be located strictly inside a single cell of the original geometry.
Error in
TwoMagnets (line 38)
gm = addCell(boxGm, sphere1Gm);
^^^^^^^^^^^^^^^^^^^^^^^^^
Is there an example script on how to simulate permanent magnets using the latest PDE Toolbox?I am trying to model and visualize two permanent magnet spheres using the below code snippet:
function twoMagnets_femodel_demo()
clc; clear; close all;
% ——————————-
% 1) Define geometry parameters
% ——————————-
% Let’s assume "1/2-inch magnets" => diameter = 0.5 inch => radius = 0.25 in
% 1 inch ~ 0.0254 m, so 0.25 in ~ 0.00635 m
R_sphere = 0.00635; % [m]
% Bounding box of ~5 cm each side
boxDim = 0.05; % [m]
% ——————————-
% 2) Build geometry parts
% ——————————-
% (a) Air region: bounding box
boxGm = multicuboid(boxDim, boxDim, boxDim);
% (b) Sphere #1
sphere1Gm = multisphere(R_sphere);
% (c) Sphere #2
sphere2Gm = multisphere(R_sphere);
% Translate spheres so that they’re centered in the box, side-by-side
% Let the box corners go from (0,0,0) to (0.05,0.05,0.05), so center ~ (0.025,0.025,0.025)
centerBox = [0.025, 0.025, 0.025];
offset1 = centerBox + [-R_sphere, 0, 0]; % shift left in x
offset2 = centerBox + [ R_sphere, 0, 0]; % shift right in x
sphere1Gm = translate(sphere1Gm, offset1);
sphere2Gm = translate(sphere2Gm, offset2);
% ——————————-
% 3) Combine geometries
% ——————————-
% addCell(…) is the new PDE Toolbox function that
% keeps each shape as a separate "cell" in the final geometry.
gm = addCell(boxGm, sphere1Gm);
gm = addCell(gm, sphere2Gm);
% Visualize & label
figure;
pdegplot(gm, FaceAlpha=0.2, CellLabels="on");
title("Bounding Box + Two Spheres");
axis equal
% At this point we should see "Cell #1 = bounding box, #2 = first sphere, #3 = second sphere"
% ——————————-
% 4) Create femodel for magnetostatics
% ——————————-
model = femodel(AnalysisType="magnetostatic", Geometry=gm);
% Specify vacuum permeability in SI (H/m)
model.VacuumPermeability = 4*pi*1e-7; % ~ 1.2566370614e-6
% ——————————-
% 5) Set material properties
% ——————————-
% We have 3 cells:
% 1 => bounding box (air)
% 2 => sphere #1
% 3 => sphere #2
%
% For air and the spheres, assume relative permeability ~ 1
model.MaterialProperties(1) = materialProperties(RelativePermeability=1);
model.MaterialProperties(2) = materialProperties(RelativePermeability=1);
model.MaterialProperties(3) = materialProperties(RelativePermeability=1);
% ——————————-
% 6) Assign magnetization
% ——————————-
% In new PDE Toolbox, you can set magnetization on the "CellLoad" for a subdomain (cell).
% Suppose sphere #1 is magnetized in +z; sphere #2 in -z.
% For a typical neodymium magnet, you might set M ~ 1e5 or 1e6 A/m.
M1 = [0; 0; 1e5];
M2 = [0; 0;-1e5];
model.CellLoad(2) = cellLoad(Magnetization=M1); % sphere1
model.CellLoad(3) = cellLoad(Magnetization=M2); % sphere2
% ——————————-
% 7) Boundary conditions
% ——————————-
% Typically, you’d set a "zero" magnetic potential (A = 0) on the outer faces
% of the bounding box to approximate an open boundary.
% Identify the face IDs of the bounding box:
% pdegplot(gm, FaceAlpha=0.1, CellLabels="on", FaceLabels="on")
% The bounding box typically has 6 faces with IDs [1..6].
% So we do:
model.FaceBC(1:6) = faceBC(MagneticPotential=[0;0;0]);
% (If you find different face numbering, adjust accordingly.)
% ——————————-
% 8) Generate mesh
% ——————————-
% We can specify "Hcell" to refine the mesh in certain cells, e.g. the magnets,
% or keep it uniform. For example, refine in cells 2 & 3 (the spheres):
internalCells = [2 3];
model = generateMesh(model, Hcell={internalCells, R_sphere/3});
% ——————————-
% 9) Solve
% ——————————-
R = solve(model);
% R is a MagnetostaticResults object with fields:
% R.MagneticPotential, R.MagneticField, R.MagneticFluxDensity, R.Mesh, etc.
% ——————————-
% 10) Postprocessing
% ——————————-
% Extract Bx, By, Bz at the mesh nodes
Bx = R.MagneticFluxDensity.Bx;
By = R.MagneticFluxDensity.By;
Bz = R.MagneticFluxDensity.Bz;
Bmag = sqrt(Bx.^2 + By.^2 + Bz.^2);
% Simple 3D plot of |B|
figure;
pdeplot3D(R.Mesh, ColorMapData=Bmag);
title("|B| flux density magnitude");
colorbar;
% Optionally, a vector plot:
figure;
pdeplot3D(R.Mesh, FlowData=[Bx By Bz]);
title("Magnetic Field Vectors");
end
While logically complete, the above code doesn’t work. The current error I am struggling with is that addCell() only works inside the original geometry. And throws the following error here:
Error using pde.DiscreteGeometry/addCell
Added cells must be located strictly inside a single cell of the original geometry.
Error in
TwoMagnets (line 38)
gm = addCell(boxGm, sphere1Gm);
^^^^^^^^^^^^^^^^^^^^^^^^^
Is there an example script on how to simulate permanent magnets using the latest PDE Toolbox? I am trying to model and visualize two permanent magnet spheres using the below code snippet:
function twoMagnets_femodel_demo()
clc; clear; close all;
% ——————————-
% 1) Define geometry parameters
% ——————————-
% Let’s assume "1/2-inch magnets" => diameter = 0.5 inch => radius = 0.25 in
% 1 inch ~ 0.0254 m, so 0.25 in ~ 0.00635 m
R_sphere = 0.00635; % [m]
% Bounding box of ~5 cm each side
boxDim = 0.05; % [m]
% ——————————-
% 2) Build geometry parts
% ——————————-
% (a) Air region: bounding box
boxGm = multicuboid(boxDim, boxDim, boxDim);
% (b) Sphere #1
sphere1Gm = multisphere(R_sphere);
% (c) Sphere #2
sphere2Gm = multisphere(R_sphere);
% Translate spheres so that they’re centered in the box, side-by-side
% Let the box corners go from (0,0,0) to (0.05,0.05,0.05), so center ~ (0.025,0.025,0.025)
centerBox = [0.025, 0.025, 0.025];
offset1 = centerBox + [-R_sphere, 0, 0]; % shift left in x
offset2 = centerBox + [ R_sphere, 0, 0]; % shift right in x
sphere1Gm = translate(sphere1Gm, offset1);
sphere2Gm = translate(sphere2Gm, offset2);
% ——————————-
% 3) Combine geometries
% ——————————-
% addCell(…) is the new PDE Toolbox function that
% keeps each shape as a separate "cell" in the final geometry.
gm = addCell(boxGm, sphere1Gm);
gm = addCell(gm, sphere2Gm);
% Visualize & label
figure;
pdegplot(gm, FaceAlpha=0.2, CellLabels="on");
title("Bounding Box + Two Spheres");
axis equal
% At this point we should see "Cell #1 = bounding box, #2 = first sphere, #3 = second sphere"
% ——————————-
% 4) Create femodel for magnetostatics
% ——————————-
model = femodel(AnalysisType="magnetostatic", Geometry=gm);
% Specify vacuum permeability in SI (H/m)
model.VacuumPermeability = 4*pi*1e-7; % ~ 1.2566370614e-6
% ——————————-
% 5) Set material properties
% ——————————-
% We have 3 cells:
% 1 => bounding box (air)
% 2 => sphere #1
% 3 => sphere #2
%
% For air and the spheres, assume relative permeability ~ 1
model.MaterialProperties(1) = materialProperties(RelativePermeability=1);
model.MaterialProperties(2) = materialProperties(RelativePermeability=1);
model.MaterialProperties(3) = materialProperties(RelativePermeability=1);
% ——————————-
% 6) Assign magnetization
% ——————————-
% In new PDE Toolbox, you can set magnetization on the "CellLoad" for a subdomain (cell).
% Suppose sphere #1 is magnetized in +z; sphere #2 in -z.
% For a typical neodymium magnet, you might set M ~ 1e5 or 1e6 A/m.
M1 = [0; 0; 1e5];
M2 = [0; 0;-1e5];
model.CellLoad(2) = cellLoad(Magnetization=M1); % sphere1
model.CellLoad(3) = cellLoad(Magnetization=M2); % sphere2
% ——————————-
% 7) Boundary conditions
% ——————————-
% Typically, you’d set a "zero" magnetic potential (A = 0) on the outer faces
% of the bounding box to approximate an open boundary.
% Identify the face IDs of the bounding box:
% pdegplot(gm, FaceAlpha=0.1, CellLabels="on", FaceLabels="on")
% The bounding box typically has 6 faces with IDs [1..6].
% So we do:
model.FaceBC(1:6) = faceBC(MagneticPotential=[0;0;0]);
% (If you find different face numbering, adjust accordingly.)
% ——————————-
% 8) Generate mesh
% ——————————-
% We can specify "Hcell" to refine the mesh in certain cells, e.g. the magnets,
% or keep it uniform. For example, refine in cells 2 & 3 (the spheres):
internalCells = [2 3];
model = generateMesh(model, Hcell={internalCells, R_sphere/3});
% ——————————-
% 9) Solve
% ——————————-
R = solve(model);
% R is a MagnetostaticResults object with fields:
% R.MagneticPotential, R.MagneticField, R.MagneticFluxDensity, R.Mesh, etc.
% ——————————-
% 10) Postprocessing
% ——————————-
% Extract Bx, By, Bz at the mesh nodes
Bx = R.MagneticFluxDensity.Bx;
By = R.MagneticFluxDensity.By;
Bz = R.MagneticFluxDensity.Bz;
Bmag = sqrt(Bx.^2 + By.^2 + Bz.^2);
% Simple 3D plot of |B|
figure;
pdeplot3D(R.Mesh, ColorMapData=Bmag);
title("|B| flux density magnitude");
colorbar;
% Optionally, a vector plot:
figure;
pdeplot3D(R.Mesh, FlowData=[Bx By Bz]);
title("Magnetic Field Vectors");
end
While logically complete, the above code doesn’t work. The current error I am struggling with is that addCell() only works inside the original geometry. And throws the following error here:
Error using pde.DiscreteGeometry/addCell
Added cells must be located strictly inside a single cell of the original geometry.
Error in
TwoMagnets (line 38)
gm = addCell(boxGm, sphere1Gm);
^^^^^^^^^^^^^^^^^^^^^^^^^
Is there an example script on how to simulate permanent magnets using the latest PDE Toolbox? pde, matlab code, mesh MATLAB Answers — New Questions