How do user-supplied b-h curves work

[Nick_99]

Hi everyone,

I'm having some conceptual issues with Elmer and user-supplied b-h curves. I don't quite understand how Elmer treats the b-h curve in different situations and I can't seem to find any information in the documents. I found the magdyn_bh example on GitHub but it has left me with more questions than answers

Let's start by using this case.sif for the discussion.

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
  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

! A cube located at the centre of a cylinder (body 2)
Body 1
  Target Bodies(1) = 1
  Name = "Body 1"
  Equation = 1
  Material = 2 ! Iron with a user-supplied b-h curve
End

! A cylinder
Body 2
  Target Bodies(1) = 2
  Name = "Body 2"
  Equation = 1
  Material = 1 ! Air
End

Solver 1
  Equation = MgDyn
  Procedure = "MagnetoDynamics" "WhitneyAVSolver"
  Exec Solver = Always
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-8
  Nonlinear System Convergence Tolerance = 1.0e-10
  Nonlinear System Max Iterations = 500
  Nonlinear System Newton After Iterations = 5
  Nonlinear System Newton After Tolerance = 1.0e-4
  Nonlinear System Relaxation Factor = 1
  Linear System Solver = Iterative
! Bicgstabl is usually performing well for AV solver
  Linear System Iterative Method = BiCGStabl
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  BiCGstabl polynomial degree = 4
! We cannot use ILU since vector potential is not uniquely defined
  Linear System Preconditioning = none
  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
  Calculate Magnetic Field Strength = True
  Procedure = "MagnetoDynamics" "MagnetoDynamicsCalcFields"
  Exec Solver = Always
  Optimize Bandwidth = True
  Steady State Convergence Tolerance = 1.0e-5
  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStab
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  Linear System Preconditioning = ILU0
  Linear System Abort Not Converged = False
  Linear System Residual Output = 10
  Linear System Precondition Recompute = 1

  Calculate Nodal Fields = False
  Calculate Elemental Fields = True

! Perform the averaging within this solver
  Average Within Materials = Logical True
End

Solver 3
  Exec Solver = after timestep
  Equation = "ResultOutput"
  Procedure = "ResultOutputSolve" "ResultOutputSolver"
  Output File Name = f
  Vtu format = True
  Discontinuous Bodies = True
  ! bloated alternative for the above maintaining all discontinuities
  ! Discontinuous Galerkin = True 
  Save Geometry Ids = True
  ! use this only in conjunction with quadratic mesh   
  Save Linear Elements = True  
  Save Bulk Only = True 
  ! Save Boundaries Only = True 

  ! Enable saving of different parts to different files
  Vtu Part Collection = False
End

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

Material 1
  Name = "Air (room temperature)"
  Relative Permeability = 1.00000037
  Heat Conductivity = 0.0257
  Heat Capacity = 1005.0
  Viscosity = 1.983e-5
  Relative Permittivity = 1.00059
  Sound speed = 343.0
  Density = 1.205
  Heat expansion Coefficient = 3.43e-3
End

Material 2
  Name = "Iron"
  Relative Permittivity = 1.0
  
  h-b curve = Variable "dummy"
    Real
      include HB
    End 
End

! BOUNDARY CONDITIONS TO BE DISCUSSED!!

This case is applied to a cylinder (which acts as the region) with a cube in the centre of the cylinder. Body 1 is the cube and body 2 is the cylinder. The sides of the cylinder are boundaries 1 and 2 and the rectangular part of the cylinder is boundary 9. The b-h curve is applied to the cube which is iron. The volume inside the cylinder is air.

A image of the setup is shown below.

setup-min.png

(26.97 KiB) Not downloaded yet

General questions
For the first part of the discussion, I'd like to ask some general questions.

1. How must I format the b-h curve file and how do I tell Elmer what format I have used? I see in the GitHub example that they specify the B-field first and then the H-field. I would've thought you'd specify H first seeing as B will depend on the applied field (H)? Please see the sample below (which is an extract from the HB file on GitHub).

Code: Select all

0 0
0.7 355
.............
.............
3.22 800000

2. What are the units for the columns in the file? Are they Tesla for the B-field and A/m for the H-field?

3. What does the variable "dummy" do (see code snippet below)? I get that its called dummy since it isn't used (?) but why are we defining a variable in the first place? Is that just the syntax for supplying an external file to the solver? The code below has been taken from lines 121 through 133 in the case.sif on GitHub.

4. What happens when I change the 'h-b' in the line h-b curve = Variable "dummy" to 'b-h'? Does that swap the way the columns are read from the datafile?

Code: Select all

Material 2
  relative Permittivity = 1.0
  Name = "Iron"
  Electric Conductivity = 0.0

!  H-B Curve(31,2) = Real
!   INCLUDE HB

  h-b curve = Variable "dummy"
    Real
      include HB
    End 
End

No externally applied field
I know this is a bit daft but let's not have any externally applied field. In other words, the boundary condition for the case.sif will be

Code: Select all

Boundary Condition 1
  ! Sides of the cylinder are 1 and 2 and the rectangular part is 9.
  Target Boundaries(3) = 1 2 9
  Name = "Farfield"
  AV {e} = Real 0
End

5. What happens when I supply a magnetisation for a material that also has a b-h curve associated with it (see example below)? How is the B-field resulting from the magnetisation computed by the solver seeing as a permeability isn't supplied? Does the solver just use the b-h curve as a lookup table?

Code: Select all

Material 2
  Relative Permittivity = 1.0
  Name = "Iron"
  Electric Conductivity = 0.0
  Magnetization 3 = 21100

  h-b curve = Variable "dummy"
    Real
      include HB
    End 
End

Externally applied field
I would like to create a uniform magnetic field that flows from the one side of the cylinder to the other side. I do this by setting the one side (i.e. one of the circles) of the cylinder to 2.2T and the opposite side of the cylinder to -2.2T. I leave the rectangular part of the cylinder with no field. The boundary conditions added to the case.sif will be:

Code: Select all

Boundary Condition 1
  Target Boundaries(1) = 1 
  Name = "Circle1"
  Magnetic Flux Density {n} = Real -2.2
End

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

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

6. Is this the correct way to supply boundary conditions for an external field or is there a better way to do it?

7. Assuming the cube has no magnetisation. With 6. in mind, how does the solver work out the net magnetic field emanating from the cube seeing as the external field is given in Tesla and not A/m. Surely this means it can't use the b-h curve to determine the permeability?

8. Slightly unrelated but are the units for the magnetic potential Dirichlet boundary condition in Weber per meter?

Apologies for the length of the post.

Thanks in advance for any help.

Kind regards,
Nick

[hielau]

I have been using BH curves like this

Code: Select all

Material 1
  Electric Conductivity = Real 0.0
  Relative Permittivity = Real 1.0
  H-B Curve = Variable coupled iter
    Real Monotone Cubic 
		0.0	0.0
		0.1	123.4
			..........
		4.0	12345678.9
    End

  ! Ch hysteris component
  Harmonic Loss Coefficient 1 = Real 123.4

  ! Ce eddy current component
  Harmonic Loss Coefficient 2 = Real 0.123
End

B (T) values in the first column and H (A/m) values in the second column.

Elmer Models manual say HB curve must be given as a cubic spline

BR.
Lauri

[Nick_99]

Hi Lauri,

Thanks very much for your input and answers to questions one and two in my post

Cheers,
Nick.

[raback]

Hi Nick,

There is some history behind the b-h curve. Initially it was given as a (x,y) table with given size and internally changed to spline. Then later we adopted more general approach where it is a function H=H(B) is given in the same style as any other real valued keyword in Elmer. The catch just is that the argument is set internally (|B| is computed on-the-fly for each integration point) and the actual argument is not used. Hence it is called "dummy".

The Real valued keyword has many different ways for evaluation and the non-constant ones are activated with the "variable" string.

The primary variable is always the vector potential and B=curl(A). As H=H(B) is usually nonlinear either fixed point or Newton linearization is used. A challenge here may be a crappy curve. Hence a monotone spline is used since it has a continuous positive derivative.

So I would recommend using the newer version. The nice thing is that the b-h curve can also be given by a user function or MATC expression etc. Probably best to go with the table still usually. However, one could imagine more complex functional depencies, for example, on temperature etc.

-Peter

[Nick_99]

Hi Peter,

Thank you very much for your reply!! I will use the newer version

Kind regards,
Nick