© Q.-P. Zhou and H. Ding, published by EDP Sciences, 2018

1 Introduction

The solving process in finite element method (FEM) includes four main steps: first step, is to discrete the physical model. Second step, is to determine the governing equations and their boundary conditions. Third step, is to give the finite element (FE) equations of discrete elements. Finally, to solve the FE equations by a joint solution. Meshing is a critical step in FEM [1–3]. It affects directly the accuracy of analysis results. The choice of mesh size in FEM is the eternal question for numerical simulation [4,5]. The accuracy of the meshing size depends on the geometries of the domain, the materials properties, the element types and the loading. Theoretically, for FEM under certain conditions, the accuracy of the mesh is more precise when it's size is quite small. However, such as small size lead to high computation cost. Concretely, one should follow the algorithm below to reach the optimal precision ε:

• Start with initial mesh size equal to h_1^;

• form K_h1 (stiffness matrix) and F_h1 (discrete loads) we solve K_h1. U_h1 = F_h1;

• Choose a new mesh size h₂ that respect h₂ < h_1^;

• Then, form K_h2 and F_h2 we solve again K_h2. U_h2 = F_h2^;

• If $\frac{\left|U_{h_2}-U_{h_1}\right|}{\left|U_{h_2}\right|}\le ϵ$ then Stop, U_h2 is the final solution, Otherwise, h₁ = h_2^;

• Return to Step 1.

However, the above approach cost a lot in term of computation time. In practical, the mesh size is often chosen by experience [6–9]. On the other hand, the material is often homogenised when performing FE analysis at a macroscopic scale. The mesh size cannot be infinitely reduced. Indeed, we can mesh the local structure, but in this case, we face the time computing and also a cross-size simulation. When we select the mesh size for FE simulation, there must be some error in the FE solution and in the solution of the differential equation. In general, this error is monotonous to the mesh size. For the FEM using the displacement method, the stiffness of the finite element model increases with the increase of the mesh size under certain conditions [10–13]. In this article, we propose a new approach to identify the mesh dependency on material properties to solve this problem.

2 Relationship between FEM accuracy and mesh size

In this section, we underline above all through an basic example of cantilever beam the relation between the meshing size and the response of the structure.

2.1 Example: the cantilever beam subjected to uniform load

As discribed in Figure 1, a representative cantilever beam subjected to uniform load is discussed in this section (Fig. 2). The beam section is I50a, where y-axial moment of inertia is I_y = 1120 cm⁴. The length of the beam is L = 2000 mm. The material properties of the beam are: Young's modulus E = 210 Gpa and the Poisson's ratio μ = 0.3. The bean is under uniform loading with the following density per length unit q = 10 N/mm.

The theoretical equation for deflection curve is $w_B=-\frac{q{x}^2}{8EI}(x^2-4lx+6l^2)$. At the end of the free edge, it is equal to $w_B=-\frac{q{l}^4}{8EI}$. Once, we consider the numerical values, W_B = 8.5034 mm.

The length of the beam, the cross-sectional dimension, the load and the material properties remain unchanged, changing the mesh size h, we obtain in Table 1 the displacement of the cantilever beam W_B(h) at different mesh size.

As shown in the above table and figure, the displacement in function of the mesh size W_B(h) converges to the theoretical solution W_B when the mesh size is small enough (h converges to zero).

Fig. 1 Cantilever beam subjected to uniform load.

Fig. 2 Displacement of the cantilever beam at different mesh size.

2.2 Relation between finite element solution and mesh size

Using the minimum potential energy principle, we obtain the following static equation $$K_h{U}_h=F_h$$(1) where, K_h is the FE stiffness matrix that dependent on Young's modulus, the Poisson's ratio. U_h is the nodal displacement and F_h is the equivalent nodal load.

The elastic deformation energy of the approximate solution of the displacement obtained by using the minimum potential energy principle is the lower bound of the exact solution deformation energy (the approximate displacement field is smaller than the exact solution). As we know the finite element solution is: $$U_h={\displaystyle {\displaystyle \sum }_{i=1}^n}{N}_i{a}_i$$(2) where, N_i and a_i are respectively the shape interpolation function and the weight parameter for each i.

When n is big, there are more parameters to be determined, the accuracy of the finite element solution is higher. When n tends to infinity, the approximate solution approaches the exact solution. So the displacement U_h converges to the exact solution U₀ when the mesh size tends to zero (see Fig. 3).

Fig. 3 Displacement of the domain at different mesh size.

3 Mesh dependent material properties

Material properties are usually obtained by physical tests. For example, the Young's modular E can be obtained from the one dimension tensile test (Fig. 4). The values will depend only on the measurement of stress and strain.

However, the measurements of stress and strain are not the objective. They depend on which kind of measurements you used: sample type, measurements of the geometry, loads and/or evaluation methods of stress and strain. Using FEM as a rule to evaluate σ and ε or U etc, to obtain E and μ. Given element size h, element type, sample geometry etc, form K_h which is depend on E and μ, solve $$K_h\left(E_h,{\mu }_h\right){{\displaystyle \overline{U}}}_h={{\displaystyle \overline{F}}}_h\text{.}$$(3)

We obtain mesh dependent evaluations: $$E_h=E(h)\text{and}{\mu }_h=\mu (h)\text{,}$$(4) where, K_h(E_h, μ_h) is the stiffness matrix when the mesh size is h, E_h is the Young's modulus and μ_h is the Poisson's radio, ${{\displaystyle \overline{U}}}_h$ is the nodal displacement when the mesh size is h; ${{\displaystyle \overline{F}}}_h$ is the load when the mesh size is h. ${{\displaystyle \overline{U}}}_h$ and ${{\displaystyle \overline{F}}}_h$ can be obtained by simulation or experiment.

Fig. 4 Stress-strain curve.

4 Numerical examples

In this section, some numerical examples are calculated to illustrate the proposed method. First, the optimal solutions to meet a certain accuracy of a typical one dimension tensile model and a shear dominant model are obtained. Then, based on the solutions, we obtain the related mesh that depends on material parameters E_h and μ_h with a certain accuracy. We consider for this example a cubic structure structure subjected to uniform loads. The mesh scale k is evaluated as: $$k=\frac{h}{l_0}\text{,}$$(5) where h is the mesh size, l_o is the minimum characteristic dimension of concern.

4.1 Example of accuracy for a cuboid under different loading case

In the following, we will present 3 cases of loading for better understanding.

4.1.1 Case 1: Tension loading

As displayed in Figure 5, the simulation model is a cubic structure 5 cm × 5 cm × 10 cm subjected to one dimension tensile load 200 Mpa on a 5 cm × 5 cm section. The material properties of the cubic structure are Young's modulus E = 210 Gpa and the Poisson's ratio μ = 0.3.

We choose the accuracy of the solution is ε = 10⁻³. Than, we reduce gradually the mesh size from 50 mm to 3.1 mm in order converge to the optimal value ε that we consider after as a reference value. The error of the solution before and after encryption meet than the parameter ε.

In Table 2 and Figure 6, we obtain the following displacement U_t^* = 9.37779 × 10⁻² mm at the smallest value of ϵ that corresponds a mesh scale equal to 0.0625.

Fig. 5 One dimension tensile model.

Fig. 6 Elongation of the cubic structure at different mesh scale.

4.1.2 Case 2: Shearing loading

As displayed in Figure 7, the simulation model is a cubic structure 5 cm × 5 cm × 10 cm subjected to shear force 100 Mpa. The material properties remain the same.

We choose the accuracy of the solution ε equal to 3 × 10⁻³. As the previous case, we reduce the mesh size from 50 mm to 1.5 mm to reach the optimal ε (that we consider as a real or optimal solution).

Table 3 and Figure 8 show that the solution related to optimal ε gives U_s^* = 8.91528 × 10⁻¹ mm with a corresponding mesh scale equal to 0.03125.

Fig. 7 Shear loading of the cubic structure.

Fig. 8 the displacement of the cubic structure at different mesh scale.

4.1.3 Case 3: Mesh dependent material properties

Based on the previous target solutions with a certain accuracy, we obtain the mesh dependent material parameters E_h and μ_h. First we identify the E_h then we identify μ_h.

4.1.3.1 Identification of E_h

As the domain, the load and the Poisson's ratio remain unchanged (see Fig. 5), we change the mesh size value h to solve the following equality: $$U_{E_h}={U_t}^{*}\text{.}$$

The material parameter E_h is identified afterward as described in Table 4 and in Figure 9.

Fig. 9 Identified Young's modulus at different mesh scale.

4.1.3.2 Identification of μ_h

As the domain, the load remain unchanged (see Fig. 7), and using the identified Young's modulus E_h we obtained above we change the mesh size h to solve the following equation to get parameter μ_h: $$U_{{\mu }_h}={U_s}^{*}\text{.}$$

4.1.3.3 Application of the mesh dependency material properties

When we use the given mesh size to solve the real application (Fig. 10), we should introduce mesh dependent material properties E_h and μ_h (Tab. 5) to form the stiffness matrix K_h(E_h, μ_h). As displayed in Figure 11, the cubic structure 5 cm × 5 cm × 10 cm is subjected to uniform loading with the following density equal to 50 Mpa.

As the domain and the load remain unchanged during EF analysis, the comparison of the cases results 1–3 with different mesh size h are shown below : Case 1 in Table 6, case 2 in Table 7 and case 3 (using the identified parameters Young's modulus E_h and the identified Poisson's ratio μ_h) in Table 8.

As we can see in Figure 12 after Young's modulus and Poisson's ratio correction, the mesh size can increase almost two levels to reach the same accuracy.

Fig. 10 Identified Poisson's ratio at different mesh scale.

Fig. 11 Cubic structure under uniform loading.

Fig. 12 Displacement of the cubic structure at different mesh scale.

5 Conclusion

In the finite element analysis, the structural response is related to the mesh size under certain conditions, the solution of FE analysis is increase uniformly to converge to the exact solution. However, with the discontinuity of the materials and non uniform shape such as composite and pore materials, the mesh size cannot be too small, otherwise, it will cause imprecise. The introduction of mesh dependent material properties correction results in more accurate results. Moreover, the Young's modulus and the Poisson's ratio correction together gave better accuracy than with only modulus correction. Finally, the correction is more efficient for coarse mesh than for fine mesh. This interval is also what we need to improve most.

References

All Tables

All Figures

	Fig. 1 Cantilever beam subjected to uniform load.
In the text

	Fig. 2 Displacement of the cantilever beam at different mesh size.
In the text

	Fig. 3 Displacement of the domain at different mesh size.
In the text

	Fig. 4 Stress-strain curve.
In the text

	Fig. 5 One dimension tensile model.
In the text

	Fig. 6 Elongation of the cubic structure at different mesh scale.
In the text

	Fig. 7 Shear loading of the cubic structure.
In the text

	Fig. 8 the displacement of the cubic structure at different mesh scale.
In the text

	Fig. 9 Identified Young's modulus at different mesh scale.
In the text

	Fig. 10 Identified Poisson's ratio at different mesh scale.
In the text

	Fig. 11 Cubic structure under uniform loading.
In the text

	Fig. 12 Displacement of the cubic structure at different mesh scale.
In the text