COMSOL 6.4 - Kirsch Infinite Plate Problem

In this example, you perform a static stress analysis to obtain the stress distribution in the vicinity of a small hole in an infinite plate. Two approximations of the infinite plate are evaluated. The first one uses a plate that is large compared to the hole while the second one employs an infinite element domain.

The problem is a classic benchmark, and the theoretical solution was derived by G. Kirsch in 1898. This implementation is based on the Kirsch plate model described on page 184 in Mechanics of Materials, D. Roylance (Ref. 1). The stress level is compared with the theoretical values.

Model Definition

Model the infinite plate in a 2D plane stress approximation as a 2 m-by-2 m plate with a hole with a radius of 0.1 m in the middle. Due to symmetry in load and geometry you need to analyze only a quarter of the plate, see Figure 1. Choose the size of the plate sufficiently large so that the stress concentration close to the hole is not affected.

Figure 1: Geometry model of the Kirsch plate with rollers defining the symmetry plane.

When modeling a plate using the infinite element domain you need to create an additional layers around the plate. Those layers simulate the part that stretches to infinity and can have an arbitrarily length along the direction that stretches to infinity, for example 0.1 m.

In the model the infinite element domain is created along the y direction only since the numerical results along x = 0 symmetry plane are compared to an analytical reference and infinite element domain in x direction only have a minor influence.

Material

Isotropic material with, E = 2.1·1011 Pa, ν = 0.3.

Load

A distributed stress of 103 Pa on the right edge pointing in the x direction.

Constraints

Symmetry planes, x = 0, y = 0.

Results and Discussion

The distribution of the normal stress in the x direction, σx, is shown in Figure 2and Figure 3. The stress contours of the finite model and the infinite model are very similar.

Figure 2: Distribution of the normal stress in the x direction for the finite model.

Figure 3: Distribution of the normal stress in the x direction for the infinite model.

According to Ref. 1 the stress σx along the vertical symmetry line can be calculated as

(1)

Figure 4 shows the stress σx obtained from the solved models, and plotted as a function of the true y-coordinate along the left symmetry edge, which are in close agreement with the theoretical value according to Equation 1.

Figure 4: Normal stress, simulated results (solid line) versus the theoretical values (dashed line).

Away from the hole, stresses from the finite model starts drifting from the theoretical values, while stresses from the infinite model matches closely with the theoretical value.

The stress error is summarized in the following table:

Table 1: Stress error relative to analytical solution.
	Finite plate	Infinite plate
Near hole	1.1%	0.2%
Away from hole	-4.0%	-0.1%

Notes About the COMSOL Implementation

The default scaling function in is rational. This type of function is well adapted to cases where the degrees of freedom vanish to zero at infinity. The present model is submitted to infinite loads at infinity, that means that constant strain and linear displacement are expected. For this type of infinite solution, polynomial functions are preferred. The relation between the stretched and geometric coordinates is

where the function f is defined by an analytic function. Here we want f as a second-order polynomial: f(ξ) = aξ2 + bξ + c. The continuity conditions at X0, f(0) = 0, and at the end of the domain, f(1) = pw, imply that the polynomial is

The infinite element domain gives the best results when meshed with rectangular elements; see Figure 5.

Figure 5: Infinite element domain modeled with rectangular elements.

Reference

1. D. Roylance, Mechanics of Materials, John Wiley & Sons, 1996.

Structural_Mechanics_Module/Verification_Examples/kirsch_plate

Modeling Instructions

From the menu, choose .

New

In the window, click

Model Wizard

In the window, click

In the tree, select > .

Click .

Click

In the tree, select > .

Click

Global Definitions

Parameters 1

In the window, under click .

In the window for , locate the section.

In the table, enter the following settings:


Name	Expression	Value	Description
pw	10[m]	10 m	Physical width of infinite element domain
deltaY	0.1[m]	0.1 m	Geometric thickness of infinite element layer

Draw a rectangle with a top layer that represents the infinite element domain.

Geometry 1

Rectangle 1 (r1)

In the toolbar, click

In the window for , locate the section.

In the text field, type 1+deltaY.

Click to expand the section. Clear the checkbox.

Select the checkbox.

In the table, enter the following settings:


Layer name	Thickness (m)
Layer 1	deltaY

Circle 1 (c1)

In the toolbar, click

In the window for , locate the section.

In the text field, type 0.1.

Click

Difference 1 (dif1)

In the toolbar, click

and choose .

Add the rectangle and remove the circle in the section.

In the window for , click

Global Definitions

First add an analytical function for stress, based on Kirsch’s theoretical solution of an infinite plate.

Analytic Stress

In the toolbar, click

and choose > .

In the window for , type Analytic Stress in the text field.

In the text field, type AnaStress.

Locate the section. In the text field, type 1000/2*(2+(0.1/y)^2+3*(0.1/y)^4).

In the text field, type y.

Locate the section. In the table, enter the following settings:


Argument	Unit
y	m

In the text field, type N/m^2.

Create an analytic polynomial function to define the scaling in the infinite element domain.

Definitions

Analytic 2 (an2)

In the toolbar, click

In the window for , locate the section.

In the text field, type (pw-10*deltaY)*x^2+10*deltaY*x.

Locate the section. In the table, enter the following settings:


Argument	Unit
x	m

In the text field, type m.

Infinite Element Domain 1 (ie1)

In the toolbar, click

Select Domain 2 only.

In the window for , locate the section.

From the list, choose .

Add a variable representing the physical y-coordinate to be used in postprocessing.

Variables 1

In the window, under right-click and choose .

In the window for , locate the section.

In the table, enter the following settings:


Name	Expression	Unit	Description
ym	if(dom==2,ie1.Ym,y)	m	Physical y-coordinate

Solid Mechanics (solid)

First set up a model without the .

In the window, under click .

In the window for , locate the section.

In the list box, select .

Click

Select Domain 1 only.

Locate the section. From the list, choose .

Locate the section. In the d text field, type 0.1.

Symmetry 1

In the toolbar, click

and choose .

Select Boundaries 1, 2, and 5 only.

Boundary Load 1

In the toolbar, click

and choose .

Select Boundaries 6 and 7 only.

In the window for , locate the section.

Specify the fA vector as


1e3	x
0	y

Now set up a model with the infinite element domain.

Add Physics

In the toolbar, click

to open the window.

Go to the window.

In the tree, select > .

Click the button in the window toolbar.

In the toolbar, click

to close the window.

Solid Mechanics 2 (solid2)

In the window for , locate the section.

From the list, choose .

Locate the section. In the d text field, type 0.1.

Symmetry 1

In the toolbar, click

and choose .

Select Boundaries 1, 2, and 5 only.

Boundary Load 1

In the toolbar, click

and choose .

Select Boundaries 6 and 7 only.

In the window for , locate the section.

Specify the fA vector as


1e3	x
0	y

Materials

Material 1 (mat1)

In the window, under right-click and choose .

In the window for , locate the section.

In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Young’s modulus	E	2.1e11	Pa	Young’s modulus and Poisson’s ratio
Poisson’s ratio	nu	0.3	1	Young’s modulus and Poisson’s ratio
Density	rho	7800	kg/m³	Basic

Mesh 1

For the finite plate selection, a customized free triangular mesh must be used for getting a better solution in the stress concentration region. A smaller element size is set at the expected location of stress concentration.

Free Triangular 1

In the toolbar, click

In the window for , locate the section.

From the list, choose .

Select Domain 1 only.

Size 1

Right-click and choose .

In the window for , locate the section.

From the list, choose .

Select Point 1 only.

Locate the section. Click the button.

Locate the section.

Select the checkbox. In the associated text field, type 0.02.

Infinite elements give better results when meshed with rectangular elements.

Mapped 1

In the toolbar, click

Distribution 1

Right-click and choose .

Select Boundary 2 only.

In the window for , locate the section.

In the text field, type 4.

Click

Study 1

In the toolbar, click

Results

To check the error in the computed results, make a point evaluation of stresses near the hole (y = 0.1) and away from the hole (y = 1) for the solution computed with and without the infinite element domain. The error can be determined by finding the difference between computed stresses and analytical stresses.

Error Evaluation

In the toolbar, click