Anisotropic non linear BH curve

[mika]

When curl-conforming basis functions are applied, there is no need to call ElmerGrid with -increase flag. Instead try with the lowest-order nodal mesh as a background mesh and add Quadratic Approximation = True to the solver section.

-- Mika

[mika]

As additional info, a second-order nodal mesh (10-node tetrahedra, 13-node pyramids, 15-node prisms, or 20/27-node bricks) should also work as a background mesh, together with "Quadratic Approximation = True" in the solver section. The second-order curl-conforming basis in principle needs also "Use Piola Transform = True" in the solver section, but this should follow automatically when "Quadratic Approximation = True".

[gforti]

i've been having convergence issues using quadratic elements. Is there a guideline?
Tested some examples from Github and convergence issues are similar.
For instance, using only quadratic approximation = true diverges. It needs a higher order element (with ElmerGrid -increase)
Also, is there difference between using "Quadratic Approximation = True" or "Element = p:2"? I though the results would be more smooth for the magnetic flux density, but still B is constant over an element.

EDIT:
Just finished a little evaluation. I got the simulation files from https://github.com/ElmerCSC/elmer-elmag ... OnIronCore and did some simplifications. I tested with and without higher order elements. The results from the solver SaveLine showed the expected improvements (for Magnetic Flux Density and Magnetic Flux Density e).
My question was with the paraview results. Using "Discontinuous Bodies" in ResultOutputSolver gives smooth B distribution for Magnetic Flux Density e. Without this parameter, the results of B seems to be constante over each element. And this is the result for Magnetic Flux Density e from the vtu file from "Post File = case_qd.vtu" in the simulation section.

Code: Select all

! post: http://www.elmerfem.org/forum/viewtopic.php?f=3&t=7010&p=22330#p22330
!
! This file has minor modification for style and economics. 

Header
  CHECK KEYWORDS Warn
  Mesh DB "." "MESH_inc"
  Results Directory "RESU"
End

Simulation
  Coordinate System = "Cartesian 3D"
  Coordinate Mapping(3) = 1 2 3

  Simulation Type = Steady
  Steady State Max Iterations = 1

  Output Intervals = 0
  Max Output Level = 20

  Post File = case_qd.vtu
End

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

! -- BODIES --
Body 1
  Name = "Inductor"
  Target Bodies(1) = 1
  Equation = 2
  Material = 1
  Body Force = 1 
End

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

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

! -- MATERIAL PROPERTIES --
Material 1
  Name = "Inductor"
  Relative Permeability = Real 1.0
  Relative Permittivity = Real 1.0
  Electric Conductivity = Real 1e+14
End

Material 2
  Name = "Air"
  Relative Permeability = Real 1.0
  Relative Permittivity = Real 1.0
End

Material 3
  Name = "Core"
  Relative Permeability = Real 10000.0
  Relative Permittivity = Real 1.0
End

! -- COMPONENTS --
Component 1
  Name = String "Inductor 1"
  Coil Type = String "test"
  Master Bodies(1) = Integer 1
  Desired Current Density = Real -1e+08
End

! -- BODY FORCES --
Body Force 1
  Name = "Source"
  Current Density 1 = Equals "CoilCurrent e 1"
  Current Density 2 = Equals "CoilCurrent e 2"
  Current Density 3 = Equals "CoilCurrent e 3"
End

!-- BOUNDARY CONDITIONS --
Boundary Condition 1
  Name = "FAC003"
  Target Boundaries(1) = 6
  AV {e} = real 0
  AV = real 0
  Jfix = real 0.0
End

!-- EQUATIONS --
Equation 1
  Name = "MGDyn for air/iron"
  Active Solvers(2) = 2 3
End

Equation 2
  Name = "MGDyn for coils"
  Active Solvers(3) = 1 2 3
End

!-- SOLVERS --
Solver 1
  Equation = "CoilSolver"
  Procedure = "CoilSolver" "CoilSolver"

  Linear System Solver = "Iterative"
  Linear System Preconditioning = ILU1
  Linear System Max Iterations = 1000
  Linear System Convergence Tolerance = 1e-10
  Linear System Iterative Method = BiCGStab
  Linear System Residual Output = 10
  Steady State Convergence Tolerance = 1e-06

  Normalize Coil Current = Logical True
  Nonlinear System Consistent Norm = Logical True
  Coil Closed = Logical True
  Narrow Interface = Logical True
  Save Coil Set = Logical True
  Save Coil Index = Logical True
  Calculate Elemental Fields = Logical True
  
  !Coil Center(3) = Real ...
  Coil Normal(3) = Real 0.0 0.0 1.0

  ! Ensure divergence-freeness of source
  ! Must be False to equal results from case.sif with Jfix
  Fix Input Current Density = Logical False
End

Solver 2
  Equation = MGDynamics
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  
  Steady State Convergence Tolerance = 1.0e-5
  !Nonlinear System Max Iterations = 200
  !Nonlinear System Convergence Tolerance = 1.0e-5
  !Nonlinear System Newton After Iterations = 3
  !Nonlinear System Newton After Tolerance = 1.0e-3
  !Nonlinear System Relaxation Factor = 1
  !Nonlinear System Consistent Norm = Logical True
  
  Linear System Solver = iterative
  Linear System Iterative Method = BiCGStabl
  Linear System GCR Restart = 100
  BiCGStabl Polynomial Degree = 6
  Idrs Parameter = Integer 6
  Linear System Robust = Logical True
  Linear System Max Iterations = 1000
  Linear System Convergence Tolerance = 1.0e-8
  Linear System Preconditioning = ILUT
  Linear System ILUT Tolerance = 1e-3
  Linear System Residual Output = 10
  Linear System Abort Not Converged = False
!  Linear System Row Equilibration = Logical True

  Fix Input Current Density = True

!   Element = p:2 ! ERROR:: MarkOuterNodes: The faces are not created!

  Quadratic Approximation = Logical True
  Use Piola Transform = Logical True
End

Solver 3
  Equation = MGDynamicsCalc
  Procedure = "MagnetoDynamics" "MagnetoDynamicsCalcFields"
  
  Nonlinear System Relaxation Factor = 1

!   Quadratic Approximation = Logical True
!   Use Piola Transform = Logical True
  
  Linear System Solver = Iterative
  Linear System Symmetric = True
  Linear System Iterative Method = CG
  Linear System Max Iterations = 5000
  Linear System Convergence Tolerance = 1.0e-8
  Linear System Preconditioning = ILU0
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 1
  Linear System Precondition Recompute = 1

  Discontinuous Bodies = True
  Calculate Nodal Fields = True 
  Calculate Elemental Fields = True

  Calculate Current Density = True
  Calculate Magnetic Field Strength = True
End

Solver 4
   Exec Solver = After Timestep
   Procedure = "SaveData" "SaveScalars"
   Filename = res_case.dat
End

Solver 5
  Exec Solver = After Timestep
  Equation = SaveLine
  Procedure = "SaveData" "SaveLine"
  Filename = BH_line_qd.dat
  Variable 1 = Magnetic Flux Density
  Polyline Coordinates(2,3) = 0 -0.03 0 0 0.03 0
  Polyline Divisions(1) = 500
End

Solver 6
  Exec Solver = "after saving"
  Equation = "ResultOutput"
  Procedure = "ResultOutputSolve" "ResultOutputSolver"
  Output File Name = case_res_qd
  Vtu format = Logical True
  Binary Output = Logical True
  Discontinuous Bodies = Logical True

  Save Elemental Fields = Logical True
  Save Geometry Ids = True
  Save Bulk Only = True
End

[gforti]

Hello to all,

I've finished another set of trials. Please refer to the git folder:
https://github.com/gforti84/ElmerFiles/ ... in/debug_4

Here I've tested several conditions, including linear and non-linear for isotropic and anisotropic silicon steel. There are a few things that I wish to point out.
With linear and anisotropic case, the simulation time is really sensible with the z-axis permeability.
Ii could not use vanka preconditioner with mpi.
Also, using MUMPS with non-linear material did not converge, even isotropic.

I have added different line plots for comparison against the commercial software. These lines are all defined in the core, two in the joints of yoke and leg and one in the middle o the core leg.
I each folder (EF_Iso_Lin, EF_Iso_Nlin, ...) there is a jupyter notebook for generating these curves from the result files. The output from the commercial software is saved in the res_AM folder.

[gforti]

Hello all,

Just for the record, the result of magnetic flux along a line, defined as illustrated in the figure below, shows good agreement with commercial software result.
Line:

Joint_line.png (90.38 KiB) Viewed 479 times

Result:

EF_Aniso_Nlin_Joint.png

(69.53 KiB) Not downloaded yet

[gforti]

Just added another case folder:
https://github.com/gforti84/ElmerFiles/ ... in/debug_5

Simulation time is similar now, and results are close. Now it is time to try this solution on a real case.
As soon as I have results for a more realistic case I can share it in elmer-elmag.

Thanks for all the help so far.