[Solved!] Mgdyn simulation problem (is it boundary conditions?)

General discussion about Elmer
kevinarden
Posts: 2411
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by kevinarden »

Try this out, the last few lines control mesh size.
sphere.geo
(2.64 KiB) Downloaded 57 times
Nick_99
Posts: 49
Joined: 12 Jul 2023, 10:07
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by Nick_99 »

Thank you so much!! That works like a treat :D
Nick_99
Posts: 49
Joined: 12 Jul 2023, 10:07
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by Nick_99 »

Ok so I do have one small problem. It seems that setting the mesh size parameter for the points constituting the magnets (lm) does not change the mesh size. Instead this mesh size is tied to the mesh size of the points in the sphere. I think it has to do with the

Code: Select all

Characteristic Length{ PointsOf{ Volume{3}; } } = lr;
Do you know of a work around?

Edit: Never mind. Helps to read the manual :lol: Thanks again for all your help (and Rich's help).


Kind regards,
Nick.
kevinarden
Posts: 2411
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by kevinarden »

My understanding is the characteristic length sets the sphere, but the transfinite per curve sets the magnets. However, that overrides the point density mesh size set when creating the points.
Nick_99
Posts: 49
Joined: 12 Jul 2023, 10:07
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by Nick_99 »

So from what I read in the manual, you are correct. I was losing my mind when changing the mesh size wasn't changing the mesh. Thanks for confirming I understand what I read :)
Nick_99
Posts: 49
Joined: 12 Jul 2023, 10:07
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by Nick_99 »

Hi everyone,

I have a new dilemma that I can't quite find a fix for.

Using the geometry that Kevin supplied (and not changing anything) along with the default solver setup for the magdyn solver, I can get a solution that converges correctly and produces results (as Kevin said in an earlier post). However, when I load the results in ParaView I get something like this:
kevin_mesh_default_case.png
(48.6 KiB) Not downloaded yet

But, I'm expecting the Elmer result to have a similar shape to this one (this is a 3D simulation of similarly sized bar magnets in Ansys Maxwell):
Magnetic flux density plot-min.png
(18.5 KiB) Not downloaded yet
The case.sif used:

Code: Select all

Header
  CHECK KEYWORDS Warn
  Mesh DB "." "."
  Include Path ""
  Results Directory ""
End

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Steady state
  Steady State Max Iterations = 1
  Output Intervals(1) = 1
  Coordinate Scaling = 1e-3
  Solver Input File = case.sif
  Post File = case.vtu
End

Constants
  Gravity(4) = 0 -1 0 9.82
  Stefan Boltzmann = 5.670374419e-08
  Permittivity of Vacuum = 8.85418781e-12
  Permeability of Vacuum = 1.25663706e-6
  Boltzmann Constant = 1.380649e-23
  Unit Charge = 1.6021766e-19
End

Body 1
  Target Bodies(1) = 1
  Name = "Body 1"
  Equation = 1
  Material = 2
End

Body 2
  Target Bodies(1) = 2
  Name = "Body 2"
  Equation = 1
  Material = 3
End

Body 3
  Target Bodies(1) = 3
  Name = "Body 3"
  Equation = 1
  Material = 1
End

Solver 1
  Equation = MgDyn
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  Exec Solver = Always
  Stabilize = True
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-5
  Nonlinear System Convergence Tolerance = 1.0e-7
  Nonlinear System Max Iterations = 20
  Nonlinear System Newton After Iterations = 3
  Nonlinear System Newton After Tolerance = 1.0e-3
  Nonlinear System Relaxation Factor = 1
  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStab
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  BiCGstabl polynomial degree = 2
  Linear System Preconditioning = ILU0
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 10
  Linear System Precondition Recompute = 1
End

Solver 2
  Equation = MgDynPost
  Procedure = "MagnetoDynamics" "MagnetoDynamicsCalcFields"
  Calculate Magnetic Field Strength = True
  Exec Solver = Always
  Stabilize = True
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-5
  Nonlinear System Convergence Tolerance = 1.0e-7
  Nonlinear System Max Iterations = 20
  Nonlinear System Newton After Iterations = 3
  Nonlinear System Newton After Tolerance = 1.0e-3
  Nonlinear System Relaxation Factor = 1
  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStab
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  BiCGstabl polynomial degree = 2
  Linear System Preconditioning = ILU0
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 10
  Linear System Precondition Recompute = 1
End

Equation 1
  Name = "Equation 1"
  Active Solvers(2) = 1 2
End

Material 1
  Name = "Vacuum"
  Relative Permeability = 1.25663706e-6
End

Material 2
  Name = "N"
  Relative Permeability = 5000
  Magnetization 1 = Real 750e3
End

Material 3
  Name = "S"
  Magnetization 1 = Real -750e3
  Relative Permeability = 5000
End

Boundary Condition 1
  Target Boundaries(1) = 3 
  Name = "BoundaryCondition 1"
  AV {e} = 0
  AV {e} 1 = 0
  AV {e} 2 = 0
  AV {e} 3 = 0
End
So I though the issue was just the mesh size. I progressively made the mesh larger and larger. However, at a certain point, the solution stopped converging. I ended up using the config Rich supplied in an earlier post but the results pretty much look the same as in the first image. I eventually made a mesh that was ~300MB and ran it on a node on the cluster I have access to. The simulation took about 30 minutes but it did produce a result that converged. However, the plot (seen below) also looks very similar to that in the first image.
kevin_mesh_rich_case.png
(45.98 KiB) Not downloaded yet
The case.sif for the above is:

Code: Select all

Header
  CHECK KEYWORDS Warn
  Mesh DB "." "."
  Include Path ""
  Results Directory "./results"
End

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Steady state
  Steady State Max Iterations = 1
  Output Intervals(1) = 1
  Coordinate Scaling = 1e-3
  Solver Input File = case.sif
  Post File = case.vtu
End

Constants
  Gravity(4) = 0 -1 0 9.82
  Stefan Boltzmann = 5.670374419e-08
  Permittivity of Vacuum = 8.85418781e-12
  Permeability of Vacuum = 1.25663706e-6
  Boltzmann Constant = 1.380649e-23
  Unit Charge = 1.6021766e-19
End

Body 1
  Target Bodies(1) = 1
  Name = "Body 1"
  Equation = 1
  Material = 2
End

Body 2
  Target Bodies(1) = 2
  Name = "Body 2"
  Equation = 1
  Material = 3
End

Body 3
  Target Bodies(1) = 3
  Name = "Body 3"
  Equation = 1
  Material = 1
End

Solver 1
  Equation = MgDyn
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  Exec Solver = Always
  Stabilize = True
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-5
  Nonlinear System Convergence Tolerance = 1.0e-7
  Nonlinear System Max Iterations = 200
  Nonlinear System Newton After Iterations = 3
  Nonlinear System Newton After Tolerance = 1.0e-3
  Nonlinear System Relaxation Factor = 1

  Use Piola Transform = Logical True
  Quadratic Approximation = True

  Linear System Solver = "Iterative"
  Linear System Preconditioning = none
  Linear System Residual Output = 20
  Linear System Max Iterations = 5000
  Linear System Iterative Method = gcr
  Linear System GCR Restart = 100
  Linear System Convergence Tolerance = 1.0e-06
  Steady State Convergence Tolerance = 1e-06
End

Solver 2
  Equation = MgDynPost
  Procedure = "MagnetoDynamics" "MagnetoDynamicsCalcFields"
  Calculate Magnetic Field Strength = True
  Exec Solver = Always
  Stabilize = True
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-5
  Nonlinear System Convergence Tolerance = 1.0e-7
  Nonlinear System Max Iterations = 200
  Nonlinear System Newton After Iterations = 3
  Nonlinear System Newton After Tolerance = 1.0e-3
  Nonlinear System Relaxation Factor = 1
  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStab
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  BiCGstabl polynomial degree = 2
  Linear System Preconditioning = ILU0
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 10
  Linear System Precondition Recompute = 1
End

Equation 1
  Name = "Equation 1"
  Active Solvers(2) = 1 2
End

Material 1
  Name = "Vacuum"
  Relative Permeability = 1.25663706e-6
End

Material 2
  Name = "N"
  Relative Permeability = 5000
  Magnetization 1 = Real 750e3
End

Material 3
  Name = "S"
  Magnetization 1 = Real -750e3
  Relative Permeability = 5000
End

Boundary Condition 1
  Target Boundaries(1) = 3 
  Name = "BoundaryCondition 1"
  AV {e} = 0
End

The geometry file used for both cases

Code: Select all

SetFactory("OpenCASCADE");

// lm corresponds to the mesh size on and inside the magnets' geometry.
// lr specifies the mesh size on and inside the region's geometry.
lm = 1;
lr = 100;

Point(1) = {-50, 7, -5, lm};
Point(2) = {-50, -7, -5, lm};
Point(3) = {-10, -7, -5, lm};
Point(4) = {-10, 7, -5, lm};
Point(5) = {-50, 7, 5, lm};
Point(6) = {-50, -7, 5, lm};
Point(7) = {-10, -7, 5, lm};
Point(8) = {-10, 7, 5, lm};

Line(1) = {5, 8};
Line(2) = {5, 6};
Line(3) = {6, 2};
Line(4) = {2, 1};
Line(5) = {1, 5};
Line(6) = {1, 4};
Line(7) = {4, 8};
Line(8) = {3, 4};
Line(9) = {3, 2};
Line(10) = {3, 7};
Line(11) = {7, 8};
Line(12) = {7, 6};

Curve Loop(1) = {6, -8, 9, 4};
Plane Surface(1) = {1};

Curve Loop(2) = {6, 7, -1, -5};
Plane Surface(2) = {2};

Curve Loop(3) = {7, -11, -10, 8};
Plane Surface(3) = {3};

Curve Loop(4) = {4, 5, 2, 3};
Plane Surface(4) = {4};

Curve Loop(5) = {9, -3, -12, -10};
Plane Surface(5) = {5};

Curve Loop(6) = {1, -11, 12, -2};
Plane Surface(6) = {6};

Surface Loop(1) = {6, 2, 1, 3, 5, 4};
Volume(1) = {1};

Point(9) = {10, 7, 5, lm};
Point(10) = {10, -7, 5, lm};
Point(11) = {10, -7, -5, lm};
Point(12) = {10, 7, -5, lm};
Point(13) = {50, 7, -5, lm};
Point(14) = {50, 7, 5, lm};
Point(15) = {50, -7, 5, lm};
Point(16) = {50, -7, -5, lm};

Line(13) = {11, 16};
Line(14) = {16, 13};
Line(15) = {13, 12};
Line(16) = {12, 11};
Line(17) = {11, 10};
Line(18) = {12, 9};
Line(19) = {13, 14};
Line(20) = {16, 15};
Line(21) = {15, 14};
Line(22) = {10, 9};
Line(23) = {9, 14};
Line(24) = {15, 10};

Curve Loop(7) = {13, 14, 15, 16};
Plane Surface(7) = {7};
Curve Loop(8) = {15, 18, 23, -19};
Plane Surface(8) = {8};
Curve Loop(9) = {19, -21, -20, 14};
Plane Surface(9) = {9};
Curve Loop(10) = {17, -24, -20, -13};
Plane Surface(10) = {10};
Curve Loop(11) = {16, 17, 22, -18};
Plane Surface(11) = {11};
Curve Loop(12) = {22, 23, -21, 24};
Plane Surface(12) = {12};

Surface Loop(2) = {8, 7, 10, 11, 12, 9};
Volume(2) = {2};

//+
Physical Volume("Left", 25) = {1};
Physical Surface("sleft", 26) = {1, 2, 3, 4, 5, 6};
//+
Physical Volume("Right", 27) = {2};
Physical Surface("sright", 28) = {7, 8, 9, 10, 11, 12};
//+
Sphere(3) = {0, 0, 0, 1000, -Pi/2, Pi/2, 2*Pi};
//+
Physical Volume("Sphere", 29) = {3};
Physical Surface("Outer", 30) = {13};
BooleanDifference{ Volume{3}; Delete; }{ Volume{2}; Volume{1}; }
//+
//+
Coherence;
//+
//+
Characteristic Length{ PointsOf{ Volume{3}; } } = lr;
//+
Transfinite Curve {1, 9, 6, 12, 13, 15, 24, 23} = 10 Using Progression 1;
//+
Transfinite Curve {3, 5, 2, 10, 7, 11, 18, 17, 20, 14, 19, 21, 22, 16, 8, 4} = 4 Using Progression 1;
I've been fiddling with other aspects of the solver's setup but have not managed to get anything resembling the output from Maxwell.

Does anyone have some advice for me? Am I maybe not using ParaView correctly? I loaded the pvtu in ParaView, used the 'Plot over line' option and positioned the line along the x-axis straight through both bar magnets (and the airgap between them).

Kind regards,
Nick
raback
Site Admin
Posts: 4849
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by raback »

Hi,

Try projecting to elemental fields. Maybe you can share the full case. I'll PM my email.

-Peter
Rich_B
Posts: 431
Joined: 24 Aug 2009, 20:18

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by Rich_B »

Hello,

Attached is another version of the sif file, this time adding SaveLine. This will save the data along the x-axis, similar to what you described for Paraview.

The elemental data in red looks much smoother than the nodal data in blue. Notice the red line inside the magnets has a value of 750,000, which matches the magnetization set in materials 2 and 3.

Rich.
FieldStrength.png
FieldStrength.png (65.61 KiB) Viewed 704 times
Attachments
case-saveline.sif
(3.26 KiB) Downloaded 41 times
fjimenez
Posts: 63
Joined: 27 Sep 2021, 23:40
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by fjimenez »

Hi,

Try modifying your materials' parameters:
- Relative Permeability of air is 1
- Relative Permeability of a permanent magnet is approx 1.05

Cheers,
Nick_99
Posts: 49
Joined: 12 Jul 2023, 10:07
Antispam: Yes

Re: Mgdyn simulation problem (is it boundary conditions?)

Post by Nick_99 »

Hi everyone,

Thank you all so much for your assistance :D

@fjimenez Thank you for pointing out my error with the permeability! I had not noticed that it was too high until you mentioned it.

@Rich_B Thank you for the new sif file. I ran your file on my computer and got the exact same results in Excel as you did. The results in ParaView are also identical. This gives me a little bit more confidence that I'm using ParaView correctly. Sadly, I'm still a bit suspicious of the "spikes" in the results. I have a working Elmer simulation for 2 bar magnets in 2D that produces this result.
2d_case_elmer.png
(38.92 KiB) Not downloaded yet
I also cross checked these results with Ansys Maxwell and got a very similar plot.
2d_case_maxwell.png
(38.65 KiB) Not downloaded yet
So, I'm expecting the 3D results to have a similar shape (i.e. a shape without "spikes") to the 2D results.

I have just sent @raback an email. If we get a solution, I'll be sure to post it here :-)
Post Reply