MagnetoDynamics; WhitneyAVSolver; Uniform B-field

Numerical methods and mathematical models of Elmer
Post Reply
peavinepeak
Posts: 8
Joined: 06 May 2020, 02:06
Antispam: Yes

MagnetoDynamics; WhitneyAVSolver; Uniform B-field

Post by peavinepeak » 19 May 2020, 03:41

Hello everyone,

I am a new user of Elmer and large numerical simulations in general. I have worked through some of the magnetic field calculation tutorials and read through relevant parts of the solver and the models manual and other slides on the Elmer website.

I work on precision physics experiments and need to look into magnetic field distribution in the vicinity of our sensing atoms.
There are some non-magnetic parts of varying permeability near the atoms and they are all immersed in a uniform magnetic field (created by a superconducting solenoid). My interest is in magnetic flux density and its gradients around the non-magnetic parts. I do not need to model the setup with the solenoid itself at the moment.

I decided to try a simple problem first. I take a cylinder filled with air and apply uniform flux density boundary conditions. For simplicity, I am looking at the X,Y, and Z components of the magnetic flux density along the cylindrical axis (the z-axis). The result is shown in the image files below. However, all the simulations I have tried give me inaccurate magnetic flux densities.

This problem is identical to the one discussed by
sarz:
viewtopic.php?f=3&t=3422&p=11290&hilit= ... 561#p11290
and
szewro:
viewtopic.php?f=3&t=6609&p=20589&hilit= ... eld#p20589

I tried 3 different ways of specifying boundary conditions, a combination from sarz's and szewro's topics.

Case 1:
SIF file:

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
  Timestepping Method = BDF
  BDF Order = 1
  Solver Input File = case.sif
  Post File = case.vtu
End

Constants
  Gravity(4) = 0 -1 0 9.82
  Stefan Boltzmann = 5.67e-08
  Permittivity of Vacuum = 8.8542e-12
  Boltzmann Constant = 1.3807e-23
  Unit Charge = 1.602e-19
End

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

Solver 1
  Equation = MgDyn
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  Exec Solver = Always
  Stabilize = True
  Bubbles = False
  Lumped Mass Matrix = False
  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"
  Exec Solver = Before Saving
  Stabilize = True
  Bubbles = False
  Lumped Mass Matrix = False
  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 = "Coupled equations"
  Active Solvers(2) = 1 2
End

Material 1
  Name = "Air"
  Relative Permeability = 1.0
End

Boundary Condition 1
  Target Boundaries(3) = 3 2 1 
  Name = "Uniform_flux_density"
  Magnetic Flux Density 2 = Real 0.0
  Magnetic Flux Density 3 = Real 1.0
  Magnetic Flux Density 1 = Real 0.0
End
Result:
Solver log:

Code: Select all

ELMER SOLVER (v 8.4) STARTED AT: 2020/05/18 15:49:37
ParCommInit:  Initialize #PEs:            1
MAIN: 
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Webpage http://www.csc.fi/elmer, Email elmeradm@csc.fi
MAIN: Version: 8.4 (Rev: unknown, Compiled: 2020-02-28)
MAIN:  Running one task without MPI parallelization.
MAIN:  Running with just one thread per task.
MAIN:  HYPRE library linked in.
MAIN:  MUMPS library linked in.
MAIN: =============================================================
MAIN: 
MAIN: 
MAIN: -------------------------------------
MAIN: Reading Model: case.sif
LoadInputFile: Scanning input file: case.sif
LoadInputFile: Loading input file: case.sif
LoadInputFile: Number of BCs: 1
LoadInputFile: Number of Body Forces: 0
LoadInputFile: Number of Initial Conditions: 0
LoadInputFile: Number of Materials: 1
LoadInputFile: Number of Equations: 1
LoadInputFile: Number of Solvers: 2
LoadInputFile: Number of Bodies: 1
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_Init0]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_Init0]
LoadMesh: Base mesh name: ./.
LoadMesh: Elapsed REAL time:     0.2679 (s)
MAIN: -------------------------------------
AddVtuOutputSolverHack: Adding ResultOutputSolver to write VTU output in file: case
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_Init]
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_bulk]
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver]
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdyn...done.
OptimizeBandwidth: Half bandwidth without optimization: 34122
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 2690
OptimizeBandwidth: ---------------------------------------------------------
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_Init]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_bulk]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields]
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdynpost...done.
OptimizeBandwidth: Half bandwidth without optimization: 5289
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 689
OptimizeBandwidth: ---------------------------------------------------------
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_Init]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_bulk]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver]
ElmerSolver: Number of timesteps to be saved: 1
MAIN: 
MAIN: -------------------------------------
MAIN:  Steady state iteration:            1
MAIN: -------------------------------------
MAIN: 
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: mgdyn
WhitneyAVSolver: Solving the AV equations with edge elements
MGDynAssembly: Elapsed REAL time:     0.5577 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.02
      10 0.2341E-02
      20 0.2018E-03
      30 0.5390E-05
      40 0.2806E-06
      50 0.6461E-08
      60 0.1097E-09
      61 0.8760E-10
ComputeChange: NS (ITER=1) (NRM,RELC): (  8.8287480      2.0000000     ) :: mgdyn
MGDynAssembly: Elapsed REAL time:     0.5669 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.02
      10 0.3128E-02
      20 0.4472E-03
      30 0.1099E-04
      40 0.1610E-06
      50 0.1511E-07
      60 0.3043E-10
      60 0.3043E-10
ComputeChange: NS (ITER=2) (NRM,RELC): (  6.6285411     0.28468211     ) :: mgdyn
MGDynAssembly: Elapsed REAL time:     0.5539 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.02
       1 0.2629E-10
ComputeChange: NS (ITER=3) (NRM,RELC): (  6.6285411     0.15173373E-11 ) :: mgdyn
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_post]
ComputeChange: SS (ITER=1) (NRM,RELC): (  6.6285411      2.0000000     ) :: mgdyn
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: mgdynpost
MagnetoDynamicsCalcFields: Computing postprocessed fields
WARNING:: GetPermittivity: Permittivity not defined in material, defaulting to that of vacuum
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:     73770
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.00
       5 0.6772E-12
ComputeChange: NS (ITER=1) (NRM,RELC): (  3.9438957      2.0000000     ) :: mgdynpost
       5 0.6344E-12
ComputeChange: NS (ITER=2) (NRM,RELC): (  3.8298566     0.29339508E-01 ) :: mgdynpost
       5 0.7986E-12
ComputeChange: NS (ITER=3) (NRM,RELC): (  3.9227380     0.23961358E-01 ) :: mgdynpost
MagnetoDynamicsCalcFields:  Eddy current power:    0.0000000000000000
MagnetoDynamicsCalcFields:  (Electro)Magnetic Field Energy:    12217953433.834257
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_post]
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: internalvtuoutputsolver
ResultOutputSolver: -------------------------------------
ResultOutputSolve: Saving with prefix: case
ResultOutputSolver: Creating list for saving - if not present
CreateListForSaving: Field Variables for Saving
ResultOutputSolver: Saving in unstructured VTK XML (.vtu) format
VtuOutputSolver: Saving results in VTK XML format with prefix: case
VtuOutputSolver: Saving number of partitions: 1
WriteVtuFile: Writing variable: magnetic flux density e
ResultOutputSolver: -------------------------------------
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_post]
ElmerSolver: *** Elmer Solver: ALL DONE ***
ElmerSolver: The end
SOLVER TOTAL TIME(CPU,REAL):         5.90        6.16
ELMER SOLVER FINISHED AT: 2020/05/18 15:49:43
Result from Paraview:
cylinder_air_case1_paraview.png
(135.53 KiB) Not downloaded yet


Case 2:
This is identical to sarz's example.
Note: Fix Input Current Density = True

SIF file:

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
  Timestepping Method = BDF
  BDF Order = 1
  Solver Input File = case.sif
  Post File = case.vtu
End

Constants
  Gravity(4) = 0 -1 0 9.82
  Stefan Boltzmann = 5.67e-08
  Permittivity of Vacuum = 8.8542e-12
  Boltzmann Constant = 1.3807e-23
  Unit Charge = 1.602e-19
End

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

Solver 1
  Equation = MgDyn
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  Fix Input Current Density = True
  Exec Solver = Always
  Stabilize = True
  Bubbles = False
  Lumped Mass Matrix = False
  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"
  Exec Solver = Before Saving
  Stabilize = True
  Bubbles = False
  Lumped Mass Matrix = False
  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 = "Coupled equations"
  Active Solvers(2) = 1 2
End

Material 1
  Name = "Air"
  Relative Permeability = 1.0
End

Boundary Condition 1
  Target Boundaries(1) = 2 
  Name = "Uniform_flux_pos_z_axis"
  Magnetic Flux Density {n} = Real 1.0
End

Boundary Condition 2
  Target Boundaries(1) = 3 
  Name = "Uniform_flux_neg_z_axis"
  Magnetic Flux Density {n} = Real -1.0
End

Boundary Condition 3
  Target Boundaries(1) = 1 
  Name = "Uniform_flux_XY_axis"
  Magnetic Flux Density {n} = Real 0.0
End

Result:

Solver log:

Code: Select all

ELMER SOLVER (v 8.4) STARTED AT: 2020/05/18 16:20:03
ParCommInit:  Initialize #PEs:            1
MAIN: 
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Webpage http://www.csc.fi/elmer, Email elmeradm@csc.fi
MAIN: Version: 8.4 (Rev: unknown, Compiled: 2020-02-28)
MAIN:  Running one task without MPI parallelization.
MAIN:  Running with just one thread per task.
MAIN:  HYPRE library linked in.
MAIN:  MUMPS library linked in.
MAIN: =============================================================
MAIN: 
MAIN: 
MAIN: -------------------------------------
MAIN: Reading Model: case.sif
LoadInputFile: Scanning input file: case.sif
LoadInputFile: Loading input file: case.sif
Model Input:  Unlisted keyword: [magnetic flux density {n}] in section: [boundary condition 1]
Model Input:  Unlisted keyword: [magnetic flux density {n}] in section: [boundary condition 2]
Model Input:  Unlisted keyword: [magnetic flux density {n}] in section: [boundary condition 3]
LoadInputFile: Number of BCs: 3
LoadInputFile: Number of Body Forces: 0
LoadInputFile: Number of Initial Conditions: 0
LoadInputFile: Number of Materials: 1
LoadInputFile: Number of Equations: 1
LoadInputFile: Number of Solvers: 2
LoadInputFile: Number of Bodies: 1
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_Init0]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_Init0]
LoadMesh: Base mesh name: ./.
LoadMesh: Elapsed REAL time:     0.2631 (s)
MAIN: -------------------------------------
AddVtuOutputSolverHack: Adding ResultOutputSolver to write VTU output in file: case
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_Init]
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_bulk]
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver]
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdyn...done.
OptimizeBandwidth: Half bandwidth without optimization: 34122
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 2690
OptimizeBandwidth: ---------------------------------------------------------
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_Init]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_bulk]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields]
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdynpost...done.
OptimizeBandwidth: Half bandwidth without optimization: 5289
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 689
OptimizeBandwidth: ---------------------------------------------------------
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_Init]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_bulk]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver]
ElmerSolver: Number of timesteps to be saved: 1
MAIN: 
MAIN: -------------------------------------
MAIN:  Steady state iteration:            1
MAIN: -------------------------------------
MAIN: 
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: mgdyn
WhitneyAVSolver: Solving the AV equations with edge elements
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdyn...done.
OptimizeBandwidth: Half bandwidth without optimization: 5289
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 689
OptimizeBandwidth: ---------------------------------------------------------
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WARNING:: JfixPotentialSolver: No Dirichlet conditions used to define Jfix level!
SolveSystem: Solution trivially zero!
MGDynAssembly: Elapsed REAL time:     0.6739 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.02
      10 0.2534E-02
      20 0.1903E-03
      30 0.2295E-05
      40 0.2687E-07
      50 0.6972E-09
      54 0.8833E-10
ComputeChange: NS (ITER=1) (NRM,RELC): (  8.8266964      2.0000000     ) :: mgdyn
MGDynAssembly: Elapsed REAL time:     0.6513 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.01
      10 0.2329E-02
      20 0.3384E-03
      30 0.6711E-04
      40 0.9014E-06
      50 0.9067E-06
      60 0.6132E-10
      60 0.6132E-10
ComputeChange: NS (ITER=2) (NRM,RELC): (  6.6300115     0.28423709     ) :: mgdyn
MGDynAssembly: Elapsed REAL time:     0.6474 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.01
       1 0.5636E-10
ComputeChange: NS (ITER=3) (NRM,RELC): (  6.6300115     0.10563946E-10 ) :: mgdyn
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_post]
ComputeChange: SS (ITER=1) (NRM,RELC): (  6.6300115      2.0000000     ) :: mgdyn
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: mgdynpost
MagnetoDynamicsCalcFields: Computing postprocessed fields
WARNING:: GetPermittivity: Permittivity not defined in material, defaulting to that of vacuum
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:     73770
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.00
       5 0.6771E-12
ComputeChange: NS (ITER=1) (NRM,RELC): (  3.9454737      2.0000000     ) :: mgdynpost
       5 0.6345E-12
ComputeChange: NS (ITER=2) (NRM,RELC): (  3.8303861     0.29601252E-01 ) :: mgdynpost
       5 0.7984E-12
ComputeChange: NS (ITER=3) (NRM,RELC): (  3.9232352     0.23949874E-01 ) :: mgdynpost
MagnetoDynamicsCalcFields:  Eddy current power:    0.0000000000000000
MagnetoDynamicsCalcFields:  (Electro)Magnetic Field Energy:    12222836596.278994
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_post]
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: internalvtuoutputsolver
ResultOutputSolver: -------------------------------------
ResultOutputSolve: Saving with prefix: case
ResultOutputSolver: Creating list for saving - if not present
CreateListForSaving: Field Variables for Saving
ResultOutputSolver: Saving in unstructured VTK XML (.vtu) format
VtuOutputSolver: Saving results in VTK XML format with prefix: case
VtuOutputSolver: Saving number of partitions: 1
WriteVtuFile: Writing variable: magnetic flux density e
ResultOutputSolver: -------------------------------------
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_post]
ElmerSolver: *** Elmer Solver: ALL DONE ***
ElmerSolver: The end
SOLVER TOTAL TIME(CPU,REAL):         6.38        6.65
ELMER SOLVER FINISHED AT: 2020/05/18 16:20:10
Result from Paraview:
cylinder_air_case2_paraview.png
(188.52 KiB) Not downloaded yet


Case 3:

SIF file:

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
  Timestepping Method = BDF
  BDF Order = 1
  Solver Input File = case.sif
  Post File = case.vtu
End

Constants
  Gravity(4) = 0 -1 0 9.82
  Stefan Boltzmann = 5.67e-08
  Permittivity of Vacuum = 8.8542e-12
  Boltzmann Constant = 1.3807e-23
  Unit Charge = 1.602e-19
End

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

Solver 2
  Equation = MgDynPost
  Procedure = "MagnetoDynamics" "MagnetoDynamicsCalcFields"
  Exec Solver = Before Saving
  Stabilize = True
  Bubbles = False
  Lumped Mass Matrix = False
  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 1
  Equation = MgDyn
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  Exec Solver = Always
  Stabilize = True
  Bubbles = False
  Lumped Mass Matrix = False
  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 = "Coupled equations"
  Active Solvers(2) = 2 1
End

Material 1
  Name = "Air"
  Relative Permeability = 1.0
End

Boundary Condition 1
  Target Boundaries(1) = 2 
  Name = "Uniform_flux_pos_z_axis"
  Magnetic Flux Density {n} = Real 1.0
End

Boundary Condition 2
  Target Boundaries(1) = 3 
  Name = "Uniform_flux_neg_z_axis"
  Magnetic Flux Density {n} = Real -1.0
End

Boundary Condition 3
  Target Boundaries(1) = 1 
  Name = "Uniform_flux_XY_axis"
  Magnetic Flux Density {n} = Real 0.0
End

Result:
Solver log:

Code: Select all

ELMER SOLVER (v 8.4) STARTED AT: 2020/05/18 16:29:43
ParCommInit:  Initialize #PEs:            1
MAIN: 
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Webpage http://www.csc.fi/elmer, Email elmeradm@csc.fi
MAIN: Version: 8.4 (Rev: unknown, Compiled: 2020-02-28)
MAIN:  Running one task without MPI parallelization.
MAIN:  Running with just one thread per task.
MAIN:  HYPRE library linked in.
MAIN:  MUMPS library linked in.
MAIN: =============================================================
MAIN: 
MAIN: 
MAIN: -------------------------------------
MAIN: Reading Model: case.sif
LoadInputFile: Scanning input file: case.sif
LoadInputFile: Loading input file: case.sif
Model Input:  Unlisted keyword: [magnetic flux density {n}] in section: [boundary condition 1]
Model Input:  Unlisted keyword: [magnetic flux density {n}] in section: [boundary condition 2]
Model Input:  Unlisted keyword: [magnetic flux density {n}] in section: [boundary condition 3]
LoadInputFile: Number of BCs: 3
LoadInputFile: Number of Body Forces: 0
LoadInputFile: Number of Initial Conditions: 0
LoadInputFile: Number of Materials: 1
LoadInputFile: Number of Equations: 1
LoadInputFile: Number of Solvers: 2
LoadInputFile: Number of Bodies: 1
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_Init0]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_Init0]
LoadMesh: Base mesh name: ./.
LoadMesh: Elapsed REAL time:     0.2603 (s)
MAIN: -------------------------------------
AddVtuOutputSolverHack: Adding ResultOutputSolver to write VTU output in file: case
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_Init]
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_bulk]
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver]
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdyn...done.
OptimizeBandwidth: Half bandwidth without optimization: 34122
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 2690
OptimizeBandwidth: ---------------------------------------------------------
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_Init]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_bulk]
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields]
OptimizeBandwidth: ---------------------------------------------------------
OptimizeBandwidth: Computing matrix structure for: mgdynpost...done.
OptimizeBandwidth: Half bandwidth without optimization: 5289
OptimizeBandwidth: 
OptimizeBandwidth: Bandwidth Optimization ...done.
OptimizeBandwidth: Half bandwidth after optimization: 689
OptimizeBandwidth: ---------------------------------------------------------
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_Init]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_bulk]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver]
ElmerSolver: Number of timesteps to be saved: 1
MAIN: 
MAIN: -------------------------------------
MAIN:  Steady state iteration:            1
MAIN: -------------------------------------
MAIN: 
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: mgdyn
WhitneyAVSolver: Solving the AV equations with edge elements
MGDynAssembly: Elapsed REAL time:     0.5689 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.01
      10 0.2534E-02
      20 0.1903E-03
      30 0.2295E-05
      40 0.2687E-07
      50 0.6972E-09
      54 0.8833E-10
ComputeChange: NS (ITER=1) (NRM,RELC): (  8.8266964      2.0000000     ) :: mgdyn
MGDynAssembly: Elapsed REAL time:     0.5700 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.01
      10 0.2329E-02
      20 0.3384E-03
      30 0.6711E-04
      40 0.9014E-06
      50 0.9067E-06
      60 0.6132E-10
      60 0.6132E-10
ComputeChange: NS (ITER=2) (NRM,RELC): (  6.6300115     0.28423709     ) :: mgdyn
MGDynAssembly: Elapsed REAL time:     0.5649 (s)
DefUtils::DefaultDirichletBCs: Setting Dirichlet boundary conditions
DefUtils::DefaultDirichletBCs: Dirichlet boundary conditions set
WhitneyAVSolver::  Boundary tree edges: 2136 of total: 6405
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:    527062
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.01
       1 0.5636E-10
ComputeChange: NS (ITER=3) (NRM,RELC): (  6.6300115     0.10563946E-10 ) :: mgdyn
Loading user function library: [MagnetoDynamics]...[WhitneyAVSolver_post]
ComputeChange: SS (ITER=1) (NRM,RELC): (  6.6300115      2.0000000     ) :: mgdyn
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: mgdynpost
MagnetoDynamicsCalcFields: Computing postprocessed fields
WARNING:: GetPermittivity: Permittivity not defined in material, defaulting to that of vacuum
CRS_IncompleteLU: ILU(0) (Real), Starting Factorization:
CRS_IncompleteLU: ILU(0) (Real), NOF nonzeros:     73770
CRS_IncompleteLU: ILU(0) (Real), filling (%) :       100
CRS_IncompleteLU: ILU(0) (Real), Factorization ready at (s):     0.00
       5 0.6771E-12
ComputeChange: NS (ITER=1) (NRM,RELC): (  3.9454737      2.0000000     ) :: mgdynpost
       5 0.6345E-12
ComputeChange: NS (ITER=2) (NRM,RELC): (  3.8303861     0.29601252E-01 ) :: mgdynpost
       5 0.7984E-12
ComputeChange: NS (ITER=3) (NRM,RELC): (  3.9232352     0.23949874E-01 ) :: mgdynpost
MagnetoDynamicsCalcFields:  Eddy current power:    0.0000000000000000
MagnetoDynamicsCalcFields:  (Electro)Magnetic Field Energy:    12222836596.278994
Loading user function library: [MagnetoDynamics]...[MagnetoDynamicsCalcFields_post]
SingleSolver: Attempting to call solver
SingleSolver: Solver Equation string is: internalvtuoutputsolver
ResultOutputSolver: -------------------------------------
ResultOutputSolve: Saving with prefix: case
ResultOutputSolver: Creating list for saving - if not present
CreateListForSaving: Field Variables for Saving
ResultOutputSolver: Saving in unstructured VTK XML (.vtu) format
VtuOutputSolver: Saving results in VTK XML format with prefix: case
VtuOutputSolver: Saving number of partitions: 1
WriteVtuFile: Writing variable: magnetic flux density e
ResultOutputSolver: -------------------------------------
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_post]
ElmerSolver: *** Elmer Solver: ALL DONE ***
ElmerSolver: The end
SOLVER TOTAL TIME(CPU,REAL):         5.73        6.02
ELMER SOLVER FINISHED AT: 2020/05/18 16:29:49
Result from Paraview:
cylinder_air_case3_paraview.png
(123.78 KiB) Not downloaded yet
I prepare the mesh file in Gmsh 4.5.6 and use Paraview 5.4.1 along with Elmer 8.4 in Ubuntu 18.04.4. elmerfem-csc-eg package was installed using apt-get. No convergence window appears while running the simulation, perhaps there is a bug.

I converted the .msh file using : ElmerGrid 14 2 mesh.msh -autoclean.

Any comment is very much appreciated.

Thank you.

kevinarden
Posts: 618
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: MagnetoDynamics; WhitneyAVSolver; Uniform B-field

Post by kevinarden » 19 May 2020, 12:33

Judging by the graph you cylinder is 12 meters long by 12 meter diameter?

kevinarden
Posts: 618
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: MagnetoDynamics; WhitneyAVSolver; Uniform B-field

Post by kevinarden » 19 May 2020, 12:51

This is what I get if it is 5 mm radius and 10 mm long
mflux.png
(102.73 KiB) Not downloaded yet
case.sif
(2.92 KiB) Downloaded 124 times

peavinepeak
Posts: 8
Joined: 06 May 2020, 02:06
Antispam: Yes

Re: MagnetoDynamics; WhitneyAVSolver; Uniform B-field

Post by peavinepeak » 19 May 2020, 19:25

Hello Kevin,

Thank you for the prompt reply.

You are using the following boundary conditions (from your SIF file):
Boundary Condition 1
Target Boundaries(1) = 3
Name = "Uniform_flux_pos_z_axis"
Magnetic Flux Density 3 = Real 1.0
End

Boundary Condition 2
Target Boundaries(1) = 1
Name = "Uniform_flux_neg_z_axis"
Magnetic Flux Density 3 = Real -1.0
End

Boundary Condition 3
Target Boundaries(1) = 2
Name = "Uniform_flux_XY_axis"
Magnetic Flux Density 1 = Real 0.0
Magnetic Flux Density 2 = Real 0.0

End

So on one end of the cylinder you are applying a field pointing in +ve z direction and another end you are applying a field pointing in the -ve direction and the field on the curved surface of the cylinder is zero for all x,y,z components. This is not a uniform field. This make sense in the result shown in your image file above. The "magnetic_flux_density_Z" is zero at the center of the cylinder and roughly symmetrical about it the center, this is an expected result.

However, I am trying to simulate a uniform field throughout , pointing along only +ve z direction (or only -ve z direction). The problem is the large deviation in the fields computed. If you look at my Case 1 and the result image file, you see that "magnetic_flux_density_Z" shows deviation of about 10% from 1 Tesla value. For a cylinder made of air, I would expect the "magnetic_flux_density_Z" should be 1 Tesla throughout with small deviation because numerical simulations are always approximates. Note that "magnetic_flux_density_x" and "magnetic_flux_density_Y" values are close to 0 Tesla--- this is again expected. The large deviation of "magnetic_flux_density_Z" from 1 Tesla and its asymmetry about the center of the cylinder is my concern. I did use rather fine mesh for the simulation.

My question is, shouldn't I be seeing better accuracy, field homogeneity, and symmetry?
Am I defining the boundary conditions accurately or adequately?

The cylinder is 5 units in radius (SI units of meters by default) and extends from -5 to +5 along the z-axis.

kevinarden
Posts: 618
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: MagnetoDynamics; WhitneyAVSolver; Uniform B-field

Post by kevinarden » 19 May 2020, 20:55

Looks like Z flux density is 1 except for the boundaries
z.png
(73.86 KiB) Not downloaded yet
case2.sif
(2.71 KiB) Downloaded 121 times

Post Reply