T H E UN I V E R S IT Y OF M I C H I G A N COLLEGE OF ENGINEERING Department of Civil Engineering STRUCTURAL DESIGN OF ASPHALT PAVEMENTS Mihai Raf iroiu Lecturer in Civil Engineering Polytechnic Institute Timosara, Romania partially supported by: Funds from the Second International Conference on the Structural Design of Asphalt Pavements Department of Civil Engineering The University of Michigan July 1974

1. Introduction This report is part of a study on aspects of Materials Resource Management in Highway Construction. In order to find a practical system (procedure) for economical de.sign and selection of materials for road structures, the structural design of asphalt pavements will be discussed. The report has 5 main parts: 2. The Problem. 3. Stresses and Strains in Homogeneous Masses. 4. Stresses and Strains in Layered Systems. 5. Deriving a New Design Method. 6. Conclusions. 2. The Problem As it was shown in Report No. 1, Chapter 2, the construction materials transportation represents about 60% of all transportation costs or 9 to 18% of the total investment in the highway. It should also be noted that the costs of the pavement amount to 24 to 50% of the total investment and no reduction of, the pavement cost is possible without a good understanding of the behavior of

-2the road body under current loads. The topic of this report deals with a semiempirical approach for a new design method based on the energy consumed during the loading of the pavement. The new design method will be used in further investigations, in order to get an appropriate procedure for the economical design of the pavement and the reduction of pavement construction costs by taking into account all future maintenance operations and costs. The new design method will be also used in order to assess the needed quantities of materials upon the highway sections before any optimization of the materials transportation is done. In this report only part 3, entitled "Stresses and Strains in Homogeneous Masses" deals with things already known. All further developments were done by the author of this report even if quotations from the technical literature are given. 3. Stresses and Strains in Homogeneous Masses The distribution of the stresses below a concentrated or circular load has been studied by many researchers. The general equations have been obtained for an ideal material which is perfectly elastic and homogeneous and obeys Hooke's law [Ref. 4, page 21]. Many experiments have demonstrated that these assumptions are generally acceptable for a low level of stresses and displacements because the measured values of the stresses and strains reasonably agree with the computed ones. The theory of elasticity [Ref. 14, page 174] demonstrates that

the stresses in axial symmetric and homogeneous masses can be computed with equations (1) where p is an arbitrary function which has to fit some boundary conditions. The general equations for the stresses are the following: Or ~z 2 zr = -r (1) a 2 a z = a [(2-1) T r= [(l1-) V24 a- 2] where all parameters are defined as shown in Fig. 1. The boundary conditions used in order to define the function ~ are generally reduced to one: 2 + - a --- + a 2 ) 2)r( A+ 2 r (2) 3r r 3r az ar r;r az usually written as follows: v2 (V 2) = 0 (2.a) The literature in this area shows that one can use many different

-4solutions for c and no restrictions are formulated to limit the use of new expressions for ~, except those in (2) or (2.a), which means that any function p used has to be at least a biharmonic function. As is well known, these types of functions are called stress functions. The general equations for the strains can be also formulated in the same way as the stresses but it was demonstrated by Gallerkine [Ref. 14, page 176] that one can get the strains with more accuracy by using equations (3), where (l, D2 and c3 are scalar biharmonic functions: aw 2 u = B * - + 1 ax aw 2 v= B * - + V 2 (3) ay aw 2 w = B ~ - + V \3 az where: asl _ 2 a3 w = --- + - ax 9y az (4) B 1 2 ~ (1- -) and 1l' m2 and e3 have to satisfy the condition in (2) or (2.a). The parameters in (3) are defined as shown in Fig. 2. It is to be added that, according to Ref. 14, pages 300-302,

-5the idea is to use the principle of the minimal energy, where the potential energy must be expressed in terms of deformation components, namely in terms of displacements at the top level of the halfspace, the integration being extended throughout all the mass of the halfspace. It is well known that Boussinesq used in his developments an Airy function for ~ and his results may be written as follows [Ref. 4, page 22]: [ (a2+ z 2)3/2]l (5) =_ _2 ( _1 + ) z z3 r= + 2. 1 2 (1 2i/2 + 2 2 3/2 (a + z') 2 (a +z ) (a + where all parameters are defined as shown in Fig. 3. One can complete the equations (5) (that are given on a vertical axis passing through the center of the circular load) with the equation necessary for shear stress computing [Ref. 14, page 108]: Trz =zr = ( - ) (6) the equation (6) being the simplest possible. Boussinesq's problem was solved for different cases of loading. Only two of them will be presented herein. The first case deals with an uniform loaded circular area as shown in Fig. 4. The deflection under the center of the loaded area (w0) is

-6given by equation (7) [Ref. 14, page 257]: 2. (1- 2 w P (7) * a E where: p is the Poisson's ratio of the loaded mass E - Young's modulus for the same mass a - radius of the circular area P - total load on the circular area; it can be computed with relation in (8): P =T D ~ p (8) where: D is the diameter of the circular area, and p - the unit load on the. circular plate. In terms of D and p, equation (7) can be presented as follows: 2 w = (1- * * D (9) E The deflection at the edge of the plate (w ) is given by a equation (10) presented also in terms of D and p [Ref. 14, page 258]: 2 (1- 2) w =- p ~ D (10) a Tr E The ratio Wo/W is: = _ 1.57 (11) W 2

-7which means that w is almost two-thirds of w a.... o It should be said that the edge of the loaded area corresponds to an inflection point on the deflection curve as shown in Fig. 4. The second case deals with a "semispherical" loaded circular area as shown in Fig. 5. The deflection under the center of the loaded area (w S) is given by [Ref. 14, page 260]: = 3 ~ff (1 - P2) WO =-. *p D (12) 8 E and the deflection under the edge of the loaded area (w a): s 3 ~ (1 -2) wa p D (13) 16 E It means that, for this particular type of loading, the ratio wo /wa will be: O a s 0 3 7' T~ 16 16 =2 (14) w s 8 3 ~ a which means that w s is one-half of w s a o This fact is very important for all further investigations done in this report. But BEZUHOV [Ref. 14, page 260] adds: "If the radius of the deformed surface at the top level of the mass is'large' in comparison with the radius of the loaded area (and usually this is the situation), equation (15) can be considered, from a practical point of view, as being the equation

-8of a certain spherical surface." Equation (15) is given below: (1- 2) 2 5) r o E D where r represents the distance from the center of the plate to a point situated inside the loaded area. The assertion is very important for the developments which will be made through this work, because the computation of the tensile stresses at the bottom of the layers will be much simplified. 4. Stresses and Strains in Layered Systems 4.1. Introduction The problem of stresses and strains in layered systems has been under study for a long time and will continue well into the future. A great many papers, articles and reports have been published on this matter and only a few will be reviewed here. The author of this report will present mainly his concepts, and only the most recent developments will be presented in more detail. 4.2. Notations D - diameter of the loaded area (assumed as a circular plate), at the top of the pavement, in units of length; a - radius of the loaded area, in units of length;

-9P - gross or total load on the circular plate (total load on single wheel), in units of force; p - unit load on the circular plate, in units of force per unit of area; N - present annual traffic, in EWL (Equivalent Wheel Load); NR - annual traffic increase rate, dimensionless; NT - total estimated traffic for the service life of the pavement, in EWL; K - wheel equivalency factor, dimensionless; Ks - season equivalency factor, dimensionless; S1 - service life of the pavement, in years; r - horizontal distance from the center of the load, in units of length; Z - depth under the center of the plate, generally expressed by the relative depth Z/D, in units of length or dimensionless; E. - Young's modulus of the material in layer i, in 1 units of force per unit of area; - Poisson's ratio of the material in layer i, dimensionless; h. - thickness of layer i, in units of length; - natural density of the material in layer i, in units of weight per unit of volume; - angle of internal friction of the material in layer i, in degrees; - angle of superficial friction between two layers i and j, in degrees; f.. - friction factor, at the interface of two layers i and and j, dimensionless, (f.i = tan kij)

-10A. - area of assumed spherical surface at the top of layer i, in units of area; a. - distribution factor of the load, at the top level 2 of layer i, dimensionless; (i. ~ D) A Tr * 1i = -4 A. - multiplier, representing the ratio between the diameter of the surface A. and D; dimensionless; - vertical stress at depth z, under the center of the plate, in units of force per unit of area; ay = - radial stress at depth z, under the center of the y r plate, in units of force per unit of area; T z - shear stress at depth z, under the center of the plate, in units of force per unit of area; Tmax - maximum tangential stress which can be mobilized z at depth z, in units of force per unit of area; ~ - unit vertical strain at depth z, under the center of the plate, dimensionless; a. - vertical compressive stress at the top of layer i, 1 under the center of the plate, in units of force per unit of area; a. - maximum bending stress in layer i, in units of force per unit of area; m T. - maximum shear stress in layer i, in units of force per unit of area; O. - maximum radial stress at the bottom of layer i (equal to radial stress at the top of layer underneath), in units of force per unit of area; of - maximum friction stress at the bottom of layer i (equal to friction stress at the top of layer underneath), in units of force per unit of area;

-11a - average compressive stress in layer i, under the center of the plate, in units of force per unit of area; gP _-,compressive stress at the bottom of layer i, due to the natural weight of the above layers; a T. - allowable vertical shear stress along the circum1 ference of the assumed loaded area, for layer i, in units of force per unit of area; ad pid - permitted pressure (compressive stress) at the top of layer i, in units of force per unit of area; MXZ - bending moment along axis x, in units of force times units of length; T - total shear force along axis x, in units of force; x - Ri - radius of curvature at the top level of layer i, in units of length; Wa - permitted deflection at the surface of the pavement in the center of the loaded area, in units of length; IW - deflection at the top level of layer i, under the edge of the surface Aj, where j is the index of another layer, in units of length; wpj - partial deflection due to the compression of layer i, under the edge of the surface Aj, in units of length; We - exterior energy transferred by the load during the deformation of the highway body, in units of force times units of length; W.int - interior energy gained by the highway body during deformation, in units of force times units of length; Wt - subgrade energy gained during deformation, in units of force times units of length;

W - pavement energy, gained during deformation, in sr units of force times units of length; Wsi - layer energy, gained in layer i during deformation, in units of force times units of length; Fd - correction factor, used in the computation of the permitted deflection of the surface, in the center of the loaded area; Fd = 1 if one uses the International Units of Measure; Fd = 0.3937 if one uses the American Units of Measure; F - deflection reduction factor, which represents the ratio.between the average deflection of the surface, inside the loaded area, and the deflection under the center of the plate; dimensionless; - correction factor, used in the computation of the deflection reduction factor, dimensionless; c - correction factor, used in stress computation formula, dimensionless; b - correction factor, used in equivalent layer thickness computation, dimensionless; FS - safety factor, dimensionless; F1 - "lack of fit" factor, dimensionless; Zf - frost penetration depth, in units of length. 4.3. Study on Stresses and Strains in Flexible Pavements Fundamental Assumptions In order to develop a particular theory concerning the distribution of stresses and strains in flexible pavements, several assumptions were made as follows: - Each layer acts as a continuous, isotropic, homogeneous,

-13linearly elastic medium infinite in horizontal extent; - The surface loading can be represented by a uniformly distributed vertical stress acting over a circular area; - The surface loading is elastic; - The interface conditions between layers can be represented as being partially rough; which means that friction stresses can occur between layers; - Inertia forces are negligible; - The principle of superposition is valid (in the pavement analysis); - The influence of the other wheels, than that which is taken into account, is negligible; - Each layer, except the subgrade, behaves as a slab of a low rigidity, for a given range of stresses and strains; - Deformations throughout the system are small; - The vertical deflection is extremely small in comparison with the thickness of the layers; - The temperature, the moisture contents and the frost penetration can be taken into consideration by their influence on the moduli of the layers and Poisson's ratios; - The failure of flexible pavements is essentially a fatigue phenomenon, though the distress mechanisms which appear to be most prevalent are: (1) facture resulting from repeated applications of vehicular loads; (2) distortion which is traffic associated, and (3) fracture resulting from non-traffic associated factors. The behavior of the highway body under the given load is imagined to be as shown in Fig. 6. Each layer is characterized by its Young's modulus, Poisson's ratio and thickness, except for the subgrade whose thickness is infinite.

-14Each layer, except the subgrade, is assumed to be stressed by both the load transmitted through layers above and by the bending of the layer resulting from its following the deflection of the subgrade. This double loading of the layer is due to the relatively low friction stresses that can be mobilized at the interface of each layer. For this reason, the most important surface which is to be characterized is the deformed surface at the top of the subgrade, which is concordance with the theory of elasticity (Ref. 14, pages 300-302, quoted above). For the m-layered system, one can present the stresses under the center of the load as shown in Fig. 7. As boundary conditions at the interface between two layers in contact one can write: ri rs r vi Ci-l = (16) afi = fs f zi 1i-1 i As a consequence of equations (16) the principal stresses under the center of the loaded plate and for the layer i can be presented as shown in Fig. 8. 4.4. The Equation of the Deformed Surface at the Top of the Subgrade Many equations have been proposed for characterizing the deformed surface at the top of the subgrade.

-15Generally, the proposed equations are of the second degree as follows: a' w = 2 (17) b' * r + c' * r + d' where a', b', c' and d' are coefficients that differ from those defined in the list of notations. After studies, Leger (Ref. 3, page 1196) proposed a simplified equation: w w = (18) r 2 a' r + b' The main difference between equations (17) and (18) is of the number of external conditions necessary to define the unknown values of the parameters a', b', c' and d'. Leger's equation needs only three external conditions as shown in Table 1. It should be noted that the conditions in Table 1 are given according to the model presented in Fig. 6 (where the load at the top of the subgrade is assumed to be distributed semispherically) and is in equation (14). It should also be noted that the shape and the general equation of the surface studied corresponds rather well to the experimental findings of Coffmann (Ref. 2, page 819-862), Brown and Pell (Ref. 2, page 487-504) and Rafiroiu (Ref. 7, page 14-16), and the procedure used for defining the values of the parameters a and b in equation (18) differs from that which was used by Leger.

-16According to the condition 1 from the Table 1, when r = 0, w must be wo. Then: w 0 W 2 (19) a (0) + b' and one finds immediately that b = 1. o D According to the second condition in Table 1, when r = 2 w w = 2~ Therefore: w w o o 0 0 (20) 2 - D 2 2 from which: 4 a' = 2 2 (20.a) = ~ a D Now one can write the complete equation for the deformed surface (the third condition being automatically satisfied): w w = (21) r 4 2 2 2r + *D o But it is known that the second derivative of the deflection gives the curvature of the surface at the point of differentiation: ~d~~~~w _ ~ ~1 ~(22) dr r

-17From equation (21): 8 wo 2 2 e *D dw o (23) -2 * r + 1)2 and: 2 2r + 282 w0 2 r ow0 dr ( r2 + 1 8 _ a2 =8 w (25) and, finally: a 2 D2 RO (26) Equation (26) corresponds with the experimental findings reported in Ref. 7, page 22. But, the validity of equation (26) can also be demonstrated by several experimental findings as follows: - Leger and Autret, in Ref. 3, page 1199 report that the value of the product R (m)' wo (i00 mm) varies, as an average, from 6000 to 7500.

-18- Scala and Dickinson, in Ref. 2, page 971 report that the product R (ft) * wo (10 o in) varies from 2250 to 12500, the average being 9000 which corresponds in Leger's units of measure to 7000. - Dehlen, in Ref. 1, page 814 reports different measured values, which, after interpretations, lead to the same values as those found out by Scala and Dickinson. - Finally, Grant and Walker, in Ref. 3, page 1161, using the same units of measure as Leger, report that the variation of the value of the product RI w~ is to be related to the ratio E1/Eo, where E1 would be the average modulus of elasticity of the pavement and to the modulus of the subgrade. Some computations were made within this work in order to verify whether equation (26) satisfies the experimental findings mentioned above. Three different types of pavements were considered: - Gravel on subgrade (two layered system), - Macadam on gravel and subgrade (three layered system), and - Asphalt concrete on Macadam, gravel and subgrade (four layered system). The thickness of the layers were maintained constant: 5 cm (2 in) for the asphalt concrete surface work, 10 cm (4 in) for the Macadam base course, 27 cm (~11 in) for the gravel subbase. The modulus of the subgrade was varied, in order to get different values for the ratio E1/E.

-19The results of the computations are given in Table 2. One can see immediately that the values given by Grant and Walker can be considered a good approach though they had enough experimental data only for E1/Eo ratios between 1 and 3. One can also see that equation (26) gives very good results, and it seems closer to experimental findings than that of Grant and Walker. It should. also be said that Leger and Autret in Ref. 3, page 1197,give the distribution of R ~ wo, which shows variation of R. w~ values related to the moisture of the subgrade (dry season —average value = 8200; wet season —average value = 6600). If the moisture influences the modulus of the subgrade, the lower it is, the lower the radius of curvature is, and this fact agrees with equation (26). Nevertheless one shall see that equation (26) gives a radius of curvature a little larger than that which can be measured at the peak of the deflection curve and a reduction coefficient will be introduced later. This is due to the fact that in reality, although the layers have almost the same radius of curvature as the top of the subgrade, there is a certain difference between radius measured at the top of the surface layer and that measured at the top of the subgrade. 4.5. Stresses in Layered Systems 4.5.1. Vertical Stress As given in equation (5), the vertical stress in a homogeneous mass under a circular plate can be computed by Boussinesq's equation, which may be written as follows (Ref. 4, page 31):

-20{= z 2 ]2 1/2 (27) (a + z (27) wherein n is called "Froehlich's coefficient" or "concentration factor." For the particular case of an homogeneous mass, n = 3. For layered systems different values for this factor have been suggested, the most used being n = 2 as shown in Ref. 1, page 538. It has been also suggested by many authors that, for layered systems, Z represents the equivalent depth of a point estimated under the center of the loaded area. It is well known that, for a two layered system (where E, is the modulus of the surface course and to the modulus of the subgrade), the equivalent thickness of the surface course is (hleg): 3 E1 hleg bl1h1. W = al h1 (28) Many values have been the correction factor bl, but the most used one is b1 = 0.9 as proposed by Odemark long ago. If one replaces the actual depth Z with the relative depth (x): X = (29) and one uses n = 2 for the concentration factor, equation (27) will become: az = p 2 1 al... (30) [ 1 -:0a25.X

-21In order to check the accuracy of the equation (30) many papers were studied and all experimental data, available in those papers, concerning this peculiar kind of problem were collected (Ref. 1, page 284; Ref. 2, pages 303, 493, 698-701, 713; Ref. 3, pages 411-420, 478-486, and 877-900). All of these experimental data are presented in Fig. 9. But, unfortunately, neither Boussinesq's equation, nor equation (30) with al = 0.9. E1/Eo, fitted the experimental data. For this reason a new equation was written as shown in (31): a2 2 x (31) 1 x 1z =p c 1 0.25 + al x and different values for al and c1 were used. In order to speed up the calculation, a computer program was written and it is presented in Appendix 1. After an accuracy check had been made, new values for al and cl were obtained: a a 31 1 a 2l = (32) a2 = 0.84; c1 = 0.97 Therefore, the empirical equation for the computation of the vertical stress is the following:

-223E 2 a2 p V E a x (33) z = p C1 2 0.25 + a2 * x a2. E x L 0 The accuracy of the results, which can be obtained by using equation (33), was checked by replacing at each particular depth, the actual data with the mean. The results are shown in Fig. 10 and they demonstrate the high accuracy of equation (33). Although one can appreciate that the differences between coefficients al and cl in (33) and those generally used (a = 0.9 /E; = 1.00) are rather small, it should be said that the latter do not fit at all the experimental data for relative depths below 2.0 where usually the subgrade-subbase interface lies and hence the vertical stress on the subgrade is not correctly estimated. 4.5.2. Radial Stress If we rewrite equation (5) by using the following replacement: a = (34) \2 +2 a +z

-23one can get: r= P [ 1 + 2 * - 2 (1 + A + A3] * (35) If we set: B = A3 (36) equation (35) will become: = 2 [ 1 + 2 * - 2 * (1 + p) A + B (37) r 2 A + B] Finally, we can write: A = on (38) wherein n - for Boussinesq's problem - is equal to 3. But, in the previous developments regarding the vertical stress, it was found: ( )2 B = c (39) 0.25 + a2 3' x) and n = 2. For this reason one can write: ( 2 312 0.25+ a2~ 3 l ~ x3 0 /

-24or: 3 Y1 a2 * E ~ x A, 3' )2 0.25 + a2'Vo* x and therefore: ~r 2- [1 + 2 ~ + 2* (1 + I) 1 B(42) r 2 L (42) The accuracy of equation (42) was checked by also using the experimental data available in the studied literature (Ref. 2, pages 495 and 712) and the results are given in Fig. 11. The computed data agree reasonably well with measurements, especially at the larger relative depths. For this reason one can consider that all developments in equations (34) to (42) are correct, at least from an empirical point of view. A computer program for the computation of the radial stress at different levels is given in Appendix 2, 4.5.3. Shear Stress The literature in the area is rather poor in information about shear stresses and no experimental data could be found. This is the reason why only a theoretical development has been done.

It is well known from the theory of elasticity (Ref. 14, page 108, and Ref. 19, page 370) that the maximum shear stress is given by the relationship: 1 Tzr = 2 (z or) (43) By replacing Or and Cr by equations in (33) and (42) z r respectively and after simplifying one gets: zr 2 l 2 + (1 + p) * A1 3' B T 2 2 1 2(44) or Tzr = [ (1 - 2 * p) + 2 * (1 + p) * A1 - 3 ~ B1 ] (45) Equation (45) was used in the computer program in Appendix 3. The results are shown in Fig. 12 for two values of Poisson's ratio, ~ = 0.35 and i = 0.50. In order to make a comparison Boussinesq's formula was also used. According to Ref. 19, page 370, Boussinesq's formula is: 1 2. z 3 (46) Tzr 2 2 2 2 2 (a2 + z2 3)2 a + z

-26It is very interesting to observe the Boussinesq's curves have the same shape, but have the maximum values at a depth two times greater than the other curves, and approximately the same values at the depths z = 0 and z = oIt should also be noted that the maximum values of the shear stress obtained using equation (45) are greater than those using equation (46), namely, in layered systems shear stresses should be greater than through an halfspace under the same loading. It is to be said that high shear stress could reduce considerably the bearing capacity of the pavement, especially when the maximum shear stress lies in a weak material as sandy gravel or sand which has a small angle of internal friction when wet. The writer would like to report his own experience with finding material displaced by the action of shear stresses, in the base course of an asphalt pavement. The distinctive action of the shear stresses was observed in a three-layered system (5 cm (2 in) asphalt surface course; 30 cm (12 in) sandy gravel base course; and silty clay as subgrade), in the Spring, after a year of heavy traffic. The equivalent load was 5000 kg per wheel, with about 30 cm (12 in) equivalent diameter of the loaded area and about 5 kg/cm2 (= 69 psi) tire pressure. At about 10 cm (4 in) under the top level of the sandy gravel base course, a region about 2 cm thick, obviously disturbed, was observed, the rest of the material remaining well compacted both above and below this region.

-27This means, in the writer's opinion, that the thin asphalt concrete surface course had not had any important influence either on the bearing capacity of the pavement or on the reduction of the shear stresses, especially because it was intensely cracked. The shear stresses displaced the material at about 0.33 h/d, which fits almost perfectly with the theoretical findings. 4.5.4 Stress Function As was shown in Part 3, the stress formulas are derived from the four basic equations (1) by the use of a stress function. The general requirement for the stress function is that it be biharmonic. But it seems that more than one requirement is necessary. As a matter of fact, at least one further requirement should be formulated: a stress function should also lead to formulas whose results must agree with experiments. Nowadays, a limited number of stress functions are used: Airy functions, Bessel functions, polynomials, and exponential functions. Each of them is arbitrarily chosen, and the usual way to derive the stress formulas is to choose one function, introduce it in equation (1) and solve it in such a way that all constants will satisfy the boundary conditions. This is usually called the "inverse procedure." The "direct procedure" is to integrate the stress formulas to get the general form of the stress function and then, determine the constants by using the boundary conditions.

-28Although the purpose of this work is not to solve this particular kind of problem, the homogeneity of the procedure used to get the stress formulas suggested to the writer that a stress function should exist even for the "empirical" stress formulas given in equations (33), (42), and (44), Unfortunately, these formulas are given only for the region right below the center of the loaded area. This is the reason why the following stress function search aims only at finding an expression which would be useful in guiding further choices. For this particular problem, r = 0 and hence: 2 ~a = 0 ar 1 *'f _ (47) r Dr Therefore: v2 (V2)= 2 = a% a)= O (48) az Dz which is the same as: = 0 (49) az for z = 0 and z = A. Although, from now on, one could use common derivatives through the theoretical developments, the partial derivatives were however maintained in order not to change the general form of the problem.

-29If one analyzes equation (33) one will see that for z = 0: a = p = constant z (50) z =0 and for z = (51) atz As a conclusion of the relationships in (50) and (51) one can state that any stress function obtained as a result of a triple integration process applied to the above stress formulas will be a biharmonic function and it will also satisfy the requirement of the experimental agreement. For a simpler development one will consider that: al = a2* E b = bi (52) c = 0.25

-30and equation (33) becomes: which will be more convenient for the integration process. As a conclusion of the relationships in (47) and (48) one can write: 22 2 V = 2 (54) ax We know that: z= — 2''2- 5 52 az a 1(2 -i)* 2 2 2 (56) which will be more convenient for the integration process. can write: 2 2 or: = ax L a2p21 (1 - ) (57)

-31Then, we can write: 3 a = (1 - ) (58) z a and: 2 2 p 1 - b 2a 2 2 (1 - ) 59) 2 2 2 3 c + a1 ~ ax The first integration can now be performed: z 2 2 2 a * X 2~ = P ~ 1 - b ~ x 2 (1 - t) 2 2 2 2)(60) Dx c + a1 X After integration one gets: 2 / 2 P ~ x+ k1b a x 2 1- 11 b a k x c 1 V- 1 tan-i-+ k (61) a1 a1 c 0

-32The constants kl and k2 are determined by boundary conditions. 2 z 0; 2 0 and k = = (62) 1 2 ax Thus: P ~ (63) __'1= (13 tan-1 c x2 1 3 c Now, one can make the second integration: (~z F c a x] = p(1- b) x -3 * tan-i1 C. dx (64) aal 0 which is: ri~ 2 = x c (65) a1 a e x c 2 2 2 1 *x tan-1 c n ( + a2 *x2) + k c 2 4 0 The constants k3 and k4 are also determined by boundary conditions: z=; = = 0; k= 0; k 2 (66) at ~3

-33And now, the last integration: =' (1 - b) *. a k tana1 vO (67) in (c + al x2) + n c2 dx 2 11 2) After integration one gets: 3 (1 - b) X c c 1 - P 2 3 5 2 a1 z a1 x c. a1 x (c + al x) tan-1 k 6 2 1 c 2 6 0 +c 1 2 2 2 + c3 b ~ a1 * x i in (c + a1' x ) - 2 ~ a1 (68) al k 2 al k + 2 ~ tan-i 1 + k 0 Finally, by boundary conditions: z = 0; = 0; k5 = 0; k6 = 0; k7 0 (69)

.34 4 The final form of the stress function, only for the region underneath the center of the loaded area, is: 22 aa 1 (1- b) x3 c * a2 + 2 * c2 - in (c2 + a2 x2. x (70) cc 1 c2 + 2 x l1 c } a2 ~ al (c+ a1 ~ x tan-1l wherein x = Z/D as presented in equation (29). If it is assumed that al = 1.0 and b = 1.0, and i represents the ratio between z and the radius of the loaded area, equation (70) becomes: = -. { 1 2 + 2 -in (1 + 2) ~I -(1 - ui) 2 L (71) + 1 (l+ 2)~ tan-1 i or 1 = - 2(1p - i- n (1- 2)1. ~ + (1 + 2) * tan-1 (72) 21- vi)

-35This simplifed form of the stress function shows that the problem of the stresses in layered systems is a very complicated one and much more investigation is necessary to solve it accurately from either a theoretical or practical point of view. 4.6. Strains in Layered Systems Although we have both stress formulas and a stress function, the problem cannot be solved throughout the whole highway body. All the formulas or equations given are only for the region underneath the center of the loaded area. This is the reason why a quite different approach for strains computation was done. Two main assumptions were made: - The surface at the top level of the subgrade is deformed as shown on the right-hand side of Fig. 6; - Each layer deforms itself under the load, and the higher the unit vertical stress is the thinner the layer becomes. If one measures the vertical deflections at the top of each layer, in comparison with the same surface before deformation and at points situated in the same vertical plane as the edges of the Ai areas, one will get the following matrix (see Fig. 6):

-36m+l m m-1 2 1 o o 0 w0 w0 w 0 w O O O O O O m+l m m-1 2 1 o W1 Wi w1iW W W1 w1 m+l m m-1 2 1 o W2 W2 Ww2........W2 W1 Wo W =.... (73) m+l m m-i 2 1 o W W W.W W. W Wm-l nm-1 Wml.2 w1 m+l m m-1 2 1 o Wm w w..W2 1 W nm- 2 1 o The computation of the elements of the matrix W can be done only after the compuation of the elements belonging to another matrix (Wp). In the latter matrix the elements, noted by wPi represent the partial strain which occurs in layer i at the projection of the edge of the surface A,. The matrix W will be: p pm+l P pm -1.p2 p 1 po w1 1 1 1 1.... w pm+l pm pm-1 p2 pl po 2 w2 w2 w2 w2 2 W = ~ ~.... (74) wPm+l WP wPm-l wp2 p1 WPO M-1 Min in1 in1 in in

-37The computation of the values of the elements in the matrices W and W is the following: p a) For the first row in W: 1- 2 ~ D wm0l = X E -E a (75a) O O where: a 2 (76a) 3 o and: mj +l F w =w 1 j (75b) j = 1,..., m where: j = Pv (76b) 3except for the surface course whose = 1.0 from geometric except for the surface course whose am = 1.0 from geometric conditions. a. being calculated by the following equation:

-38a2 h r-i P 1 - C 2 m 0.25 + 2 hl D =j0 which is derived from equation (33). It should be noted that equation (75a) is derived from equation (12) which is obtained from the theory of elasticity, and equation (75b) is derived from equation (15) which is also obtained from the theory of elasticity. Computations made with the computer-program, Appendix 5, demonstrated that the deflection calculated by using equation (75b) are very close to those calculated by using equation which is based on the theory of elasticity. Actually, the differences between the results of the two equations are in the range of 0.2 to 0.06 mm (0.003 to 0.009 in.) hence are very small. b. For each row in W: p If ~ is the strain at the level z under the center of the z loaded area, according to the theory of elasticity: E [iz - i' (ax + oy)] (78) z -' x v

-39In our particular case a = a = a. Then: x y r = Cy[a-2e i Y (79) z E v If med and med vIf and ri are respectively: v v med a-1 +1 med i+ + a Vi 2 one can write: wpm+ hi ed med wi E.= - ri (81) In equation (80), the radial stress ai can be computed from the following equation: r. = 1 + 2 1 - 2 * (1 + v) 1 + i i (82) 21 2 l0.25 + di 0.25 + d. where: ( m 2 di hl l-=i~ ~I~d"~eE -) (83)

-40Because of the similitude of many terms in the i's and (r,'s formulas, after computations one may write: 1 h r avi 21 med wpm+l = i'(l + P) -med 2 p * 1 -. (84) i... E. vi p 1 L If one accepts that the partial deflection of a layer is constant within the loaded area Ai, one may write: pj = wpm+l (j =..m) (85) 1 1 The main problem now is to define the partial deflection of a layer outside of the surface Ai. In order to solve this problem two equations were used. The first equation is suggested by both equations (15) and (21): wp wm+l i J (j 1,.., (86) The second equation is suggested by the developments of M. M. Filonenko-Borodits (Ref. 14, page 474): v 1 1 D c /M=C' * t r

-41Computations made as part of the computer-program in Appendix 5 demonstrated that there are no real differences between the results of equations (86) and (87), at least for the computation of the radii of curvature. c. For the rest of m rows in W: wj = w + wj (i = 1,... m) (88) 1= -'1 (88) j = 1,..., m+l The procedure described above and which is based upon the developments concerning the stresses in layered systems allows the computation of the deflections at a net of points situated in the most stressed part of the pavement, namely inside a cylinder whose radius is equal to o ~ D/2. Outside of this cylinder, the computation of the strains would be very difficult and, in fact, the computation of the deflections of the top surface of the pavement by using equation (21) is sufficient for practice. The computer-program, given in Appendix 5, was used in order to check the validity of equations (73) through (86). The basic data for the computations were obtained from Referenee 1 (pages 511-573; 655; 404), Reference 2 (pages 303-305, 677-683; 698-701; 872-874; 1084-1085) and Reference 3 (pages 411-420; 478-486; 789-794; 796-797; 877-900; 1133; 590-603). These data were completed with all necessary information from Ref. 5, Ref. 7, Ref. 9, Ref. 10 and Ref. 11.

-42It should be said that the majority of data given in the references noted above were not ready to be used in the present developments and a process of evaluating the data had to be done. Therefore, a stochastic approach was used especially as only average values were used for the Young's moduli. The conclusions resulting from the computer processing of these data are the following: - The vertical compressive stresses under the center of the loaded area were predicted correctly by equation (77) (which is the same as equation (33)) as one can see in Fig. 13. The horizontal tensile stresses at the bottom of the layers are about 15/4 of the measured stresses, as shown in Fig. 14. The deflections at the surface of the pavement under the center of the loaded area are about 14/5 of the measured deflections as shown in Fig. 15. - If one takes into account the energy necessary for the change of the form of the layers (procedure according to Ref. 7), the deflections computed by using the procedure given in equations (75) to (88) represent about 14/4.8 of the new deflections as shown in Fig. 16. For this reason the deflections calculated by taking into account the energy necessary to change the form of the layers are expected to be: 5 1, therefore from a statistical 5 14 point of view, equal to those measured. - If real deflections are 5/14 of those calculated by equations (75) through (88), the curvature of the layers and implicitly the calculated tensile stresses at the bottom of the layers are expected to be 5/14 of those calculated by the computer

-43program in Appendix 5. Then, the calculated tensile stresses 5 15 would represent about: 15 4 1.32 of the measured ones. It should be noted that Nijboer (in Ref. 2, pages 699-701) has shown that the tensile stresses at the bottom of the layers, calculated by any method resulting from the theory of elasticity, are almost always higher than those measured (deviations up to 35%). This means that any theoretical approach cannot evaluate correctly both the influence of the friction between layers and the real tensile stresses at these levels. Anyway, it seems to be better to overestimate the tensile stresses at the bottom of the layers, because the higher the tensile stresses are the more cracks at the surface of the pavement occur. It should be noted that all these results may be incorrect as to their exact values. The purpose of their presentation is to put in evidence a qualitative process whose consequences will be studied in the next part of this report. The average ratio E1/Eo, where E1 is the average modulus of the pavement and to the modulus of the subgrade, for the set of given data, is 16; the average number of layers is 5; and the 1 average value of the product R * w (cm ~ TiO mm) is 18000, the range of the latter being between 16000 and 20000. If the results given in Table 2 are interpreted graphically, and if one accepts a possible extrapolation with the increasing number of the layers and the increasing of the E1/Eo ratio, one gets the results shown in Fig. 17. This means to the writer that the products R. w should be a function not only of the E1/E ratios as suggested by Grant and Walker

-44in Ref. 3, page 1161, but also of the number of the layers of the pavement. If this is true one can consider that any approach using an "equivalent halfspace" is wrong. It will remain for further developments to demonstrate the validity of the above assertion. - The last conclusion is that the tensile stresses at the bottom of the layers measured under the center of the load put in evidence that the real radius of curvature of the layers is the same as the radius of curvature of the top surface of the subgrade, under the center of the loaded area. 4.7. Energy Approach 4.7.1. Introduction During the deflection of the highway body under loading, two types of deformations occur: - Deformations due to the change of the volume; - Deformations due to the change of the form. The stresses and strains studied in the last three sections are connected with the change of the volume. The stresses and strains due to the change of the form are studied next. During the deflection of the highway body, an energy phenomenon occurs: the exterior load transfers energy to the highway body and the highway body gains energy. When the exterior load is removed, the highway body (if loaded within elastic limits) loses the energy by recovering its previous form.

-45It is obvious that the phenomenon described implies that no plastic or viscoelastic strains occur in the highway body during its loading. For this is the reason that additional stresses occur in the pavement during deflection and these also should be taken into account. 4.7.2. Bending Stresses As it was demonstrated in some previous works (Ref. 7; Ref. 8; Ref. 9; Ref. 11) the highway layers, except the subgrade, usually behave as slabs of a low rigidity. The explanation of this fact lies in the low friction stresses which can be mobilized at the interfaces of the layers and the prestressed state (compressive stresses) which occurs even in the granular or noncohesive materials during compaction and which disappears only when the pavement has completely failed. Experiments reported in Ref. 11 have demonstrated that the friction factor at the interface of two layers (fi) lies in the range of 0.3 to 0.665, the lower values (0.3-0.4) being reported for rolled granular materials, the medium values (0.4-0.5) for bituminous mixtures placed upon any kind of compacted material, and the higher values (0.5-0.665) for crushed stone, Macadam, or materials which penetrate into each another during compaction; for example, crushed stone on clay or clayey materials, bituminous mixtures spread upon an open textured crushed stone, etc.

-46Consider the region under the center of the loaded area. If one makes the assumption that, in this particular region, one can add radial stresses with shear stresses, and if one considers an average value for fi, for example 0.5, one will get —for the particular depths where the interfaces of the layers lies (see Fig. 9) —the values shown in Table 3. Certainly, the above computation does not observe the rules of the theory of elasticity and it was made only to get a better idea about the general state of stresses in that region. It can be observed immediately that the friction stresses which can be mobilized at the interfaces of the layers are less than the global horizontal stresses which occur at these levels. This is the reason why slip is possible between the layers. But there is neither perfect slip nor rough contact there but an intermediate situation which can be taken into account within computations if the friction factor is known. The experiments reported in Ref. 11 have shown that the friction factor at the interface of two layers is equal to the * smallest value of the tan.i of the materials the layers are built with, except when materials penetrate into each another. In this case fi would be the average of the two values of tan f%. A more interesting problem is that involving the prestressed state of the layers. The prestressed state has been taken into account by many authors. Among them are: Burmister (Ref. 1, page 444), Bonitzer * See list of notations

-47and Leger (Ref. 2, page 783) and Livneh and Shklarsky (Ref. 1, page 345-353) with the mention that the latter used also cohesion characteristics in their developments. The compressive stresses which occur at the bottom of the granular layers and which remain approximately at the same level even after the exterior loading ceases its action were studied in depth by Vasile Izdraila*. He has shown that a prestressed situation can occur even if the material of the layer is composed of steel balls of the same diameter. Also, that the remaining compressive stresses increase with the internal angle of friction and the nonuniformity of the grading of the material. These are the reasons why even the granular layers can behave as slabs of low rigidity. Nevertheless it should be taken into account that the prestressed condition of the layers cannot be correctly estimated within a computation because the remaining compressive stresses at the bottom of the layers are functions of the characteristics of the means of compaction and the number of coverages, and not only of the characteristics of the materials in layers. If the layers of the pavement behave as slabs, one can compute the stresses under the center of the load. Associate Professor at the Polytechnic Institute from Timisoara, Romania, Doctoral Thesis.

-48According to the theory of elasticity (Ref. 14, page 313) one can write: E 2w a 2 1- s ax ay 2+ 2 aY 1 - ay x2 2 E a w 2 2 aI W~ = -; 2 =(90) z 1x ay y Because of the axial symmetry of the load one can write: a2w a2w 1 (91) X2. ay2 R If one replaces the partial derivatives from (89) by the values in (91) and one considers: z -2 (92)

-49one can write:* E. * h. m 1 1 - 2 * R. (1 93) 1 i E. * h. m 1 1 T. - 94) 1 2 Ri (1 + iv) The bending moment Mi will be (Ref. 11, page 440): E. ~ h. (95) 1 12 * Ri (1 - pi) This bending moment will occur as a result of each layer following the subgrade during deformation. But two other moments will occur in the layer: - The moment due to the nonuniform distribution of the radial stress within the thickness of the layer; - The moment due to the friction stresses which occur at the interfaces of the layers. The moment due to the nonuniform distribution of the radial stress can be approximate as suggested in Fig. 18: 2 r r r 1i M = (U - a ) (96) i i+l i 12.0 See the list of notations.

-50The moment due to the friction stresses which occur at the interfaces of the layers can also be approximate as suggested in Fig. 19: f h i f f / M. l,- + CF. ~ (97) Mi ai+l 2 _i 2 i+l 2 The total moment would be: T Mm + Mr + Mf (98) i 1 1 i If one analyzes the variations of the radial and friction stresses within the thickness of the pavement, it will be observed that the bending moments due to the friction stresses always reduce the bending moment due to the following of the subgrade by each layer while the moment due to the radial stresses increases the bending moment. In the writer's opinion the moment due to the nonuniformity of the distribution of the radial stresses should not be necessarily taken into account because the radial stresses and their influence were taken into consideration when the deformations due to the change of volume were studied. For this reason the moment as given in (96) will not be considered further. Previous studies made in this particular area (Ref. 11, page 440) have shown that the smaller the radius of curvature of the surface the smaller the influence of radial and friction stresses. The relationship in (98) agrees with the previous results because if M and M. are independent of the radius of curvature of the 1 1

layer, the bending moment Mm. is directly related to the radius of curvature. Thus, one can check immediately the relationship in (95) and see that the smaller the radiu's of curvature, the greater the bending moment. A logical inference of the fact mentioned above is that the greater the bending moment, the greater the horizontal stresses in the layers. Therefore, if the layers have to support higher stresses than was assumed when only the deformations due to the change of the volume were taken into account, it may be hypothesized that the participation of the layers in the load bearing energy is greater. For this reason, the deflection at the top level of the subgrade would be less than that computed before. The above hypothesis agrees with the conclusions of the previous section and this shall be made use of in the next chapter. 4.8. Shear Stresses on the Contour of the Load In Ref. 3, pages 823 to 843 it was demonstrated that important shear stresses can occur along the contour of the load. This fact fits one of the writer's previous ideas, that, for subgrades with a very low bearing capacity, the shear stresses along the contour of the loaded area, at each level of the pavement, should be taken into account. According to these facts a model for computing these stresses was assumed as shown in Fig. 20. The magnitude of these shear stresses should be less than that which can be supported by the layer (Tai).

-52For this reason we can write (for a layer i): (1/ D) v 1 Ta - r * 4 *CT < - Ti D h (99) V V ad ad If o.. and Cr are replaced respectively by pi and Pi.1 and one considers hi as function of the other parameters, one can write: ad ad. ~ D ~(Ppi -Pi-) h..T1 (100) The relationship in (100) is obviously an approximation but it can help us with computing at least the degree of magnitude of these stresses (or the requiredthickness of the layer to support them). But, the sum of all shear stresses along contours should be equal to the exterior load (as a limit). So one can write: 2 m P < hk a D (101) k=2 and after simplifications: m p D < / aT. h *(102) 4 - E k k k (102) k=l

-53Equation (102) is to be used only to verify, at the level of the first layer (layer number 1), what the remaining hi could be: m 4 k k k k=2 hi > ak (103) 101 If hi calculated with the relationship in (103) is less than h1 calculated with equation (100), it means that the pavement will not be distressed under the load. If the assertion above is not true, one: must increase the thickness of the layer number 1 at least up to the value computed from equation (103). One must also take into account that all these computations must be made for the heaviest load which is expected to use the highway. 5. Deriving a New Design Method 5.1. Design Criteria As a consequence of the developments presented in the previous part, the following criteria were considered: 1. The energy gained by the pavement system including subgrade must be equal to the exterior energy, transferred by the load during deflection. 2. Vertical stresses in the top of the layer (subgrade included) must be less than the allowable stresses.

-543. Tensile stresses at the bottom of the layer must be lower than the allowable stresses. 4. Shear stresses (horizontal) in the center plane of the layer must be lower than the allowable stresses. 5. The maximum vertical stresses should be located in the surface course. 6. The location of maximum shear stresses should be in the base-course. 7. Vertical shear stresses along the circumference of the load must be lower than the allowable stresses. 8. The depth of frost penetration into the pavement must be less than the thickness of the pavement. 9. A minimum practical thickness must be assigned to each layer. One could consider this list of design criteria as being too long. In fact, previous studies (Ref. 7) have demonstrated that, for most subgrade, if the first design factor has been accounted for, the following three ones are also satisfied automatically. But all the design criteria are necessary for determining all the constraints which are to be observed for a structural or economical (optimal) design of the pavement. 5.2. Design Parameters Many parameters are to be taken into account in the design procedure in order to solve all the problems raised by the design criteria.

-55Within the further developments the following design parameters will be taken into consideration (for each layer i): E. - Young's modulus, units of force/unit of area; Pi. - Poisson's ratio, dimensionless; hi - thickness of the layer, units of length; - angle of internal friction, degrees; 1 av a. - allowable vertical compressive stress at the top ad of the layer, units of force/unit of area (= pi ); am aam - allowable horizontal tensile stress at the bottom 1 of the layer, units of force/unit of area; am Ta - permitted shear stress (horizontal), units of force/ units of area; Ta - permitted vertical shear stress along the circumference 1 of the loaded area, units of force/unit of area; fr Cfr - coefficient of resistance to the frost penetration. fr 1 Cfr = (104) where Si is the thermo conductivity of the material in layers i. Many other parameters will be used for the exterior load characterization (environment included): D - diameter, units of length; P - total load, units of force; p - unit load on the loaded area, in units of force/unit of area;

N - present annual traffic, in EWL; NR - annual traffic growth rate, dimensionless; N - total estimated traffic, in EWL; T K - wheel equivalence factor, dimensionless; K - season equivalence factor, dimensionless; S~ - service life of the pavement, years; w - permitted deflection at the surface of the pavement, units of length; F - correction factor, as defined in the list of d notations part 4; Zf - frost penetration depth, in units of length. The main problem is to define the parameters characterizing the behavior of the layers. Some of them are given in the Table 4, and they are a result of several statistical manipulations presented in Ref. 3, 8, 10, 11 and 20. It should be noted that the Young's moduli are defined as being those found out within the plate test, at the second loading of the plate, on the "return" curve. 5.3. The Development of the New Method 5.3.1. Limit Thickness of the Layers As a consequence of the results obtained within the previous part one can define both the minimum and maximum thickness of the layers. The limit thicknesses of the layers can be determined by the aid of the design criteria, namely numbers 2, 3, 4, 5, 6, 7

-57and 9. Having been used once they will not be used again in this development. It should be noted that number 3, 4 and 9 are common for all the layers, except the subgrade. In order to make this development more systematic the equations will be determined for each layer separately. It will be assumed that the pavement is composed from four main layers: surface, base; subbase and drainage courses. 5.3.1.1. Surface Course 5.3.1.1.1. Minimum Thickness a. According to the design criterium number 2 and equation (33) and after computations one can write: Pi h = D (105) min 3. 1 a i \E b. According to the design criterium number 5 the maximum vertical stresses has to lie in the surface course. The maximum vertical stresses are defined herein as being those whose range lies between 0.9 and 1.0 p. For the lower limit, according to the curve shown in Fig. 9, the corresponding value Z/D is 0.3. So one can write: min. = 0 3 a * D (106) 0

or, after computations: h = 0.252. D / (107) mini E0 c. According to the design criterium number 7 and equation (100) we have: D(pad pad ) h (108) min a d. According to the design criterium number 9, one can write: hmin = Constant (109) where Constant i is a value chosen arbitrarily by the designing agency according to its experience. In Romania this value, for the surface course, is 5.0 cm (2.0 in). 5.3.1.1.2. Maximum Thickness a. If one wants to limit the tensile stress at the bottom of the layer, according to the design criteria number 3 and equation (93), one can write: 2 * am (1 - i ) R. hi 1 1 (110) E. 1

-59But, at the beginning of the computations, one does not know what the value of R. is. If one takes into account the fact that R. is a function of the E1/E ratio, the number of the layers and the value of the deflection (which is this case would be the permitted deflection) one can write: R = f (E1/Eo; Nk; wa) (111) After a stochastic approach the function f in (111) was defined as follows: R. 500.0 1(1 Ri w E (112) a wa if Ri and wa are given in the metric system of measure, and 1=200.0 3 500.0 1 R. wa E Fd w (112.a) I W o d Wa Eo if Ri and wa are given in the American system of measure. If one computes the value of hi by equation (110), the 1 thicknesses found would be too conservative, namely too small. For this reason the value of hi calculated is to be multiplied by the safety coefficient, F, used in the design method. As a first approximation one could use the value of 17 which was obtained after more than 2000 experiments reported in Ref. 7, this value being only a provisional one. Now one can write:

-60am 2 * a. * (- i.) R. max. E. Fs (113) 1 1 b. In order to limit the horizontal shear stresses in the center plane of the layer one must observe the design criteria number 4 and use equation (94). If one uses the development presented above one obtains the following equation: am 2 T am (1 + pi) R. h = F (114) max. E. s am But the value of T. in (114) is not a constant one and 1 it can be calculated by Coulomb's equation: am v T. = Co. + a. * tan i (115) 1 where Co. is the cohesion of the material i. am One can see immediately that the value of T. is a function 1 of the vertical stress o. which cannot be calculated until one knows the thickness of the layer. The most conservative approach is to consider that all layers above the layer i have their maximum thickness. Then, the vertical stresses on the layer beneath will be am ~am the least ones and hence the value of T. will be the least Therefore the values of Ti will be a little higher than the least ones.

-61For the surface course one will have: Tam = Co + tan i (116) where the cohesion Co can be measured in the laboratory before the design procedure is started. By using equation (105, (107), (108) and (109) one gets four values for h. It is obvious that the greatest one will min. be used in the computations. Equations (113) and (114) give the maximum values for h.. It is also obvious that the smallest one will be used in the further computation. Because most of the developments are the same, in the presentation below only those which show some differences will be presented in detail. 5.3.1.2. Base Course 5.3.1.2.1. Minimum Thickness a. Design criterium No. 2. ad- ai P1l h = D. (105.a) min. E 1 \IE 0/E

b. Design criterium No. 6. It requires that the maximum horizontal shear stresses lie in the base course. This means that the maximum shear stress should lie in the middle third of the thickness of the base course, in order to better use its resistance to shear failure. According to Figure 12 and equation (46) the maximum shear stresses under the center of the load will occur at a relative depth 0.35 The relative depth of the surface course being: h3 Em h = alDm 3 (117) egm o one can determine immediately the minimum thickness of the base course: hmi 32D al0.35 alm 3 1 (118) After computation and reduction of the constants, equation (118) becomes: h = 0.63 D - 1.5 h 3 (119) min. m E c. Design criterium No. 7. (ad ad h = D ((108.a) min a T.a~ 1

d. Design criterium No. 9. hmin. = constant i (109.a) In Romania the minimum thickness of the base course is 8.0 cm (~3.0 in). 5.3.1.2.2. Maximum Thickness The equations are the same as those in (113) and (114) and nothing is to be gained by repeating them again. 5.3.1.3. Subbase Course 5.3.1.3.1. Minimum Thickness The same equations are used as for the previous layers except for (106) or (119). Therefore, one uses equations (105), (108) and (109) only. In Romania the minimum thickness for the subbase is 15 cm (6.0 in). 5.3.1.3.2. Maximum Thickness One uses equations (113) and (114).

5.3.1.4. Drainage Course The drainage course has only minimum thickness and no maximum thickness. Previous experiments and theoretical developments report in Ref. 7 and 9 revealed the fact that if the thickness of a layer is greater than 2.5* D, then it behaves as having an infinite thickness. It can be considered as an actual subgrade if the real subgrade is not too weak and can support the highway body without appreciable deformations (E > 2 100.0 kgf/cm or E > 1380 psi). It is not the purpose of this work to study this phenomenon but the previous findings reported in Ref. 7 and 9 will be used in the optimal design of the pavements. The minimum thickness can be either a calculated or an evaluated one. The calculated minimum thickness is related to the frost penetration. In Romania it is required that the thickness of the drainage course be such that the thermo-conductivity of the pavement will be the same as that of the natural subgrade with a depth that is equal to the frost penetration depth. It should be noted that this requirement is applicable only to subgrades which are frost susceptible. If the coefficients of resistance to the frost penetration, fr Cf, are known, one can write the equation which will allow the determination of the computed minimum thickness of the drainage course:

m fr Zf - i h. ~ C. hmin. Cfr (120) min. fr o 1 C1 One can observe immediately that the thickness of the drainage course is related to the thicknesses of the other layers. If one accepts the bearing capacity of this layer is not taken into consideration, one can compute its thickness after the thicknesses of the other layers are found. If not, an alternative calculation will be necessary. According to Ref. 7 and 9, if the thickness of the drainage course is less than 0.4 ~ D it is not worth taking it into account. If the thickness of this layer lies between 0.4 ~ D and 0.8 ~ D it is taken into account as a common layer. Between 0.8 ~ D and 2.5 ~ D a correction function % is used: = 0.51/(hi/ai. D)3 (121) Finally, if hi > 2.5 D it is considered as being an actual subgrade and all computations are done over again. The evaluated minimum thickness of the drainage course varies from one design agency to another. In Romania this value is 5.0 cm (2.0 in) but it will not remain at this value as many proposals have been made to increase it up to 15.0 cm (6.0 in).

5.3.2. Actual Development of the New Method The fundamental principle of the new method is the principle of the conservation of the energy. This requires that the exterior energy, transferred by the load during the deflection of the pavement, will be equal to the interior energy gained by the highway body during its deflection. This principle can be presented as follows: We = Wi (122) The interior energy is formed by the subgrade energy Wt and the pavement energy Wsr defined as presented in the list, of notations at the beginning of the fourth part. So, one can write: We = Wt + Wsr (123) Now, the main problem is to find the equations whereby these energies can be computed. Each energy will be considered separately. 5.3.2.1. Exterior Energy At first sight this energy seems to be easy to compute: one multiplies the load by the deflection and is done.

For rigid pavements or for flexible pavements with high bearing capacity, the error could be less than 3%. But for flexible pavements with medium or low bearing capacity the error can increase to 10 and 30% respectively, which cannot be accepted in a design method. It was demonstrated in Ref. 7 that the load being flexible follows the pavement in its deformation and the average displacement of the load is less than the deflection under the center of the load. This is shown in Fig. 21. The problem is to calculate the value of the deflection reduction factor F r In Ref. 7 a geometric approach is presented and the following formula is given: F = 1.0 - 1-.0 (124) r 2 3 * a 0 In order to have amore accurate definition of the value of Fr a new development was done as follows. The deflection at the distance r from the center of the load is given by equation (21): w w (21) r 2 1+ a r where: a = 2 (20.a) ~ 2 D

-68For more convenience in further compuations one will consider: 2.0 a D (125) o and: w Wr 2 (126) r 2 2 + a ~ r The load pressure is assumed to be constant throughout the loaded area and equal to p. Now, one can compute the exterior energy as follows: D/2 w 1 + a2 r 0 Equation (127) can be presented as follows: D/2 W = 2 ~ p I W a r~ d (128) o 0 2 2 r then: D/2 D/2 1 + a2. r2 dr | 2 in (1 + a ) + C (129) L~~~~~a2 r2 in(

-69By boundary conditions C1 = 0. Then, after replacing "a" with its value from (125) one gets: 2 D 1 (1.0 0 ) (130) e o 8 n2 (130) The factor F is then: r e 2 1.0 Fr - n 1.0 + 2(131) =io 2 2 2 H * D w p o w+ a Equation (131) is very different from (124). In order to verify the validity of equations (124) and (131) their results were compared with measured values of Fr as shown in Fig. 22. Equation (124) better fits the data than equation (131). Nevertheless, neither of them is satisfactory and an experimental correction factor is added to get a better agreement between calculated and measured values of F r The reason of this lack of fit may be explained by the non-uniform pressure of the tire and by the fact that the real F factors were measured under a dual wheel area whereas the r calculated values of F were computed under a unique circular area.

-70So, for the further developments one may use one of the following equations: F = 1.0 0.- (132) r 2+ 0.03 2 = ~1 1 0 +1"0 Fr Cn (1.0 + 2) + 0.046 (133) and then: W =F ~P ~w (134) e r s In order to get a more consistent theory the relationship in (131) will be used. 5.3.2.2. Subgrade Energy There are two ways to compute the subgrade energy: - By using equation (15) given by the theory of elasticity (Bezuhov variant) or - By using equation (21) which is an empirical approach (Empirical variant) Both of them will be discussed. Bezuhov Variant The equation which gives the distribution of the load on

-71the subgrade (and which is not given in Bezuhov's developments) is: r 1.0 [ ] (135) 0 where: a D b = 2 a (135.a) One can demonstrate easily that if equation (135) is integrated throughout the loaded area one obtains the actual exterior load, so there exists a real equality between actions and reactions. In order to simplify the further computations equation (135) will be presented as follows: 2 2 r = C1 (b - r) (136) where C P (137) 1 2 2 b *a o Equation (15) can be presented as follows if required computations are made: 1' 2 W=4 2 E b r) (138) Wr 4 2 E cx o

-72or, in a more convenient form: WrC= (2 b r2) (139) where: 2 -P. (140) 2 4 2 E (140) The subgrade energy is: b Wt = 2 T H b * r' w ~ dr (141) 0 or: b W= 2 b * Wr dr (142) r r r 0 and, after replacing rr and wr with their respective expressions from (136) and 139) one gets: b P 2 2 2 2 Wt= 2 * ~* b C1 C2 (b2 - r (2 *b2 - r) dr (143) 0 Equation (143) can be written as follows: Wt = C3 f (2 * b4 - 3. b2. r2 + 4) d

-73where: C3 = 2 b ~ C1 * C2 (145) After integration one obtains: b 3 5 4 2 r r Wt C 2* b4 *r -3.. b2. r + +C (146) 3t C3 5 4 0 By boundary conditions C4 = 0 and after computations: 6 5 W = C 6 b (147) If one replaces C3 by its expression in (145) and all required computations are made -one obtains: 2 1 2 2 3 18 — (148) t 160 E. (148) o o or, in a more compact form: 1 - B 2 Wt= 0.45 E * (149) o o Empirical Variant The equation for 5r is the same as in Bezuhov variant,

-74The equation of the subgrade deflection is the same as in (126). Using the same approach as for the Bezuhov variant, b W= 2 1 H * b W dr (150) t r r 0 b 2 2 Wt = 2 ~ b w C1 (b - r ) 2 dr (151) 1+ a r 0 C2 = 2. 7r * b * w0 C1 (152) b 2 b r t = C2 + ~ dr (153) 1 + a r + a r O After integration one obtains: Wt = Cctan(a r)- a + L arc tan(a * r)+ C3 (154) 0 After replacing C2, a and b by their expressions and after making all required computations one gets: 3 = 2 (1 - 2) 2 p2 Wr-3( 32 E* (155) o o

-75and, in a more compact form: Wt = 0.43 E D (156) 0 0 o o The difference between equations (149) and (156) is about 5% which means that the empirical approach used to define the equation of the deformed surface in (21) is very good. 5.3.2.3. Pavement Energy According to Ref. 7, the energy of a layer can be computed by using the following equation: W 1 le E h W +f w + s 22 2 2'12 (1 - 12) J ax y A (157) 2w 2w 1 w a2 a2 + 2 (1 - a) ( )J d ~ d 2 2 ax y x y where the domain A is the area of the deformed surface. Equation (157) can be simplified by assuming that the loaded area and hence the deformed area are circular. Recognizing that: x = y= r (158) equation (157) becomes:

-763 2 r2 / s 2 12 (1-[2 +2 r2 A ( )2] d2 (159) + 2 (1 - r) r dr But, since the deformed area is a circular one and it is assumed that, at the same distance r from the center of the loaded area, the same deflection w occurs, the last term under the integral is equal to zero (Ref. 7, page 20). Therefore equation (159) can be simplified as follows: Ws 1 E; - 2.J w. 2 + dr (160) s1 2 12 ) + \2 A and finally E2 (; 2) -~ dr (161) 12Although equation (161) seems to be easy to integrate, in Although equation (161) seems to be easy to integrate, in reality it is very difficult. For this reason some simplifications were made as follows. First of all it is known that the curvature of the surface at a given point can be computed by equation (90). This is the reason why one can write: s 312 ( ff) ( Rn - dr *~ dr (162) Ws~E = (1'h

-77where R is the radius of curvature of the deformed surface at the distance r from the center of the load. It it can be assumed that the radius of curvature is constant along the deformed surface, the integration of equation (162) will no longer be difficult. Experiments made by Brown and Pell (Ref. 2, pages 487-504) and which were analyzed in Ref. 7 (page 14) have demonstrated that there is not too great a difference between an "ideal" radius of curvature (constant) and the real radius of curvature, as shown in Fig. 23. For this reason Rr is replaced by the radius of curvature at the top surface of the subgrade, calculated by equation (26). Then one can write: E h ff We = 12 - ~ dr (163) A But: r ~ dr = 0 (164 Jy n (ao D * 2) (164) 4 A because the deformed surface is a circle with a diameter equal to 2.06? D. =3 l (2 ~ D)2 E O (165) s 12 (1 -p ) R 4 or 2 2 W= D..h (166) s 12 R 1

-78If W is the energy accumulated by a layer during deflection, the energy of the pavement would be: m Wsr W=i (167) sr - Si i=l and, after replacing Wsi by equation (166): 2 2 m 3 ~ a ~' D E. h. W = (168) sr 12 ~ R 1 - i i=l 1 But W does not represent the whole energy of the pavement sr because the effect of the friction stresses is to be added. If one compares equations (95) and (166) it will be immediately observed that the energy of a layer is a direct function of the bending moment: w.]. a 2 ~ D2 si o M.5~ -_~-~ 0~ ~(169) M. R 1 The moment due to the friction stresses is given by equation (97). So-,one can write: W' o 2 h (170). 0 _ f f 1 W + + 170) sr R (1 i+l 2 1

-79and: 2 2 2 ( 2-ci 2D Ei h h' i f f i W = (aR. Z+ll-l) (i il) Cr (171) Ir o R e2 (1 2 i+l 2l in 2 Fr P (w + Z wpi, (j=m+l) = 0.43 ~_ P r 01i'' E a D o o i=l (172) 2 3E. * h H a0 D E hi f f hi 16 w1 (1 +8w + (173) 3 R D21- 2 i i+l0 2 i=l z [' + i+')] 11 Fr o + 3 2 2' ~'

-80If one considers that: 2 w 1.5 *D (174) o E a D O o which agrees with the theory of elasticity (Ref. 14, page 260) and one replaces in (173): 2 -,O p w (175) E' D 1.5 0o it follows that: Fr * P(wo + wpiJ(i=m+l))= 0.286 * wO * i=l m 2 3 (176) i=l i hi. (aE h++l After grouping terms one can write: 16 11 Ei hi 2 2 0o i=l 1. (crf Cy f B = -F * P + 0.286 * P + 4 * 1 hi (a. + a. i i=i m i=l 3 ~

-81and finally, 2 A w + B * w + C = 0 (178) o o In order to verify equation (178) a computer program shown in Appendix 6 was written. The results are presented in Fig. 24 and 25. In Fig. 24 a comparison between calculated and measured deflections (regardless of the layer on which they were measured) was done. Several observations are to be made. The correlation between computed and measured deflections is much better using the energy approach. The computed deflections show a definite trend to be higher than those measured (- 34%) which is explained by the following reasons: - As was demonstrated* several years ago, for the same pavement structures the deflections measured under real wheels are about 8% higher than those measured under flexible plates. In the plot only measures under plates are given. Thus, the real difference between calculated and measured deflections is about 26%. It should be noted here that in the definition of F a r very flexible load was taken into account, so the difference mentioned above is justified. * Between 1968-1971 a research group formed by: Assist. Prof. Horia Zarojanu; B. Cososchi and N. Vlad from the Polytechnic Institute of Iassy Romania, lead by Emm. Prof. D. Atanasiu, worked in this area.

-82- Two important simplifications were made for these computations. First, load is distributed upon a limited surface at the top level of the subgrade (which leads to a concentration of the load in the central part of the loaded area and hence a higher computed deflection). Second, the radius of curvature of the layers was considered as being a little greater than the real one as one can see in Fig. 23 (and hence the energy of the layers was underestimated). It is hoped that further studies will allow the discarding of these simplifications, after these studies have been done, better correlations can be expected. In Fig. 25 a comparison between computed and measured tensile stresses at the bottom of the layers is presented. It can be seen immediately that actual tensile stresses show a definite tendency to be greater than the computed ones. The average ratio between calculated and measured stresses is 0.79. It should be pointed out that the average values mentioned above (either for deflections or for tensile stresses) were determined statistically and they differ from the means (arithmetic means) which were computed within the computer-program in Appendix 6. The latter are greater than the statistical values used. One can assume that the magnitude of tensile stresses is proportional to the energy accumulated in the layers. If this energy is underestimated the tensile stresses are also underestimated.

-83If there exists some proportionality, the product of deflection ratio (1.26) and the stress ratio (0.79) should be close to 1. In fact, 1.26 * 0.79 = 1.0. This means that further studies have to start with a better approach to evaluate the energy of the layers. At this time a coefficient of lack of fit will be introduced. This coefficient can be calculated as follows: If the tensile stresses computed at the bottom of the layers represent, as an average, 0.79 of the actual tensile stresses, the average computed radius of curvature is about 1/0.79 greater than the actual one. As a consequence the formula in (26) becomes: 2 2 ~ D R = 0.79 8 (179) 8 w Then, by replacing R from equation (172) by equation (179) one obtains: m r o 1 o Fr P ( w0 + E wP.J (j=m+l)) = 0.286 w0 ~ P i=l 16.0 Wo 2 L * h3 1.84 2 2 1 ~0 i=l +6 5 * n 2hi (of6.+ o. (180) i=l

-84Then one gets: m 3 A = 16.0 11 E i 1.84 ~ t ~ D i=l m B = - P (F - 0.286) + 6.5 7 h ( + +l) (181) m C = - F P w (j=m+l) r i ( i=l But now a new problem occurs. As can be seen in Fig. 24 and 25 the ratios between calculated and measured values have a stochastic distribution which is natural. This stochastic distribution is caused by the fact that all values chosen for the parameters characterizing the quality of the materials in the layers are not constants, but they have normal distributions as shown in Ref. 10. For this reason a safety coefficient is to be introduced. The computation of the safety factor is a stochastic problem and will be amply discussed in the next report. However some considerations will be made in this report. It can be required that the actual deflections measured upon a new road be less than or equal to a given value, for instance

-85the permitted or the computed deflection. Although it is possible to have this requirement for all the surface of the road, it would be too expensive to build such strong pavement structures. For this reason, in Romania for instance, only 85% of the surface of the road has to have deflections less than or equal to the permitted deflection. Using this rule, which was provisionally accepted in this work, the value of the safety factor is determined in Fig. 24. This value is 1.5 which means that the permitted deflection is to be reduced in computations by 1.5. It is interesting to observe that the new value obtained for the safety coefficient is less than that mentioned before. But, as shall be seen in the next report, the value 1.5 for the safety coefficient "covers" only the variations due to the nonuniformity of the materials characteristics. The influence of the variation of the thicknesses of the layers will be also taken into account. And it is to be expected that the latter will lead to an increase of the value of the safety coefficient. The calculation of the permitted deflection may be done by different means. One of them is the formula offered in Ref. 7 page 26: Wa = Fd (0.280 - 0.031 log Nt) (182) where Fd is a correction factor defined in the list of notations Chapter 4. The result of the compuation by using formula in (182) is given in cm, if Fd = 1 or in inches, if Fd = 0.3937.

-86Nt is, as one knows, the total estimated traffic for the service life of the pavement and it can be' easily computed by means of the present annual traffic (N), the annual traffic increase rate (NR), and the service life of the pavement (S9): S Nt N (1 + NR) (183) i=l But, the total estimated traffic has to be affected by a coefficient taking into account the seasonal influence (Ks). This coefficient may be calculated, and in Ref. 7 page 33 a way to do this computation is presented but this computation is a little too difficult to do through a structural design computation. For this reason a more simple and convenient formula can be used: K 1 + 0.30 F (184) s 1000 where F represents days x ~C below 0~ C during one year. If F calculated in degrees (below +320 F), equation (184) becomes: F K 0.968 + 0.166 1000 (185) s 1000 and finally, S Nt = Ks ~ N lNR (1 + NR86) i=l

-87It should be said that the formula in (184) fits both with data in Ref. 7 and data reported several years ago by Kucera (Czechoslovakia). Another problem is to derive a formula for the equivalence wheel load factor computation. Many methods are given in the literature, the most of them resulting from correlations based on experience. The writer would like to present here the results reported in Ref. 7, because these results were based on a theoretical development which fitted rather well the AASHO road test findings. The formula presented in Ref. 7, page 30 is the following" log N = (9 - log NB) (n - 1) + n ~ log N. (187) e 1 where: NB - the basic traffic evaluated for the service life of the pavement. It can be evaluated at first by the equivalent wheel load traffic expected over the highway and then it can be corrected by an iterative procedure. Finally NB is to be equal to Nt; N. - number of load repetitions over a lane estimated for the vehicle type i; Ne - number of load repetitions over a lane calculated for the equivalent vehicle; - factor, dimensionless, that is computed by using the formula in (188) or in (188.a):

-88I e (188) 4i 41P. * P...e (188.a) Ve Pi where all teyms are known. The equivalence load factor will be: N. K = (189) w N e 6 For NB = 10 (which corresponds to the traffic applied in AASHO Road Test) and for three different types of equivalent vehicles, comparisons between the results given by the above procedure and AASHO test results are given in Fig. 26, 27 and 28. If the permitted deflection under the standard wheel is known (wa reduced by the safety coefficient) one can develop the design formula. In order to make the appropriate development to obtain the design formula one has to observe that: m wa wO + E wpj7(j=m+l) (190) i=l

-8 9and then to go back to equation (173) which is written in another form: 2 1 - o p2 F P ~ w 0.43 r a E [ w D o O m Fr 16 P wa8i=6 l [aEi hi (191) 2 2 1 -. O D i=l 1 f f I wa - wPi'(j=m+l ( i i+l i=l i=l Now introduce the lack of fit factor (Fi = 0.79) and equation (191) becomes: m Fr P w = 0.286 w - wp (j=m+l) P m m + * a -E i hi.i i+l +i= 1 i =l 1 - i

-90It can immediately be observed that there is a common term on the right side of equation (192). This will allow us to simplify the form of the design formula as follows: m F a P w w -~ IP(j=m+l) i0.286 P r a a m E. h.3 + 8.7 ~ 1 T 1 1 + 6.5 7r(193) i=l m, hi 1(of ++fi+l ] i=l The ultimate form of the design formula will be: m Fr *Pwa [wa w J(j=m+l) 286 i=l + 6.5 * II [134 1 + (194) i=l 1 h. (f + of i+l)

-91A discussion of this equation follows: - There is no direct way to determine the thicknesses of the layers from the formula even if all exterior data (exterior load and its characteristics; subgrade and its characteristics; Young's moduli of the layers and their Poisson's ratios) are given. - The stress and strains state within the pavement structure varies with the absolute and relative thickness of the layers, their position in the pavement structure, and their relative rigidity. - The design process must be an iterative one. A set of thicknesses is to be defined by the designer and then the iterative computation can be continued until limit differences between "initial" and "final" thicknesses are obtained. These limit differences should differ from layer to layer according to some practical rules. In Romania, for instance, the thickness of the surface course varies by 1/2 cm and the thickness of the other layers by 1 cm. - Any set of thicknesses must observe the constraints concerning the minimum and maximum thicknesses as it was discussed in the previous section of this report. For all these reasons no "manual" procedure can be used for the structural design of asphalt pavements and only computers are able to do this job with a real efficiency. If only computers can be used it is natural to try to optimize the design, namely the minimization of the construction costs, all technical requirements or constraints being observed. This will be the subject of the next report.

-926. Conclusions Although partial conclusions were presented throughout the previous parts, several main conclusions are to be summarized here. a) There exists empirical formulas, developed through a consistent theory, which can predict with a high accuracy (under the center of the load): - the vertical stresses at the top of the layers (equation 33), - the radial stresses at the bottom of the layers (equation 42), and - the shear stresses at the bottom of the layers (equation 45). b) The integration of these equations has demonstrated that the stress function which could solve the equations of the theory of elasticity (equation 1), if only one,could be a transcendental function as shown in equation (72). This means that the functions which are usually introduced in equation (1), for instance Airy, Bessel or polynomials types, are suspect and could lead to wrong conclusions. c) The basic equation proposed by Leger, presented in equation (18), describes rather well the top surface of the pavement when distorted. This equation may also be used to predict the strains (deflections) at the top surface of each layer of the pavement. d) The study of the distorted state of the highway body, under load, taking into account only the energy necessary for

-93the change of volume is not sufficient for describing the real state of the stresses and strains. This is because the layers slip with respect to each other and the energy for the change of form must be taken into account. e) The layers can actually be considered as slabs of a low rigidity even if they are built with noncohesive materials. f) In order to eliminate any factor which could lead to the distress of the pavement, the following design criteria are to be taken into account: - The energy gained by the highway body during deflection must be equal to the exterior energy, transferred by the load. - Vertical stress at the top of each layer must be lower than the allowable stress. - Tensile stress at the bottom of each layer, except subgrade, must be lower than the allowable stress. - Shear stress in the center plane of each layer, except subgrade, must be lower than the allowable stress. - The maximum vertical stresses should be located in the surface course. - The location of maximum shear stress should be in the base-course. - Vertical shear stresses along the circumference of the loaded area at the top of each layer, except subgrade, must be lower than the allowable stress. - The depth of frost penetration must be less than the thickness of the pavement.

-94g) The design process is an iterative calculation under technical constraints and cannot be done except with the use of a computer with a large memory (3 million bits or more) and with a high speed of execution. This is because of the large quantity of information that is to be processed and the length of the computer-programs. h) The method described in the previous part is subject to important improvements that will result from replacing a number of simplifying assumptions with more realistic ones. Nevertheless this method predicts correctly (at least from a stochastic point of view) the vertical stresses at the top of the layers, deflections at the top of the layers and tensile stresses at the bottom of the layers. i) Although the deflections permitted at the top surface of the pavement appear to be the principal design criterium, it should be pointed out that this is as a consequence of using the radius of curvature as a parameter (see equations 171 to 176).

-95Acknowledgments The material in this paper was developed in a research project sponsored by the National Academy of Sciences of the U.S.A., The Ministry of Education of Romania and The University of Michigan, Department of Civil Engineering. The opinions, findings and conclusions expressed herein are those of the author and not necessarily those of the sponsors. The author would like to express his appreciation to Professor Robert Goetz for his marvelous lectures, full of common sense, which helped the author with his not forgetting the true reality of the structural design of asphalt pavements. The author would also like to express his particular acknowledgments to his adviser, Professor Egons Tons, who helped him in many different ways with his improving the knowledge in the area.

References 1. International Conference on the Structural Design of Asphalt Pavements, Proceedings, Ed., Braun-Brumfield Inc., Ann Arbor, Michigan, 1963. 2. Second International Conference on the Structural Design of Asphalt Pavements, Proceedings, Ed., Braun-Brumfield Inc., Ann Arbor, Michigan, 1968. 3. Third International Conference on the Structural Design of Asphalt Pavements,Vol. I Proceedings, Ed., CushingMalloy Inc., Ann Arbor, Michigan, 1972. 4. Yoder, E. J., Principles of Pavement Design, Ed., John Wiley & Sons, Inc., New York, 1959. 5. Carnahan, B., Wilkes, O. J., Digital Computing, FORTRAN IV, WATFIV, and MTS, Ann Arbor, Michigan, 1972. 6. McLeod, N.W., The Asphalt Institute's Layer Equivalency Program, RS-15, March, 1967. 7. Rafiroiu, M., Contributii la studiul problemei dimensionarii sistemelor rutiere nerigide (Resum6 of the doctoral thesis), Iasi, 1971. 8. Rafiroiu, M., Une nouvelle methode pour le dimensionnement des chaussves souples, Revue Generale des Routes et des Afrodromec Paris, Nr. 431, Apr., 1968, page 49-74. 9. Rafiroiu, M., Quelques considerations sur les chaussees bicouches, Revue Generale des Routes et des Agrodromes, Paris, Nr. 452, March, 1970, page 85-92. 10. Rafiroiu, M., Essai d*etude statisti"ue des constantes elastiques des materiaux routiers, Revue Generale des Routes et des Aerodromes, Paris, Nr. 467, July-Aug. 1971, page 51-56. 11. Rafiroiu M., Sur le frottement a linterface des couches des chausszes, Strasse und Verkehr, Zurich, Nr. 9, Sept., 1971, page 438-440.

12. Rafiroiu M., Le dimensionnement aux bords des chaussees souples, Strasse und Verkehr, Zurich, Nr. 2, Feb. 1971, page 61-65. 13. Dimensionierung der Strassenoberbaues, Strasse und Verkehr, Zurich, Nr. 8, Aug., 1972, page 373-415. 14. Bezuhov, N. I., Teoria elasticitatii si plasticitatii, Transl. from Russian, Ed. Tehnica, Bucuresti, 1957. 15. Hilton, P. J. Partial Derivatives, Dover Publications Inc., New York, 1965. 16. Smoleanski, M. L., Tabele de integrale nedefinite, Ed. Technica, Bucuresti, 1972. 17. Ahlvin, R. G., et al., The Principle of Superposition in Pavement Analysis, HRB, 52nd Annual Meeting, Report No. 5 in Session 32, Washington, Jan. 22 - 26, 1973. 18. Barksdale, R. D. and Hicks, R. G., Material Characterization and Layered Theory for Use in Fatigue Analyses, HRB, 52nd Annual Meeting, Report No. 2 in Session 22, Washington, Jan. 22-26, 1973. 19. Timoshenko, S. and Goodier, Y. N., Theory of Elasticity, Ed. McGraw-Hill Book Company, Inc., New York, 1951. 20. Rafiroiu, M., O INCERCARE DE REZOLVARE A PROBLEMEI DIMENSIONARII ECONOMICE A SISTEMELOR RUTIERE NERIGIDE, Transporturi, Bucuresti, 2 (19), Nr. 5, Mai 1972, page 259-264.

List of Figures Figure 1. Element of stress. 2. Element of strain. 3. Stresses acting on an element. 4. Uniform loaded plate. 5. Semispherical loaded plate. 6. Model for flexible pavements analysis. 7. Stresses under the center of the loaded area for an m-layered system. 8. Stresses in the i-th layer. 9. Vertical stresses under a uniform loaded flexible circular plate. 10. Average vertical stresses under a uniformly loaded circular plate. 11. Radial stresses under a uniformly loaded circular plate. 12. Shear stresses under a uniformly loaded circular plate. 13. Comparison between measured and calculated vertical compressive stresses under the center of the loaded area. 14. Comparison between measured and computed tensile stresses at the bottom of the layers under the center of the loader area. 15. Comparison between measured and computed deflections at the surface of the pavement under the center of the loaded area. 16. Comparison between deflections calculated without taking into account the energy necessary to change the form of the layers (w, computed) and those which take it into account (w, LMD).

17. R * w values for various structures of asphalt pavements. 18. Bending moment due to the radial stresses. 19. Bending moment due to the friction stresses. 20. Shear stresses along the contour of the loaded area, at different levels. 21. Average deflection under load. 22. Computed versus measured values for the factor Fr 23. Theoretical versus actual deflection profile. 24. Computed versus measured deflections. 25. Computed versus measured tensile stresses at the bottom of the layers. 26. Theoretical versus AASHO Road Test results for the wheel equivalency factor, for EWL = 8.1 metric tons (18000 pounds) standard axle. 27. Theoretical versus AASHO Road Test results for 10.4 metric tons (22500 pounds) standard axle. 28. Theoretical versus AASHO Road Test results for 13.6 metric tons (30000 pounds) standard axle.

y dz _L_ Figure 1. Element of stress. Y Ax y dy, ---— 7 dz I // Fg.l/ Figure 2. Element of strain.

I.. P I-MI "'-,, i N t Figure 3. Stresses acting on an element. a P Inflection Point Figure 4. Uniform loaded plate.

~~Eo |*,~~~~(m a D/2 -m 6Jm >'/ >...... Eigr /6. Moeolxbepvmn/sn2 M Ei rm p/a Wm-I, ai *D/2 m ~M-l hi Ei ri z p/a.2 ~m-I Ri Wm-1 ai-D/2 —~` "..iD/2 El "Semisphericol".,M ht Distribution 0 ai'D/2... E L ao' D/2 M'w. /lom m-I Ro >Ri >-".>Ri >...>Rma r 6 MoDe/2. Figure 6. Model for flexible pavements analysis.

D/2 D/V f ri rn I rn fi _, rfs r rs or m m rs fs m-I m-I m-i m-l m-1 m-I hm-c I m-2 I fi ri m fi IMm- Im- m- mn cm - ~- -, - -. o_ -m. i_+._. _I o'_ h, M 1.' — - - |rfs rs'-m 0. " rs O'fs m-2 m-2 ~''2 v m-2 m-2 m-2 fi ri m ri fi -i I i I i+ I I i+l i +I Os ___ _i r rs fs ri m 0-n ri fi'1i Of'-~ C ci cr 0-fs rs m Om r s fs're O. O, 0-rf m- e'fs ~m' o m rs ~fs hi (fi g ri m m nr fo yT fs rs m rs fst Mo Oa 0 Figur 7. Stessudrtecnef th lode arafra h~~nlaee sytm

Axis of Load |civ Ei+I |!f er m m r f hi f r m rn r r f I'rI ii..+ i — i i i i i |I Ei Figure 8. Stresses in the i-th layer.

Oz/p 0 0.2 0.4 0.6 0.8 1.0 0.2- 0 0.4/ 0.60.8- o Bousinesq e0 z/D Range for I 0 -.~o,,~ Surface / Base Interface z/D + 1.6 + z/D Range for /.8 0 Eq uation Base/Subbase x + 2.0) 2.2- O+ 2.42.4 1 Author Ref. Page 2.6 - z/D ~ Vesic 1 284;, Range for o Brown & Pell 2 493 2.8 I Subbase/ r Brown& Pell 2 493 Subgrade + Waterhouse 2 303 32.80 Interface |o0 Nijboer 2 698-701' 3.2 ~ Morgan & Holden 2 713 a Hicks & Monismith 3 411-20 3.4t6 x Miura 3 478-86 3 Ledbetter, et al. 3 877-900 3.6 Figure 9. Vertical stresses under a uniform loaded flexible circular plate.

o-z/p 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 0.4 /I I 0 0.2 0.4 0.6 0.8 ~~~1.0~ z/D 1.2 1.4 - 1.6 0 1.8 2.0 2.2 2.4 Figure 10. Average vertical stresses under a uniformly loaded circular plate.

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 = L0.35.. 0.2 Boussinesq o 0 =0.50.6//. Equation (42) 0.8 1.0 5' z/D 1.2. / 1.4 1.6 Author Ref. Page 1.8 0 Morgan 8 Holder 2 712 2.0 ha x Brown 8 Pell 2 495 2.2 2.4' Figure 11. Radial stresses under a uniformly loaded circular plate.

0 0.1 0.2 0.3 0.4 0.5 0.6. \x =0.3 z/D 0.2 t max=0.3 z/D max =0.35 z/D 0.4. 0.4 Empirical Curves 0.6 /.,=.50 maxBOSsl NESQ =0.70 z/D 0.8 i =0.35 /i:. BouBou ssinesq, =0.35 1.8 I: Fu I2 Sa sBoussinesq,u = 0.50 1p2e 1.2 I I I 1.6 -I 1.6 Figure 12. Shear stresses under a uniformly loaded circular plate.

1.6 1.4 0. ~** 0 7l 1.21 V (:3 * 0* 0.8 0 0 r ~~0.6. ~0 00 O0 0.4 % 0 0 0 0.2 t* / 0.2 0.4 0.6 0.8 1.0 1.2 1.4 fv, COMPUTED, KGF/CM2 Figure 13. Comparison between measured and calculated vertical compressive stresses under the center of the loaded area.

om, MEASURED, KGF/CM2 0 2 4 6 8 10 12 K 1 I I I I I I I I 2I.'~ 4- \' - 12* Ideal b 6 t. Correlation U..0n 8_ w 2E Actual b.'Correlation(4/15) 16 \ 18 Figure 14. Comparison between measured and computed tensile stresses at the bottom of the layers under the center of the loaded area.

1.4 1.2 c 1.0 - E < 0.6 o 0.4. *. 0 0.2. ~ * 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 w, COMPUTED, mm Figure 15. Comparison between measured and computed deflections at the surface of the pavement under the center of the loaded area.

0.24 / S / 4 0.22 / 0.20 0.18 * / * ~* 0 lf E 0. * / ~ / 0.0814 0 h / I 0.061 0 0. / t 0.04 /t 0.02 t 0 1 1 1 1 1 I I I I I... I I I I I I 0.02 0.04 0.06 0.08 0.10 0.12 0.14 0.16 0.18 w, LMD, mm Figure 16. Comparison between deflections calculated without taking into account the energy necessary to change the form of the layers (w, computed) and those which take it into account (w, LMD).

20 Range of R, w Values 18 in Appendix 5 1 16 ~4 |_/// // r14//// / 12 / A12 /' /5 E1/ / / / E wt8 lo-,// 8/ X I/ / / o 82 4 0//E/E.5 2 3 4 5 NUMBER OF LAYERS Figure 17. R- w values for various structures of asphalt pavements.

J I Ior *f I/3h hi =..... __ h1~~~~~ ~1/2 h -6 I/2h Figure 18. Bending moment due to the radial stresses. h ef Iff igure 1. einmfic oi Figure 19. Bending moment due to the friction stress.

i-' oad #'I Figure 20. Shear stresses along the contour of the loaded area, at different levels. ao.D/2-j aD/2 Tr Wo Wo Figure 21. Average deflection under load.

*Measures rra2 Tr=a P2 -nt nI- — 0.046 1.00'fr 1 Tr r T r +0.03' T r 3o 0.90 s'Tr n ) (1+ 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 3.0 3.2 ao0 Figure 22. Computed versus measured values for the factorTr.

20 /,/ o40 zar Experimental 60j _ Curves: a _ Z. —— a =15 cm. -. 0...a 22 cm 80- Theoretical Curve' F -— a =30 cm 100 0 0.5 1.0 1.5 2.0 r/ao Figure 23. Theoretical versus actual deflection profile.

w, MEASURED, mm 0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 I'..1 I I I 1 1111 1 I I x x~ xx 0.2 -~\ >~v o = w/!O o* x = w'10 0.4' 0,* 0 0.6 X \ E C,0.8 0 - 1 ~' O~0, C/ 0 1.6 18gr 4 Figure 24. Computed versus measured deflections.

12 / 10 U.!P -'" D 8 ~~6c d. E *. 0 b~~~~ 0 2 6 8 IC 12 g'm,, COMPUTED, KGF/CM 2 Figure 25. Computed versus measured tensile stresses at the bottom of the layers.

T K' 13.6 30.0 AASHa X I 10.2 22.5...... 90 19. 19.8 J 8.2 18.0 7.0 15.4 E =8.2 T = 18.0 K............. | 5.4 12.0 12 10 8 6 4 2 0 EQUIVALENCY FACTOR Figure 26. Theoretical versus AASHO road test results for the wheel equivalency factor, for EWL 8. 1 metric tons (18, 000 pounds) standard axle.

T K -\AA - 13.6 30.0 AASHO LMD \ 10.2 22.5 9.0 19.8 w -J 8.2 18.0 < 7.0 15.4 E, = 10.2 T =22.5K 5.4 42.0 5 4 3 2 1 0 EQUIVALENCY FACTOR Figure 27. Theoretical versus AASHO road test results for 10.4 metric tons (22, 000 pounds) standard axle.

- __ ___ ____ T K 13.6 30.0 AASHO i 8.2 18.0 < LMD`\\ 9.0 19.8 w 7.0 15.4 E, =13.6T -30.0K ]5.4 12.0 1.0 0.8 0.6 0.4 0.2 0 EQUIVALENCY FACTOR Figure 28. Theoretical versus AASHO road test results for 13.6 metric tons (30, 000 pounds) standard axle.

List of Tables Table 1. External conditions to define Leger's equation. 2. R * w values. 3. Stresses at the interfaces of the layers. 4. Design parameters for some materials.

Table 1 External Conditions to Define Leger' s Equation Condition m r w 1 0 w 0 * D w 2 0 0 2 2 3 co 0

Table 2 R * w0 Values Type of the Pavement After E1/EQ Two- Three- Four- Average Grant and WalIker layered layered layered 0.5 1,540 4,330 7,500 4,450 3,670 1 3,250 5,370 9,100 5,900 5,200 2 4,560 6,280 10,100 7,000 7,350 4 6,100 7,600 11,100 8,300 10,400 8 7,300 9,600 11,400 9,500 14,700 Average'4,550 -"6,650 -"9,650 - 7, 000 -8,450

Table 3 Stresses at the Interfaces of the Layers' = 0.35; f. = 0.5 Z/D av z r Tzr r zr r fi 0.0 1.000 0.850 0.075 0.925 0.500 1.0 0.271 0.063 0.106 0.169 0.136 1.7 0.133 0.027 0.055 0.082 0.066 2.6 0.079 0.016 0.035 0.051 0.040

TABLE 4 DESIGN PARAMETERS FOR SOME MATERIALS Temperature Young's Modulus fA Material U T Cfr Observations Coefficient ~C OF kgf/cm2 p.s.i. of Variation Degrees kgf2 p.s.i. % cm Loess - - 300 4200 26.5 0.4 20 - - 1.5 Loess = 40% Clay - - 240 3300 40.5 0.4 20 - - 0.8 porosity. Silt - - 320 4400 28.5. 0.4 20 - - 0.65 All materials Silty Clay - - 230 3200 38.0 0.4 20 - - 0.7 at the optimum Sandy Clay - - 410 5700 48.5 0.4 21 - - 0.8 moisture conClayey Sand - - 340 4700 42.5 0.35 20 - - 0.9 tent and 100% Silty Sand - - 360 5000 41.0 0.35 20 - - 0.8 standard AASHO Fine Sand - - 740 10200 42.1 0.30 22 - - 1.0 compaction. Medium Sand - - 830 11400 34.4 0.30 24 - - 1.0 Coarse Sand - - 930 12800 26.9 0.30 25 - - 1.0 Clayey Gravel - - 1450 20000 28.8 0.30 26 - - 0.8 Sandy Gravel - - 1340 18500 19.8 0.30 32 - - 1.0 Common Gravel - - 2200 30500 12.1 0.30 30 - - 1.2 Well Graded Gravel - - 3650 50000 34.8 0.30 52 - - 1.45 Soft and Medium Rocks - - 2040 28200 33.1 0.25 28 - - 0.75 Stone Tough Rocks - - 4950 68500 20.0 0.25 36 - - 0.75 Medium Rocks - - 4000 55300 - 0.25 40 4.65 64.0 0.75 Macadam Tough Rocks - - 4680 64800 11.1 0.25 35 - 0.75 0 32 66000 915000 33.3 0.35 84 - - 0.50 All bituminous 15 59 32000 445000 2.6 0.45 61 7.40 102.0 0.50 materials 30 86 10900 150000 2.6 0.49 45 - - 0.50 according to 40 104 7450 103000 25.5 0.50 40 - 0.50 the Romanian Standards 0 32 15000 207000 - 0.35 49 - 0.50 which are alBituminous Concrete 20 68 9250 128000 16.4 0.48 42 - 0.50 most the sail,,e 40 104 5000 69000 - 0.50 36 - 0.50 as AASHO Standards. 0 32 19000 263000 - 0.35 53 - - 0.80 Sand Asphalt 20 68 10600 146000 4.3 0.48 45 - - 0.80 40 104 8000 110000 6.7 0.50 41 - - 0.80 0 32 15000 207000 - 0.35 49 - - 0.60 Bituminous Binder 20 68 5200 72000 - 0.48 36 - - 0.60 40 104 4150 57000 - 0.50 34 - - 0.60 Mastic Bituminous Concrete 10 50 20800 288000 28.5 0.42 55 - - 0.85 Bituminous Penetrated Macadam 20 68 7250 100000 31.0 0.48 40 - - 0.75 Bitumen Stabilized Gravel 20 68 8920 124000 35.0 0.48 44 - - 0.85 4.5 - 5.5% Bitumen Content Bitumen Stabilized Sand 20 68 1200 16600 42.3 0.48 25 - 0.65 - Cement Stabilized Gravel - - 12500 172000 - 0.16 47 - - 0.65 4 - 6% Portland Cement Cement Stabilized Sand - - 5000 69000 - 0.20 36 - - 0.55 Content Rectangular(1 Stone Blocks - - 4780 66000 30.0 0.25 - - - 0.55 pavements with Irregular Stones(2 - - 3350 46500 31.8 0.25 - - 0.55 Cobblestones(3 - - 3200 44300 33.2 0.25 - - - 0.55 - 1)Cubes 15 x 15 x 15 cm3 (6 x 6 x 6 in.3); 2)Stones obtained by crushing large blocks in quarries. Irregular blocks with average dimensions between 12 and 20 cm (5 to 8 in.); 3 Cobblestones size range 6 to 10 in.

List of Appendices Appendix 1. Computer program for the values of coefficients a, and bl (in equation 31) determination. 2. Computer program for the computation of the radial stresses under a uniform loaded circular plate. 3. Computer program for the computation of the shear stresses under a uniform loaded circular plate (multilayered systems). 4. Computer program for the computation of the shear stresses in homogeneous mass under a uniform loaded circular plate. 5. Computer program for equations in (75) through (88) checking. 6. Computer program for equation in (178) checking.

APPENDIX 1 M ICHkGANi TrPil JAL SYSTEF1 C34TA G{t 4i3ib) SI VEP 03-15-73 15:14.57 PAGE POOl OuOl WR I TE( 6, 10) 0002 1 0 Ft< M 4T ( m 1' t 0DEMARK' ) 0003 B = 0. 94 0004 7 A = 0.84 0005 5 EXPR = 0.0 000o 3 C = ( A**2 *( tXPR**2) 0007 T = 1.U-B*C/(0.25-sC) oood WRITE(6,1) A,B,EXPRt, 0009 1 FORM!J4AT','A=',F5*2,p5X, B=t',F5.2, 15Xt'EXPR=' F5.3,5X,' r=', F5.3 0010 E.XPR = EXPR+0.2 0011 IF(EXP,.GT.3.61GO Tr 2 0012 GO T I 3 0013 2 A = A+0.01 0014 I F(A.GT.0.-94 ) GO TO 4 0015 GO TO 5 0016b 4 B = 8+0..01 0017 -Ft{'B..GT'r-*.l00o GO TO 6 0018.GO TO' 7 0019 6 WRITE(16tli) 0020 1'1 FORMAT('I',' JUSSINESQ' 0021 D = 96 0022 17 A = 0.88 0023 15 EXPR = 0.0 0024 13 8 = (DEFXPR)**2 Q0025 C = ( A* 3 *( t XPR **3) 0026 c= (-A**2 )g( EXPRZ*2 0027 T = 1.0-B.*C/( 0.5+C )**3/2 JO26d RITE(6,21J) A,8,EXPRT 0029 21 FORYMT(''t A=,F5.Z 5X, t d=,F5.2, 15X,' i (P= F 5Xt r='tF.3) J0 3 )'XP L- +3. 0031 IF(.XPR.GT. 3. 0) GO TO 12 0032 GO T'l 1003 3 12 A = O. O 10 00 3' IF(A.GT.O. 92} GO TO 14 00)OV G LO Ti 15 0036 14 D = O+J.01 0031 IF().GT. 0.92).GO TO 16 003:3 G T I 17 0.:334' 16 STOP 00'+3 IEND

APPENDIX 2 r r.. r)A1 ^-5-lS73;)8:18.51 PAGE POOl ~ 000%,~,t =. "f)008. i! to'TP: - - &Lr:. TIa n ":!3 -,'' s' -- =".-!~ ~ f) -+T i~ )"T ~. _,,.,'1)8 - "! rtf- _ V.* " 3,3) 1 ) ^r _ T r 1 t. /! / (. ) 5 +. L ( * / ) i tR I Tr- {, S1 r r 3,! t t r V) - 3 S 1 P? T I k i).:) l 7 2 ")' r v = d "? I' +t o |57 ) t 1r [ T j ^- (.~j. —~'X Tn -,:),'" Tn 4.5nrlc~ G A T9 5 00"0' 4 Sr4o'79? 1. Ct.19: "'~nT fr" -cC eTM Tn =. 7"r " i': " " t L7''", ITNt A ".._ L.* T' r.1 T ~ ~ 7: T T r. T. "i T "J"" "'t''` r':'".- 1 Th;..... ~ T.,~- r-....:', — ~-,,. ~T F

APPENDIX 3 N1TC-jIf;\: TEK NNAL SYSTEM FURTR41 G(4i3.t) TANGiEN 02-15-73 08: 19.0b PAGE POl 000]. COF14 = 0.88 0032 C0F23 = 0.96 0003 AMIU = 0.35 0004 5 PARI = 1.0-2.0*AMIUt 0005 PAR2 = 2.0*( O+AMIU) 0006 EXPR = 0.3 0007 3 ALFA = COFLA*EXPR 0008 ALFPAT = ALFA**2 0009 BETA = (CJF2B**(1.0/2.J))*ALFA/(i.25+MLFPAT)**(1.0/2.0) 0010 GAMA = BETA**2 0011 TAU = 11.0/4.0)*(PARI+PAR2*BETA-3.0*GAMA) 0012 WRITE(6,1) AMI UEXPRTAU 0013 1 FORMAT(' s,'MItJ=',F7.45X,'H/D='FT.4,5XsTAU=,F7.4) 0014 EXPR = EXPR*0.05 0015 IF(EXPR.GT.2e5) GO TO 2 0016 GO TO 3 0017 2 AMIU = AMIU+0.15 001$ IF(AMIU.GT.O.50) GO TO 4 0019 GO TJ 5 0023 4 STOP 0021 END *OPTI.JNS IN EFFECT* I;r),EBCDIC, SSOORCENOLI ST NOOECKLOADNOMAP *OPTIONS IN EFFECT* NAME = TANGEN, LINECNT = 57 *STATISTICS* SOURCE STATEMENTS = 21,PROGRAM SIZE = 782 *STATISTICS* NO DIAG;4OSTICS GENERATED NO ERRORS IN TANGEN NO STATEMENTS FLAGGED IN THE AiOVE (-0MPILtTlIjNS.

APPENDIX 4 AICHIGAl' T@t AINAL SYSTEM FORTA,,I S(4i'j+i3) tJSTAN 02-15-73 08:18.35 PAGE POOl uO)L AMIUJ = 0.35.C'X:2 6 WRITJ(6,1)`4.'AIJ 0003 I LURM.AT(' t 1 M =', l F10. ) 00 04 A = ( 1.-2.JAM IU) /2. J 000 D = 1i.0+A'11 0006 X = 0.0 0007 4 5 = SQRT T(1.J+X2-Z) 000:~ C = X/B 0009 TAU = 0.5*(+DO*C-i.5*C**3) 0010 WRITE(6t2) X,TAU 0011 2 FORMAT('',' ECH=',FS5.2,5Xl,'ITAU=',FS.3 001 t X = X+0.05 0013 IF(X.GT.3.0) GO TO 3 0014 GO TO 4 0015 3 AMIU = AMIU+,.15 0015 IF(AlIUu.GT.J.5J) GO TOI 0017 GO TO 6 0013 5 STOP 0019 END *OP TIJNS IN EFFECOT* IO,Ei CDIC,SOURCE,NOLIST r NDECK(LUADNOMAP *OPTION4S IN EFFECT* NAME = BUSTAN, LI. LENT = 57 *STATISTICS* SOURCE STTEMENTS = 19,PROGRAM SIZE = 662 *STATISTICS * NO DIAGNOSTICS GENERATED NO ERRJRS IN BUSTAN

APPENDIX 5 MI!CHIGANJ TEI4MNAL SYSTEM F3)TRA4 7(41336) VERTER 03-36-73 18:08.56 PAGE P001 C SIGVFP(I) SITGFR(IT)SIGRA4tl) ARE GIVEN 14 KGF/:m2 A4J THEY ARE MEASURE) AT THE TIP (7 )F THE LAYER. C nEFt(IJ) AN) 0rEPAR(IJ) ARE GIVEN C IN t-M AND TfEY ARE MEASURED AT T4E T3P C OF TIF LAYFR. C C R(I) IS GIVEN IN CM AND IT IS MEASURED C AT THE TIP )F TIE LAYER. C -~ ST7M7M(I) IS GIWEN TN KGFICM2 A4Nt IT IS C MEASURE) WITI SIGN - AT THE T3P OF C THE LAYEP AND WITH SIGN + Ar THE BOTTOM C F THF LAYFR. C C SIGSUP(I) IS THE HORIZONTAL STIESS AT THE'~ T3P 1;: THc LAYER, SIGINF(I) C, IS THF VpnIZINTAL STRESS AT THF 30TTIM!F THE LAYER. C BITH )= THF4 ARF MEASURE) IN <(P/CM2. C C EXPR(I),AI),ALF PAT(I),ZECH(I),3( I) AIE CAL':ULATED AT C T-IE TOP )F THF LAYFRII) C C Y(LT)tPREVE(LIh)?1IZ(LI) ARE GIVEN AT THE T3P 9F THE LAYER II). C F(), IS GIVEN AT THF 8OTTV1 ]F LAYER (I). C OEPL(L,I) IS THr CrMPUTIOD OFFLLCTIN C AT THF TIP 3F T-4E LAYER I p'N?) IT IS C MEASURE?) IN MM. 0101 1 I AFNS TI N iO ir3o, oro,~I2o o 0302 yMFNSI:NI Y(? Il'I),PREVE(?20J9,1)),IIIZ'?0 Ill) 01303 n I m r his I nV WL ( 209 ) r )( 00) I OM 310) I ZFCH( 1 0) 0004 I NOMFNST i T 7 f 1 ) )N FX R (13), RI Vl),S I3VER( LI),ALFP4T( 10) 0)06?)IMYNSI1N RU?)~(1),SIG ~D(3flMF(131)),I0d(?O~OSGIEo,1 000 C6 r Y SI3Y9 I m i= NS Tpq, SS ~D 0 ri0 S IGGI NFF( 330,1 I 0307 DIMFNSIIN 71PL('00, 10),(?),l( l1(?0,1),LF(l0) 13 08 "I Tm 1'N SI I 1,?AZ(10),SgST-,:ZAISII(11),SISJCI(0)SfUP(t) 0009 DIMENSION DEF(10,1j)),SIMAR1(IO),SI4AR2(10) 1010 DIM!:NSIJN SjNFl(10),SINE?( 113) 1 11 OI MF NS I J" H3(l0) 01,r3H3(0L ok r I 13, 0012 TITMENSTIN Q(10),DEF9'R(l3,O1),EF33(l0,l0)) 0013 L = 1 0014 6 REAnI5,1) WL(V ) r)(Lk, 1(L),9K.Y 0015 1 F1PAAT(3F8.?,T6) 0316 IF(KAY.7T.0) S' TI 8 0017 5 RFAnoI5,?) (LI),I=1,7) 0020 REAn(5t,3!) H(II),I=)T=1,7) 1019 PPAD(5,2) (c:(L, I),T=1,7) 0020 REAn(5,3)1~r)~=t)l~

A1CHIGAN TE~kM1NAL SYSTFM Fl'TPA4' G( 41136 VERTEZ 03-06-73 18:03.56 PAGE P)02 0021 2 F3RM4AT(7F8.21 0022 3 F3RMAT(6F8.2,1121 033!F(N1N.FQ.1) GO TO 13 0024 IF(INOI.FQ.2) GO TO'3 0325 TI( INn1.,FQ.) GO TO 33 0026 IF([nIl.EQ.1?) GO TO 120 007 7 TF(FNOI.FOI.1 3 0 TO 130 0 23 IF( 14D.EO.23) 10 TO 233 0029 (FtINOL.EQ.123) GO TO 123 0030 10 REAf(5,140 (Y(t.,j1t=1,7)I1NO12 0031 4 FIR"AT(7F. 2,r4) 0032 7 L = 1+1 0033 1 f( I Jl?.F 0. 0) GO TO 1 5 0034 G7 TO 6 0335 15 L = L-1 0336 FORTA WL(L) 0037 On =) 0(1) 0038 PREST PUtI) 0039 L * L*i 0040 Wt(L) = Fr-RTA 0341 n(L) =')TA 0342 PU(L) = DREST 0043 GO Tn 5 0044 z0a- REAn( 5, 4) ( PIEV4L, f I ),!=1: 71, TN12 0045 GO, T0 7 0046 30 REAO(5,4) (01IZtf:trT-T=1,7),4ID? 0047 GI Tn 7 0048 120 READ(5,4) (YULt),1=1,7) 0049 AO(') (RFVr (L,T),I =1,7),yN3? 9)050 GO Tn 7 0051 130 READ(9,0) (Y(LTiT=1,7) 0052 REA0(5,f) (O# 1Z( 1,1)!=1,7,P40 0053 GO TO 7 0054 233 READ(5,4) (PREVE(L*HI)=1tr) 0055 RFAD(5,4) (311Z(LI),1=l,7),1\fl? 0056 GO TO 7 0057 1 3 READ(5,4) (Y(LL),r=l=7) 0058 R'Ac t5,4) (PRFVr(LiT)T=1,7) 0059 RFA0(5,'1) (IO TZ( LT 1) =1,7), P141 0063 50 TO 7 3061 8 LIM1TT L-1 00621 CrrFlA = 0.88 0063 0IF?B = 0.96 0064 00 1300 1=1,LIMTTA 0065 WIITF(6#3031) L 9066 3001 FORMmAT(1j1,5XpL==,y31 0067 W ITF(6,800l) 0308 8001 FIR'AT('',5~,f0'NTRIL WRITTIS'S/,5"/6x,********X******.) 0069 WP I TE(6, 8002) WLRL)fn( L) PM L) 0070 8002 F3RMAT(' 0071 00 1001!=1.,10 007?2 T(r(L.T).r1..n..o) rT 100o 00373 T 0074 1Qo0f)1 —.'N T I TN 9'

MTCHIGAN TFRMINAL SYSTEM F TPR^N G(611336) VERTER 03-36-73 18:)8.56 PAGE P003 0075 1002 f = 1 0076 ARFC = F(L,) 0077 I = M 0078 FXPRit) = 0.3 0079 A I = 3.0 0080 1004 1 = 1-1 00R1 IFT.I. T.) 30 T1 1333 0082 K = 1+1 0083 EXPR(I) = (I(L,K)*(E(LLK)IA3E),E **(1.0/3.0)I**2 0084 G3 Tn 1004 0085 100 I = M 0086 SPA = 0.0 0087 5001 T = 1-1 0088 IF(.LT.T1} 3I T) 1006 0089 ZEC )l r = S.'4+CF1lA*(EXPR( I)**(1.0/2.0 ) ) O /D( L ) 0090 SPA = ZFCH I) 0091 GO T! 5301 0092 1006 1 = M 0393 B(I = 0.0 0094 1106 1 = 1-1 0095 IF(.tT.1) 30 TO 1005 0396 A( }) = ZFCHt)Il **2 0097 BI ) ='l=F28*A(l )/(0.25+A(I)) 0098 GO TO 1106 0099 1005 1 = 1 0100 1008 SIGVEP({) = PUIL*Ilo(1-B(II} 0101 ALFPATII) = PUIL)/SIGVER(I) 0102 AL F I) = ALFPATf I )**(1.1O2.3) 0103 1 = 1+1 0104 IF(I. T.M) -,3 TO 6100 0105 GI Tl 1308 0106.5100 1 = M 0107 S!GF( T) - 3.0 0108 SIGRANT(I = ( (P'(L)/2.30)*i.O+2.0*AMIUiL, I-)) O109 N = Mt0110 W I TC(.6,300?) N, SIlrVr{ I ),N, ALFPAT( I),N,SIGFRE I),' t S I 3' I{ I, N, ZFC:, ( T ), N, L f I) 0111 6200 I = 1-1 0112 K = T+l 0113 T'( I.tT.),n Tr 6300 0114 S*3FRt(I) = SIG;EP(IIT*FLK) 0115 SIGRA( I ) = PUiL /2.O*[ 1.0+?.3*4AIU(L, )-2.0*( 1 1.O+/'lfLJ, T )e( T )**(1.0/ l.o)+B( 0 16 4 = T-1 0117 WRITEIS,O?) ) N,13VFQ( I I),!,ALFPTAT( T ),N,$SI RE(I ), J, SGAO)( 1) I N,ZF-4(I t),,,L F (I) 0118 3002 FOPMAT("',SIm VFP (l, ll,)= r6.tX,AlFP TlrrI.J=m', I F5.~, S X,' S S F Q E I i t I I,, ) =', F 5. R, 5 <, S I G R ~ Q( 1 r f l, ( ) =', F 5. 3, 2X,'ZFC'(',tll,')=' FS..2SX,'5 LF (', 1,' =',F5.2 ) 0119'3 T 5O0 0120 53030 = 1 0121.1 = +l 1L 22 F r I ( r = (.. *.141'3 / 5, 9 )' | r t.-A, I.! * ) / F (L, I ) ) 1 ~1 O(I) $*[ (L/ -L (T I ) r1723 F. c. = IFF (I,J)

MICHIGAN TE~MINAL SYSTEM F3RTRAJ G(41336) VERTER 03-06-73 18:08.56 PAGE P004 012.4 J = I 0125 DtL(!.J1 = ALFA/2.0 0126 K = M+1l 0127 00 1009 J=?,,M 0123 OEFL(IvJl = 3FFL(1,K)/(AIFPA.T(J,)/ALFPAT(!*1.0) 0129 DESE(I,J) = ALFA*(1.0-0.5*(AtFU(J)/ALFI[1)**2) 0130 1009:3NrINUE 0131 J = m+1 0132 00 1012 T=?tl 0133 K = 1-1 0134 SI T (Sr'VF(1)*SrrFR(K))/?.0 01 35:81 B (( I)+B( K) /2.0 0136 -F''PA4(fJ) = H(L,T )*( 1.0*AM4T1(LT I)F(L.11* 0137 1012 CINTINJE 01.38 MM = m+l 0139 103 1010 122,4 0140 001011 J = t'm 0141 rWFPAR(fJ) = DFFPARUMM) 0142 loil cnNTINIJE 0143 1010 CINTINJE 0144 K = M+1 0145 KKK = M-L 0146 DO) 101,3, 1349(K-K, 0147 J = 1 0148 DPFPAR(IJ) = 0.0 0149 4 = 1-i 01 5 0 0 1014 J=2pJ 0151 DEFPAR(I,J) = oFPPARI,K)*( (ALF( I )-ALF(J ))/ALF( I) )**2 0152 DFfWlR( I,) =')FFPAU1 T MM) *F-XO-( ~SIGVEP )f9( L))fALF(J)/ALF([))I 0153 1014 CINTINUF 0154 1.013 CONTINJE 0155 Ti 1015 1=2,4 0156 K=44-1 0157 N = I-I 0158 J = 1 01 59 1.01,6 DEFL(IJ) 9 OEFL,(NJ)+OEFPARUJ) 0160, IFB0(I,J) = rFFL(NJOEF,3-lPtiJ) 0161 1 = J+l+ 016' IF(-J.,(;T.) K 3 TI 1.15 0163 GI Tr7 1016 0164 1015 C'NTININW 0165 K =M+ 0166 M4 = m41 0167 WRITr(6,3500) (ALF(I),T1=,M4) 0163 3500 FORMAT1'' I 5X,'ALFlSX,8(F13.4,-X)) 01.69 1 01170 WR I TP(6,'A751 PR"' I, J), = 1 v MI 0171 3751 FORMATlI', IX tD9F'I5X, 8( Ft.,3,?X) 0172 W3 I TE(6,3 700) 0 17 3 373')O l rP3PM TI' I,5X tr X, Lt.) ) 0174 11'A3701. T=1,4 01.75 N = 1-1 0176 W:Z1Tr3(6,370?) NOEFL(IJ),J=I,M!) 0177 373? P 3m TA'',5,K l,7X, 3(f10.1,' )

4!rHTGAN TFIMINAL SYSTrm =P7TPAN G( 41336) VFRT 7 03-06-73 18:08.56 PAGE P005 0178 3701 C3NTNLJF 0179 1 = M 0180 1018 R(U) =&LFPAT{ 1 r)*(L)**2)/(8.0*( OFFL(IKH)) 0181 RAZAM) = D(Li**2/(4.0*(0EFL(IK) )3EFL(!,j))) 0182 DERO = r)EOF3( 1, K )-3FFBO( f #M) 0183 IF(DFRI.FQ.3.0) GO TO 401 0184 IRAFM( I 0(L)* 2/(4.O*(flFF8OU,<Q 0EF80(1,N))) 0185 S3 TI 432 0186 401 C)NTINUE 0187 KKK = 1-1 0188 RAFS(1) = RAFB((KK) 0189 402 CJNTINUE 0190 NJ = 1-1 0191 WRITF6,P4001) NRUI),NRAZA(T) NRAFB(1) 019' 1031 FIRMArTU','R('11,'P=', FIO.?,9X,'Ir 1 9X,'RAFB(*,I1,')=',F1O.2) 0193 1 = 1-1 0194 F(I.LT.1) fGO TO 2501 0195'33 T9 1018 0196 2501 1 = M 0197 2500 STGMOM(1) = E(1,[)*H(Ltl)/(2. t*(l.0-AMIJ(L,T))*R ()) M198 S1S&ZM1) = F(Lt, I *4(Lt,!/(2.O*1.-AM1J(L,!))*RAZA(T)) 0199 1 = 1-1 0?00 IFUT.tT.?) 0 T" 1017 2201 "O To 2530 0?02 1017') 1019 I=2,M 0203 N = T-1 0204 S1SUP(LO) = - SIG'~omr)-sIGTRI(I)+SIGFRF(11 0?05 STGINF(Lt1) = SGM'IM(T)-STGI )(N)-SIG FRF(N) 0206 1019',3NTINUE 0207 K = M+l 0208 00 102? I=1,4 0209 DEPILtT) = 10.3*nEFLI(K) 0212 1 070 C3N TT'JJ 17 0211 WRITF(S,1021 0921? 1021 r IRMAT('?',53K,' PE SJ1t',f,59,'*******,/,5X,'F,12X,'4 1,8X, OImlJI t 9X, IFI,1X, IYI, 17 I SIVI,17XI I rSGSUP I, I t xII' S I G N, / X, I < -, F (-'A' 9 y, I - M', 3 f, X, 1 Y Y v, 1 6 X, v rKGGF / r %4 i i,?7 X, K 1 F / r "A?1 3,/.48X, I - SI PR-,X,'CnMP!ITCf),1Xr I 4SU C' r3X,', MDUTFnl 1,P3Y F A ) I JR F I, 3 X f r nmpSIIT EPl) 3.. X,U, SI F! I 1213 -)I'177 I r, 0 21.1 Iol I 1 T T NI f 01 217 1 =7 2218 PA Z I= ) 0219'Il' = R47A(1I 020 I =7 0221 4300 S MI1M1(I) I STG MOr I(*I(R)I 4Z 0?22 ST MA7I(T, I= GA Z( I P)AA T ) 7 2 I'Pl (I) = -SF 4API ( I 10 )(1)+SVFR EE(1) 0224 0225 Si v1(1) = 6l A l(T)-Sr oJ- fQ Vu 0276 SVJc P'( T) = - ST Y R?( I)-S I' + F l

MICHIGAN TE),MINAL SYSTEM rlRTRAN G(41336) VERTER 03-06-73 18:08.56 PAGE P006 02?27 SINF2{ [} = STMAR?(f)-SIGRAD(N)-SIGFRE(N) 0228 I = 1+1 0229 IF(I.GT.M) GO T! 4400 0230 GO Tr 4300 0231 4400 WRITE(5, 4500) 0232 qt503 FIRMAT('',/4X,5X4X'S ISUPl, t7X' S NFl' 13 X L,'SISUP2',.17X SINF2tt/,b6Xt'R',23Xt R',22X 2,'RAZA',19Xt' RAZA',/t 7X,4('MEASIIEO t.3X I COMPUTED' 3,3X)) 0233 T) 4600 I=2,4 0234 N - t-1 0235 WRITEi(6,4701) ORIZtLI,$JE$tJP(LI),DRIZ(LI,SNFkL(t) 1,tORIZIL, I)SISUP2[1),tORIZIL, 1),SINF2([) 0236 4700 FOR'ATt'',5X,8(FF9.4.2X) 0237 4600:INTINUc 0238 I M 0239 Pt001 DEPLittLt )-*RAZL 0240 PROD2 = DEPL(L.i:k*RAZ2 0241 WRf TE6lt6,4800!.OOW PRO2 0242 4180 F3RFMAT(' *S,t5X,'PRTOl.RI'10.2,5X,'PR002tR ZA)'El —'E0.2) 0243 l M 0244 PA003 = R(I-)*DFPL(L,I) 0245 WITF( 6,4853 ) P 0r)3 0246 4 850 FORMAT(0' - X,' PR)3 ( I ) =' F02 ) 02 47 VAL = 0.0 0248 I = 20249 301 Qti =- Ft2 Lt, (H(L,I )**3)/(I.0-AIUILtI)) 0256 VAL = V4L+QOi ) 0251 I I1+1 0252 IF I.GT.MP) GO T3 300 0253 G3O r3'301 02546 300 r = I 0255 ANUMER = 0.35*(1.O-AMIIJ(L [ )**2)* 1 )L(L**2*4LF ( I) *( L ) 0256 ANUMIT = F(L, I)*VAL 0257 DEFI... = SQ3TlANUMER/ANJ'4IT.*1.? 0258 WR ITE(6,302) OF tI Mr 0259 302 Fn1RMAT('',/ 5x,'FP 4LMn=,F1o.6) 0760 DFP.. = OEPL L,1 026t t = M 0262 DEPLPl = DEPL(L, T) 0263 DFFCOR = DEFLMn0*-EPl?R-QErLAS 0264 WRIT7I5,301) nr:CR 0265 333 FJRMATI'',SX,'DFFCOR=%rQI0.6) 0266 SIM.F = 3.0 0267 S'M.H = 0 0 0268 0O 90003 =20, 0269 H3( 1 = H(LI)**3 0270 SIJMH = M4+3 ( I ) 0271 F3H3(I) = E(L,I)*H3( ) 027? SJMF = S:JME+-F3H3(T) 0?73 900 CrNTT NJF

MIC-4I';AVi TrmTNSAL SYSTr- A FP1TRAV i G(41336) VF T 03-b-73 18:)8.'6 PAGF D0')7 0275 RAPMOD = EMED/A3EC 0276 WRITTE(5#909) EEDIRAPMOO 0277 9001 FJRMAT(' I,/5X,'EWFD0'IF12.2,5X,'R&PMCJDRFF1..3) 0218 DI 1100 1=11, 0279 A(Ml = 0.0 0280 EXPR(I) = 0.2 0281 ZECH(I) = 0.2 0282 A(M) = 2.0 0283 SIGVFR(I) = 2.0 0284 ALFPAT(T) = 0.0 0285 ALFUl) = 0.0 0286 STGFREf1) = 2.0 0287 STGRAD(I) = 2.0 0288 Z(T) a 0.0 0289 Si~u~y~)T5 ) = 2.) 0290 H3(I) = 0.0 0291 F3H2(I) = 0.2 0292 00 1101 J=197 0293 )-FttI,.J) =.0 02946 FFPARUIJ) = 0.0 0295 1101 l~NT1NIJE 0296 110)''N TrrNJE 0297 1000 2NTT'4UF 0298 STF)P 0299 Fan *OPTIONS T4 FFFFFCT* T), F R tI,S3URCF, LTX I ST, 43 ECKL0OAD9NflMA P *OPTIONS IN EFFFCT* NAME = VEV1R, LINFCNT = 57 *STATISTICS* S3URrE STATEMENTS = 299,PR)GRAM SIZE 111452 *STfATISTTCS* WI1 2)1&NOS7TCS 3ENERATFD 40 ERPR~RS IN VFRTEQ j3 STATEMFNJT FL&'rr' 1' rT4 PAPOVF -IMO1TIArThS.

APPNDIX 6 MICHIGAN TERMINAL SYSTEM FORTRAN G(41336) VELMED 03-23-73 09:19.56 PAGE POOl C C SIGVEh((I),SIGFRE(I),SIGRADf I),ARE GIVEN C IN KGF/CM2 AND THEY ARE MEASURED AT THE TOP C OF THE LAYER.. C C DEFL(,J ) AND DEFPAR(IJ) ARE GIVEN C IN CM AND THEY ARE MEASUJED AT THE TOP C OF THE LAYER. C C R(I) IS GIVEN IN CM AND IT IS MEASUREL' C AT THE TOP OF THE LAYER. C C SIGMOM(I) IS GIVEN IN KGF/CM2 AiD IT IS C MEASUJRED 4ITH SIGN - AT TtIt TLP OF C THE LAYER AND WITH SIGN + AT THE BOTTOM C OF THE LAYER. C C SIGSJP(I) IS THE HORIZONTAL STRESS AT THE C TOP OF TH E'LAYER, SIGINF(I) C IS THE HORII]NTAL STRESS AT THE _BTTOM OF THE LAYER. C'B[TH OF THEM ARE MEASURED IN KGF/C'42. C C EXPR{I),A(I),ALFPAr(I),iEC1 (r(I),B(I) ARE CALCULATED AT C THE TOP OiF THE LAYiERti I) C C Y(L,I),PREVE(L,I),O'RIL(L,l) ARE GIVEN AT THE C TOtP O' TH l LAYER ( I.. C F(I), IS'IVEN AT THE OUtTTO OF LAYER (I). C C DEPL(L,I) IS THE COMPUTED DOFLECTIONj C AT THE TOP OF THE LAYER I, AND IT IS C MEASURED[ IN IM. C IME NSI)N E'(Z,I0lO),AMIUCZR20 10) tF2OL O10) I,H2 200,1I) DIMENS I N Y( 200,O), PRE E(20 0,I),ORlZ ( 200 10) DIMENSION WL(200),!)(200),P(2OO),LZECH(1O) )IMMENS IUN Al O)1 EXPR(.U),O( l (/),SIGVERIlO),ALFPAT(1lO) I M E NSION ) SIGFRE(1 J),SIGRA)D(10),)EFL(IU,10),DEFPAR( 10,10) DIMENSION R(10),SIA;UOM(10),SIGSJP(200,10),SIGINF(2D0,10) DIMENSION:)EPL(200,D)',P(200,I0,),0(200 t10)tALF(10) DIMENS ION AMR(10),AMF(I0),AMSUPL(1l0) DIMENSION ANTPAT(1O),ANT( i ),RAPDEF(200,10) DIMENSION RAPREV 200,10),RAPORS(200,k A PORI(200) OIME -SION LIST 200) C L = 1 6 REAC(5#1) WL(L1,D(L)IPUL),IKAY 1 FORMAT (3 F8.2, b16 IF(KAY.GT.O) GO TO 8 5 READ(5,2) (E. L, I), I=1,7) READ(5,2) (AMIU(LI),I=l,7) READ(5,2):(FL', I)',:I=1,7) kREAD(5,3) (H(L, I I)=2.,?),INDI

MICC-11;4'Ai ri fMINLNAL sYTrt, -' I T;..;,J:-Lt ) 3 3-23-73 09:19.56 PAGE POOZ 03J2 3 FORM Tr(6F J. 2, [12) 0022 IF(I )1.tE~.LI GO TJ 10 0023 ftF W(O1.EQ.2 t GO TO 20 0024 tF(I Nt4.E.. 3) GO TJ 30 0325 IF( [!).tOI.12 ) GO3 ra 12 0026 IfF(INoL. i.EQ. 13 GO To 130 0027 IF( ti,)1.EQ.23 GO TO 230 0o28 iF(iN)i1.EQ.123) GO TJ 123 0,i j34 15 - = t-1 0335 FiRTA = WL+ L) O-3O8 Lt = L+L G039 Ok.40-t L Fkf r A'~'.;+ Pt P L) = PREST uJ41 T; 0~4, 2,0.rAnAD 5.4) lPREVEIL, I'-'-'-i,I7iO2 0145 0 u 4 zt I 1- (1, 0,46..Gd- tO Qo04-'.1; Zb:':"~ AC,~- tx.t 1 o f"L:,I=',,?.A 0048,RE-A D ( 5, E.. V E 4.L i } I I. NN2 00.9 GI TO 7.,C.,0 t,'a.Ao ( s,.4J -AtD 5 4': -- fZ: 0031 i.F- EA (~,z,.j) U 141t QIRZ..t: -, [},ii.,I, DIh 0052 GrO t 7'-:' 0C 5 3 -j230 READ,4 t (P.RE VE( L, I,I = 1 7 oos04 a E Ai`( S E, 4 ~: ZiZ ( 1, II, 1 7 t.N...2 0055 GO -7''}5. 1' 3' t*-' -,"t) t('_,4,tV (Lt =i,?);.'",; o>L~ e) L.AI 1 1VM I L I = i=, 1 O5-i - R'EAO[5,~4) ftRIZ(L,'l}), I~,,'-Fi:'tNh Z2 O:,;,GO T 1 7 *-.1 i S~Lifr = L-1 C = b )52 _J Of' = 0. rt Cr,,00 [li),O L-/,jLI*TTA v>' ="* r, ~I TiE(,, OOi3 l ) L JO5s 33J0j L-t,;MAT. { m' ATi2 -I'...'L', I3) CC'~,,r/':,1T E t(-,,c k dOOU { b0?7. O.)I FOF MAT(i',5X,'CONTRUL R I T I S',/,-X,'*******1: *******'*")' J(T)h~ y 9 W T' (6 /i3o,8002 L (L 0I,. l L,P (L ) c1. i~3 -0JZ) FutMAr(''SX,'WL=' tFd./-,t)A,'o=,Fs.ZbX,'UJ=F o F5,) J3C70 0 1'01 I=lt l, "'71 IFi(LttI ).r.i.')., ilb? J07!' 4 = I,' 7:.) -'! I'-!.JC,- i I: i

MICHIGAN TERMINAL SYSTEA FORTRAN G(41336) VELMED 03-23-73 09:19.56 PAGE P003 0015 LIST(L) = M 0076 ABEC H E(LI) 0077 1M 0078 EXPR(I) = 0.0 0079 A(I) = 0.0 0080 1004 I-= 1-1 0081 IF(I.LT.1) GO TO 1003 0082 K = 1+1 0083 EXP(I) = (H(LK)*(E(LKJ/ABEC)**( 1.0/3.0) )**2 0084 GO TO 1004 0085 1003 I = M CQ86 SPA = 0.0 0C87 5001 I =,1-1 0088 IF(I.LT.1) GO TO 100o C089 LECHMI) = SPA*COFIA*(EXPR(I)**(1.0/2.0))/O(1) 0090 SPA ZECH(I) 0091 GO TO 5001 0092 1006 I = M 0093 B(I.) = 0.0 0094 ALF(I1 = 1.0 0095 1106 I 1-1 0096 IF(I.LT.1) GO TO 1005 C097 AM) = ZECH(I)**2 0098 B(I) = COF23*A(I)/(0.25+A(I)) 0099 GO TO 1106 CLOO 1005 1 = 1 0101 1038 SIGVER(I) = PU(L)*(1.0-B(lJ) 0102 ALF(I) = SQRT(PU(L)/(SIGVE{( I)*2.0/3.0)) 0103. ALFPAT(I) = ALF(I)**2 0104 1 = I+1 0105 IF(I.GTAM) GO TO 6100 0106 GO TO 1008 0107 6100 I = M 0108 SIGFRE(I) = 0.0 0109 SIGRADlI) = (PU(L)/2.0i*(1.0,2.0*AMIU(LI)) 0110 N = M-1, 0111 6200 1 = 1-1 0112 K = 1+1 0113 IF(I. LT.1) GO-TO 6300 0114 SIGFRE(I) = SIGVER(I)*F(LK) 0115 SIGRADO(I) = PU(L)/2.0*(l.0-2.*AMIU(LI)-2.0*( 1 1.O+AMIU(LI )) *B(I) **.(1.012. 0) 8(fI) 0116 N=I-I 0117 GO TO 6200 0118 6300 J = M+l 0119 SUMOEF = Jou 0120 O 1012 I=2,t 3121 K = I -l='A 0122 SI =, (SIGVER(I)+SIGVER(K) )2.0 0123 81 = (B(I)+BIK))/2.0 0124 DEFPA I,(19 J) = H((LI) *(1.044MIU(Li))/E(LI)* 1 (SI-2.0*AMIU(L,1I)*PU(L)*LL.0-Bl**(1.0/2.0))) 0125 SUMOEF = SUMOEF+OEFPARlIJ) 012'& 1012 CONTINUE C

M[C-4Tl lu I c313 i'4E3MibL) D 3-23-7i 39: 14. 5b PAGE P004 C A E TriJ 1 1 INSISTS JF C XI = SQRT((j*XG-C)/A) C -9?H4tY) 2 LfJiSISTS UF ~ -/t(A*Xo-b) C METHOD 3 I0 T1SS JF C Xl = 4ABsX0*YL+C)/ C M-FTHQOr 4 ~GC)NS~J:VS, 0F C XI = (B*X)-C)/(A*XU) C 0127 SUMP 3. Cf IL jDC 13.) >I2,1'4 OI2c- SUML &l3ai3 1134-0. 473TIN4JE 0131 I = 0152 A1r. 4.03 *j.1416/L/3.J*ALFPAr ()*i( L)**2 ) *SUM:L"4. I~~~f:3 ~ ~ gbat,a ~~:F P~ATH At +1 WAd ICL+ Lii A FPATI~rf 0135 00 1200 l=2,' 1] 3' K 1-1 t L31 aMR ~d A) C PAi I I i -1(L i /Z f J S1(FiGF E IU*SJGF R Ei(AK) I ~LVL3 4M*SUDLL1) =hi()~ c144~f 1201) ~~t4T14UE 01 ti (0 D Ow Cl' StH4ON~~~~~~~~~~~I I.2 I)C 12.1 >2Z,1i 0144 SUMC61 = S J J#M+A A UP L I) 0145 1221 CO NTINU l CL46 =5 I C 0147 1471 C-, L) F* S UMDEF CL>-i 3fE =. I>"AI T = 0 -3 IT;'.~ 6tlZ57)~7.3 1> ~1 i 7 FLur:AAT(',3X, 1A I-: EGATCivF') i2J2 1) L A S ) 1203 12)2 U).LLr z 14AI 113U~ 5157 i209 It:L'L T-uEr-Lli.3 LT..J1Jj),cLFu 12 4 46IlF( LA T( )1 0U LEC).GT. 1 0'i. 0) OJ i 120 6 0i 1 59 1F (AfS( (DEFL JT) LT.0*00JOi-) Cu — -T L23 015~81 T 6 12'i. -1 127 ( )',5X,',iw~ ~ r\JO CJuufCU UjTlZtJFLcU-'TlJd\4 l: TUO) ~oAAL Li I164 I, —)5121 O[FLL;C J -rLAr (1 I 2 C) (6, 1207) ~;\7 1 3' 1 P -', I [ T 1: J _rL! L

I4'CH[GAN TERMINAL SYSTEM FORTRAN G(41336)'IELMED 03-23-73 09:19.56 PAGE PO05 0169 1203 DEFLAT = -C/(A1*DEFLEC-BIC) 0170. JALE = JALE+1 0171 IF(JALE.GT.JO) GO TO 1290 0172 GO TO 1291 0173 1290 WRITE(6,1292) 0174 1292 FORMAT('',5X,'WARNING METHOD NO 2 IS NOT CONVERGENT') 01.75 GO TO "1293 0176 1291 MIMI = 2 0177 IF(ABS(DEFLAt-DEFLEC).LT.O.0001) GO TO 1204 0178 I.F(A3S(DEFLAT-DEFLEC).GT.100.O) GO TO 1206 0179 DEFLEC = DEFLAT 0180 GO TO 1203 0181 1293 JALE = 0 0182 1295 D'EFLAT = (Al*DEFLEC**2+C-)/tiC 0183 MIMI = 3 0.184 IF(ASS(DEFLAT-DEFLEC).LT.O.0OO1) GO TO 1204 0185 IF(ABS (OEFLAT-DEFLEC )',GT.100.0 GO TO 1206 0186 DEFLEC =.DEFLAT 0187 JALE = JALE+1 0188 IFIJALE.GT.301 GO TO 1294 0189 GO TO 1295 0o190 1294 WRITE(6, 1296) 0191 1296 FORMAT(''iSX,'WARNING METHOD NO 3 IS NOT CONVERGENT') 0192 1360 DEFLEC = 1.0 0193 JALE = 0 0194 1351 DEFLAT = (8IC*DEFLEC-C)/(Al*DEFLEC) 0195 MIMI = 4 0196 IF(-A3S(DEFLAT-DEFLEC).LT.0.0001) GO TO 1204 0197 IF(ABS(OEFLAT-DEFLEC).GT.1.0.00) GO TO- 1206 0198 OEFLEC = DEFLAT C199 JALE = JALE+1 0200 IF(JALE.GT.30) GO TO 1350 0201 GO TO 1351 0202 1350 WRITE(6, 1352) 0203 1352 FORMAT(',5X,'WARNING METHOD NO 4 [S NOT CONVERGENT') 0204 GO TO 1000 0205 120-4 K = mi+1 0206 I = 1 0207 DEPLCL,I)'= DEFLAT 0208 DO 1501 I=2,M 0209 N = I-1 0210 DEPL(L,I) = DEPL(L,N) + DEFPAR(I,K) 0211 1501 CONTINUE 0212 I = 1 0213 R(I) = (ALFPAT(I)*O(L)**2)/(8.0*DEFLAT) 0214 PROD = DEFLAT*R(I)*10.O 0215 WRITE(6,4800) PROD,R(l) 0216 4800 FORMAT(''/,5X,' PRD.=' F10.25X,' R=',F10.2,'CM' ) 0217- = M 0218 K = M+1 0219 DEPLAS = DEPL(LtI)*10.0 0220 0221 K =1 0222 2500 SIGMOM(i } = ELI) *H LI)/ (2.O*(1.0'-AMIU(L, 1)*RIK') 0223 I = I-1

* f ~ i I ~~ J ( i~ 1 z ) ~: i J..i~:_ - J ~- L -3-3 09: 1 9. 5 6 PAGE P006 4- 1 /+ I " (!.L.2); Oi r I i1 7 1017 0U Lil. 1=2, 02~' 4 1Jc -1 3 - i j 5 I IG, tP (L,1) = -S I GM ( I)1Gi)(L)+S[.,F E I) Si61!14-'( l ) = I (;; IJMLI)-SIoiPA t) 4-SIGFRE(N) 1)9 CJNT J - 1'4 =`2733 MfPL(L,1) 1U..O: EPLlLI) 0234 12.i2 C NTTI WJE. 1021 FlO&MAT(' E,58X9 EIRES UtSe, 59X ru1 L**** / t5Xv*E,1ZX-t FI 8 X r i1 lu I i 9X v F vlW:1 v A Y' L 7 X i S IGCllER # I? T-XCt I S I GSU-P #-I~ Zj~ IJG f'~i~s~, 2,2X';F/k7 i',ix1AC -,`CIM-, OA I 44i4',1 I oX, K / C M21 t22X,*)GF / CM2' - 6.Xi #Ati-A-iQ-KfE f/4 #3A,*a`-:; P:T,3 0,3X,4 pMEASJRE0,3XVOCMPUJT E.D* 1~-,.3 3 2 I t, I L t,tt L. L 1, - t 4 *: L e )rP I-1 F-f- r I I y i c OEP-P ( L I I i #PR 7 L I V E r trt U L v I S S Jr' i L t I J S I fH ~ij 0239 4W 23 FG KKAT(' tF I 0. Z iq,1 ZF9. 0 I40 212 2 coh[NffTrE it-(Y(l * J.Lc'.U.UJ 03 To 9 %-jO Q2t4 IF t A4) )L —'L(L,I..1A)jU-J GJ Tij 9-iJ ~-I L:4L ~ ~'1 j LI 0O Trj'95i1 02. 505J2;(4P iL i ) =.0 0 ~i c'''ih1TI~ (o,- b t,d't~i) (KRAPUL)LFt (L, Ioi I=1 A) vdF i ( 6, 9Y-j) C k4 PR E Vt, 1),rI =1,t4) C, X,' "K Cii),rP I 7( F. )A)i 4' I,A~t zIL) = ~'" r 0~ )' I/ L t'A.EI..P.. 02:,~~~ ~~aPI~'~ ~) 7, 32 T ~. ii i L r.).E~ ~,. ~, r I, 5U.- p,' C ~ - 02ti 970C 4pf! ( 6,r 9fSL -iA) P lORfLkR3i4L) C, ~ ~ 5,09 ~M Al~ ASO 2- <SIF. I5x, 2-'21', F-. 3) $702 DO 11.73 11,A 0271 ACI) = 2? EA V"

MICHIGAN TERMINAL SYSTEM FORTRAN G(41336) VELMED 03-23-73 09: 19.56 0273 ZECH(I) = 0.0 0274 B(1) = 0.0 0275 SIGVER(I) = 0.0 0276 ALFPAT(I) = 0.0 0277 ALF(I) q.0 0278 SIGFRE(I) = 0.0 0279 SIGRAO(I) = 0.0 0280 RM) = J.0 0281 SIGMOM(I) = 0.0 0282 DO 1101 J=1,7 0283 DEFL(I,J)= 0. 0 0284 DEFPAR(IvJ) =0~.0 C285 1101 CONTINUE 0286 1100 CONTINUE 0287 1000 CONTINU4 0288 NUMAR1 = 0 0289 NUMAR2 = U 0290 NUMAR3 = 0 0291 NUMAR4 = J 0292 SUMRAP = 3.0 0293.SUR4PR = 3J. C294 SURAOS = 3.3 0295 SURAJI = 3.0 0296 00 9531 L=1,LIMITA 0`97 M = LIST(L) 0298 D0 9630 I=19.M 0299 [F(R.APOEF(Lfl.EQ.U.0) GO TO 9632 0300 NUMAR1 4U.IARL+1 0301 SUMR.P = SUMKAP+RAP0EF(Lti) 0302 9602 IF(R4PKEV(L,,I).EQ.0.0) G;U Tn 9630 0303 NOMAR2 = NUiAR2+1 0304 SURAP3 = S-URAPk+RAPREV(LI) 0305 9600 CCNTINUE 0306 IF(RAPORS(L).EQ.o.0) GO TO 9604 0307 INUMAR3 =;4JOMIIAKR+1 0308 SURAJS = SUiAOStRAPJRS(L) 0309 9604 IF(RAPOi)I(L).EQ.0.0) GO rO 9601 0313 NUMAF4 = NUMAR4+1 0311 SURAII = SWUAI4I+RAPORH(L) 0312 9601 CONTINUE 0313 WRITE(6,26U.O) OUM4ARiNOMAR2,NUOMAR3,. NO MAR4 0314 2600 FORMAT($ J',X,4(I3,p5X)) 0315 IF(HNUMARI.&.0.0) GO TC 1900 0316 AMDEP = SUMRAP/NUMARI 0317 GO TO 1901 0318 1930 AMOEF = 0.0 0319 1901 IF(NUMAR2.EQ.0.0) GO TO 1902 0320 AMPREV = SURAPR/NUMAR2 0321 GO TJ 1903 0322 1902 AMPREV = 3.0O 0323 1903 IF(NUMAR3.EQ.U.0) GO TO 19,04 0324 AMFORS = SURAOS/NUMAR3 0325 GO TO 1905 0326 1904 AMFORS = 0.0 0327 1905 IF(NJMAR4.E6.0.0) GO TO 1906

.11 -,'' L' ~ J (L13bL D 23- 73 ij' g.?5) P'GE PO08 O32,'':= jRt4U1tNiMAR4 9....'."'. F -," (i t -)F- r{5.',:' 4,,,,'' Pi J.EV=, F5,2 t5X,' I,- ), s= t,. t..-,,F L i _.., A fFOK I =.Z) i,'' ~R..". 1 = I ~ 1 A?", ~: i ='! 0. 03 t' = LIST(LJ o- 7 3ij ojg6 i i.M O -.'' { i;r'n' C 1 L, [) I E' J:. O': 0 T'9 6.53: 0339 itl~' 1 - SJJUQl tP1Eftd -AYUEF e 0340 g65 i ( 1; \ "-'V IL I).Ef ).3.0:~-: TT 9i) -62 0 4, - + R ~ E T t. i L ) -I'' R EV l*.2 vv43 I r J R; 1R S (L 3 ^ 1;j.b ) J J, 9 u 0Y3-44 SuR'i SU S0RA:'.+: R Ai,'R S' FORS )'*2,. i;''-654 Ifr l.:Oil LI.E 3..) AU i13 9'65I L;34t 4Uf. S-+ = 4UF 4+ -, PJi]<I (L)-AAF3~,[)i*~2 c5;' 7 -J64t' C or u T I -UF: - C-~'ts::' It-(',),l01. -:L.:4o.O) IA; TC 19[50 el3,'.9,,.', — = AsFt j(:U, AI/ NU4MAR L 0 350 3 rJ L?qi C351l. 4_Mc 1 nj'Z Ji \-::i-._'- = 0.0 0352 Ll F1 1FNu,'4AR,2.E.O') SGu TU 13D2 i353;-,^F.Pl: =.jFR, IR Ir URA I/i'UAR2- t u:': C 3,- 3 T'.I. 9...53 I:.:'5 -19~z la~EOD E _: - J, 3 b,, l 9 53 I Ft' N'JJAR.3 E1:.'~{iO' GO T -, 9 >4 ) ) ta 7 AE F 1;3:} i OG T l':~:5'D.-:, t? - - 4' 15! 4 X,-i 5 = t-) (-...,.) i 7')' IF (, Ut,'J-.:/,,,E~.u.;)t Gu r" i956,:.oi A'.'!i: - I = SJ:Iw SJ:~~t':WJMAH4 0 30 q 1': \': L- ~:[ U 0:i -:<"......, ~..'. "~:,.' I'T"{,::) 6' —' M,E tEDEfF,:IE,.,,i FA. }S i AAE' I 031~~~L 195~~~: ~ (I T L r7:,',',:. )tJo~ r:J K;-,( A A:`A**M: E F-' F 3S,4.',, 4lr_.-PRE=#9Ft-.3,A 1' AMEFG =',:i.. 3, 5At M F 1= A, F 5. J I J * J1 J *c:X * ~ ilS I -i -.F.-riT:' [ )fF, ui,L;', id;, "L'..;L.,,L.,JAl),NUMAP *tOPTIONS IN FF-FF,: T*?AME- = ~ELMED, LSNECGiF -? 7'TATISTIC ",';;:'E.TAT4E~NT -: 36o,P:':,-, -' Sl: = l)Cg7 J -,,T:';rI". A' j )ST i,.,,-;l.I" ['k NO ERR.', IhN VEt:j ) "\i',.j:- 4E:, r S F L.A E JE I,, T Jt LA., l:I;4 P ILAT I l.,S.

U NIVERSITY OF M ICHINA li 3 9015 03695 1922