Use of computers in component design.

Rubber components have been used widely for many years and the exact geometry of the component has been found to have a dramatic effect on its overall performance and life expectancy. A great deal of attention is given by the rubber industry to the composition and chemistry of the rubber material used in a component, but nearly as significant is the effect that design has on overall component performance. Adjustments to chemistry, polymer composition and additives can have a large effect on the durability of a component and its life expectancy. But the design of geometry has the primary impact on the mechanical performance of the component. Advances in the use of computer design technology have enhanced the ability of scientists to predict the mechanical characteristics of any given design geometry using basic elastomer properties. The intent of this article is to give a general overview of the process of using computer technology to develop and evaluate rubber component design.

The authors have focused their design efforts to a newer group of elastomers in the class of thermoplastic elastomers (TPEs) called elastomeric alloys (EAs) (refs. 1-3). Component design for these elastomeric alloys should be of general interest to anyone developing lower cost components with improved part consistency and uniformity. The elastomeric alloys are described as TPEs generated from a chemical combination of two or more polymers to give an alloy having better elastomeric properties than those of the corresponding blend. They are available commercially in a wide variety of durometers and polymer chemistries which are suitable for general purpose up to high fluid, high temperature environments. The most widely used EAs are the EPDM/polypropylene alloys products such as Santoprene thermoplastic rubber, which the authors have used primarily during the design work.

Design of rubber components using elastomeric alloys is very much like that for any other rubber. The process, properties and evaluations of elastomeric alloy component design is the same for virtually any rubber chosen. So the reader may generally follow the procedures for designing and analyzing the performance of any general elastomer part. In the case of elatomeric alloys designs are prepared that also account for the thermoplastic process by which the part will be manufactured. These processes allow for wider latitude in design similar to wide design flexibility seen for thermoplastic parts. Designs also benefit from the tighter tolerance and greater consistency achievable with an elastomeric alloy. Many thermoset rubber parts are manufactured by compression molding or transfer molding and have wider part variation including flash that must be removed. These issues lead to tolerance variations for which the designs must account.

A simple computer based technique for designing rubber parts is computer aided drafting (CAD). A very wide variety of computers and software companies have very useful products in this area. But this article will focus on the uses of the more advanced part design capabilities possible with the finite element analysis (FEA) technique for elastomeric alloys. Less computer software is available for FEA because of its greater complexity and less common use. However, its use for elastomer part design is growing as computer power increases and pricing decreses. These changes have been a driving force for the development of better software that continues to fuel further interest yet in its application to elastomer part design. The rest of the article will contain a look at what is required in FEA for elastomer parts. This includes a brief recounting of the equations to describe the behavior of elastomers, a description of how to determine the properties for FEA, an overview of the finite element method and a few examples of elastomeric alloy part design using computer based finite element analysis.

Theory

One of the key elements to using a computer to model the behavior of an elastomer in a component is the use of a mathematical equation that properly describes the actual behavior. It is beyond the scope of this article to fully address these theories. But a few brief comments are needed since they have a large effect on FEA. There are several models that can be used to varying degrees of accuracy. Unfortunately the more sophisticated models that have improved accuracy also use more coefficient terms in the equations. This increased number of terms requires increased computer power to use when analyzing actual elastomer part geometries. And it also requires more testing to develop the data necessary about the elastomer behavior and to allow a proper fit of coefficients over a practically useful range of stress and strain values.

For elastomer part analysis three equations of state are commonly used: Linear elastic, Mooney-Rivlin and Ogden equations. The equation for linear elastic theory is Hooke's Law (ref. 4) used for a linear spring where the strees, [stega], is defined as a function of the spring constant (or three dimensional stiffness matrix, K, and the local strain in the part, [sigma]. For a unidirectional deformation in the 1 direction the stretch ratio, [[lambda].sub.1] = [sigma] + 1:

[stega] = K [not equeal to] [sigma] (1)

So for unidirectional deformation one would have, [stega] = K ([[lambda].sub.1 ] - 1), in terms of a stretch ratio, [[lambda].sub.1].

This linear model is not very valid for elastomers because they exhibit highly non-linear stress vs. strain behavior, especially at high strain values. But if elastomer components have a small amount of deflection during use, then they may be modeled with only a small amount of error with this simple material law.

A widely accepted law that describes the behavior of elastomers for higher strains is the Mooney-Rivlin equation (refs. 5 and 6). Abaqus software uses the strain energy (U) function form of the equation which can be represented in terms of the stretch ratio, [lambda]i, in the three orthogonal directions as:

[MATHEMATICAL EXPRESSION OMITTED]

The equation models the behavior of the elastomer well through much of its stress vs. strain range. But at high strain values an elastomer usually begins to exhibit strain stiffening behavior that is not modeled by this simple form of the Mooney-Rivlin equation. For the softer elastomeric alloys the behavior is modeled well through 25 to 50% strain with this simple form of the equation.

An improved fit is obtained by using a second order version of the Mooney-Rivlin equation (ref. 6), which is given by:

[MATHEMATICAL EXPRESSION OMITTED]

This higher order equation provides a mathematical model that exhibits an inflection in the stress vs. strain behavior. Unfortunately the fit is rather poor in many instances at the more meaningful low strain values when the full range of elastomer properties is used to fit the equation. If some of the data at high extensions or compressions is not used and the coefficients fit to the lower strains, then the equation predicts wildly varying values for the high strains. So a better fit at high strains can be obtained with this higher order equation, but at the penalty of a relatively poor fit at low and moderate strain values.

The better equation to use for the full range of strain behavior is the Ogden equation (refs. 7 and 8) which also is higher order and has an inflection that matches an elastomer across the whole extension to compression strain range. The Ogden equation for the strain energy function, U, in terms of the stretch ratio, [lambda], the coefficient, [micro], and the exponent, [alpha], used for polynomial expansions for N = 3 to 6 are (ref. 6):

[MATHEMATICAL EXPRESSION OMITTED]

Across a wide range of strain the Ogden equation is a better description for the behavior of soft elastomers and elastomeric alloys. The equation uses several coefficients and the higher order fit costs more in terms of computer time and power required for analysis. But no trade off is made on low strain behavior to be able to model the high strain accurately.

A comparison of the relative fit of each of these models is shown in figure 1 for a simple cube shape in tension along one axis (uniaxial tension). The plot shows the correlation of nominal stress versus nominal strain, that is these values are compared to the original sample surface area, not the true area as the sample is deformed. Here it can be seen that the linear equation follows the elastomeric alloy behavior well for low strains. Even though the true [alpha] varies linearly versus the true [sigma] as described by equation 1, the values measured in tensile test are the nominal stress which is based on the original cross-sectional area, and not true values which account for the decreasing area and the sample draws down during extension. The Mooney Rivlin equation is a good fit up to moderately high strain values. But the best match with the true behavior over the full extension range is seen with the Ogden equation.

[CHART OMITTED]

Biaxial extension of an elastomer can be modeled by compression in one uniaxial direction. Compression for incompressible materials, like elastomers, causes an equal proportionate extension in both off-axis directions. Compression of a specimen is compared for hese models in figure 2. The linear equation shows an acceptable fit only for very low strain values. Mooney-Rivlin equations yield a very good fit to compressions of over 50%. But the Ogden equation has the best fit to compression behavior matching a 55 Shore A elastomeric alloy well up to 70% compression.

[CHART OMITTED]

The shear behavior of a planar specimen in extension is compared in figure 3 for a 55 Shore A elastomeric alloy. The very low strain shear is actually best modeled by the linear equation. But for the overall behavior to higher shear strain is better modeled by the Mooney-Rivlin and Ogden equations.

[CHART OMITTED]

Experimental

Property measurements

The test procedures to determine the proper coefficient values for the mathematical mdoel equation are determined from physical property tests and are documented by Hibbitt, Karlsson and Sorensen (ref. 6). Also a method is discussed by Finney and Kumar (ref. 9) for determining the coefficients from test data. The cofficients are determined by measuring three basic modes of deformation for an elastomer:

* Uniaxial tension where the material is extended in just one direction;

* Biaxial tension where there is a strain in two directions simultaneously; and

* Shear where there is a non-uniform strain across the material, such as when the two opposite sides move in opposite directions.

The uniaxial data are collected using a dumbbell shaped specimen as in ASTM D412 die C pulled in a tensometer at 50 cm/min The data collected must be reduced to the proper format for the curve fit program being used. Some programs will use the measured nominal stress and nominal strain values, while others wil require the input of true stress and strain determined from the data. These uniaxial data are all that is required to use the linear elastic model.

The biaxial data are most simply collected for incompressible elastomers using a uniaxial compression of a pellet shaped specimen. Because an elastomer is incompressible a compression in the one direction causes an equal, biaxial extension in the two directions perpendicular to the compressio direction. This simple technique replaces the more complex apparatus required to pull a sheet shaped specimen in two directions simultaneously. Although this biaxial tension test would provide more accurate measurements at very high strain values. But for most practical situations the simpler pellet compression test is adequate through 10 to 25% biaxial extension (25 to 60% uniaxial compression).

Shear strain behavior can be evaluated by extension of a planar sheet specimen of wide aspect ratio. The wide sample is physically constrained to not pull in at the edges by the clamps, so the sample will thin out but not draw in during the extension test. This mode places the specimen in a pure shear strain mode, at least until the extension becomes large and the sheet sample begins to draw in.

The computer software, Abaqus by Hibbitt, Karlsson and Sorensen, will automatically calculate the proper coefficients for the material model desired using the property measurements. The software will even provide feedback on the range of strain where the coefficients have a valid fit. Coefficients determined for a 55 Shore A elastomeric alloy are shown in table 1 as calculated from fit data shown in the previous figures for the various state equations discussed above.

Table 1 - coefficients for elastomer state equations (55 shore A

elastomeric alloy)

State equation             Coefficient             Value
Linear                         K                    5.66  MPa
                               v                   0.495
Mooney Rivlin, N=1            [C.sub.10]         0.5779   MPa
                              [C.sub.01]        -0.0861   MPa
Mooney Rivlin, N=2            [C.sub.10]         0.9660   MPa
                              [C.sub.01]        -0.3354   MPa
                              [C.sub.20]        -0.1091   MPa
                              [C.sub.11]         0.1149   MPa
                              [C.sub.02]        -0.0199   MPa
Ogden, N=3 coefficients       [[micro].sub.1]    3.0125   MPa
                              [[micro].sub.2]   -0.1282   MPa
                              [[micro].sub.3]   -1.8279   MPa
Exponents                     [[alpha].sub.1]     1.232
                              [[alpha].sub.2]     1.799
                              [[alpha].sub.3]   -0.0855
Ogden, N=5 coefficients       [[micro].sub.1]    1.6845   MPa
                              [[micro].sub.2]   -0.1432   MPa
                              [[micro].sub.3]    0.0183   MPa
                              [[micro].sub.4]   -0.5806   MPa
                              [[micro].sub.5]    0.0707   MPa
Exponents                     [[alpha].sub.1]     2.191
                              [[alpha].sub.2]     4.662
                              [[alpha].sub.3]     6.005
                              [[alpha].sub.4]    -1.187
                              [[alpha].sub.5]     1.872

Analysis

Finite element analysis method

The analysis of a component by computer hinges upon the ability to describe the part mathematically and detailed discussions of this are available (refs. 10 and 11). There are relatively few geometries where there is an explicit analytical expression to use for predictng the component behavior. Therefore an approximate method must be used to describe the behavior. If a component geometry is broken down into small discrete pieces of a geometry that is possible to be described explicitly with an analytical expression. Then if an analysis is performed in each individual piece and the effect each has on those pieces nearest to it one would have an adequate description of the overall components behavior mathematically. This is precisely what the finite element method does. The components must be broken into a fairly large number of pieces (elements) and all the appropriate equations solved simultaneously. The calculations required are quite immense and can only be practically done using a computer.

The part geometry is meshed into elements of type selected to be appropriate for the type of analysis desired. The representation of part can be broken down to three levels of increasing sophistication, and complexity. The first is to represent the part by a one dimensional element with an assumed cross-section profile, especially long, uniform crosssection parts. One dimensional elements have been developed that have mathematical behavior suitable to describe a wide variety of cross-section profiles. These elements all can be broken down into two general classes:

* First order elements with just two nodes on each end (so stresses and strains can only be described to vary linearly between them); and

* Second order elements that have three nodes with one at each end and another in the middle (so stresses and strain can be described by a second order parabolic equation).

Two dimensional elements can be used to describe another class of parts that are very thin and can be simulated as a shell. Again the elements will come in first and second order types with nodes at just the corner interesections or also including a mid-side node, respectively. In parts symmetrical about an axis, a two dimensional element with axial symmetry is used for the model. These two dimensional elements can have three or four sides as triangles or quadrilateral shapes. Or the elements can have an additional mid-side node, which can allow for curvature along a parabolic shape to be modeled. The higher order elements requre more computation to solve but give higher accuracy and they avoid the quirks of elements locking together and giving artificially high stiffness seen in simple linear elements, especially in a bending mode. Similarly parts with axial symmetry are well modeled by elements on a single plane through the model with two dimensional elements.

The most extensive type elements are the three dimensional elements that are used to represent full non-symmetrical three dimensional parts. The two general shapes used for these are quadrilateral shapes and tetraheadral. They are also used in first and second order types. These elements use the most nodes to describe the part since they have no advantage of reducing the part description due to symmetry and also due to the larger number of nodes required per element.

Some key element types for use in rubber part analysis include a linear type, a hyperelastic and reduced integration point versions of these. The mathematical formulation of these elements is optimized for certain types of problems. The linear type is formulated for use with just linear material models, the hyperelastic is able to use the Mooney-Rivlin and Ogden equations for rubber problems. The reduced integration elements make use of fewer internal points instead of all the nodes of the element, so the solution requires fewer calculations. Element mathematical formulations also include a number of other features that space does not allow for discussion. But they are suited for certain geometrical symmetries, and a variety of modeling situations. Proper element selection is critical to the accurate representation of an elastomer part and reducing the computation hardware size and time required.

Computer requirements

The calculation of finite element models requires fairly high powered computers. Simpler linear material models with two dimensional elements can be calculated in a few minutes or hours on a small workstation consisting of an RISC central processor unit chip (cpu) rated at 25 to 75 mips, 16 MB of main memory and 500 MB to 1 GB of disk space. Larger problems that require the use of the nonlinear material models for rubber on smaller two dimensional element models can still be run on workstation class computers although the calculations generally are stretched out to a few hours. These problems benefit greatly from increased main memory and disk space as well as faster or multiple cpu's. The more ultimate problems using three dimensional elements and the nonlinear rubber equations consume even more computer power. A three dimensional problem benefits greatly from >100 mips cpu's with disk space of over 5 GB being necessary for modes sized problms and main memory of over 100 MB being most useful. Still calculation times will stretch into multiple hours and possibly days on workstation class computers. A super computer/super workstation class computer is needed to reduce these calculation times to shorter times.

Representation of the results of the calculations requires a computer or workstation with high speed graphics capabilities. Results are displayed with 1280 x 1024 bits in full color on graphics workstations. Colored contours of stress, strain, strain energy, etc. can be displayed through the cross-section or on the surfaces of even three dimensional models. These high speed graphics capabilities also allow the representation of the part as three dimensional, shaded images. Calculations of the mass, center of gravity and simple geometrical factors can be displayed and documented visually. The resulting display of deformed geometries can locate high stress and strain areas. The effects of various load conditions and use parameters can be displayed in color shaded contours to allow the prediction and identification of problem areas. Altternate designs can be compared for predicted suitability match with acceptable design and material parameters.

Fast color graphics displays allow the animation of the loading and deformation of parts. The animations can be watched on the workstation screen and studied for acceptability per design criteria. Computer workstation visual display can also be captured in playback files on the computer or on video tape for use in presentations to document the results of the analysis. Also color printers available do an excellent job of capturing the colors and high resolution available on the computer workstations. The color prints provide a very good tool to document the predicted performance and comparative benefits of designs or elastomeric alloy selections.

Design examples

Perhaps the clearest way to show how the part design process works is to discuss an example of a typical analysis for an elastomeric alloy. Examples can be found for a variety of rubber applications (refs. 12-15). A first example here will deal with the approach typically taken when trying to design a 55 Shore A elastomeric alloy weather seal to meet a specific force deflection curve. Weather seals are most often specified by a geometric size or geometric window they must operate within, and by a specific force deflection curve. This curve defines the amount of force needed to create a life long seal and/or the maximum force needed/allowed to close the door.

In this case, a door will be closed against a lip seal. The door is represented by the vertical line in the contour plots shown in figure 4.

[CHART OMITTED]

The overall movement of the lip is a maximum of 4.5 mm. Because this movement is very small in relationship to the arc movement of the door, a simple linear horizontal movement of the line can be used. As this line is moved to the left, it deflects the seal. This deflection can be defined in a force/deflection curve based on the movement of the line and the reaction force (seal force) it sees from the lip seal.

The seal is defined by a geometry defined by eight noded quadrilateral elements. These second order elements have nodes at their corners and the midpoints of their sides. Because the seal has a consistent geometry down its length, only a cross section of it must be modeled and meshed with plain strain elements. These elements, by definition, take into account the geometric consistency of an extruded profile without compromise in the mathematical answers.

If a working seal exists in a similar geometry and function, it is wise to model it to check the validity of the problem set-up. This includes elements chosen, material model to be used, physical constraints, and mesh densities. In this example, an existing 55 Shore A elastomeric lip seal exists and its force deflection curve was measured. This seal and its cure will be defined as the actual profile or lab.

This actual profile's geometry was built and meshed as seen in the contour plots. This model contains 197 elements and 674 nodes. Several types of elements and material models were then used to try and match the lab force deflection curve. All the computer run force deflectin curves ignore friction where it was seen in the lab sample. The base element used was a CPE8 (eight noded second order plane strain quadrilateral). Variations of this element used were CPE8H (taking into account the hyperelasticity of elastomeric alloys) and CPE8RH (reduced integration version of CPE8H). The reduced integration used only four internal (to the elements) calculation points as opposed to the eight defining nodes. This, in theory, leads to faster run times and less stiff elements in large deflections.

In this experiment, we experience an approximate 40% savings in computer CPU time using the reduced integration elements over their full integration equivalents. As can be seen in table 2 the stresses and strains solved were only slightly lower. A graphic represntation of the total reaction force vs. contact surface displacement is shown in figure 5 using the hyperelastic elements CPE8H. The force deflection curves of the reduced integration elements, CPE8RH, in figure 6 are essentially identical to their full integration versions in figure 5 using the same element mesh. This comparison shows that for this type of deflection and model size the reduced elements are as good of a choice as full integration equivalents.

[CHART OMITTED]

Table 2 - weather seal profile calculated results 55 Shore A elastomeric

alloy

Elements    Material   Mises strain    X strain      Y strain
             model        MPa         max.   min.   max.    min.
CPE8H      MR N = 1      .8069        +8.4   -15.6  +15.6   -8.4
CPE8RH     MR N = 1      .7172        +8.3   -16.7  +16.7   -8.3
CPE8H      MR N = 2      .9172        +8.4   -15.6  +15.6   -8.4
CPE8RH     MR N = 2      .9310        +8.3   -16.7  +16.7   -8.3
CPE8H      OG N = 5      .7655        +8.4   -15.6  +15.6   -8.4
CPE8RH     OG N = 5      .7724        +8.3   -16.7  +16.7   -8.3
CPE8H      Linear       1.3586        +8.3   -17.2  +16.0   -8.5
CPE8       Linear       1.2621        +8.2   -15.0  +15.0   -8.5

A check was also made comparing material model choices. This included a linear, N = 1, N = 2, N = 3, and N = 5 (see figure 5). Because the strains are low, "8-16%, a linear material model appears to give the best fit in this case. This is not typically the case. This is probably true in this case because the majority of the part is in less than 3% strain or essentially neutral. Most materials are quite linear in this strain range. There is probably some inaccuracy at this very low strain with the higher order curve fit for material models. In this case, there was no computer run time savings for the linear approach compared to the other material models. As can be see in figure 7, there is no real difference in the force deflection curves of various elements with linear data. The use of CPE8 elements did take slightly longer to solve, but it is thought this difference may be atypical.

[CHART OMITTED]

When designing an elastomeric part, geometric changes should be made to minimize stress concentrations. In figure 4a the maximum Von Mises Stress is 1.36 MPa. This maximum stress can be due to mesh density problems or to part design. In this case, it is due to a coarse mesh in the curved region of the seal. The stress around this peak is much lower and probably more realistic at 1.0 MPa. If this were not the case, the part should be redesigned with a larger radius curve, different angled lip, or a change in deflection distance as allowed by the definition of the seal's function.

With elastomeric alloys, the strain levels are more important than the absolute stress levels. For a very small increase in strain, a large increase in stress can be seen. This makes it difficult to base part function on stress levels (although high stress concentrations should be minimized). If strain levels are analyzed in this lip seal, it is seen that the peaks range from 8-16%. This is not a problem for a 55 Shore A elastomeric alloy.

Other examples

FEA analysis of elastomeric alloys is not restrained to long symmetrical parts. Parts of all shapes and sizes can be analyzed. The more complex a part is, the more computer intensive the problem becomes. The more the problem can be simplified through symmetry, as the weather seal was, the faster an analysis can be completed.

Another typical type of problem would be a round bumper, mount or spacer. This type of part can be modeled as a section through the part, such as the bumper variations shown in figure 8. The section is then meshed as an axisymmetric part. These axisymmetric elements are formulated to account for the symmetry of the round part. In this case, a force deflection curve may need to be matched, the deformed geometry analyzed, a stress and strain minimized through design changes as was discussed in the first example. The match of part design, grade and hardness of an elastomeric alloy with the performance criteria can be evaluated without making prototypes. So fewer prototypes will be needed to develop a part design and select the proper material.

[CHART OMITTED]

A third typical type of problem deals with an elastomeric convoluted bellows. This type of part an be produced by blow molding the elastomeric alloy. This process allows a designer the opportunity to design the shape and function of the bellows without concern for part ejection as in injection molding. The part can be modeled as a cross section symmetric down the center line. Typically, the part is analyzed for its desired fit or movement in the system. Any interferences, excessive collapse forces or high stresses can be corrected in this design phase.

Conclusions

The use of computer design tools such as geometrical modeling software, and finite element analysis software provide powerful tools in the development of better elastomeric alloy part designs. Part designs can be compared visually for suitability. The effects of loads and contact on the part's ultimate deformed geometry, the internal stress and strain can be visualized. The computer tool allows the faster evaluation of elastomeric alloy part designs and the relative suitability of the part can be evaluated. The ultimate effect is to allow faster iterations through the design-development-testing cycle. Generally fewer prototype parts are required and a faster iteration to the most appropriate design is achieved. The computer analysis currently available does not preclude the use of prototype parts. One must still use the appropriate testing on actual parts prior to approval for production and commercial use. The computer models only account for a limited number of factors and cannot account for all potential environments, aging effects and overall system interactions. But elastomeric alloys have been predicted to be suitable replacements for typical thermoset rubbers based on computerized design analysis. The results on prototype tests have confirmed these results in a variety of commercial applications for rubbers. The selection of the appropriate grade of elastomeric alloy is made more clearly for more complex parts using computer analysis. Fine tune geometry adjustments can be done to a part with computer analysis to render the part design more appropriate for an elastomeric alloy.

References

(1.)Rader, C.P., Rubber & Plastics News, October 28, 1992, p 63.

(2.)O'Connor, G.E., Fath, M.A., Rubber World, 184 (12), 1981; ibid. 185 (1), (1982).

(3.)Walker, B., Rader, C.P., Eds. "Handbook of thermoplastic elastomers," Van Nostrand Reinhold, New York, 1988.

(4.)Tapley, B.D. ed., "Eshbach's handbook of engineering fundamentals," John Wiley & Sons, New York, 1990, ch. 4, p. 6.

(5.)Mooney, M., J. Appl. Phys., 11, 582 (1940).

(6.)"Abaqus theory manual," Hibbit, Karlsson & Sorensen, Pawtucket, RlI, 1992.

(7.)Ogden, R.W., in "Elastic deformations of rubberlike solids," Hopkins, H.G., Sewell, M.J. Eds., Pergamon Press, Oxford, UK, 1982, p. 499.

(8.)Ogden, R.W., Rubber Chem. Technol., 59, 361 (1986)

(9.)Finney, R.H.; Kumar, A.; Rubber Chem. Technol.; 61; 879 (1988).

(10.)Kardestuncer, H., Ed. "Finite element handbook," McGraw-Hill, New York, 1987.

(11.)Bathe, K.J., "Finite element procedures in engineering analysis," Prentice-Hall, Englewood Cliffs, NJ, 1982.

(12.)Morman, K.N., Rubber Chem. Technol., 61, 503 (1988).

(13)Seki, W., Fukabori, Y., Iseda, Y., Matsunaga, T., Rubber Chem. Technol., 60, 856 (1987).

(14.)Schrung, D.J., Clark, J.D., Rubber Chem. Technol., 61, 669 (1988).

(15.)Lau, J.H., Jeans, A.H., J. Appl. Mech., 56, 751 (1989).