Stress Relaxation of a Viscoelastic Tube

This example studies the temperature effects on the viscoelastic stress relaxation in a generalized Maxwell material with four branches.

Assume that the viscous part of the deformation is incompressible, so that the volume change is purely elastic. The relaxation shear modulus function is approximated by a Prony series as

(1)

where Gm represents the stiffness of the spring in the mth Maxwell branch, and τm is the relaxation time of the spring-dashpot pair in the same branch.

The instantaneous shear modulus is defined as

The alternative form of Equation 1 is

where the constants μm = Gm/G0 are such that

For many materials the viscoelastic properties have a strong dependence on the temperature. A common assumption is that the material is thermorheologically simple (TRS). In a material of this class, a change in temperature can be directly transformed into a change in the time scale. Thus, the relaxation time for a TRS material is modified to aT(T)τm, where aT(T) is a shift function. One of the most commonly used shift functions is defined by the Williams-Landel-Ferry (WLF) equation:

where a base-10 logarithm is assumed. Usually T0 is the glass transition temperature of the material. Note that aT(T0) = 1. If the temperature drops below T0 − C2, the WLF equation is no longer valid. The constants C1 and C2 are material dependent.

Model Definition

A long thick-walled cylinder has an inner radius of 5 mm and an outer radius of 10 mm. The inner surface is subjected to a prescribed radial displacement of 1 μm. A four-branch generalized Maxwell model represents the viscoelastic material.

In this example, you study the decay of the stresses during a period of two hours under the influence of the temperature field. You model a quarter of the cylinder and use the Solid Mechanics interface to compute the displacements in the cylinder cross section. The geometry is shown in Figure 1.

Figure 1: Geometry and mesh.

Study the following two cases:

•

Apply a stationary temperature field, causing the problem to lose its axisymmetry.

•

Make the temperature field transient.

Material Parameters

•

Elastic data: instantaneous shear modulus G0 = 27.46 GPa and long-term shear modulus of G = μ0G0, bulk modulus K = 39.88 GPa.

•

Viscoelastic data: four-branch Generalized Maxwell material with the shear modulus per branch defined from Gm = μmG0

μ1 = 0.04, τ1 = 30 s

μ2 = 0.08, τ2 = 300 s

μ3 = 0.09, τ3 = 3000 s

μ4 = 0.25, τ4 = 12000 s

•

Thermal properties: A WLF shift function is used with C1 = 17.44 and C2 = 51.6. These values are reasonable approximations for many polymers.

•

The reference temperature is 500 K.

•

Heat transfer properties: ρ = 7850 kg/m3, Cp = 2100 J/(kg·K), k = 6·10−2 W/(m·K).

Constraints

•

The circumferential displacements are constrained on the radial edges.

•

The inner edge of the hole is constrained to a radial displacement of 1 μm.

Heat Transfer Boundary Conditions

•

Stationary analysis: The inner and outer circular edges both have a temperature distribution varying linearly with the y-coordinate from 500 K at the y = 0 symmetry section to 506 K at the x = 0 symmetry section.

•

Transient analysis: The temperature distribution obtained from the stationary analysis is used as initial condition, while all edges are insulated.

Results and Discussion

In the first case, the temperature field is stationary and is shown in Figure 2.

Figure 2: Stationary temperature distribution.

The stress relaxation is faster where the temperature is higher. In Figure 3, the plots of radial stresses at radius 7.5 mm display this effect. The solid line is taken at y = 0 (cold) and the dashed line at x = 0 (warm).

Figure 3: Stress relaxation under a stationary temperature distribution.

In the second case, the temperature initially has the same distribution as shown in Figure 2, but it is allowed to settle in time to a final homogeneous value. Again, compare the radial stresses at radius 7.5 mm. The initial behavior is similar to the previous case, but as the difference in material properties decreases, the curves approach each other (Figure 4). The strain rate in the initially warm point decreases, while it increases in the initially cold point.

Figure 4: Stress relaxation under temperature settling conditions.

Notes About the COMSOL Implementation

When you solve a problem where a load is applied instantaneously at the beginning of a transient analysis, you can choose between two approaches to compute the initial stresses: either apply the load over a short period of time at the beginning of the time-dependent study. An alternative approach is to set the viscoelastic material to use an instantaneous static stiffness, and to add a separate stationary step before the time-dependent step. In the present example, you use the later method.

The results computed during the stationary step are stored and used as the initial value for the consequent transient analysis. COMSOL Multiphysics handles this automatically if you use a single study with one solver sequence. In this example, however, you study two different relaxation histories with different thermal boundary conditions. In the second case, the thermal boundary conditions differ between the computation of the initial state and the transient analysis. For this reason you use two separate studies for the stationary and the time-dependent studies.

Structural_Mechanics_Module/Material_Models/viscoelastic_tube

Modeling Instructions

From the menu, choose .

New

In the window, click

Model Wizard

In the window, click

In the tree, select .

Click .

In the tree, select .

In the window, under click .

In the window for , locate the section.

In the table, enter the following settings:


Name	Expression	Value	Description
G_inst	27.46[GPa]	2.746E10 Pa	Instantaneous shear modulus

Geometry 1

In the window, under click .

In the window for , locate the section.

From the list, choose .

Circle 1 (c1)

In the toolbar, click

In the window for , locate the section.

In the text field, type 10.

In the text field, type 90.

Click to expand the section. In the table, enter the following settings:


Layer name	Thickness (mm)
Layer 1	5

Delete Entities 1 (del1)

In the window, right-click and choose .

In the window for , locate the section.

From the list, choose .

On the object , select Domain 1 only.

Click

The equations only need to include first-order time derivatives since the inertia effects can be neglected for this class of problems.

Solid Mechanics (solid)

In the window, under click .

In the window for , locate the section.

From the list, choose .

Linear Elastic Material 1

In the window, under click .

In the window for , locate the section.

From the list, choose .

Viscoelasticity 1

In the toolbar, click

and choose .

In the window for , locate the section.

In the table, enter the following settings:


Branch	Shear modulus (Pa)	Relaxation time (s)
1	0.04*G_inst	30

Click

In the table, enter the following settings:


Branch	Shear modulus (Pa)	Relaxation time (s)
2	0.08*G_inst	300

Click

In the table, enter the following settings:


Branch	Shear modulus (Pa)	Relaxation time (s)
3	0.09*G_inst	3000

Click

In the table, enter the following settings:


Branch	Shear modulus (Pa)	Relaxation time (s)
4	0.25*G_inst	12000

Locate the section. From the list, choose .

In the Tref text field, type 500[K].

Locate the section. From the list, choose .

Symmetry 1

In the toolbar, click

and choose .

Select Boundaries 1 and 2 only.

Prescribed Displacement 1

In the toolbar, click

and choose .

Select Boundary 3 only.

In the window for , locate the section.

From the list, choose .

Locate the section. Select the check box.

In the u 0 n text field, type -1[um].

Heat Transfer in Solids (ht)

In the window, under click .

Temperature 1

In the toolbar, click

and choose .

Select Boundaries 3 and 4 only.

In the window for , locate the section.

In the T0 text field, type (500+6*y/sqrt(x^2+y^2))[K].

Materials

Material 1 (mat1)

In the window, under right-click and choose .

Select Domain 1 only.

In the window for , locate the section.

In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Bulk modulus	K	39.88[GPa]	N/m²	Bulk modulus and shear modulus
Shear modulus	G	0.54*G_inst	N/m²	Bulk modulus and shear modulus
Density	rho	7850[kg/m^3]	kg/m³	Basic
Thermal conductivity	k_iso ; kii = k_iso, kij = 0	0.06[N/(s*K)]	W/(m·K)	Basic
Heat capacity at constant pressure	Cp	2100[J/(kg*K)]	J/(kg·K)	Basic