Problem understanding boundary conditions

Numerical methods and mathematical models of Elmer
Post Reply
hauptmech
Posts: 8
Joined: 06 Nov 2016, 07:17
Antispam: Yes

Problem understanding boundary conditions

Post by hauptmech » 08 Nov 2016, 01:54

Hi,

I'm trying to get valid results from the linear elastic solver and failing. I think I'm misunderstanding the boundary conditions.

I have a 100mmx20mmx4mm aluminum beam, fixed at one end, with a load of 100N applied in the thin direction on the other end. I expect about 4.3mm of deflection.

I'm using ElmerGUI (nglib) to load a STEP file and create a mesh. I have 2 different meshes (2 different projects), one of 264 elements and one of 933 elements.

I'm using Coordinate Scaling = Real 0.001 as a Simulation parameter to deal with the STEP import scaling problem.

After some initial confusion, I expect that the boundary condition "Force" should be Desired_Force/Boundary_Area or 100N/8e-5m^2, or "Displacement i Load" could be Desired_Force/Nodes_in_boundary, when using Boundaries as a Target.

However, the results are not near to what I expect. Moreover, the two different meshes produce different results. I'm getting 0.7mm deflection for the 264 element simulation and 0.99mm for the 933 element mesh.

Does anyone know what boundary conditions I should be using to get more accurate results? Why are the 2 simulations giving such different results?

An example 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.ep
Coordinate Scaling = Real 0.001
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 = Linear elasticity
  Procedure = "StressSolve" "StressSolver"
  Variable = -dofs 3 Displacement

  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 = Diagonal
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 1
  Linear System Precondition Recompute = 1
End

Solver 1
  Equation = Result Output
  Procedure = "ResultOutputSolve" "ResultOutputSolver"
  Output File Name = case
  Output Format = Vtu
  Exec Solver = After Timestep
End

Equation 1
  Name = "strain"
  Calculate Stresses = True
  Active Solvers(1) = 2
End

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

Material 1
  Name = "Aluminium (generic)"
  Heat expansion Coefficient = 23.1e-6
  Heat Conductivity = 237.0
  Sound speed = 5000.0
  Heat Capacity = 897.0
  Mesh Poisson ratio = 0.35
  Density = 2700.0
  Poisson ratio = 0.33
  Youngs modulus = 72.0e9
End

Boundary Condition 1
  Target Boundaries(1) = 1 
  Name = "Fixed"
  Displacement 3 = 0
  Displacement 2 = 0
  Displacement 1 = 0
End

Boundary Condition 2
 Target Boundaries(1) = 3 
  Name = "Force"
  Force 2 = $100/8e-5
 !Displacement 2 Load = Real $100/8

End
Any ideas or help is appreciated.

Thanks in advance

raback
Site Admin
Posts: 3568
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Problem understanding boundary conditions

Post by raback » 08 Nov 2016, 23:08

Hi

The 1st way of giving the force should be ok. I would not use the 2nd method as you don't know the representative area of each boundary element.

For 3D computation your meshes seem very small. The difference in results is probably due to not having proper mesh convergence.

You could try to use "Element = p:2" in Solver section for better accuracy.

-Peter

hauptmech
Posts: 8
Joined: 06 Nov 2016, 07:17
Antispam: Yes

Re: Problem understanding boundary conditions

Post by hauptmech » 09 Nov 2016, 04:30

Hi peter,

Thanks, the results are a lot closer in agreement and both close enough to the expected solution with that change. (For others hitting this issue, the docs are Appendix E and later of the ElmerSolverManual).

The note about small number of meshes is helpful too. I'm new to FEM and am still learning about element type and density.

I'm setting up user functions to speed up my workflow such that I can apply forces (simplified engineering style) and torques to boundaries. I figured out functions to calculate the number of nodes and the area of each boundary. For setting forces on a boundary I have:

Code: Select all

Boundary Condition 2
  Target Boundaries(1) = 3 
  Name = "Force"
  Total Load 2 = Real $100
  Force 2 =  Variable Time  
  Real Procedure "simpleforces" "BoundaryForce2"
End
In which the user function calculates the boundary area and divides the load by it.

For torque I'm not sure what is the approach I want. However your note on my second method suggests that the displacement load is not a per-node kind of thing? Is it a function of the area of the elements that share the node? I'm reading through the code in SolverUtils but I don't have a background in FEM math.

raback
Site Admin
Posts: 3568
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Problem understanding boundary conditions

Post by raback » 09 Nov 2016, 15:37

Hi

Good that you've found better agreement.

There is a poorly advertised feature that any Real keyword on boundary can be normalized by area with

Code: Select all

Somevalue = 1.23
Somevalue Normalize by Area = Logical True
This keywords are at the start of simulation then checked whether they need to be divided by the area where they live. If so then the expression is always multiplied with the scaling factor. No changes needed in the code.

-Peter

hauptmech
Posts: 8
Joined: 06 Nov 2016, 07:17
Antispam: Yes

Re: Problem understanding boundary conditions

Post by hauptmech » 10 Nov 2016, 01:43

Hi,

This sounds like just the thing I needed, however, it does not seem to work.

I have added the keyword to my local SOLVER.KEYWORDS as it fails without it.

Code: Select all

bc:logical:          'Force 2 Normalize by Area'
When I add the keyword to my sif file the result is different (more than an order of magnitude) than when I normalize manually, or when I normalize via my own area calculation. Any ideas why this might be?

Code: Select all

Boundary Condition 2
  Target Boundaries(1) = 3 
  Name = "Force"

  !!Manual normalization
  !Force 2  = $100/8e-5

  !!My normalization user function
  !Total Load 2 = Real $100
  !Force 2 =  Variable Time  
  !Real Procedure "torque" "BoundaryForce2"

  !!Built in normalization 
  Force 2 = Real $100
  Force 2 Normalize by Area = Logical True
End

raback
Site Admin
Posts: 3568
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Problem understanding boundary conditions

Post by raback » 10 Nov 2016, 20:24

Hi

You should be able to see from the ElmerSolver output the normalization factor used by "ListSetCoefficients". Is that not correct? Area you running in serial or parallel?

-Peter

hauptmech
Posts: 8
Joined: 06 Nov 2016, 07:17
Antispam: Yes

Re: Problem understanding boundary conditions

Post by hauptmech » 11 Nov 2016, 01:53

I am running serial. The relevant output is: (I just now added some dummy boundary conditions so that all boundary weights are shown.)

Code: Select all

 
BC weight:           1   2.0906241574201690E-003
 BC weight:           2   2.0861969993537025E-003
 BC weight:           3   1.4964948152328917E-002
 BC weight:           4   3.6268986341234169E-003
 BC weight:           5   3.2263823175843562E-003
 BC weight:           6   1.6604450039076200E-002
 BF weight:           1   0.0000000000000000
 Body weight:        1   8.0000000000000047E-006
 Mat weight:           1   8.0000000000000047E-006
ListSetCoefficients: Normalizing > force 2 < by  4.79341E+02
The expected area for boundaries 1,2 is 8e-5 (Which is what my user function, using ElementArea on the boundary elements, calculates).



The full sif, in case it's helpful, is:

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.ep
Coordinate Scaling = Real 0.001
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 = Linear elasticity
  Procedure = "StressSolve" "StressSolver"
  Variable = -dofs 3 Displacement
Calculate Strains = True
Calculate Principal = True
Calculate PAngle = True
  Calculate Stresses = True
  Element = p:2
  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 = Diagonal
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 1
  Linear System Precondition Recompute = 1
End

Solver 1
  Equation = Result Output
  Procedure = "ResultOutputSolve" "ResultOutputSolver"
  Output File Name = case
  Output Format = Vtu
  Exec Solver = After Timestep
End

Equation 1
  Name = "strain"
  Calculate Stresses = True
  Active Solvers(1) = 2
End

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

Material 1
  Name = "Aluminium (generic)"
  Heat expansion Coefficient = 23.1e-6
  Heat Conductivity = 237.0
  Sound speed = 5000.0
  Heat Capacity = 897.0
  Mesh Poisson ratio = 0.35
  Density = 2700.0
  Poisson ratio = 0.33
  Youngs modulus = 72.0e9
End

Boundary Condition 1
  Target Boundaries(1) = 1 
  Name = "Fixed"
  Displacement 3 = 0
  Displacement 2 = 0
  Displacement 1 = 0
End

Boundary Condition 2
  Target Boundaries(1) = 3 
  !Total Load 2 = Real $100
  !Force 2 =  Variable Time 
  !Real Procedure "torque" "BoundaryForce2"
  Name = "Force"
  Force 2 = Real $100
  Force 2 Normalize by Area = Logical True
End

Post Reply