Local Adaptive Meshing

Chapter 37: Local Adaptive Meshing

37

Local Adaptive Meshing



Summary



Introduction



Modeling Details



Mesh Refinement Process Definition



Material Modeling



Loading and Boundary Conditions



Solution Procedure



Results



Modeling Tips



Input File(s)



Video

554 555 556

559

564

563 564

557

558

559

558

554 MD Demonstration Problems CHAPTER 37

Summary Title

Chapter 37: Local Adaptive Meshing

Features

• 2-D structure mesh refinement • Region to be refined defined by property identifier • Mesh adaptivity criterion based on error indicator • Free-Free structure

Geometry H = 0.4 m d = 0.2 m s = 0.02 m F = 280 N

σmax F

d

F

H

σnominal

s

Material properties

E = 69GPa ,  = 0.33 ,  = 3200kg/m

Analysis characteristics

Linear static analysis using local adaptive meshing functionality

Boundary conditions

Automatic inertia relief (INREL = -2)

Applied loads

Tensile axial loading acting on the shortest edge of the plate (F = 280 N)

Element type

4-node, 3-node 2-D elements (QUAD4/TRIA3)

FE results

1. Stress at the most critical point for each refinement cycle

3

2. Stress concentration factor to compare with theoretical results Stress Concentration 2.5 Kt=

σmax σnominal

Theoretical (2.157) Numerical

2.0

Refinement Cycle

1.5

0

1

2

3

4

5

6

7

CHAPTER 37 555 Local Adaptive Meshing

Introduction This problem demonstrates the ability of the mesh refinement capability to converge on the correct solution in terms of stress distribution. A very simple structure has been considered to enable a comparison between theoretical and numerical results. F Theory, based on net section, states that if the nominal stress  nom = --------------------is defined as the stress acting on the net  H – d s

section (defined as the section that results from the difference between the width of the plate and the diameter of the   nom

max hole), then the stress concentration factor due to the presence of the hole is k t = ----------- , where  max is the actual stress

at the critical point. The stress concentration factor can be calculated from the empirical relationship shown in Figure 37-1.

Kt

3.0

2

d ⎞ d ⎞ d ⎞ ⎛ ⎛ ⎛ kt = 2 + 0.284 ⋅ ⎜1 − ⎟ − 0.600 ⋅ ⎜1 − ⎟ + 1.320 ⋅ ⎜1 − ⎟ H H H ⎝ ⎠ ⎝ ⎠ ⎝ ⎠

3

2.5

2.0

1.5 d H 1.0 0.00 Figure 37-1

0.25

0.50

0.75

1.00

Graphical Representation of Stress Concentration Factor versus d/H

The theoretical results and input data are shown in Table 37-1. Table 37-1

Input Data and Theoretical Results

Applied Stress Load Concentration (N) Factor Geometrical Data (m)

Nominal Stress (Pa)

Maximum Stress (Pa)

F

H

d

s

kt

 nom

 max

280.0

0.4

0.2

0.02

2.157

70000.0

150990.0

556 MD Demonstration Problems CHAPTER 37

Modeling Details A numerical solution has been obtained with MD Nastran’s solution sequence 101 for a 2-D representation of a freefree plate with a hole in its central region. The details of finite element model, contact simulation, material, load, boundary conditions and solution procedure are discussed below.

Figure 37-2

Initial Finite Element Model with Zoom in the Critical Region

The initial finite element model consists of: M O D E L NUMBER OF GRID NUMBER OF CQUAD4 NUMBER OF CTRIA3

S U M M A R Y POINTS = ELEMENTS = ELEMENTS =

964 872 12

The case control section of the input contains the typical entries for a linear static analysis. The only command that has been added to activate the mesh refinement is HADAPT. This, in turn, specifies the use of the bulk data entries, HADAPTL and HADACRI that control all the refinement process (see next section for details): ECHO = NONE HADAPT = 1 PARAM POST 0 PARAM,INREL,-2 SUBCASE 1 TITLE=First mesh LOAD = 2 DISPLACEMENT(PLOT,SORT1,REAL)=ALL STRESS(PLOT,SORT1,REAL,VONMISES,BILIN)=ALL Furthermore, the INREL parameter has been included with a value of ‘-2’ to activate the automatic inertia relief process. It is needed (automatic or manual) because the structure is in free-free conditions (unrestrained). The output request for displacement has been considered only to check the congruency of the deformation while the stress output is what we really need for comparison with the theoretical results. The Bulk Data Section contains the standard options for a linear static analysis plus the specific option for mesh refinement control.

CHAPTER 37 557 Local Adaptive Meshing

Mesh Refinement Process Definition The following options have been added to the standard linear static analysis Bulk Data section to define the mesh refinement process: $---------------------------------------------------------------------$ C A R D S F O R R E F I N E M E N T $---1---$---2---$---3---$---4---$---5---$---6---$---7---$---8---$---9--HADAPTL 1 7 1 PROP 2 1 + HADACRI 1 1 0.9 $----------------------------------------------------------------------

The HADAPTL option specifies the local adaptive mesh refinement control parameters. In particular, referring to the specific name associated to each field in the MD Nastran Quick Reference Guide, the process has been defined in this way: • REPEAT = 7 (5th field): maximum number of refinement cycles executed before the process is stopped • CRITID = 1 (6th field): associated HADACRI option identifier • WHEREMET = PROP (7th field): method used to specify the mesh refinement region subjected to the adaptivity criteria referenced in the associated HADACRI. PROP means that all the elements associated to a specific property option are considered by the refinement process • WHEREID = 2 (8th field): Identifier of the mesh refinement region subjected to the selected adaptivity criteria. Considering the WHEREMET value and the elements used in the finite element model, all the elements associated to PSHELL which identifier is 2 will be involved in the refinement process • SNAPMETH =1 (9th field): Method to project, snap, or relax new grid points on mid-edge or mid-face during the refinement process. The selected value allows the projection onto a smooth approximation of the analysis domain boundary interpolated from the mesh boundary. • MAXLEVEL = default (2nd field in the second physical option for HADAPTL): Maximum refinement level allowed for each individual element in the mesh. No elements will be refined to a level higher than the specified value. The default value is equal to that one defined in the REPEAT field. The HADACRI option specifies the mesh adaptivity criterion and the corresponding parameters. In this case, the method based on a scalar error indicator has been chosen (TYPE = 1 in the 3rd field). According to this criterion a scalar error indicator Ee is computed in the finite element mesh and an element ‘e’ will be refined if: 2

Ee --------------------------------1 N 1 2 F 1  ----  E f N f= 1

where N is the total number of elements in the element set to which it belongs and F 1 is the value specified in the 4th field of this option (in the specific case F 1 = 0.9 ). Note that the elemental error indicator is computed using the grid point stresses following the procedure utilized by the ELSDCON Case Control command.

558 MD Demonstration Problems CHAPTER 37

Material Modeling Isotropic elastic material properties of the deformable body are defined using the following MAT1 option as follows: MAT1

1

6.9+10

.33

3200.

The Young’s modulus is taken to be 6.9 GPa with a Poisson’s ratio of 0.33. Mass density ( = 3200 Kg/mm3) has also been specified to support the inertia relief process. Note that the results are not affected by the value included in this field.

Loading and Boundary Conditions No boundary conditions have been defined because the structure is not constrained. The loading involves application of concentrated forces (to simulate uniform load distribution) at the nodes located on the shortest edges of the plate: $ Loads for Load Case : Case_2 LOAD 2 1. 1. 1 1. 3 1. 5 $ Nodal Forces of Load Set : TENSILE_FORCE_Central_NEG FORCE 1 252 0 20. -1. 0. FORCE 1 253 0 20. -1. 0. FORCE 1 254 0 20. -1. 0. FORCE 1 255 0 20. -1. 0. FORCE 1 256 0 20. -1. 0. FORCE 1 257 0 20. -1. 0. FORCE 1 258 0 20. -1. 0. FORCE 1 735 0 20. -1. 0. FORCE 1 736 0 20. -1. 0. FORCE 1 737 0 20. -1. 0. FORCE 1 738 0 20. -1. 0. FORCE 1 739 0 20. -1. 0. FORCE 1 740 0 20. -1. 0. $ Nodal Forces of Load Set : TENSILE_FORCE_Central_POS FORCE 3 504 0 20. 1. 0. FORCE 3 505 0 20. 1. 0. FORCE 3 506 0 20. 1. 0. FORCE 3 507 0 20. 1. 0. FORCE 3 508 0 20. 1. 0. FORCE 3 509 0 20. 1. 0. FORCE 3 510 0 20. 1. 0. FORCE 3 959 0 20. 1. 0. FORCE 3 960 0 20. 1. 0. FORCE 3 961 0 20. 1. 0. FORCE 3 962 0 20. 1. 0. FORCE 3 963 0 20. 1. 0. FORCE 3 964 0 20. 1. 0. $ Nodal Forces of Load Set : TENSILE_FORCE_Ends_POS FORCE 4 503 0 10. 1. 0. FORCE 4 958 0 10. 1. 0. $ Nodal Forces of Load Set : TENSILE_FORCE_Ends_NEG FORCE 5 251 0 10. -1. 0. FORCE 5 734 0 10. -1. 0.

1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.

4

CHAPTER 37 559 Local Adaptive Meshing

Solution Procedure According to the HADAPTL and HADACRI control options, the refinement process starts with a preliminary calculation (CYCLE = 0 or ANALYSIS number 1) using the initial finite element model. Then, the refinement process starts and continues up to a number of cycles equal to REPEAT (3rd field in HADAPTL). During these cycles, each element involved will be refined up to MAXLEVEL value (2nd field in the second physical option in HADAPTL). As result of each refinement cycle the following files will be generated (xxxx.bdf is the generic name of the input file and): • xxxx.n.bdf – It contains the grid points, the elements and the MPC options related to the refined mesh created at the specific refinement cycle • xxxx.n.xdb – It contains the model and the results for the specific refinement cycle. where n is the number of the generic refinement cycle. Furthermore, the standard files xxxx.log, xxxx.f04, xxxx.f06 are generated. In the last one, it is possible to read some information about the refinement process show in the example below: ^^^-----------------------------------------------------^^^GLOBAL NUMBER OF ELEMENTS: 1096 ^^^AVERAGE ERROR INDICATOR: 1.766260E+03 ^^^CHANGE IN AVERAGE ERROR INDICATOR: 5.402161E-01 % ^^^-----------------------------------------------------^^^* * * E N D O F A N A L Y S I S #: 2 * * * ^^^-----------------------------------------------------by which it is possible to verify how it is proceeding and when the specific cycle is finished.

Results The first result to analyze is the way in which the finite element mesh is changed during the refinement cycles. In the figure below all the refined models are summarized. Note that the MPC relationships used to establish the congruency between regions with different meshes are not displayed to make the images clearer.

560 MD Demonstration Problems CHAPTER 37

Refinement Cycle 1

Refinement Cycle 2

Refinement Cycle 3

Refinement Cycle 4

Refinement Cycle 5

Refinement Cycle 6

Refinement Cycle 7

Figure 37-3

Refined Finite Element Models

CHAPTER 37 561 Local Adaptive Meshing

Looking at the refined meshes obtained in the subsequent cycles, it can be seen how important it is to activate a projection onto a smooth approximation of the analysis domain boundary from the mesh boundary (SNAPMTHD field in HADAPTL option). In fact, it avoids the creation of kinks that create two problems: • Driving the refinement process around the geometrical singularities • Generating stress concentration in the singular regions Displacement output has been required only to verify the correctness of the solution in terms of deformed structure. The use of PARAM.INREL,-2 enables a meaningful deformed structure in the case of free-free boundary conditions (Figure 37-4). The deformation seems to be congruent with the applied loads.

Figure 37-4

Deformed Structure

Also, the stress distribution is as we expect.

Figure 37-5

von Mises Stress Distribution Relative to the Last Refinement Cycle

The stress level in the critical point is compared with the theoretical one and the relative stress concentration factor is calculated. The resulting data are summarized in the following table together with other general information related to refinement effects on mesh size and error indicator. The error percentages are calculated according to the following relationship: Calculated value – Theoretical value Error% = ----------------------------------------------------------------------------------------------------------  100 Theoretical value

while, as already mentioned, the stress concentration factor is calculated as:  max k t = ---------- nom

562 MD Demonstration Problems CHAPTER 37

Referring to Table 37-1, the theoretical values for maximum stress and stress concentration factor are: Maximum Stress

= 150990 N/m2

Nominal Stress

= 70000 N/m2

Stress Concentration Factor = 2.157

Table 37-2

Results Comparison

Refinement Cycle

Global Number of elements

Average Error Indicator

Maximum von Mises Stress (N/mm2)

Stress Concentration Factor

Stress Concentration Factor Error (%)

0 1 2 3 4 5 6 7

880 1100 1568 2900 6476 16326 41060 98996

1756.77 1766.26 1582.74 1211.11 838.27 530.66 338.77

115461.46 130292.07 137223.50 141407.04 145797-84 148861.35 150535.57 150545.23

1.649 1.861 1.960 2.020 2.083 2.127 2.151 2.151

-23.53 -13.71 -9.12 -6.35 -3.44 -1.41 -0.30 -0.29

Two important considerations can be seen in Table 37-2: • The evaluated stress concentration factor is close to the theoretical one. The relative error is about 0.3%. • The differences between two subsequent maximum stresses decreases increasing the refinement level (Figure 37-6).

CHAPTER 37 563 Local Adaptive Meshing

Stress Concentration 2.5 Kt=

σmax σnominal

Theoretical (2.157) Numerical

2.0

Refinement Cycle

1.5

0

Figure 37-6

1

2

3

4

5

6

7

Theoretical/Numerical Stress Concentration Factor Comparison

Modeling Tips Some suggestions can be helpful to define the best refinement process: The refinement can be limited using the field MAXLEVEL in the HADAPTL option. None of the elements in the mesh will be refined to a level larger than MAXLEVEL. Limiting this process is necessary to avoid run-away refinement. In this example, the default value (MAXLEVEL = REPEAT) has been used not only to test the right convergence towards the theoretical stress but also the limited improvement introduced in the latest refinement cycles. Kinks (e.g., sharp internal corners that lack C 1 continuity) should be avoided in order to limit their influence: • on the refinement process (if they exist, the refinement is concentrated around the geometrical singularities) • on results (avoiding kinks prevents fictitious stress singularities) Kinks can be controlled defining SNAPMTHD = 1 in the HADAPTL option. In this example, the relaxation/projection method has been activated for the grid points created by the procedure; to verify its positive effect, change SNAPMTHD from 1 to 0 and see how the refinement process behaves. The refined meshes are concentrated along the geometrical singularities (sharp corners or kinks of a polygonal hole) and the results (the maximum value always increases) will continue to subdivide elements near the kinks.

564 MD Demonstration Problems CHAPTER 37

Setting SNAPMTHD = 1 ensures the geometry of the hole is correctly represented during the refinement process. By creating a cylindrical coordinate system at the center of the hole, all the grid points that have been generated on the boundary are all at R = 0.1 m, exactly the radius of the circle (the error is on the fifth decimal digit). It confirms the need to use the SNAPMTHD = 1 relaxation/projection procedure.

Input File(s) File

Description

nug_37.dat

MD Nastran input file for local adaptive meshing example

nug_37.bdf

MD Nastran input file for local adaptive meshing example used in video

nug37.SimXpert

MD Nastran input file for local adaptive meshing example

Video Click on the image or caption below to view a streaming video of this problem; it lasts approximately 30 minutes and explains how the steps are performed.

H = 0.4 m d = 0.2 m s = 0.02 m F = 280 N

σmax F

d

σnominal

s

Figure 37-7

Video of the Above Steps

H

F