2 Verification Problems
2.3 Stresses in long pressurized pipe with butt-welded supporting trunnion pipes |
2.1 Cantilever beam with end shear load
2.1.1 Objective
This verification problem demonstrates that higher-order elements can capture the exact solution of a simple cantilever beam problem with coarse meshes. In this problem we analyse the normal and shear stress profiles of a cantilever beam with a square cross-section that is subjected to a tip load using quadratic and cubic basis functions. The results show that both choices of basis functions exactly capture the axial stress and that cubic elements capture both the axial and shear stresses exactly on any mesh.
Although real world engineering problems will not in general have solutions that can be captured exactly with quadratic or cubic basis functions, these higher-order functions will typically be more accurate for a given number of degrees of freedom. This allows coarser meshes to be used in simulations, which can lead to reduced computational cost in many problems.
2.1.2 Geometry Description
2.1.3 External Interactions
A transverse shear traction is applied at the free end of the beam with an average shear stress of in the direction. As the cross-sectional area is , the total applied force is . The shear traction applied to the tip is distributed according to the parabolic solution from Euler–Bernoulli beam theory given in equation , to facilitate precise quantitative accuracy testing. The cantilever support at the fixed end is modeled by fixing displacement in both transverse directions ( and ), and applying a linearly-varying pressure distribution, coming from the Euler–Bernoulli axial stress solution given in equation . The remaining rigid modes are eliminated by constraining axial displacement at only the corners of the cantilever support cross-section. Although the Euler–Bernoulli solution is formally derived through an assumption of cross-sections remaining planar (alongside other approximations), that does not hold exactly for the full three-dimensional displacement solution, and the support cross-section needs to remain free to warp in the axial direction. The remaining rigid modes are eliminated by constraining axial displacement at only the corners of the cantilever support cross-section.
Due to St. Venant’s principle, qualitatively-similar results could be obtained using a fixed boundary condition for the cantilever support and a uniformly-distributed shear traction for the tip load. However, the cubic solution would no longer be exact and significant deviations from the beam-theory solution would occur near the ends of the beam.
2.1.4 Material Description
2.1.5 Method
2.1.6 Solution Characteristics
Since the cross-section of the beam is a square and is composed of a single homogenous material the neutral axis is located at the center of the beam. Thus the normal stress profile can be computed as where is the moment at the axial location being probed and is the area moment of inertia of the cross-section, which is given by for a rectangle. The shear stress profile can be calculated as In the present problem, the shear stress calculation can be simplified to which we note is a parabolic distribution with bounds
2.1.7 Results
2.2 Infinite plate with a hole
2.2.1 Objective
This verification problem shows that the Flex method is able to achieve optimal theoretical mesh convergence rates on the well-known infinite plate with a hole benchmark problem.
This problem is chosen as a verification problem not only because the exact solution is known, but in particular because the geometry includes a hole. The presence of a hole in the geometry presents two problem features that are critical for verification of flex-method stress analysis. First, the hole will produce a stress concentration. Second, the hole will require elements of the background spline to be trimmed in the area of the stress concentration. Accurate representation of stresses near features that produce stress concentrations is critical for stress analysis, and such features will almost always produce trimmed elements. The results presented for this problem show that we are able to achieve theoretical mesh convergence with trimmed elements. We also show that stress gradients can be captured to a reasonable level of accuracy on relatively coarse meshes. These results show that the flex method is a poweful tool for stress analyis of complex parts.
2.2.2 Geometry Description
The in-plane geometry for this problem is a quarter-symmetry representation of a square plate with a hole in the center, shown in figure 5. For this problem, and . The out-of-plane thickness is selected to be the mesh element size in each computation, so that there is always only one element in the thickness direction.
Figure 5: A quarter-symmetry representation of a infinite plate with a hole.
2.2.3 External Interactions
The model is constrained using symmetry boundary conditions on the left and bottom surfaces, and plane strain deformations are enforced by constraining out-of-plane displacements on both the front and back surfaces. The exact stress is given by in a polar coordinate system with the origin at the center of the hole, where is the traction applied in the direction at infinity. This stress is used to apply the exact tractions on the top and right surfaces, as shown in figure 5.
2.2.4 Material Description
2.2.5 Method
Figure 6: The background mesh for mesh #3. There are approximately eight elements along the long edges of the part.
2.2.6 Results
2.3 Stresses in long pressurized pipe with butt-welded supporting trunnion pipes
2.3.1 Objective
In this verification problem, the maximum von Mises stress is computed for a supported, pressurized pipe under a tension load. The pipe is supported by two trunnions that have been butt-welded to the pipe as shown in figure 10. The results are compared with the reference solution published in the NAFEMS report [4].
The Infinite plate with a hole verification problem showed that the flex method is able to accurately represent stresses in trimmed elements for a simple geometry where an exact stress field is known. In this problem, the geometry is more representative of a real-world part wherein an exact stress is not known because of the geometric complexity. In general, as the complexity of CAD parts increases, the complexity of the trimmed elements will also increase. However, it is important to keep in mind that the trimmed elements are only used to define the region in which evaluations of the background spline will be performed. The computations themselves will always be performed on the well-behaved functions of the background spline, which provides accurate results in trimmed elements. The presented results show that the flex method is able to accurately predict the maximum von Mises stress for this problem using a relatively coarse background spline.
A comparison between a conformal body-fitted approach and the fully-immersed flex approach is also shown in the results. The body-fitted approach required additional work to decompose the geometry into hex-meshable sections. The results show that both approaches produce excellent stress approximations for relatively coarse meshes. The flex approach, however, requires much less work because there is no need to decompose the geometry.
2.3.2 Geometry Description
2.3.3 External Interactions
2.3.4 Material Description
2.3.5 Method
Figure 11: The 5mm body-fitted mesh at top, followed by the untrimmed background spline and trimmed U-spline for the 5mm immersed case.
2.3.6 Solution Characteristics
Mesh Convergence Iteration
Mesh Size (mm)
Maximum Von Mises Stress (MPa)
Table 8: Reference solution convergence results. The mesh size column reports the mesh edge length in the area where the maximum stress occurs.
2.3.7 Results
Mesh scheme
Mesh size (mm)
Max. von Mises stress (MPa)
Percent difference
Body fit coarse
10
509
4.8
Body fit medium
5
516
3.4
Body fit fine
2.5
526
1.5
Immersed coarse
10
523
2.1
Immersed medium
5
533
0.1
Immersed fine
2.5
534
0.1
Table 9: Convergence results for the body-fitted and immersed simulations.
Figure 12: Comparison of immersed versus body-fitted von Mises stress and displacement for the fine mesh. The von Mises stress is shown in units of megapascals (MPa) and displacements are given in units of millimeters (mm). The minimum value for von Mises stress contours is set at 400 MPa to better show contours around the maximum value.
2.4 Composite cantilever beam with end shear load
2.4.1 Objective
This verification problem demonstrates the accuracy of tied constraints between adjacent parts by computing the stress profile through a composite cantiliver beam composed of a stiff and flexible section. The normal and shear stress are compared to an exact reference solution.
Tied constraints play a critical role when computing results for assemblies using the flex method. In the typical body-fitted approach for assemblies, adjacent parts are often tied together by creating a conformal mesh between the two parts and specifying that coincident nodes at the interface map to the same degrees of freedom. For the flex approach, adjacent parts will each have a different background spline, which in general will not align or share coincident points at the interface. Because of this, tied constraints are used to tie adjacent parts together. In this simple problem, two adjacent cantilever beams are tied together and then subjected to an end load. Results are computed for two cases. In the first case, the background mesh for each section has the same element size. For the second case, the element sizes are different. The results for both cases show that these interfaces can be accurately modeled with the flex method’s tied constraint approach.
2.4.2 Geometry Description
Figure 13: A sketch of the composite cantilever beam showing the stiff and flexible sections that comprise the beam.
Dimension
Value
2.4.3 External Interactions
Displacements in the transverse ( and ) directions are contrained on the fixed end of the cantilever beam. The axial displacement is only constrained at the corners of the cantilever support face of the beam’s lower section, to remove rigid modes, while the analytical distribution of axial stress from Euler–Bernoulli beam theory given in equation is applied on the fixed end instead of constraining axial displacement. This is because the analytical solution depends on allowing cross-sections to warp slightly, and we want to be able to perform precise testing of the method’s accuracy. A transverse shear traction is applied at the free end of the beam with a non-dimensional average shear stress of , as shown in figure 13. The beam has a cross-sectional area of , which results in a total applied force of . The traction is distributed according to the analtyical beam-theory solution given in equation , to facilitate accuracy testing. The interface between the two sections is modeled as a material interface, using a dimensionless penalty scale factor of . This is selected higher than the default value recommended for practical analyses, to enable precise comparisons with the analytical solution.
2.4.4 Material Description
Property
Flexible
Stiff
Young's Modulus
Poisson's Ratio
Table 11: The non-dimensional material properties for the stiff and flexible layers of the cantilever beam.
2.4.5 Method
2.4.6 Solution Characteristics
A reference solution for the composite beam can be calculated using the transformed section method, wherein the geometry of the beam is transformed into an imaginary section of uniform material with equivalent elastic properties. We assume here that the notional uniform beam material is equivalent to the material found in the composite beam’s flexible section, so we transform only the geometry of the stiffer material.
The normal stress profile can then be calculated as
where is the distance to the neutral axis of the section and is radius of curvature. Intermediate calculations can be computed via the following expressions:
The shear stress profile can then be calculated as
where intermediate calculations can be computed via the following expressions:
2.4.7 Results
The normal stress is plotted along a transverse line at a point midway along the length of the beam in figure 15. The reference solution given in equation is also plotted. All plots are coincident because both quadratic and cubic basis functions are able to represent the solution exactly.
Figure 15: Normal stress profiles. The computed normal stress is coincident with the reference solution in all cases.
Figure 16: Shear stress profiles. The shear stress computed with the cubic basis is coincident with the reference solution in all cases.
In all cases, the interface between the stiff and flexible sections of the beam was modeled as a tied constraint. These results demonstrate that we are able to accurately represent material interfaces with tied constraints.
2.5 Thin structures: Scordelis-Lo roof
2.5.1 Objective
The Scordelis-Lo roof is a benchmark problem from the well-known shell obstacle course [3]. It is often used to evaluate shell formulations because it contains both membrane and bending modes. Here we use this problem to evaluate the effectiveness of using the flex method to model thin structures. A commonly-encountered challenge when thin structures is having sufficient resolution in the background spline to capture bending and shear deformations through the thickness of the structure while remaining computationally tractable. This is particularly difficult for non-planar structures in cases where (i) the background spline does not align with the structure, and (ii) a uniform element size is used in all directions. The results presented for this problem show that higher-degree splines are able to accurately model thin structures with very coarse background splines. In particular, we show that when using a degree-four spline, we are able to get accurate results for both displacements and stress with an element size that is eight times larger than the thickness of the structure. These results are in contrast to the typical practice of choosing element sizes such that there are multiple elements through the thickness of a structure.
2.5.2 Geometry Description
2.5.3 External Interactions
The roof is supported at each end by a rigid diaphragm, and a non-dimensional body load of 360 per unit volume is applied to represent the weight of the roof. The body load for this problem is typically given as 90 per unit area when modeling with shells and is converted to a volumetric load here by dividing by the thickness of the shell. The rigid diaphragm support is modeled by constraining the and displacements at the end of the roof to be zero.
2.5.4 Material Description
2.5.5 Method
Taking advantage of symmetry, we model only one-fourth of the roof. The roof is immersed in an axis-aligned background spline with uniform size in each direction. The background spline, with an element size of , is shown in figure 18. We compute solutions on a series of refined background splines for degrees two, three, and four until the center displacement at the midpoint of the unsupported edge converges to a reference displacement of -0.30148. The reference displacement was computed using a refined quartic background spline.
2.5.6 Results
The quantities of interest for this problem, shown in figure 19, are (i) the displacement at the midpoint of the unsupported edge, and (ii) the shear and membrane stress along a curve on the midsurface at the midpoint of the roof. The displacement-versus-element count along the unsuported edge of the quarter-symmetry model is shown in figure 20. The coarsest background spline has an element size of 2.0, which is eight times larger than the thickness of the roof. For this mesh size, all displacements are within 1% of the reference solution of -0.30148, and the cubic and quartic results are within 0.1% of the reference solution.
Figure 19: A sketch of the Scordelis-Lo roof showing the probes where the displacement and stress are evaluated. The displacement is evaluated at the midpoint of the unsupported edge. The shear and membrane stresses are evaluated along a curve on the midsurface, at the midpoint of the roof.
Figure 20: The displacement at the midpoint of the free edge of the roof.
Shear stress results are shown in figure 21. Here we see that the quadratic results do not capture the shear stresses accurately and exhibit symptoms of locking for large element sizes. This behavior is reduced when using degree three splines, and almost completely eliminated when using degree four splines, where only the coarsest background spline shows noticeable error. Similar results for membrane stress can be seen in figure 22.
2.6 Frequency Extraction of Free Circular Plate
2.6.1 Objective
This verification problem shows that the Flex method is able to calculate several of the lowest fundamental frequencies of an unconstrained circular plate.
This problem is chosen as a verification problem because while exact solutions known they are quite challenging to derive and thus this problem represents, perhaps, the simplest structural modal problem for which known solutions exist yet it is still preferred, in practice, to produce a numerical approximation.
2.6.2 Geometry Description
The geometry for this problem is a circular plate shown in figure 23. The parameters used to describe this geometry are provided in table 13.
Parameter
Value
t
0.01
r
1
2.6.3 External Interactions
There are no boundary conditions applied to the circular plate (i.e., “free vibration”).
2.6.4 Material Description
Dimension
Value
Young's Modulus
Poisson's Ratio
0.3
Density
1
2.6.5 Method
Results for this problem were computed using two approaches, each using Coreform IGA. One approach consisted of a traditional bodyfit-meshing approach while the other utilizes an immersed meshing approach. In addition, three levels of refinement were performed using each approach. These refinement levels are listed in table 15. A quadratic basis with maximal continuity was used in each case.
Refinement
Quadrature points
DOF
Immersed
Bodyfit
Immersed
Bodyfit
Coarse
243
324
225
360
Medium
864
1296
540
864
Fine
2943
5616
1413
2700
Results from both bodyfit and immersed meshes were compared against the analytic solution to the first fundamental frequency of a free vibrating circular plate. For a circular plate with radius , thickness , density , Young’s modulus , and Poisson’s ratio , the analytic solution for the natural frequencies in a vacuum are where is the non-dimensional frequency parameter and is the flexural rigidity. For a circular plate, the flexural rigidity is The non-dimensional frequency parameter is tabulated in the literature for values of and , where is the number of nodal circles and is the number of nodal diameters (see figure 25). A shortened list of calculated non-dimensional frequency parameters is shown in table 16.
Figure 25: Examples of nodal circles () and nodal circles (). The bolded lines represent the nodal lines of the associated eigenmode.
0
1
2
3.00052
0
2
5
6.20025
0
3
8
9.36751
0
4
11
12.5227
0
5
14
15.6727
1
1
3
4.52488
1
2
6
7.7338
1
3
9
10.9068
1
4
12
14.0667
1
5
15
17.2203
2
0
1
2.31481
2
1
4
5.93802
2
2
7
9.18511
2
3
10
12.3817
2
4
13
15.5575
2
5
16
18.722
2.6.6 Results
Table 17 shows the calculated first frequency of the finest mesh for both approaches, as well as the relative error between each solution and the analytic solution.
Approach
Calculated Frequency ()
Relative Error (%)
Immersed
0.8200
0.47 %
Bodyfit
0.8266
1.29 %
Table 17: Relative error in fine mesh solutions for both approaches
Table 18 shows the solution characteristics of the finest mesh for both approaches.
Approach
Initialization Time (s)
Linear Equation Time (s)
Nonlinear Equation Time (s)
Solver Time (s)
immersed
1.66e-04
0.22
0.29
2.14
bodyfit
2.66e-04
0.19
0.66
132.75
Table 18: Solution characteristics for fine mesh solutions for both approaches
2.7 Neo-Hookean cylindrical pressure vessel
2.7.1 Objective
The main purpose of this test is to demonstrate the accuracy of the neo-Hookean hyperelastic material model in a finite-strain problem with a known non-trivial analytical solution. In particular, we solve the problem of a plane-strain section of a cylinder subjected to internal pressure. This problem was used to verify an academic implementation in Section 3.2 of [1]. A secondary purpose of this test is to demonstrate the effectiveness of higher-order immersed analysis in a situation where classical low-order finite element analysis would be affected by severe "locking" due to a large bulk modulus.
2.7.2 Geometry Description
2.7.3 External Interactions
To model a full cylinder subject to plane-strain kinematics, we apply the symmetry boundary conditions shown in figure 28, and additional symmetry boundary conditions on the faces. The inner surface of the cylinder is subject to a pressure follower load of magnitude . This is treated as given data in the problem setup, but chosen based on an analytical solution to result in a known radial displacement of the inner surface, via equation . The outer surface of the cylinder is traction-free. The symmetry boundary conditions use a dimensionless penalty scale factor of 1000 to perform a controlled study on the accuracy of the constitutive model and potential for locking in high-order splines.
2.7.4 Material Description
2.7.5 Method
2.7.6 Solution Characteristics
An analytical solution is available from the literature [1] for the fully-incompressible case, which we are approximating via . For the mesh sizes considered here, the error due to this approximate incompressibility is small relative to the approximation error inherent in using polynomial splines to approximate a transcendental exact solution. The exact solution is most conveniently parameterized in closed form by considering a target radial displacement of the inner surface, , and then calculating the corresponding pressure as where This solution is valid for arbitrarily-large deformations, and we select for calculating the results presented here, corresponding to an applied follower pressure of psi.
