maggfifig. a. _ . 2a.... .....r..r..: .. .. ... .. . ., , _ . .. .. . .. . . . . : . , .. 1...... . . :5... z,“ .. . .. .71.... . r . .. . . . . . . . : . . . . . .. ,. .. . . . xc.r..ewnu . v [If . 2 J 11 .... i .. v .r- 4,. . . . . . . . . . . . . .. .. r m.«‘!r..1.!v 01.9. .1 a}... .n; 3. 7 . whh..\....$. . .. . .. ... . .. . . . . . . . .._?..~ 4.2:... . lama»? é... .. 21.. .. 1m , . , .. . .. .. . : . . D . .1. .. _ . 1 a . .1. .. t . . . .. I. . . . . . . . ..u . : .4 ”NET 33”.... a2“. 33.. ‘ H kmBfiM «Wumm. .. .. . x: «fin-WW... . . r . _ .. . .. . ‘ .. , l . _. . ... . . .. . . J .. .. . . . . . .. . . .. .. . .; . .. . . , . . 3. . . .. . .. . : .. . .. . . ... . . m. . L . .. . . .5} z . . . s .. a. . . . ~ .R. .f ‘ .19 .rmwfi .. 5.1%.... .3»... .s. . . o~l.7¢. .. . . . . .rr . . t . .2...» t . . . . .. .. . . . .. .. . .. . ., “Cu. $2..-.Rttuzf law. 2T . . .. i. . . .. . .. . . . . . . . . sky}. wry!!! v»... - .. . .. v.9 . . _ .. . . . . .1 . . . u. -. (I: ... 6 .. ... , . ... ._ . 1%... r... .4..olk11(u..il..u.fll\. 1. .L . .. . E . . . . . .u v I! I. . . .. . .. .. a. . . . . . . 4. . w; a. .1 .2 .1. .... .: . . .. . . . . . .wunuh.i..¢,a. it... .t . . .5“ . 3.... . f7. Mn. . r. . - . . . . . . .. . Elvin...” .gfiu ! 1‘ 171:\ . 3.. .\2 L3 .WQ 15% .. . I.1. .I in} fig. . n a! ,. .. I.l1\ on .. -hwvximfiasu. ...HL......... . .. 112.3... h....fi...x¥ ... . mm... A... . an... . «$3.13.? . _ ”X5 fig . . . . . . . .. . . ) .. :1 3 . . . ). . .. o. . . .. n. 31:41 I . . . . v. n‘!!!.!!v¢n _ l O! . . . I 1 . . '1‘ .'G .‘tll'i'i \v '4 v t . v _. .h s A .‘ . c .D .. J. n, '2}... Ill... . . . .. - ._ . . .. . . . . . .. 3.}. 131m... 5: .....:i....!.k.2i1.flz 1. vi...l!!...1?$}§o¥. .. y . . . .. .. _ .. . . . . .. 3.0.... . F11? .8..-ir’in ill}; £512.? . so...) . . r .1\.fl3.. u .. , . . .. x .. . .r‘! . :1 . 505‘. seizure. 13...... 1.5.9 . . .. . \Tt..\é\ _ .1 a... .w‘ . V..m...-.bf.l .. 45..-»..3. . 1.63.“. I)! .. s. a X Wuvhfl!ka3nb .3 . . 2 a»... .. .I.‘ 483:}... 31x13 ... . '1’; {5.1 . ,. r £1? . u. S... .f. .1 . . . .. . . :fiuwhunu.“§li .594... .iszafizm umwbwgs .amWI... - . . a... .33., buuvthz 9...: . . . . .u. .5013). 112133."? 5.: Eu u :d:!.a.r .. . P{P§.1.hw§.:v) . . . . a . .. 1L8. 4.9.3.5129 .1. 1.5.... . ......§.3........8,..? . 2.... . .. .. . . .. 1 "m... . «123...... $525.31 . 1.24.18... 1.. .. "a . . . 4 . . A .n _ v .0 5 ’I {ivriufflfiifl .. . . . #333523315iuk1. a. Lqu 74%.. . . . n .ii. u... . . .. . . . . . . .. .3... . . .fiyh. . . . 1...». .. fig... .. :31; :2 1:8. . r. . . .. . . . . . ..... . z w... alum-L1 3n... ifitgr. 33.3.}...m.‘ . taunfionryiv.‘ » .gflnflunfli 0T“ . V\»\Nu3x..a__r§\ 7. . . u . J)‘ runwxl . ”V4.35“? gran? 3». r . . 3L. . .3 .. ‘ .~ .\ 3.9L»: Klilil? 31.. he? . 5.3.5... . . . figg fl . . . .. .3 . ... am, .3. Mann ._ ,. . . 3:38.551 . .. r. . .. 59.3. .. {(51) I}. .. . . .. . 5%... 8.3»... r . fl... . . . "t uh Q . IJ’RFtLqu-fl. It. . . . line... .7- . JIIIutf... . grflfirafiwm .z1.fim.4&zi}.uhu3§i{; . . . . . .4 t .v . M . 2...). :‘Q 18'! . . . an}... . .3 .2633. . F... . ...... . . . f0. 13:. p.71... kn . . .n . . .. .. . : ,. . . Mauziflkmm. 1n.r..........p...... .. . . \...flnvtn stick. ‘5 . ., . . . ... . . . ... x... .r. .49....) .i. . . 17L... 2.1.. . . . n. . . . _. twafimzétx .wguhfifihufiwf i . {573.11.1410 .{fixflnz’fi . . . . . . “EWQQ‘I ‘ 9,311 . . .. . J . I. . ‘ k . I" . ”an-c . it . . .. . tutti). lint vhu 7: . » s. . fl :1... . .. . . 19'... .\.I L . . . .. .. . ....t..§nu..o...t....k. .. . .. . 1.2.5.... .1 2 . . . . . . 1t............. 1...". 4.3m. . II. . . . . . . . I . (rt... (.9 . .. . . . a. . .nnnfifili...w...£i.ft , .. u . .. . . .. .. .. . 3h . . . . . . . 3.3... 3%.}...2122. s.-. . . . .. , . .. . . .1 z... w. a. k .3. ti. .3“... . .. . . . .. . . .1... nuanfiiwfit i......t...... r. f. . . . . . I E), .. rhuu...<:. .. 5 .3 gui . . .. . . . .iurh..¢..t . 0 )ii. » . 1 . 21.55.5453,}... .5. - .. .iaoflrafihi..fif..hnsaunilnutl:é....n . .. . . . - . . .....-I€....-.l!w.anuz.l.ll.tH.. .. . . ‘x .i: . . .N: . \‘ . 4.... . ..... .2. . . . all... Va: 31?...2‘...) x.’ git:~1§. .. $1.. . - 1.11!!-§~1-....E§E}N.(n.rr 1..... .73... .5I...Ic.....f-¢.n.§1{. 08...... 2.... {dwflmwlu “5... (Jun :it.:1.u..fi.-.rxiv.i.3 : . Hu..ru$(s1...9h. .H. Wag LIBRARY MichiganScnc , U" .'ty This is to certify that the thesis entitled FINITE ELEMENT FORMULATION OF THE ’IHERMO—HYDRO STRESS PROBLEM IN SOYBEANS presented by Kamyar Haghighi has been accepted towards fulfillment of the requirements for Ph . D. degree in Wing and Metallurgy , Nbchanics and Material Science Major professor DateMarCh 8, 1979 0.7 639 lllllllllllfllflllllllllll Illlllllfllllllllll 082 OVERDUE FINES ARE 25¢ PER DAY PER ITEM Return to book drop to remove this checkout from your record. ‘ zoo A324 ‘x W‘- . / 300 A265 FINITE ELEMENT FORMULATION OF THE THERMO-HYDRO STRESS PROBLEM IN SOYBEANS by Kamyar Haghighi A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Agricultural Engineering Department of Metallurgy, Mechanics and Material Science 1979 ABSTRACT FINITE ELEMENT FORMULATION OF THE THERMO-HYDRO STRESS PROBLEM IN SOYBEANS By Kamyar Haghighi The objective of this study was the formulation of a finite element model that could be used to analyze the stress crack formation in soybean kernels resulting from temperature and moisture gradients during the drying pro- cess. The soybean kernel was taken as an isotropic sphere, symmetrical with respect to its center and a two dimen- sional axisymmetric grid was used to model the temperature and moisture gradients. The finite element method was used to obtain numerical solutions to the simultaneous moisture and heat diffusion equations describing moisture removal and heat intake pro- cess for the soybean model. All material parameters of the soybean were assumed to be independent of temperature and moisture content. The distribution and gradients of tem- perature and moisture deve10ped inside the kernel during the drying process were established. The simulated drying curve for the soybean model was obtained and compared favorably with the experimental results reported in the literature. The calculated temperature and moisture gradients were used in a finite element analysis of the thermo-hydro viscoelastic boundary value problem to simulate the stresses in soybeans under the thermo-hydro loads. The kernel was assumed to have no initial stress and free ex- pansion and contraction of the kernel was allowed. Due to lack of information on viscoelastic prOperties of the soy- bean, the simulated results reduced to an elastic stress analysis. Tangential stress was shown to change from com- pressive to tensile stress as it approaches the surface. It reaches its peak value at the surface in one hour and then decays slowly. The effect of different drying con- ditions on the stresses was studied and the results were discussed. It was found that the magnitude of peak tan- gential stress is directly proportional to the step in— creases in temperature and the time to reach its peak value is independent of the change in drying temperature. This study concludes with suggestions regarding the need for additional investigations on viscoelastic pro- perties of the soybean and the master curves to represent these properties over the entire drying period. This information will improve the soybean model and allow the L‘Mih _‘ _ ____,_ , .._._ establishment of specific guidelines for crack free drying of soybeans. Approved: Major Professor Department Chairman To Iranian Workers And Peasants ii ACKNOWLEDGMENTS The author sincerely appreciates the guidance and close c00peration of his major professor, Dr. Larry J. Segerlind. Equally, his gratitude is extended to Dr. W.A. Bradley, Dr. G.E. Mase, Dr. G.E. Merva, and Dr. C. Cress for their time and counseling. Thanks are due to his wife Atossa and his parents for their patience and encouragement throughout this study. The author is particularly indebted to the Iranian masses from which, through a scholarship, this graduate program was made possible. iii TABLE OF CONTENTS LIST OF FIGURES LIST OF TABLES l. 2. INTRODUCTION REVIEW OF LITERATURE 2.1 2.2 Development of Stress Cracks Process of Moisture Movement Properties of the Soybean Kernel Structure of the Soybean Kernel 2.4.1 Description 2.4.2 Morphology and Anatomy Viscoelastic Characterization as Affected by Temperature and Moisture Gradients . SIMULTANEOUS MOISTURE AND HEAT DIFFUSION 3.1 3.2 3.3 3.4 Governing Equations Initial and Boundary Conditions Finite Element Formulation of Coupled Heat and Mass Diffusion Finite Difference Solution in the Time Domain Computer Implementation and the Soy- bean Model 3.5.1 The Soybean Model 3.5.2 Mass Average Temperature and Moisture . . iv Page vii 13 16 16 18 20 .24 .24 .25 .26 .37 4o 41 44 fiH_Al_- p - ---.~..__ ~._ - 3.6 Distribution of Moisture and Tem- perature and Drying Curves 4. SIMULATION OF THERMO-HYDRO STRESSES 4.1 Equations of the Linear Theory of Viscoelasticity . . . . . 4.1.1 Viscoelastic Constitutive Equations 4.1.2 Laplace Transform of Viscoelastic Equations 4.1.3 Viscoelastic Boundary Value Problem 4.2 Thermo-Hydro Viscoelastic Constitutive Equations . . . . . . . . . . . .. 4.2.1 Thermo- Hydro Rheologically Simple Materials and Shift Factors 4.2.2 Constant Temperature and Moisture States 4.2.3 Non-Constant Temperature and Moisture States 4.3 Finite Element Formulation in Thermo- Hydro Viscoelasticity 4.3.1 A Variational Theorem . 4.3.2 Radially Symmetric Viscoelastic Solids . . . . . . . . .. 4.3.3 Discretization of a Region. A. Nodal Displacements B. Element Stresses 5. RESULTS AND DISCUSSION 5.1 Discussion on the Failure Criteria 5.2 Analysis of the Results .46 .51 51 51 .54 .56 .57 .58 61 68 73 .73 .74 .79 .86 .88 .90 93 6. SUMMARY AND CONCLUSIONS 7. SUGGESTIONS FOR FUTURE RESEARCH LIST OF REFERENCES APPENDIX A vi Figure 3. 3. l 2 LIST OF FIGURES Page Simplex Triangular Axisymmetric Elements 31 Numerical Determination of the First Derivatives of {M} and {T}. . . . . . . 38 The Two-Dimensional Axisymmetric Finite Element Grid with Simples Triangular Elements . . . . . . . . . . . . . . . . . 42 Drying Curve for the Soybean Model Compared with Experimental Results . 47 Simulated Mass Average Temperature of the Soybean Model vs. Time . . . . . . . . . . 48 Moisture Content Distribution of the Soy- bean Model at Different Radii for the Whole Drying Period (2= 0) 49 Distribution of the Soybean Moisture Con- tent at Different Radii for the First Four Hours of Drying . . . . . . . . . . . . .50 Temperature Dependence of Relaxation Function . . . . . . . . . . . . 59 Stresses in Radially Symmetric Problems . .75 Radially Symmetric Line Elements and Corresponding Concentric Shells . . . . . .80 Simulated Maximum Shear Stress at Different Radii Inside the Soybean Model . 94 Simulated Maximum Shear Stress at the Surface . . . . . . . . . . . . . . . . . 95 Simulated Radial and Tangential Stresses at the Time Tangential Stress Reaches its Maximum (One Hour) . . . . . . . . . . . . 97 Figure Page 5.4 Simulated Tangential Stress at the Surface vs. Time . . . . . . . . . . . . . 98 5.5 Simulated Tangential Stress at the Surface as Affected by Drying Temperature 100 viii LIST OF TABLES Page Properties of the soybean kernel . . . . . 17 Values of the different material pro- perties of the soybean used in modeling the kernel . . . . . . . . . . . . . . . . 44 Mean values of temperature and moisture shift factors for soybean (Herum, et al., 1973) . . . . . . . . . . . . . . . . . . 66 Rheological properties of the soybean used in stress analysis . . . . . . . . . . . . 90 ix 1. INTRODUCTION The field harvesting of soybeans at a higher than Optimal moisture content necessitates the artificial drying of this product. During drying, temperature and moisture gradients are superimposed within the soybean kernel, both resulting in volumetric change. Expansion occurs in heat- ing and shrinkage or contraction occurs because of mois— ture removal. The stress states arising throughout the thermal and hydro-loaded material often results in the initiation of cracks. Stress crack phenomena is one of the most important problems in grain drying. The cracks lead to splits, decreased storage life and lower germination. Thompson and Foster (1963) reported that stress cracks develop on heating or cooling of corn kernel when the thermal and hydro stresses exceed the failure strength of the material. Overhults et a1. (1973) noted severe physical damage to soybeans during their drying tests. They believed these cracks occured because the seed coat shrank before the endosperm had begun to dry. Milner and Shellenberger(1953) has studied the formation of stress cracks in wheat through the use of radiography. High moisture content or high drying temperature gave visible fissures. Soybeans initially damaged during harvest are likely to be more susceptible to damage from subsequent handling 1 and processing. Mechanically damaged beans, including those which are split or have seed coat checks, are of less commercial value because of reduced canning and cook- ing quality. Mechanical damage not only affects market value of the bean but also impairs germination vigor. Increased concern has been shown by bean growers, processors, and shippers, to minimize the mechanical dam- age to beans caused by artificial drying and impact load- ing during harvesting and handling. The stress crack phenomenon is a function of rheoloi- cal prOperties of the material like bulk, shear and tension moduli, Poisson's ratio, failure stress and time-tempera- ture and time-moisture shift factors. In order to apply the theories of engineering mechanics to the solution of stress field equations, it is necessary to know the rheological properties of the stressed material as affected by time, temperature and moisture content. It is apparent that physical properties and more impor- tantly viscoelastic properties of the soybean kernel are required with respect to the important variables such as temperature and moisture content. There is a physical basis for assuming that the soybean kernel behaves in a similar manner as other viscoelastic materials and that similar mathematical techniques may be employed to define the be- havior of the soybean kernel under the influence of external loads. The basic unit of a kernel is the living cell. The cell consists of wall and protoplasm. The cell walls are composed of cellulose microfibrils embedded in an amorphous matrix. The cell walls exhibit a high degree of elasticity. Inside the cell, protoplasm consists of cytOplasm, nucleus and vacuoles. The cytoplasm shows both elastic and vis- cous properties. Vacuoles are made up of droplets of solution called cell sap. While the elasticity of the cell walls is the main factor responsible for the elasticity of the kernel, the cell sap is responsible for exerting a hydrostatic pressure called turgor pressure on the cell walls. The combined effects with the elastic cell walls determine the viscoelastic properties of soybean kernel. Based on experimental evidences due to Zoerb and Hall (1960), Mohsenin (1968), Timbers (1964) and Morrow (1965), agricultural products are viscoelastic. The soybean is processed in large quantities. Besides its oil which is used for different kinds of oil products and protein rich livestock feed, the soybean has been used in the production of synthetic food in recent years. Stress cracks, in fact, account for increased break- age during storage, handling and processing. Such damage contributes to susceptibility to molds and insect damage and to the production of low quality and high cost soy- bean grains and soybean products. It is hardly possible to give a dollar value to the damage resulting from thermo-hydro stress cracks. Martin and Stephens (1976) reported that corn initially consisting of 2% breakage had 15.7% breakage after 21 handling Opera- tions. If the stress state arising from both temperature and moisture gradients could be determined, and if the nec- essary material prOperties are known or could be evaluated, a failure criterion could be established which would de- fine the critical temperature and moisture gradients. This would enable drying facilities to adjust ambient tempera- tures and relative humidities to soybean moisture content in order to minimize stress crack formation in soybean kernel. The objective of this work was to develop a numerical technique that could be used to analyze the stress crack formation in the soybean kernel resulting from temperature and moisture gradients during the drying process. The specific objectives were: (1) To numerically study the phenomenon of coupled moisture and heat diffusion within the soybean kernel by use of finite element method. (2) To determine the stress state due to the tempera- ture and moisture gradients inside the soybean kernel, using computer simulation. The basic assumption of isotropy, continuity and homo- geneity which is a good macrosc0pic approximation was made in order to solve the theoretical stress analysis problem. 2. REVIEW OF LITERATURE The cracking phenomena is of great concern to soybean industry. Cracked kernels are objectionable because they are quite susceptible to breakage during handling and cause problems in storage, manufacturing processes, grading and shipping. Over much of the moisture content range encountered in grain drying, a kernel will shrink as heat is applied and water is removed. In this process, a moisture gradient and a temperature gradient are established with the outer- most layers of the kernel. The outer layer thus tends to skrink more causing tensile stresses tangent to the surface. These stresses in the outer layers increase as the moisture and temperature gradients increase. Kernel breakage re- sults if the magnitude of the moisture and temperature gradients are such that the outer layers are stressed beyound their ultimate strength (Hammerle, 1970). 2.1 Deve10pment of Stress Cracks Under favorable weather conditions, soybeans dry in the field to moisture levels safe for storage. Soybeans require lower moisture content than corn or other cereal grains to store safely under similar conditions (Hall, 1961). However, the soybean stalks and pods become brittle at a moisture content of 11 per cent, resulting in high harvest losses. Byg and Johnson (1970) reported 8 per cent shat- tering losses when soybeans were harvestedenza moisture content of 10 percent as compared with only a 2 percent? loss when harvested at 17 percent. Harvesting soybeans at moisture contents much above 11 percent would require conditioning the beans to a lower moisture for safe stor- age (Alam and Shove, 1972). The actual mechanisms for mechanical breakdown of seed coat and loss of seed viability during drying are com- plicated and frequently involve the previous history of the seed including any mechanical damage and disease infection. Seed coat cracking induced by drying is associated with the drying rate imposed upon the seed and with the corresponding rate of shrinkage as well as the level of moisture content. Thus the high percentage of cracking associated with the conventional drying method is essentially due to the high percentage of the overdried seed. Sabbah, et a1. (1976) achieved a great reduction in percentage of the overdried seeds by reversing the air flow. The use of the reversed-direction drying method also resulted in more uniform final bed moisture contents, with less seed coat cracking and germination loss than with the conventional (one- direction-air-flow) method. Overhults et al. (1973) determined thin-layer drying 7 characteristics of soybeans as affected by harvest moisture content and drying air conditions. In a study by Ekstrom (1965) certain physical pro- perties of corn kernels related to stress cracking were investigated and an attempt was made to determine whether or not stress cracks can be caused by temperature gradients alone. Zoerb (1958) concluded that moisture content has the greatest influence on the mechanical prOperties of grain and that all of the strength properties generally decrease in magnitude as the moisture increases. He further con- cluded that initial rate of deformation has more effect on the rate of stress relaxation than moisture content or the initial amount of deformation. Goncharova (1962) studied the dependence of structural strength characteristics of grain on humidity, temperature, and the speed of loading. His results were in agreement with the conclusions from the the work of Zoerb (1958) and Chizhikov (1960). Extensive studies of the causes of stress cracks and breakage in artifically dried corn were made by Thompson and Foster (1963), who found that drying speed, expressed in percentage points of moisture removed per hour, was the most significant factor in stress crack development. The amount of drying, as well as the rate, also appeared to affect stress crack development. According to Chizhikov (1960) stress crack formation in corn can be blamed on the compact construction of the shells of the kernels. This structure makes evaporation difficult since the main channel through which moisture can leave during the drying process is located around the germ. He further observed that when the drying temperature was increased, the structural strength decreased. Henderson (1954) observed that the cracking in rice standing in the field began at the center of the kernel and progressed toward the minor circumference. He concluded that stress cracks resulted from an increase in either temperature or moisture. Kunze and Hall (1965) investigated the effects of relative humidity and moisture gradient on stress cracking in brown rice. They found that a thermal gradient of 17°C did not produce fissures in rice as long as the grain were maintained at a constant moisture content. In a study of the penetration of wheat grains by water, Grosh and Milner (1959) concluded that the forces which produced cracking as a result of tempering might have two causes: a residual stress set up within the wheat endosperm during the maturing stage of kernel development, and a gradient of swelling forces produced when the moisture is absorbed into the wheat kernel. Ekstrom et al. (1966) and Mannapperuma (1975) con- cluded that stress cracks in shelled corn and brown rice, respectively, probably are not caused exclusively by tem- perature gradients alone in the kernel. Therefore moisture gradients or a combination of moisture gradients and temp- erature gradients appear to produce greater stresses than those caused by thermal gradients alone. Arora et al. (1973) also found the temperature gradient in the rice kernel is dependent upon the moisture gradient. Milner and Shellenberger (1953) detected an increased number of fissures in weathered wheat under increased initial moisture content and elevated temperature. The cracking was primarily due to stresses arising from unequal moisture gradients. Wang (1956) used intermittent application of dry air in a test to dry pea beans. He was interested in seed coat cracking as well as splitting. Perry (1959) noted in examining the checked beans that the cracks ”seemed to radiate from the hilum" creating a common check pattern. The analysis of navy bean seed coat strength and maxi- mum shear stress acting on the bean was studied by Hoki (1973) based upon the thin ring theory and Hertz's contact theory, respectively. He concluded that Young's moduli and ultimate stresses increased with decreasing moisture content for both the seed coat and cotyledon. 2.2 Process of Moisture Movement Sarvacos and Charm (1962) and Chirife (1971) have studies the diffusion process in fruits and vegetables and 10 have indicated that the linear diffusion equation does not apply over the entire period of drying. Chen and Johnson (1969) gave a mathematical analysis of moisture diffusion assuming thermal diffusion to be negligible. The removal of moisture from grains has been the sub- ject of considerable research. However, most of this re- search has been devoted to a study of bulk or batch drying in which only the average effect on a relatively large quantity of material has been considered. In order to understand better the effects of certain drying and curing practices, it is important that more emphasis be placed on moisture movement within individual kernels. There- fore, the drying rate of the whole bed, as well as that of each layer, depends upon the accurate description of mois- ture movement from a single kernel which is fully exposed to drying air of a changing temperature and relative humidity. Moisture movement has been described reasonably well with the diffusion equation when the kernel temperature does not change substantially. Becker and Sallans (1955), Hustrulid and Flikke (1959), Henderson and Perry (1967), and others treated a grain kernel as an isotrOpic Sphere, symmetrical with respect to its center and used the diffu- sion equation to describe moisture movement, without regard to the kernel temperature. Chittenden and Hustrulid (1966) and many othersfbund that the diffusion coefficient varied 11 with the initial moisture content of corn and concluded that diffusion coefficient must depend on moisture con- centration. Whitaker, et al., (1969) solved the equation for a non-homogeneous sphere with a moisture and position dependent diffusion coefficient and a non-uniform initial moisture concentration drying under changing air conditions. They concluded that the diffusion coefficient was dependent more on temperature than on moisture and that the moisture diffusion equation alone is inadequate for describing the process completely. Wang and Hall (1961) used the diffusion equation to characterize heat and moisture transfer in a corn kernel assuming constant boundary conditions. They pointed out that the effect of temperature changes due to moisture vaporization within the kernel is important and has a pro- found influence on the rate of moisture diffusion. Whitaker and Young (1972) evaluated the diffusion equation, characterizing moisture movement in a homogeneous body, as a model to describe the drying rate of peanut kernel and further studies are needed to evaluate the diffusion equation characterizing moisture movement in a composite body as a method of predicting the drying rate of peanut pods. Young and Whitaker (1971) listed the solutions of the diffusion equation characterizing moisture transfer in plane sheets, infinite cylinders, finite cylinders and spheres, for the constant boundary conditions. 12 The influence of a temperature gradient as a moisture- driving potential was investigated by Whitney and Porterfield (1968) as a means of controlling the moisture gradient during the drying process. The investigations revealed a net moisture transfer in the direction of de- creasing temperature. They further concluded that after a given amount of water has been removed, the moisture gradients are similar whether or not a temperature gradient (resulting from internal heating) exists in the same direc- tion as the moisture gradient and this can be used to de- crease the total drying time, without increasing the shrinkage stresses in the solid due to the moisture grad- ient. Absorption or evolution of water by solid results in evolution or absorption of heat, respectively. This heat diffuses through the solid, causing changes in temperature, which affects the ability of solid to absorb or evolve water (Whitney and Porterfield, 1968). Thus the transfer of moisture and heat are coupled together and in general, should be considered simultaneously. Young (1969) described a mathematical model for a drying porous sphere using the diffusion equation for both moisture and heat transfer assuming that moisture diffu- sivity is linearly function of moisture. He defined a modified Lewis number and suggested that the moisture diffu- sion equation alone is sufficient if the number is greater than 60 (negligible temperature gradient). 13 Singh, et a1. (1972) solved the simultaneous heat and moisture diffusion equations under continuously changing boundary conditions. This study improved the kernel model by using a more accurate heat conduction equation, but did not account for the temperature dependence of the moisture diffusion coefficient. 2.3 Properties of the Soybean Kernel Swelling and shrinkage are natural occurring phenomena as biomaterials are heated, cooled, or absorb and deabsorb moisture. This movement of the material results in internal strains causing formation of stress cracks which reduce the quality of the product, permit excessive loss of moisture and provide accessability for disease. The differential expansion of a composite body ex- periencing a thermal gradient has been known for a long time. The stresses which result from a thermal gradient are called thermal stresses and their calculation is discussed by several authors including Timoshenko and Goodier (1970). The mass transfer due to moisture loss (or addition) also produces stresses in composite bodies. These stresses are called hydro stresses. The combined existence of thermal and hydro (thermo-hydro) stresses requires a modification of the theory for thermal stresses as outlined by Timoshenko and Goodier (1970). This modification was performed by Hammerle (1972) who expanded the basic theory to include 14 the viscoelastic nature of agricultural products and the expansion or contraction resulting from a change in mois- ture content. The solution of this theory will give the stresses in the endosperm and seed coat or hull. An analysis of stresses would indicate what temperature and moisture conditions give rise to stresses that are likely to crack the seed coat or hull. A knowledge of failure stresses of the hull or seed coat must also be known in order to make this decision. An important contribution of Hammerle's work is the indication of which material constants must be determined before an analysis can be performed. Several are required and they include the viscoelastic modulus and the coeffi- cients of thermal and hydro expansion of each material com- ponent. If there are temperature and moisture gradients within the kernel, then the theories governing the tem- perature and moisture distributions must also be considered. These theories require a knowledge of the thermal conduct- ivity and specific heat of each material as well as the moisture diffusion coefficient. Several of the required constants have been determined but some such as the moisture diffusion coefficient and thermal pr0perties of the individual material coponents have not been determined and no specific standard or pro- cedure exists for their determination. Alam and Shove (1972) obtained a linear relationship between specific heat and moisture content of soybean which 15 is in agreement with the results obtained by Bellinger (1964) and Watts and Bilanski (1970). Jasansky and Bilanski (1973) found a linear relationship between soy- bean thermal conductivity and kernel temperature at 11.2%. moisture content (wet bulb) using transient heat flow measurement method. They further indicated that the ther- mal conductivity of soybean is dependent upon its particle size, moisture content and temperature. Their results for whole soybean thermal conductivity compare favorably with those of Watts (1967), even though his conductivity evalua- tion was calculated from measurement of density, specific heat and thermal diffusivity. Using the method of transient heat conduction, Watts and Bilanski (1973) found an average value for soybean thermal diffusivity. A small inter- action between temperature and moisture content was noted but was not considered. The coefficient of linear thermal expansion of soybean was found by Bacchus (1971). He further studied the effect of moisture on this property and concluded that this property is linearly related to kernel moisture content. Paulsen and Brusewitz (1976) determined the coefficient of thermal expansion for spanish peanut kernels (dk) and skin (as) and as drying occurs, as de- creases with moisture loss more than does ak. As further drying continues there would be a greater tendency for the kernel to expand more than the skin causing an increased stress on the skin. It also appeared from their work that peanuts and possibly other oilseed grains would have 16 expansion coefficients larger than those of grains con- taining relatively low amount of oil. The values of the different properties of soybean which are available in the literature are tabulated in Table 2.1. 2.4 Structure of the Soybean Kernel 2.4.1 Description According to the botanical classification, soybean is a member of the family Leguminosae, subfamily Papilion- oideae and its genus is Glycine (L.). Wilson (1955) men- tions the correct botanical name as Glycine Max (L.) Merrill. Classification relating to the life of the plant and its response to its environment i.e. on the basis of its growth habit describe the plant as a summer annual with a maturity requirement of approximately 75 days for the early varities to 200 days or more for the very late varities. According to the market classification, there are five market classes of soybeans based on the color of seed: yellow, green, brown, black and mixed soybeans. The chemical composition of the soybean kernel is, protein 33.98 percent, fat 16.85 percent, nitrogen free extract 28.89 percent, fiber 4.79 percent, ash 4.69 per- cent and water 10.80 percent (Saxena, 1972). 17 chauwpogaop comm esp moumuflwcw flew Afiamav manuuam - om-w.e ma.eouxww.mH-AOHxNe.om Q 938: mass» Haemav acme: - AH-m.e ma.woaxm.N- mafixmk.eN flenmfiv acme: AHAGHS msauumm - - e. fl.o ouuam m.:ommuoa meamfic acmuz” - - fl.n.e.uea-ev\e -ouXameN. Ame acumeagxm cues: co ucouuucmmou flfiamfic mseuoam oa-m~ AH-A u.\ -onflw.wm-w.mmv Axe coumaaaxm HmEhozh mo ucowofimwoou meamfiv .Ha um seesaw A com-w AN-ON .uom\Ns . .Hn-ofixwm.a-NH-Omem.m fine peouuficmoou acumsmwum Amamflc uxmcafium eca meow: - - .umm\Neu.meooo.o A.fle breasts €62: Amnmfiv uxmaaflum was xxmcamme m.me-oH N.HH eo-e\muua2 .mMH.-Amo. A.ev xufl>fluuswcou Hmsuonh naemflc mega: - mm-e me\mx.omNH flfiamfic msgouam - NH-A me\wx.eONH mac spumaon fleoafiv cmmeufifiom - - ma\mx.owHH fleemfiv Homqufifiom ww-mm HN u..mg\e .m.umo~ floemfiv axmqafium eea mega: w.ANH-N.om 4.5 eo.wx\e .o.ewwH flue new: eacuuoam AHAQHS msauuam - aa-m.e sfim suuuucezam ea an... eezwwwm oocohomom 3320?an 83 we 02 633/ 5.89on .Hochox Gdonxom exp we mofluquOpm .H.N ofinme 18 2.4.2 Morphology and Anatomy Most mature seeds are made up of three basic parts: the seed coat, the embryo, and one or more food storage structures. The soybean seed has only two parts. Its two food storage structures, the cotyledons, are part of the embryo. That is, the embryo of the soybean seed con- sists of two thick, fleshy cotyledons which account for most of the bulk and weight of the seed (about 90 percent). They contain nearly all the oil and protein in the bean. The seed coat is relatively thin and protects the embryo from fungi and bacteria before and after planting. If this protective coat is cracked, the seed has little chance of developing into a healthy seedling. A variety of distinct layers can be recognized in the seed coat and cotyledons. The outermost layer, the epi- dermis, remains uniseriate and develops into the palisade layer characteristic of leguminous seeds. The palisade cells have a definite cuticular layer (Smith and Circle 1976). Two palisade layers and a compact group of trach- eids of unknown role occur in the hilum region. The palisade layer has attracted much attention because its structure in certain hard legume seeds such as soybean is assumed to be connected with their high degree of imper- meability and thus germinability (Esau, 1960). The column cells are hourglass or I-shaped, similar to those of other beans, but in soybeans theyare thicker and longer. They 19 vary from 30 to 70 microns in length and from 16 to 36 microns in width. The spongy parenchyma consist of 6 to 8 layers of thin-walled, somewhat compressed, boxlike, empty cells ranging in width from 20 to 120 microns but mostly 40 to 60 microns. One of the most characteristic features of the soybean seed is the presence of the endo- sperm. This is represented by a layer of aleurone cells and several layers of thin-walled cells which have been crushed by compression. The aleurone cells are rather thick-walled and are filled with dense protein. The sur- face of the cotyledons are covered with a typical epidermis made up of small cubical cells filled with grains of aleu- rone. The remainder of the cotyledon is largely made up of layers of elongated palisadelike cells with thin walls and filled with aleurone and oil. The hard legume seeds such as soybean achieve and maintain a very low percent moisture which is not affected by fluctuations in moisture content of the surrounding air. The attainment of this high degree of desiccation is ascribed to a combination of intensely impremiable testa with valvular action of the hilum (Hyde 1954). The hilum is said to act like a hygrosc0pic valve. A fissure occurs along the groove of the hilum. This fissure opens when the seed is surrounded by dry air and closes when the outside air is moist. Thus the entry of moisture is prevented but loss of moisture is permitted. The occurrence of highly impermeable seed coats is one of the important factors in 20 delayed germination of seeds in Leguminosae. In this connection the occurrence of cuticular layers in seeds is of particular interest. Some good pictures of the soybean seed structure are shown by Esau (1960) and Smith and Circle (1976). 2.5 Viscoelastic Characterization as Affected by Tempera- ture and Moisture Gradients In recent years the subject of Viscoelasticity has received considerable attention from both theoreticians and experimentalists, which is evident in the reviews of Hilton (1954), Lee (1964) and Morland and Lee (1960). Within the engineering requirements of accuracy many material have been found to satisfy the assumption of thermorheologically simple behavior. Mérland and Lee (1960) have shown how this assumption can be used to derive iso- trOpic constitutive equations with transient temperatures. Hilton (1954) has investigated the thermal stresses in thick walled incompressible cylinders exhibiting temperature dependent viscoelastic prOperties of the Kelvin type. Sharma (1964) and Rosen (1964) have discussed in detail the time—temperature dependence of linear viscoelastic materials and the Boltzman superpostion theory. Verifica- tion of these principles by experimentation has been done by Ferry (1961), where the transient tests described are basically creep and relaxation for various shaped specimen. 21 The experimental determination of stress relaxation moduli from tensile stress-strain curves at various temperatures has been discussed by Sharma (1964). Hammerle and Mohsenin (1970) determined the stress relaxation modulus of corn at various temperature and moisture levels. Muki and Sternberg (1961) have worked on the theoretical analysis of thermal stress distribution within infinite and finite viscoelastic plates, and viscoelastic spheres. The time dependent yielding for the above problems was investigated by Landau et al. (1960). Bland (1960), Perry (1961) and Sharma (1964) have arrived at the relaxation integral laws where deviatoric and isotropic stresses are given as a function of strain-time loading. Utilizing the constitutive laws of a continuous Maxwell model and the superposition method, Rao (1971) arrived with the solution of the stress state in a viscoelastic cylinder subjected to transient temperature and moisture gradients. Later on Rao et al. (1975) solved the same problem for a viscoelastic sphere subjected to radially symmetric temperature and moisture gradients. They further concluded that this approach is applicable to approximately spherical seeds or formed foods. Studying the physical prOperties of small grains by a triaxial testing method, Stewart (1964) found that an interrelationship exists between their moisture content and their viscoelastic prOperties. Bacchus (1971) studied the mechanical properties of 22 soybean by assuming Hookean behavior using the Hertz con- tact theory for convex bodies. He concluded that yield point, Young's modulus and maximum compressive strength decreased with an increase in moisture content and initial rate of deformation had more effect on the rate of stress relaxation than moisture level or initial amount of de- formation. He also found that the coefficient of linear thermal expansion has a linear relation with the kernel moisture content and the assumption that soybean is a sphereholds good for all practical purposes. Silberstein and Rao (1976) evaluated the uniaxial modulus of elasticity at constant loading rates for Florruhher peanut variety. The time-temperature and time moisture shift factors were evaluated and a three element Maxwell model was found to adequately apprximate the pea- nut behavior. This information is vital for stress analy- sis of peanut kernels to determine the stress cracks and shrinkage. Herum et a1. (1973) evaluated the time depen- dent uniaxial moduli of intact soybeans, determined by re- laxation tests in parallel plate compression, at four tem- peratures and four levels of moisture content. They con- cluded that the time scale shift factors were not constant but were sufficiently similar to permit description of intact soybeans as thermo-rheologically and hydro—rheolo= gically simple. This means that viscoelastic properties such as uniaxial modulus function of the soybean can be described by one master curve. This curve allows the 23 calculation of the viscoelastic modulus over the range of temperatures and moisture contents which would be en- countered during the drying process. Further, they found that an increase or decrease of 1% in moisture content, computed on the dry basis, has approximately thirteen times as much effect on the uniaxial modulus in compression as does a 1° Fahrenheit increase or decrease in temperature. Saxena (1972) reported that the soybean kernel behaves as a viscoelastic material in response to the applied forces and is characterized by a generalized Maxwell model having as few as three Maxwell elements and one spring in parallel. He further observed that the orientation of stress cracks depends on the test position of the kernel and also there was no way to ascertain whether seed coat or the cotyledons cracked first. 3. SIMULTANEOUS MOISTURE AND HEAT DIFFUSION 3.1 Governing Equations Since the temperature distribution within the soybean kernel may not be uniform during the initial stages of drying, simultaneous equations of moisture and heat diffu- sion are needed to describe moisture movement within the soybean kernel. Crank (1964) gave the mathematical model characterizing moisture and heat diffusion, in a sphere (in cylinderical coordinates) as follows: ,a_1\_1_ 1 a 23M 32M 3t — ¥ 5? (rD 5?) + D.__? (3.1) 82 2 3T=l 3 .3T .8T M DC at F ar (TK 8r) + K 822 + LO 3 (3‘2) where O M = mositure concentration, 6 dry basis T = temperature, °C r = radial coordinate, z = vertical coordinate t = time, sec. D = moisture diffusivity, mZ/sec. K'= thermal conductivity, w/m°C or J/m sec. °C c = specific heat, J/kg°C p = density, kg/m3 L = latent heat of vaporization of water, J/kg 24 25 Parameters D, K', p and c are in general functions of local temperature and moisture values inside the kernel. No analytical solution to this complex problem has been found. Computer solutions would also be extremely hard to obtain, especially when the nature of some of these pro- perties such as D is not known. It is felt that consider— able information on the model about temperature and mois- ture distribution within the kernel could be gained by treating p, c, D and K' as constants in a finite element analysis. Therefore this study assumes that p, c, D and K' do not change substantially in the temperature and moisture ranges under consideration. Equations (3.1) and (3.2) now reduce to 3M 6...”! =13 3121.2“:“4. D ——7: (3.3) 3T _ , 2T K' 3T , 232T 3M 4 OCa—t’K 2"?5?*K 2*L05‘t (3°) 3r az 3.2 Initial and Boundary Conditions (1) Kernel is at uniform temperature and moisture initially (2) (3) 26 for t>o and at the surface (r=R) 3M 3M _ = D(§? + EZ)| _ + hm (Msurf. Meq.) 0 r—R where hm = mass transfer coefficient, m/sec. for t>o and at the surface (r=R) VET fl .. = K (3? + 82)] =R + h(Tsurf. Teq.) 0 where h = heat transfer coefficient, W/m2°C or J/mzsec.°C Finite Element Formulation of Coupled Heat and Mass Diffusion In most problems the integral of a functional is minimized. X1 = This functional possesses the property that be written as l any function which makes it a minimum also satisfies the governing differential equation and the boundary conditions. As for this case, the functional form of the governing equations (3.3) and (3.4) by use of variational calculus may 82 h + fig (Ms - MMm)2dS (3.5) S ,- 27 x2 = If [K't3333 2 + K'(33§i)3 — 2Lo (33M}){T} + 2pcc3333){r}]dv + x% (TS - Tm)2ds (3.6) S where M00 and T00 denote ambient moisture content and tem- perature respectively. V is the total volume of the body and S is the boundary. The solution of the boundary value problem stated in (3.3) and (3.4) can be obtained by finding the stationary value of the functionals x1 and x2 in (3.5) and (3.6) with respect to the set of nodal values {M} and {T}. The finite element method can be used as a numercial technique based on the minimization of these functionals. Define four matrices: T 3M 3M {g1} = [332' 9 a—Z (3.7a) {g2}T = [3r , 33] (3.7b) [D] E)(j (3.8a) o D [H 0] (3.8b) 0 K' Equations (3.5) and (3.6) can be written as [k] x1 = I— 2[{g1}T [D]{g1}+ 2(33“}) {M}]dv 3511 . ;_7 (M: - ZMSMW + M:)ds (3.9) S 28 x2 = i %[{g2}3[k1{g2} - 2Lp(3§¥l){r} + 2pc(3%%l){T}]dv mm: .2 (15 + f - zrsrm + T:)ds (3 10) S The volume and surface integrals in (3.9) and (3.10) can be expressed as a sum of integrals over a set of sub- regions or elements. E 3MOB) at {M(e)}dv = f 1{ 33)}3 D33) { 33)}d f X1 8:1 V(e)2 g1 [ ] 81 V + V(e) (e) + f E? (M§G)M(e) - 2M333M3e) + M(:)M£e))ds (e) S S w s (3.11) E . (e) (e) x2 = 2 f %{g§e)}T[k(e)]{g§e)}dV - Lpf 3§¥——{T }dv 8:1 V(e) V(e) (e) (6) 3T h ( ) (e) (e) ‘3t--{T(e)}dv + f(e)—7"-ITse Ts ' S +pcf v 2T333T333 + Tie3T£e))ds (3.12) 29 where the superscripts e indicate the elements or sub- regions and E is the total number of elements. These two equations can also be written as E (e) _ (1 2 E - Xl _ Xl ) + Xi ) + ,,, + XE ) — g=lxl (3.13) E x2 = x31) + x32) + ,,, XEE) = 2 x39) (3.14) e=1 (e) where X1 is the contribution of a single element to X1 . . (e) and s1m11arly for X2 . The minimization of X1 and X2 occurs when 8x B E 3 __l_= __a_. 2 x333 = z ——l+—= 0 (3.15) a{M} 3{M} e=1 e=1 3{M} E E (e) 3X2 3x = _2_.. 2 x33) — z 2 = o (3.16) a{T} 8{T} e=1 e=1 a{T} . (e) . C ) . . 9X . oX2 The der1vat1ves -—————- in (3.15) and can not be a{M} 8{T} evaluated until the integrals in (3.11) and (3.12) have been written in terms of nodal moisture contents and tem- peratures, {M} and {T}. The moistures and temperatures in each subregion or element are approximated by algebraic poly- nomials relating them to moisutres and temperature of nodal points of that element (Zienkiewicz, 1971) or 30 M333 = [N363] {M} (3.17) T(e) = [N(e)] {T} (3.18) [N383] is the matrix of shape functions relating moisturesand temperatures in the element, M38) and T36), to the nodal moistures and temperatures, {M} and {T}. The similarity between the axisymmetric and two dimensional problem makes the solution of the axisymmetric problem quite straight forward. An example of an axisym- metric element is given in Figure 3.1. Equations (3.17) and (3.18) can be written as (e) _ , . M - NiMi + Mij + NkMk (3.19) T33) = N T. + N T + N T (3 20) i 1 j j k k ' where N. = l— (a. + b r + c z) 1 2A 1 i i _ 1 . N. - 2A (aj + bjr + cjz) _ 1 . Nk ' 2K 33k + bkr 3 Ckz) A is total area of the triangular element and constants a, b and c are defined as Figure 3.1. 31 3 (”fir-IF ~.\ \\ \..J>‘l - \~._ ,.” J G#_"’,AV' Simplex Triangular Axisymmetric Elements 1 j k X J j k 1 1 b1 = zj - zk bj = 2k - z1 ci = rk - rj cj = r1 - rk aL - rin - r 21 bk = 21 - zj ck = rj - rl We can now evaluate (3.7a) and (3.7b) which along with (3.19) and (3.20) can be substituted in (3.11) and (3.12). aM3e) ”33133) 3N3e) angen Mi {338)} = 3r = Br Br Br M. J aM3e) 8N3€3 3N3e) 8N£e) M 32 L 32 az ‘33 . k = [333)] {M} (3.21) and similarly {gée3} = [B363] {T} (3.22) 33 Element integrals in (3.11) and (3.12) now become X1Ce) = I %{M}T[B(e)]T[D(e)][B(e)]{M}dV + I [N(e)]{M} V(e) V(e) h [N333133%3 dv + I —%IM}T[N33)]T[N33)]{M}ds - S(e) I hme[N(e)]{M}ds (3.23) S(e) (e) (e) X2(e) = I %{T}T[B(e)]T[k(e)][B3e)]{T}dV - Lpf [N383] V V {T}[N3e)]3ifl3 dV + ocf [N33)]{T}[N3e)]3il}dv + at v36) at I §{T}T[N3331T[N(e)]{r}ds - I hTw[N(e)]{T}ds S(e) 5(6) (3.24) where [B383]T is the transpose of matrix [B] and element superscripts are going to be deleted for clarity from now on. The order of integration can be changed. Performing differentiation in (3.15) and (3.16) for the purpose of minimizing x1 and x2 in (3.23) and (3.24) yields 34 3X1 -§—— = (IINJTINidv133¥}+ (lelTlDllBldV + fhm[N]T[N]ds) {M} v v s {M} - fhme[N]Tds = o (3.25) S 3X2 a{M} ——— = (oc IINJTINidV)3§§}- (Lof[N]T[N]dV)—§f 8{T} v v + (f[B]T[k][B]dV + fh[N]T[N]ds){T} - thm[N]Tds = 0 V S S (3.26) Evaluating the integrals in (3.25) and (3.26) (Segerlind, 1976) yield " 2 2 bi+ci sym. 3 ZWA 3{M} ZETD 2 2 {—60 [Cij])_§f' + (_4A_' bibj+cicj bj+cj 2 2 .bibk+cick bjbk+cjck bk+CkJ . 0 0 o + ZnthJk O 3r.+rk r +rk ){M} - Zflthjka 2 o r.+rk r +3rk 6 o o 0 r1 0 2 1 r. = 0 (3.27) o l 2 ri 217A 3{T}_ 211A 3{M} ank' (DC 60 [51)” t (LOTO‘ 3C1)“ at ‘3 (7T r N bz+cg 1 1 ZnhL. o o o b.b.+c.c. bg+cg + ————l£ 0 3rj+r rj+r§ 1 j 1 j j j 12 o rj rk rj+3 k 2 2 ubibk+cickbjbk+cjck bk+CkJ ZnhL. T o o o r. {T} - ———L3‘—2 2 1 r3 = 0 (3.28) o l 2 rk _ ri+r.+rk where r = ————l——— 3 L. jk = length between nodes j and k of the triangular element where convection takes place. and [Cij] = (l + éij)(3r + ri + rj) (BrOCC1, 1969) where aij = Kronecker delta (3.27) and (3.28) can be written in condensed form [CHE-33$} + [K11]{M} - {fl} - 0 (3.29) M}. [c211 t [62213—333 + [KZZHTl - {f2} = 0 (3.30) and in matrix form C O M K11 O M f 11 = o (3.31) 36 where [all]. [c211, [c221, [K11], [K22],{f1}and (£2) are~appropriatecoefficients in (3.27) and (3.28) and - _ 3{M} ' = 31M} {M} - —§f , {T} at Equation (3.31) can also be written as M M [C] . + [K] - {F} = 0 (3-32) T T Note that this equation is for one element and for the whole body a summation over all the elements must be carried out. We now solve (3.29) for §%%} or 3333= T—l—T ({f }-[K 1(M)) at C11 1 ll inserting this in (3.30) yields . _ 1 . [C22]{T}+ [K22]{T}-{f2}-[C21][TCIIT({fl}-[K11]{M})] = 0 After simplification and substitution of (3.29) the equations become [c111{N}+[K111{M}—{f1} = 0 (3.29) [CZZ]{T}+[K22]{T}-{f2} + Lp([Kll]{M}-{f1}) = o (3.33) Equation (3.33) indicates that the moisture values {M} are needed for every time step in order to obtain the temperature values for the same time step 37 3.4 Finite Difference Solution in the Time Domain To solve the linear differential equations (3.29) and (3.33) to obtain the values of {M} and {T} at each point in time, a central finite difference technique to approximate the time derivatives is used (Segerlind, 1976). Given two points on a curve (Figure 3.2) an approxi- mation for the first derivative at the midpoint of the increment is {M}: d{M}= 3M3i+1 ' {M31 (3.34a) at At {T}~ - {T}. ' _ d{T}_ 1+1 1 {T}- —dt - At (3.34b) Since {M} and {T} are evaluated at the midpoint of the time increment, we should also evaluate {M}, {T}, {f1} and {£2} at this point. So we have {N}* = 3 ({M}i+1 +{M}i) (3.35a) * 1 {T} = 7 ({T}1+1 +{T}i) (3.35b) it _ 1 - {f1} - 7 ({f1}i+1 + (£111) - {£1} (3.36a) * 1 .- {f2} = 7 ({fz}i+1 + {f2}i) - {£2} (3.36b) 38 MM. Ti+1 M, T = f(tl .— At.._+ )- t ti ti+1 Figure 3.2. Numerical Determination of the First Derivatives of {M} and {T} 39 since {fl} and {f2} are constant and do not vary with time. Substitution of (3.34a), (3.35a) and (3.36a) into the differential equation (3.29) gives ([Clll + A%[K111){M}i.1 = ([c111-451K11])(M}i+ At {f1} (3.37) Given the nodal moisture values at time t, equation (3.37) can be solved to yield the nodal moisture values at the time t+At. The column vector {fl} consists of known para- meters; hence, its evaluation at t and t+At can be carried out before (3.37) is solved. Let us rewrite equation (3.33) as [c221{i} + [K22]{T} - {£5} = 0 (3.38) where {fi} = {f2} - Lp([Klll{M} - {f1}) (3.39) using (3.34b), (3.35b) and (3.36b) equation (3.38) may be written as ([sz] + A% [K221){T}i+1 = ([sz] - A§3K22]{T}i + At 332}* (3.40) where f' * - 3 ( f' + f' ) (3 41) 3 2} ' 2 3 23i+1 3 231 . {32}i+1 ‘ {£2} - Lo([K11]{M}i,1 - (£11) (3.42a) 40 and 3f231 = {f2} ' Lo([K11]{M}i - {f1}) (3.42b) substituting (3.41), (3.42a) and (3.42b) in (3.40) yield ([szl + A; [K22]){T}i+1 =333322] ' A”323K22333T3i 3 A‘EIHfZ} - Lp([K11]{M}i+1 ‘ 3313)) 3 3332} ’ 143333311] {M}i - {fl}))] Simplifications gives ([CZZ] + 4%[K22]){T}i+1= ([szl - 9%IKZZ]){T}1 + - %9[K111({M}i+1 + {M}i)] (3.43) AtIIfZ} + LDIfl} This equation along with equation (3.37) defines the process of coupled heat and mass diffusion inside the soy- bean kernel. 3.5 Computer Implementation and the Soybean Model A finite element method for the solution of coupled moisture and heat diffusion within a spherical body was presented in the previous sections. A computer program for two dimensional transient field problem such as the one described by equation (3.37) was written by Segerlind (1976). 41 This program was modified to incorporate the coupling effect of moisture content and temperature for each time step. Modification also included a change to axially symmertric triangular elements. The modified program first solves equation (3.37) for a given initial nodal values. For every time step At and given {M}i values, a set of nodal moisture values {M}i+1 will be obtained and stored on a tape. Then the last term of equation (3.43) was evaluated and stored. Fin- ally equation (3.43) was solved for a given initial values of {T}i to yield the new nodal temperature values {T}i+l° Ninety six iterations were calculated using a time step At, of 10 minutes. Total drying time was 16 hours. Nodal moisture and temperature values were printed for one hour intervals for 16 hours. 3.5.1 The Soybean Model The following assumptions for the soybean model were made: (1) Soybean kernel was taken as an isotrpic sphere symmetrical with respect to its center. A two dimensional axisymmetric finite element grid as is shown in Figure 3.3 was used. This was due to the fact that no knowledge of temperature varia- tion within the soybean kernel was available. 42 FINITE ELEMENT GRID 255 ELEMENTS 152 NODES 8 ch .00 o . 08 or. 16 03.24 0232 0140 R-HXIS CM Figure 3.3 The Two-Dimensional Axisymmetric Finite Element Grid with Simplex Triangular Eleme nts (Z) (3) 43 Some researchers such as Whitaker et a1. (1969), Misra and Young (1978) assumed a uniform tempera- ture variation, solved only the moisture diffu- sion equation and used a one dimensional model. The finite element grid that was used in this model consists of 255 elements and 152 nodes. The thickness of the soybean skin is incorporated in the grid and the computer program has the option of having different material prOperties for the skin and the cotyledones. Moisture removal and heat intake take place by convection process through the whole surface area of the kernel skin as was shown in section 3.2. It is noteworthy to mention that some unsuccess- ful attempts were made to incorporate the effect of the hilum in the process of moisture removal. Failure was due to the fact that there was no way to find what percentage of moisture was re- moved from the kernel through the hilum region even though some preliminary tests proved that more moisture was removed through the hilum than through the skin. All material parameters of the soybean kernel were assumed to be constant for two reasons: (1) it was felt that considerable information on the model about temperature and moisture 44 distribution within the kernel could be gained. (2) There was a lack of information on the material properties of the soybean as most of the properties which where needed for this model are functions of temperature, moisture, time, etc. Some of these properties as are available in the literature are tabulated in Table 2.1. The values of the different material properties which were used in this model are presented in Table 3.1. Table 3.1 Values of the different material properties of the soybean used in modeling the kernel Material PrOperty Value Dimension Diffusion coefficient (D) 7.0)(10'11 mZ/sec. Thermal conductivity (K') 0.1 w/m-°C Specific heat (c) 2000.0 J/kg-°C Density (p) 1200.0 kg/m3 Heat convection coefficient (h) 60.0 w/m2-°C , Mass convection coefficient (hm) 0.05 m/sec. 3.5.2 Mass Average Moisture and Temperature In order to have a measure of the kernel moisture content and temperature as a whole, rather than the nodal 45 moisture and temperature values, the concept of mass average moisture and mass average temperature was used. Mass average moisture and temperature of a body are de- fined by (Segerlind, 1976) f M(r,z)dm v for every time step (3.44) fidm fT(r,z)dm v T = for every time step (3.45) fdm v where dm is an element of mass. Using (3.19) and (3.20) and after some simplification (3.44) and (3.45) yield E (e) NA , 2:1 ——6_—'332ri+rj+rk3Mi+3ri+2rj+rk3Mj+3ri+rj+zrk3Mk3 N: E (e) ZnA Z (r.+r.+r ) e=1 3 1 J k (3.46) and similarly E nA(e) 2 g=l_—6__332ri+rj+rk3Ti+3ri+2rj+rk3Tj+3ri+rj+ rk)Tk] T: E (e) )3 -2—TL131———(r.+r.+rk) e=1 1 3 (3.47) 46 Both equations are applicable only for the axisymmetric triangular simplex element. 3.6 Distribution of Moisture and Temperature within the "Kernel and Drying Curves Based on the results for M, a drying curve was obtained and compared with the experimental results from two diffe- rent sources, as is shown in Figure 3.4. The simulated drying curve is between the two experimental drying curves and agress very closely with the results given by Overhults et a1. (1973). Figure 3.5 shows the variation of the ker- nel mass average temperature vs. drying time. It is obvious from the figure that there is very sharp increase in simulated mass average temperature of the model in the first hour of drying. Mass average temperature reaches the equilibrium temperature (drying air temperature) in about one hour and then stays constant for the rest of the drying period. Distribution of the nodal moisture content of the soybean model at different radii for the whole drying per- iod is shown in Figure 3.6. It is apparent from the figure that nodes which are located on the skin reach the quili- brium moisture content much faster than the internal nodes. It was also felt that there is a need to monitor the nodal moisture distribution for the first few hours of drying since according to several researchers such as Misra (1977) this is the period in which stress cracking occurs. This is shown in Figure 3.7. MASS AVERAGE MOISTURE CONTENT (”lo D. B.) 47 34 l-SABBAH et al., 1976 2-OVERHULTS et al., 1973 3-Fl N ITE ELEMENT SOLUTI ON 30 26- 22‘ 18' 14‘ 10 .4 1 03 .l I l I I fi I o 2 4 6 811) 121416 DRYING TlME (HOURS) Figure 3.4. Drying Curve for the Soybean Model Compared with Experimental Results MASS AVERAGE TEMPERATURE (°C) 48 30 " 25 ‘ 20 0 I l r I T T I T 0 2 4 6 8 10 12 14 16 Figure 3.5. DRYING TIME (HOURS) Simulated Mass Average Temperature of the Soybean Model vs. Time MOISTURE CONTENT ("lo D. B.) 49 33 30. a-RADIUS OF SOYBEAN 27 ? R") R: 24 alZ 21 3 18 ‘ a BMW 15 " 12 - R=a 0 T T 1 T T I r T 024'6810121416 DRYING TIME (HOURS) Figure 3.6. Moisture Content Distribution of the Soybean Model at Different Radii for the Whole Drying Period (z=0) MOISTURE CONTENT (% D. B.) 50 a-RADIUS OF SOYBEAN 33 30 ‘ R-O R=a/4 26 - R=a/2 22 . 18 . R-3I4Ia) 14 ~ R=a 10 3 0 W T T F 0 1 2 3 4 DRYING TIME (HOURS) Figure 3.7. Distribution of the Soybean Moisture Con- tent at Different Radii for the First Four Hours of Drying 4. SIMULATION OF THERMO-HYDRO STRESSES It is assumed that the total stress at a point is some combination of the thermal and hydro stresses superimposed upon whatever stress may exist due to mechanical loading. The thermal and hydro stresses may arise from either tension or compression, and it might be expected that since both high temperature and high moisture content cause swelling or expansion of the soybean kernel and since tem- perature and moisture gradients are of opposite sign in a heated material undergoing drying, a very complicated tension-compression state may occure within the kernel. This complicated stress state, when it exceeds some critical state at a point, results in internal failure. 4.1 Equations of the Linear Theory of Viscoelasticity 4.1.1 Viscoelastic Constitutive Equations The theory of Viscoelasticity is adequately described by several authors such as Flfigge (1975), Christensen (1971), and Ferry (1961). There are several equivalent forms of the constitutive equations of viscoelastic materials: hereditary integral forms, differential operator forms, and 51 52 complex modulus form. The integral form of the constitu- tive equations can be written as (Christensen, 1971)1 t de.. - 1 (T) Sij - I” G1(t-T) _—a%—_— dr (4.1) t dc okk - I G2(t-t) k3333d (4.2) where G1(t) is a relaxation function appropriate to states of shear and Gz(t) is a relaxation function appropriate to states of dilatation. The deviatoric stress and strain tensors are n— - 1 _ Sij ‘ Oij 3 3ij3kk ’ Sii ‘ 0 (4'3) and - - l - respectively with Oij Eij = strain tensor 6.. = Kronecker delta, zero for ifj and 13 = = 3ii 31133223333 3 Gkk = first invariant of the stress tensor stress tensor Ekk = first invariant of the strain tensor In order to use the same notation as in elasticity, the relaxation functions in simple shear and dilatation are taken as 1The standard indicial system for a rectangular Cartesian re- ference frame is employed whenever applicable: Repeating the subscripts i, j, k or 1 implies summation, Kronecker's delta is denoted by 8 differentiation with respect to space is indicated by1335ubscripts preceded by a comma. 53 Glct) G(t) = 2 (4.5) czct) K(t) = T (4.6) These relaxation functions are equivalent to the elastic shear and bulk moduli, respectively. An alternate form of the stress-strain relation is obtained by using creep functions to represent the current strain as determined by current value and past history of stress (Christensen, 1971) CD I I J (353.33) d 4 7 1].- 1(t—t) —3TL—— T (. ) dokk(t) t Ckk = f JZCt'T) T dT (4.8) -oo where J1(t) and J2(t) are creep functions for states of shear and dilatation. They can be related to the relaxation functions by use of Laplace Transforms or other interconver- sion techniques. Shear and bulk modulus are two indepen- dent constants characterizing a homogeneous elastic solid, and the relations between these and more commonly used engineering parameters like Young's modulus and Poisson's ratio have been established. Similar relations exist between the Laplace Transforms of the viscoelastic 54 relaxation functions (Christensen, 1971). A widely used form of the relaxation function is an exponential series representation known as the generalized Maxwell Model G.e't/Tj (4.9) n G(t) = z J j: 0 where Gj and Tj are the shear moduli and relaxation times for each element in the model. In a similar manner the time dependent bulk modulus K(t) could be defined for the dilatational response (or volume change). Creep functions are sometimes approximated using elastic elements and viscous elements in parallel, as represented by the series (Gradowczyk and Moavenzadeh, 1969) m J(t) = z Jj(1-e‘t/'j) (4.10) 1 4.1.2 Laplace Transform of Viscoelastic Equations The convolution integral form of the viscoelastic constitutive equations was given in the previous section, equations (4.l)-(4.8). The relationship between the different relaxation and creep functions can be established by introducing the Laplace transformation. Let f(t) be a continuous function over o0 (4.29) f(T) is measured relative to some arbitrary temperature To' Since the rate of change is increased with increase of temperature, the modulus curve will shift towards shorter times with increase of temperature, as illustrated in Figure 4.1. We introduce a change of variable in shift function by setting f(T) = in aT(T) (4.30) or in another way the "temperature shift factor" aT(T) is 63 aT(T) = eprfITII (4.31) Relation (4.29) now implies that d T aT(TO)=l , aT(T)> 0 , _EI£_1 > 0 (4.32) dT The temperature shift factor is therefore a positive, monotonically increasing function of T throughout the range of validity of equation (4.28) (Muki and Sternberg, 1961). Substituting equation (4.30) into (4.28), the relaxa- tion modulus can be expressed as ga(t,T) La[2nt + 2n aT(T)] La2n[t-aT(T)] (4.33) Now, recalling (4.27), (4.33) is written as ga(t,T) = G t°aT(T)] = Ga(€) (4.34) o3 provided the ”reduced time" is defined by a: t.aT(T) (4.35) Thus, the relaxation function ga(t,T) at any temperature can be obtained directly from the relaxation function Ga3t) at base temperature TO by replacing t with g from (4.35) or once the temperature shift function f(T) is 64 known for the temperature range considered. The tempera- ture shift function f(T), and hence the temperature shift factor aT(T), represents an inherent property of the visco- elastic material. Christensen (1971) describes the experi- mental technique for determining the shift functions and hence, the shift factors. The viscoelastic mechanical property when plotted versus the logarithm of time can be superimposed to form a single curve merely by shifting the various curves at different temperatures along the logarithm of time axis. If the curves do coincide within experimental error, the basic postulate of thermorheologically simple material is verified. He further discusses that there are no general inclusive guidelines that can be given to answer the question as to whether a given material can be expected to exhibit the thermo-rheologically simple type of be- havior. The only safe and certain answer lies in experi- mentally verifying or invalidating the shifting procedure for every material of interest. Based on temperature shift factor, an analogy is es- tablished for the moisture shift factor, aM. It is also defined as aM(M) = exp[f(M)] (4.36) where f(M) = moisture shift function 65 From the development of a time-temperature shift factor the relaxation modulus can now be expressed by 0, (t.T,M) = f[Ga(t). aT(T). aM(M)1 (4.37) where the concept of a time-moisture shift factor is in— cluded. Equation (4.37) actually means that, by knowning any material property (here relaxation modulus) at the base temperature and moisture and also knowning temperature and moisture shift factors, one can get the relaxation modulus at any constant temperature and moisture content. The thermo-hydrorheologically simple postulate is also sometimes referred to as the time-temperature-mois— ture superpostiion principle, or the method or reduced variables. Saxena (1972) experimentally verified that the soy— bean grain is both thermo and hydrorheologically simple. Herum, et a1. (1973) found the temperature and moisture shift factors as 2n aT(T) -.049(T-32.2) (4.38) 2n aM(M) 1.61-.63(M-l6.5%) (4.39) where base temperature and moisture content were assumed to be 32.2°C and 16.5% (dry basis) respectively. His results for the mean values of aT and aM are tabulated in Table 4.1. . 1:321 .- -, 66 Table 4.1 Mean values of temperature and moisture shift factors for soybean (Herum, et al., 1973). Temperature Mean aT Moisture content (%d.b.) Mean aM O 21.1 88 11 500 32.2 1.0 16.5 1.0 43.3 .38 21.3 .045 54.4 .14 30.2 .0017 Comparing equations (4.30) and (4.38) we see that for soy— bean f(T) = -2n aT(T) OI‘ aT(T) = exp—[f(T)] this would change the equation (4.33) to ga(t.T) = t La[£nt - 2n aT(T)] = La2n[§ETT)] and relation (4.34) becomes ga(t,T) = GaigifT)1 = Gq(€) 67 where now g: t aTiT) The "shifted time" (which relates the physical time and shift factors) incorporting and summing the time-tempera- ture and time-moisture shift factors (aT and aM) which are superimposed can be written as (Rao, et al., 1975) g: g (4.40) Equations (4.37) and (4.40) enable us to pass from equations (4.1) and (4.2), which hold at the base tempera- ture and moisture to the corresponding relaxation integral law applicable at any constant temperature and moisture. This transition is evidently effected by replacing Ga(t‘T) in (4.1) and (4.2) with Ga(g-gf),where g is given by g=§ and gkg. provided the body (in the absence of loads) is considered to be in the unstrained state at the uniform tem- perature T. The preceeding section described the equivalence re- lation for the body held at different uniform temperatures and moistures, and this must now be extended to cover the case of a general temperature and moisture fields (also 68 dependent on space coordinates), T(r,t) and M(r,t). First consider the temperature and moisture dependent only on the space coordinates, T(r) and M(r), a steady state temperature and moisture fields, then assuming that the shift law applies to each particle independently, the re- lation (4.35) holds for each particle with a pseudo-time (reduced time) which depends on its position. The stress- strain law again has previous form (Ga(t—T) replaced by Ga(g-g')in equation (4.1) and (4.2) where now g=§). It should be noted that gkk(r,g) and similarly for gkk(r,g), means the mapping of the function Okk(r’t) from (r,t) plane into the (r,€) plane, and not the same function with 5 replacing t (Morland and Lee, 1960). 4.2.3 Nonconstant Temperature and Moisture States It is of interest to extend the results just described for dependence of mechanical properties upon constant tem- perature and moisture states, to model material behavior under nonconstant, nonuniform temperature and moisture states. The reason for considering this extension would be as an application in solving thermo-hydro-viscoelastic boundary value problems. The effects to be studied here are outside the scope of the first order linear theory, and consequently, a coupled thermo-hydro-viscoelastic theory which includes the nonconstant temperature and 69 moisture dependence of mechanical prOperties is necessarily nonlinear (Christensen, 1971). Starting with the general nonlinear theory Morland and Lee (1960) arrive at the stress—strain relations for nonconstant temperature and moisture states with two basic modifications. First to allow for the temperature and mois- ture dependence of the response functions in the presence of a time—dependent temperature and moisture distributions, the definition (4.35) of the reduced time g must be generalized consistant with the assumed time-temperature- moisture equivalence. Second, since constitutive law for bulk response governs the dilatational response (or volume change), it is necessary to incorporate the temperature and moisture induced expansion. The generalized relaxation constitutive equations are S _ t r aeij(r,r) ij(r’t) - fG1(€'€ ) 3r dr (4.41) t okk(r,t) = [62(g-g')%; [ekk(r,r)‘3aAT(T,T)'3BAM(T,T)]dT -oo (4.42) where a = coefficient of thermal expansion, constant B = coefficient of hydro expansion, constant AT(r,T) AM(raT) temperature deviation from the base temperature To moisture deviation from the base moisture MO 70 and according to Hammerle (1972) t t — dt" dt" €(r’t) ‘ g aTrT(r.t")] * é aM[M(r,t”)] (4'43) and ' _ T dtH T dtH E (I‘VE) - f aTTT(I‘,t")] + f aM[M(r,t”)j (4'44) 0 0 Note that a and B were taken to be constant throughout this study. For a thermo-hydrorheological1y simple material (aT and aM= constant), relation (4.43) reduces to t t t £(t) = — + — = — (4.45) a T a M a where a _ aTaM aT+aM Relation (4.45) is comparable to relation (4.40). Since both t and g are involved in equations (4.41) and (4.42), a reduced form for these stress-strain relations involving only 5 is obtained. Relation (4.43) may be inverted with respect to t, so that t = g(r,€) (4.46) Also by (4.43) and (4.46), using only a temperature gradient, 71 3 _ 1 27% " aT[T(r,t)] ’ 2112—.- 8_£.-1 3E ( ) (4.47) Suppose F(r,t) is any function of position and time. Then, to avoid ambiguity, we shall consistently adopt the notation F(r,t) = F(r,g(r.g)) = ficr.g) (4.48) Using (4.46), a different stress, strain, temperature, and moisture functions are defined as cij(r,t)‘=Efij(r,g(r,€y)==gij(r,€) (4.49) oij(r,t) = oij(r,g(r,€)) = Sier,e) (4.50) T(r,t) = T(r,g(r,g)) = i(r.g) (4.51) and M(r,t) = M(r,g(r,g)) = M(r,a) (4 52) Through the use of (4.49)-(4.52) and taking into account the relations (4.5) and (4.6), (4.41) and (4.42) become E A I 23 C(E’E')aflié215_ng' (4,53) Sij(r,g) A E A A okk(r,£) = 3éK(g‘£')§%T[€kk(r,€')'3aAT(T,€') -35AM(r,g')]dg' (4 54) 72 where strain history before time zero was neglected for simplicity. Relations (4.53) and (4.54) now involve con- volution integrals, and because of this, it would seem that the Laplace transform might conveniently be used to solve boundary value problems involving thermo-hydrorheological1y simple materials. Performing the Laplace transforms we have Sij(r,s) ZsG(s)eij(r,s) (4.55a) okk(r,5) = 38K(s)[ekk(r,s)-3aAT(r,S)-3BAM(r,S)] (4.55b) or more simply éij = zséeij (4.56a) Skk = 35x(8kk-3aai-ssrn) (4.56b) In general, by incorporating the temperature and mois- ture effects, constitutive equatoin (4.19) becomes "'1 O _ 2 _ ' 3 .. I _ t _ t t (4.57) which by using convolution notation, may be written as = k _ - Oij Gijkl (Ekl aATkl BAMkl) (4.58) 73 4.3 Finite Element Formulation in Thermo-Hydro Visco- elasticity Constitutive relations governing a viscoelastic boun— dary value problem were presented in the previous sections. These equations are essential for the analysis of the be- havior of these materials under different loading con- ditions. A finite element method for solving thermo-hydro loaded viscoelastic boundary value problems in now pre- sented. Finite element methods have already been used by several authors in solving these kind of problems. The following derivations are similar to those by Taylor et al., (1970), Heer and Chen (1969) and De Baerdemaeker (1975). 4.3.1 A Variational Theorem Let V be a functional defined as (Christensen, 1971) 7': 7: . _ * * .. ijkl Eij £k1 Gijkl “ATij €k1 - * ._ o * Fi Ui]dV f (“i Ui)da (4.59) * 6 "k1 80 where it is assumed that the diSplacement boundary conditions 74 are satisfied. V is the total volume of the solid. It can be shown that the first variation vanishes when the equilibrium equations and boundary conditions are satis- fied. In other words, the solution of the boundary value problem stated in (4.18), (4.19), (4.22) and (4.23) can be obtained by finding the stationary value of the func— tional V in (4.59). The finite element method can be used as a numerical technique based on the minimization of this functional. Note that in equation (4.59) aT and AM are prescribed functions of time and place associated with a solution of the simultaneous heat and mass diffusion boundary value problem for the body. 4.3.2 Radially Symmetric Viscoelastic Solids In this section the previous results are specialized for a particular class of problems for isotropic solids. Consider a solid sphere subjected to an arbitrary radially symmetric temperature and moisture distribution. The sphere is assumed to be istropic, continuous and homo- genous, and let the origin be at the center of the sphere. Because of symmetry, Figure 4.2, using spherical coordinates, the displacements are Ur = Ur(r,t) ; U = U¢ = O (4.60) 75 Stresses in Radially Symmetric Problems Figure 4.2. 76 The strains therefore, are .; M392. _3__ (4,61) rr ar ’ 96 r ¢¢ ere = 56¢ = €¢r = 0 (4.62) The dilatational strain is e = 1(6 +5 +3 ) (4.63) 3 rr ee ¢¢ The deviatoric strains are err = Err-E = g: %?(:£) (4.64) eee = e¢¢ = %£ %F(g£) = '%err (4.65) Again, because of symmetry there will be three non- zero components of the stress tensor, the radial component 0 and two tangential components 0 rr and 0¢¢ such that 89 066 = O¢¢ (4°66) while Ore =06¢ = O¢r = 0 (4.67) Now, decomposing the stress tensor into a uniform nor- mal stress term (spherical or isotropic part) and a pure shear term (deviatoric part), one finds the isotrOpic stress is 77 _ 1 _ 1 o — 3(Orr+oee+0¢¢) - 3(orr+2066) (4.68) The deviatoric components of the stresses are _ _ _ 2 - Srr - Orr o — 3(Orr 066) (4.69) s = s = -1s (4 70) 00 ¢¢ 2 rr ‘ and sre = 59¢ = s¢r = 0 (4.71) Combining (4.57) and (4.20) the constitutive equations are a de _ l _ d _ rr _ d Orr ‘ é[s(62 Gl)dET(€rr+€ee+€¢¢)+ Gl‘HET GZHET (aAT+BAM)]dg' (4.72) or = 1 — * l _ 1‘ _ * Orr 3[(G2 Gl)+SGl] Err+3(G2 GI) (Eee+€¢¢) GZ aAT - GZ*BAM (4.73) 78 Equation (4.73) can be rewritten as = 4 9c _2 96,1 -2 7t _. * - Orr (K+-3—G) Err+(K ‘s-G) ”66+(K 73-6) €¢¢ 3K QAT 3K*BAM (4.74) Similar expressions can be derived for the other stress com— ponents 066 and o¢¢. The expression for the functional (4.59) becomes = 1 i: k _ * if - 3% 9': _ 9': v £[7Aij ei ej 3K aATi ei 3K BAMi ei F U]dV -f (s*U)da (4.75) B o where i,j=l,2,3 Aij is a 3x3 symmetric array whose components are (4.76) 79 4.3.3 Discretization of a Region A. Nodal Displacements The volume and surface integrals in (4.75) can be expressed as a sum of integrals over a set of subregions or elements E E 1 f 8*A..* 8 dV— 2 8* K* A 8 V é=l V(e)(El 13 EJ) e=1 £(e)(€1 3 a T1)dV E e e E e e E _ _ * _ V(e) v 9* e fB(e) (U S )da (4.77) o where the superscript e indicates the elementand E is the total number of elements in the body. The displacement in each element are approximated by algebraic polynomials relating them to displacements of nodal points of that element (Zienkiewicz, 1971) or {Ue} = [Ne]{U} (4.78) An example of radially symmetric element is given in Figure 4.3, where each line element represents a con- centric shell. The shape functions for the one dimensional element, written in terms of the radial coordinate r are (Segerlind, 1976) r. '1‘ T'I‘. Ni = ?l??‘ , N. = ___£ (4.79) j i ' L'ft‘auu—ov- 3‘ .7 80 ulunllulslillcll?la 'l I! -l' ||'|'II|IAVO| ill-I 'IIIIIIIIUIIIIII-I“A . n Radially Symmetric Line Elements and Corres- ponding Concentric Shells Figure 4.3. 81 The strain vector is obtained by apprOpriate Space differentiation of (4.78) aNi 3N.“ _1 _1 q 5? r FL L {89} = ii— fi—j— Ui N—i N—j— Ui = [139 ){u} r r Uj r r Uj N- N. N. N. _l _l 1 _l “1‘ r J hr 1‘ d (4.80) where L=rj-ri is the length between nodes i and j of the line element. The evaluation of the [8] matrix is no longer a simple procedure since [E] matrix contains terms that are functions of the coordinate r. One procedure to evaluate the [8] matrix is by using the r values at the center of the element. Also note that in (4.77) terms AT and AM are AT [T(r,t)-T(r,o)] [T(r,t)-TO] {e} and (4.81) AM [M(r,t)-M(r,0)] [M(r,t)-MO] {m} So {9} and {m} are prescribed functions of place and time associated with a solution of the simultaneous heat and mass diffusion boundary value problem for the body, which are already obtained from Chapter 3. Recalling (3.17) and (3.18), we may write 82 e _ {9 }"' NiGi + N-O' e J J [N ]{e) (4.82) {me}: Nimi + N.m. e 3 J [N 1(m) where the matrix of shape functions [Ne], relates tempera- tures and moistures in the element to the nodal temperature and moisture values. Substitution of (4.78), (4.80), (4.81) and (4.82) into (4.77) and the use of matrix notation instead of indicial notation and also taking {I}T = [l l l] , yields _ _ E _ % 7((U}T[B]T*[A]*[B1(U})dv- z 7((U}T[B]T* =1 v e=1 v E _ 3K*a{1}[N]{9})dV - z 1 f({U}T[B]T*3K*B{I}[N]{m})dV 8:. V E f({U}T[N]T*{F})dV - z 7 ((U)T[N]T*{S})da v e=1 BC I (DMtTl =1 (4.83) where the element superscripts are deleted for simplicity. The order of integration can be changed. The integration over space is performed first, then the convolution. Hence, ..——.- ‘RL‘. ' ' m. V‘ 83 E — _ E v = z %{U}T*[7[B1T[A][B1dVJ*{U}- z {U}T* [sax e=1 v e=1 v - T E T ' T K[B] {I}[N]{9}dV] - 2 {U} *[3BIK[B] {I}[N]{m}dV] e: V 1 (4.84) where the last two integerals were drOpped since it is our assumption that the body is only under thermal and hydro loads and there are no body and external forces involved. This functional can be written in a simpler form as v = %{U}T*[K]*(U} - (U}T*{R} (4.85) where the stiffness matrix [K], is [A][8e]dV (4.86) and the force vector is E e E ‘e T e {R} = z {r } = z (3afeK[B I {I}[N ]{@}dV + 3B e=1 e=1 V feK[Be]T{I}[Ne]{m}dV) (4'87) V It should be noted that the displacement vector {U} now contains all the nodal displacements of a region and {R} 84 is the vector of the nodal thermo and hydro forces. Taking the first variation of (4.85) and setting it equal to zero yields 6V = 6{U}T*[K]*{U} - 5{U}T*{R} = 0 OT [K]*{U} - {R} = 0 (4.88) Equation (4.88) is very similar to the finite element equations of elasticity. However, the explicit form of (4.88) contains a time integral' t f [K(E-E')]d{U(T)} = {R(t)} (4 89) T=0 where, K is an assemblage of element stiffness relaxation functions for the body, and R is an assemblage of thermo- hydro loads. The solution of equation (4.89) yields the nodal point displacements. The strains and stresses can then be computed from equations (4.61) and (4.74). The integral equations in (4.89) can be solved numeri- cally by the use of time increments (Gupta and Heer, 1974). Rewriting the integral in (4.89) as a summation over time steps, results in M23 [K(En-Em){AU(tm)} = {R(tn)} (4.90) =1 3 85 where {AU(tm)} is a vector of displacement increments from time tm to t approximated as a displacement at the m+l’ beginning of the time increment. The last displacement increment can now be found from the previous displacement increments and a possible initial step displacement {U0} at time t=0. This is done by re- arranging the terms in the summation in (4.90) n-l [K(En-Em)]{AU(tn)} = {Rctn)) - z [K(an-am)1(AU(tm)} m=1 -[K(an)](Uo) (4.91) Therefore the equation for the first few time steps are [K(o)](UO) = {R(to)) [K(0)]{AU(t1)} {R(t1)} - [K(a1)]{uo} (4.92) [K(O)]{AU(t2)} {R(t2)} - [K(E1)]{AU(t1)} - [K(EZ)]{UO} A disadvantage of this method is that the number of terms on the right-hand side in (4.92) increases with an increasing number of time steps. If the stiffness matrix for each time interval is stored, an enormous amount of computer storage space is required. An alternate method is to rebuild these stiffness matrices each time they are used. This procedure, however, rapidly increases the 86 required computer time, and, therefore, the cost. The latter method was used in the computer programs written for the study of stresses in viscoelastic materials with constnat temperature (De Baerdemaeker, 1975). The same method with some modification in the computer program to incorporate the temperature and moisture effects on the force vector R, was used in this study. Another approach is the use of logarithmic time increments which can be organized such that only one stiffness matrix has to be stored at each moment (Gupta, 1974). Using Zak's method, White (1968) arrived at a relation where only the immediate past two solutions were required to be saved to account for memory effect. More simplifications could be made in the case of exponential series representation of relaxation functions (Taylor et al., 1970; Heer and Chen, 1969). B. Element Stresses In order to derive the stresses in each element, first by combining (4.74), (4.76), and (4.81), we may write {08} = [41*{29} - 3K*(eee} - 3K*{Bme} (4.93) Substitution of (4.80) yields (Ge) = [A]*[Be]{ue} - 3K*{eee} - 3K*{Bme} 87 01‘ {0e} = [A][Ee]*{Ue} - 3K*{aee} - 3K*{5me} (4.94) Writing (4.94) as a convolution integral and dropping the superscript notation gives t _ t {0(t)} = f[A(€‘€')][B]d{U(T)} ’ 3afK(5'5'){9(r)} o o t '38 fK(€-E'){m(r)} (4.95) O Replaceing the integral by a summation over discrete time steps yields the stress as a function of the temperature, moisture and displacement increments _ n - {0(tn)} = [4(an)1[B]{UO} + z [4(an-em)nB](AU(tm)) m=1 n n '3q 2 [K(gn—am)re(tn)} - 3B 2 [K(en-am)]{m(tn)} m=1 m=1 (4.96) Note that the calculation performed in (4.96) has to be re- peated for each element. Maximum shear stress for every element at every time step was also calculated from the relation Tmax.(r’t) = %Iorr(r,t) ' 0¢¢(T,t)| (4.97) 5. RESULTS AND DISCUSSION Finite element formulation as it applies to the solution of thermo-hydro viscoelastic boundary value pro- blems was presented in the previous chapter. De Baerde- maeker (1975) modified an existing finite element program for two-dimensional elasticity to accomodate the uniform temperature and moisture viscoelastic problems related to the loading of apples. The modification included a change to axially symmetric triangular elements, allocation of computer storage for the displacement increments calculated during each time step and the recalculatinn of the force vector on the right hand side of the equation (4.92) after each time step. For stress analysis of soybean, however, due to some existing difficulties and exceptions, the above mentioned computer progrom had to be modified to accomodate the special case of the soybean model. These difficulties and exceptions were: (1) even though soybean is reported to behave viscoelastically (Saxena, 1972), due to lack of in- formation on most of its material properties such as visco- elastic shear and bulk modulus as they vary with time, the computer program reduces to the one for elastic stress analysis. (2) from the reuslts of chapter 3, since the distribution of moisture and temperature was assumed radially 88 89 symmetric with respect to the geometrical center of the soybean model, every radius which extends from origin to the surface is representative of the whole two-dimensional model and so the formulation reduced to the one-dimen- sional case and this was also incorporated in the com- puter program. (3) the existing program was also modified to accomodate the effect of varying temperature and mois- ture on the nodal displacements, resulting force vector and the final outcome of the stresses. The resulting computer program that was used to analyze the stresses in the soybean model (STRESS) as well as the program to determine the moisture and temperature distributions (SIMHMTR) are listed in Appendix A. The program STRESS has the option of having different material properties for the skin than the kernel itself. The stresses were simulated only for the first four hours of drying since it has been reported by several researchers such as Misra (1977) that most of the cracking occurs during the first few hours of drying. The kernel was assumed to have no initial stress and free expansion and contraction (skrinkage) of the kernel was allowed. Therefore all stresses within the kernel were due only to the temperature and moisture gradients occuring within the kernel. The values of the material prOperties of the soybean Ivhich were used in analyzing the stresses are tabulated in 'Iable 5.1. 90 Table 5.1 Rheological properties of the soybean used in stress analysis Material Property Value Dimension Poisson's Ratio (H) .4 - Elastic Modulus (E) 100 MPa Skin Elastic Modulus (ES) 800 MPa Elastic Bulk Modulus (K) 166.7 MPa Elastic Shear Modulus (G) 35.7 MPa 5.1 Discussion on the Failure Criteria Because of the unknown nature of cracking process and the small size and cellular structure of the soybean kernel, the conditions under which cracking occurs were studied rather than the mechanics of the cracking itself. Failure is characterized by exceeding a critical state of stress and since the stresses for the soybean model are available, it is necessary to apply the knowledge of the entire stress distribution to prevent the critical state of stress at which failure will occur. Basically there are four types of failure as determined from: normal stresses, shear stresses, energies, and normal strains (Dal Fabbro, 1979) which may exist for uniaxial, biaxial, or triaxial states of stress or strain. 91 It is reported by Miles and Rehkugler (1973) that shear stress is the most significant failure parameter for apples. Fridley et a1. (1968)in conducting compression tests on peaches and pears, found that bruises occured in the area of maximum shear stress. Rao et a1. (1975) re- ported that, although the numerical values of radial and tangential stresses will be greater than the shear stresses, the shear stress failure theory is probably suitable for many biological materials. In the lateststudy on failure of soybeans, Misra (1977) also assumed the maximum shear stress to be the criteria for failure. Drying is a very different phenomenon than bruising, and grain materials may fail in tension. Studying the failure in a viscoelastic slab subjected to temperature and moisture gradients, Hammerle (1972) reporteddthat since no shear exists it seemed reasonable to apply a tensile failure criterion instead of a shear criteria. Misra (1977) also reported that most stress cracking is expected at the surface of the soybean where tangential stress is twice that of the maximum shear stress if the soybean is a perfect sphere. Perry (1959) noted in examining the navy beans that the cracks "seemed to radiate from the hilum" creating a common check pattern. He further observed that under end loading conditions of the beans, the two cotyledones tended to separate, subjecting therseed coat to a tensile stress. Ultimately this led to failure 92 of the coat in the form of cracking. Because of greater strength of the seed coat in the hilar region of the beans, the cracks tended to turn away from this region, producing the various check patterns often observed on damaged beans. Doing compression tests on soybeans, Saxena (1972) re- ported that the orientation of cracks depends on the test position of kernels. Soybeans in their most stable position (long axis horizontal and hilum on the side) did not split into two cotyledons, but the crack started from the sur- face at the hilum and oriented perpendicularly to the pre-existing failure plane of the specimen. From the works of Hammerle (1972), Perry (1959), Saxena (1972), Misra (1977), and the results of this study (which will be discussed in the next section) it was con— cluded that, in the drying processes, tensile stress is the most significant factor to be studied in relation to the failure of grain materials. In this study, tangential stress was assumed to be the criteria for failure during drying of the soybean kernel. Principal stresses for this problem can easily be shown to be 01 = o¢¢(r,t) = of(t) 011 = 086(r,t) = O¢¢(rst) = Of(t) (5‘1) 0111 = Orr(r’t) 93 where of(t) is the time dependent failure stress. Because of the temperature and moisture loading con- ditions, the effects of the time-temperature and time- moisture shift factors must be included modifying equation (5.1) to O¢¢(r,t:T,M) = Gf(r,t:aT(T):aM(M)) (5-2) 5.2 Analysis of the Results Simulated maximum shear stresses at different radii in- side the kernel and at the skin are plotted in Figures 5.1 and 5.2, for the first 4 hours of drying. As we see in Figure 5.1, maximum shear stress is always zero at the center of the kernel and increases as we approach the sur- face. In comparing Figures 5.1 and 5.2, a sudden jump in the magnitude of the maximum shear stress at the surface is very noticable. This is due to stiffer properties at the skin in compare to the rest of the kernel. At all different radii, maximum shear stress reaches a peak value in about an hour and then decays slowly. With the elastic solid assump- tion, the dried soybean will always have some stress, as was discussed by Misra (1977). With the viscoelastic assumption, the stresses reach a maximum and then relax slowly, approaching stress free conditions after a long period of time. This type of behavior is typical for a Maxwell model. The 94 D 5‘ RzRRDIUS 0F SOYBERN [I] R:R/4 o 0 R=H/2 :3 A R=3R/4 20‘ 0. Z: 038 U) '- mo 0: p. 00 €33- ch> I 03 3:48 2:3" 8 .. .. _. “0.00 413.00 86.00 {202.00 160.00 TIME (SECONDS) 3810 Figure 5.1. Simulated Maximum Shear Stress at Different Radii Inside the Soybean Model 1.20 1.00 95 SKIN ELEMENT c6.00 40.00 80.00 180.00 TIME (SECONDS) x102 Figure 5.2. Simulated Maximum Shear Stress at the Surface 96 maximum shear stress in the skin element reached its peak value of 1.06 MPa in one hour. The simulated radial and tangential stresses along the radius at the time the tangential stress reaches its maximum are plotted in Figure 5.3. The radius of the soybean was divided into 16 elements and all the stresses were evaluated at the center of the elements. It is clear from the figure that radial stress was negative or com- pressive all along the radius. It reached a maximum at about the mid-radius and approached zero at the skin element. Tangential stress on the other hand, was compressive in the first thirteen interior elements and became tensile be- tween elements 13 and 14 and then experienced a sudden jump at the most outer element. It reached its peak value of 2.056 MPa at the skin element. The change of tangential stress from compressive stress to tensile stress indicates that the 3 outer layers were shrinking faster than the 13 inner layers as expected. The same trend of variation in radial and tangential stresses along the radius of soybean has also been reported by Misra (1977). Since it was assumed that failure is due to tensile stresses and tangential stress is the only tensile stress which occured in the 3 outer layers of the soybean model, so it became necessary to monitor the variation of the tan- gential stress as a criteria for failure. Figure 5.4 shows the variation of the simulated tangential stress with time for the skin element where the big jump takes place in the -—- 0‘.“ --.~ ‘ o -3. .5. “.5 “v... . 97 I J Eh.00 0.10 0.20 0. 0:40 Figure 5.3. Simulated Radial and Tangential Stresses at the Time Tangential Stress Reaches its Maximum (One Hour) -.-~-.‘_—~-*¢..e‘m- n... , . 98 SKIN ELEMENT D D o c0.00 40.00 80.00 180,00 180.00 TIME (SECONDS) x10 Figure 5.4. Simulated Tangential Stress at the Sur- face vs. Time 99 value of o¢¢. The shape of the curve is very similar to the one for maximum shear stress. Tangential stress is reaching its peak value of 2.056 MPa is one hour and then decays slowly. Again with the viscoelastic assumption, this stress should approach zero after a long period of time in order to make the dried soybean to reach a stress free condition. However, due to a lack of visco- elastic material properties, existing elastic properties were used and the solution actually reduced to an elastic one and, therefore, it is expected that the dried soybean will always have some stress. Drying air temperature, relative humidity of the surrounding air and the equilibrium moisture content of the soybean used in the simulation of the stresses which are plotted in Figures 5.1 through 5.4 were 35°C, 65% (dry basis) and 11% (dry basis) respectively. The effect of different drying conditions on simulated tangential stress at the surface was also studied and is shown in Figure 5.5. Two observations about the effect of increase in drying air temperature can be made from the model results. First the magnitude of peak tangential stress and the maximum temperature and moisture differnetials across the kernel are both directly proportional to the step increases in temperature. Second, the time to reach the peak tangential stress value is independent of the change in temperature. Gustafson et a1. (1977) observed similar effects of temperature in studying the stresses in the corn kernel. It can be observed in Figure 5.5 that 222 C363 —.. F1 3 '0 II II II \1010) 010101 DOC) STRESS 3.00 J THNGENTIRL 2.00 l l .00 l .00 40.00 00.00 {20.00 100.00 TIME (SECONDS) 8102 cp.00. Figure 5.5. Simulated Tangential Stress at the Surface as Affected by Drying Temperature 101 the slope of the tangential stress curve increases with an increase in temperature. Based on his experimental results, Misra (1977) related this increase in the slope of the stress curve to the higher percentage of cracked soybeans. That is, there was less cracking at lower temperatures. Soybeans fail (crack) if stresses due to thermal and hydro loading exceed the ultimate strength of the kernel. Even though the ultimate strength of the soybean is not available, for a given value of this property, it is apparent from Figure 5.5 that for higher drying temperatures the material will fail in shorter periods of time, if it fails in tension. 6. SUMMARY AND CONCLUSIONS A numerical technique for the analysis of the stress crack formation in the soybean kernel resulting from tem- perature and moisture gradients during the drying process was deve10ped by using finite element method. This tech- nique first solves the simultaneous moisture and heat diffusion equations for the soybean model to obtain the distribution and gradients of moisture and temperature deve10ped inside the kernel during the drying process. Then, the temperature and moisture gradients are used in a thermo-hydro viscoelastic finite element analysis to simulate the stresses in soybeans under thermo-hydro loads. The following conclusions can be drawn from this study: (1) The simulated drying curve agrees closely with the experimental results in the literature. (2) Mass average temperature of the kernel reaches equilibrium in about one hour and then stays con- stant for the rest of the drying period. (3) External layers reach the equilibrium moisture con- tent much faster than the internal layers. (4) Within the assumed model, the distribution of temperature and moisture are symmetric with respect to the center. 102 (5) (6) (7) (8) (9) (10) (11) 103 The maximum shear stress is always zero at the center and reaches its peak value of 1.06 MPa at the surface in one hour. Tangential stress reaches its peak value of 2.056 MPa at the surface in one hour and then decays slowly. Radial stress is negative or compressive through- out the kernel and approaches zero at the surface after reaching its maximum value at about mid- radius. Tangential stress changes from compressive to tensile stress as it approaches the surface which means the outer layers are shrinking faster than the inner layers. The magnitude of peak tangential stress and the maximum temperature and moisture differentials across the kernel are directly prOportional to the step increases in temperature. The time to reach the peak tangential stress value is independent of the change in temperature. For a given ultimate strength of the soybean, the material will fail in shorter periods of time for higher drying temperature if it fails in tension. “a“ in. A ~.n-....u .; -.’~ -_ - 7. SUGGESTIONS FOR FUTURE RESEARCH A specific recommendation for crack free drying can not be made from the results of this study. This was primarily due the lack of information on most of the material properties of the soybean kernel. In general this dissertation was able to introduce a new concept for attacking the stress cracking problem in grain drying. This study formulated the development of temperature and moisture gradients and the state of stress arising from them. With the required material properties in hand, some specific temperature and relative humidity guidelines could have been developed for crack free drying. This research has built a sound theoretical basis for such a study. Some of the important areas and material para— meters related to stress cracking of the soybean which must be investigated are as follows: (1) The dependence of parameters such as c, p, K' and D on temperature and moisture content was discussed in chapter two. In the drying process, any specific time corresponds to some specific levels of moisture con- tent and temperature inside the kernel. For a finite element analysis of coupled heat and mass diffusion, the values of these parameters are needed at some specific increments of time. Therefore these 104 (Z) (3) 105 parameters should be represented either as functions of local temperature and moisture content or time. Amoung these parameters, the diffusion coefficient (D) needs additional study because this investigation showed that D is the most dominant factor in the diffusion process. Additional study is needed to determine whether a significant amount of moisture diffuses through the hilum region of the seed coat. This will help to im- prove the model and may eventually result in a nonsym- metric distribution of moisture within the kernel. Due to lack of information on viscoelastic prOperties of soybean, stress solutions were reduced to the one for elastic analysis. Two specific properties of soybean that are needed to obtain viscoelastic solutions are viscoelastic shear and bulk moduli. Similar to the well-known relationships between the various moduli for isotrOpic elastic materials, there exist corresponding relationships between the various moduli for isotrOpic linear viscoelastic materials. That is to say, if any two of the seven characteristic functions of isotropic materials (Heer and Chen, 1969) are known from experiment, the others can be determined analytically. The interrelationships among four primary mechanical properties, bulk modulus (K), extension or tensile modulus (E), shear modulus (G) and Poisson's ratio (u) are given by (Muki and 106 Sternberg, 1961) _ 96(t)K(t _ E(t) E“) " W) ’ K(t) ‘ 3[1-2LI(t)] _ 3‘K(t)'-2G(t _ E(t) Mt) - WK t +6 in , G(t) - 7W7] Note that in the above relations, if the value for the required property is desired at a specific time, the two properties used in the calculation must be the values at that specified time. Due to temperature and moisture content dependency of these viscoelastic properties, time-temperature and time-moisture shift factors, aT and aM, are needed to relatetthe behavior of these properties over the range of temperatures and moisture contents which would be encountered during the drying process. This is re- presented by master curves for different viscoelastic properties and its concept was discussed in chapter four. For soybean, the parameters aT and aM are avail- able (Herum, et a1. 1973 ) but no attempt has been made for the determination of the four primary pro— perties (E, n, G and K). Saxena (1972) determined a "modified uniaxial compression modulus” (compressive force divided by constant deformation, kg/cm) for soybeans and presented the master curve for this ”property", which has no use in this study. (4) (5) 107 It is clear that experimental determination of any two of the four primary properties is very essential for a successful study of this kind. With required properties in hand, their master curves could be obtained by using aT and aM. From these curves, spec- ific values of the viscoelastic shear and bulk moduli corresponding to specific time increments could be cal- culated and used in the finite element analysis and the computer program STRESS to give the resulting state of stress. Because of a different cellular structure, soybean skin has different material properties than the cotyledons. The soybean model and the computer programs SIMHMTR and STRESS used in this study have the option of in- corporating this factor, but again due to lack of information on the properties of skin, the same mate- rial prOperties were used as for cotyledons. More research is needed to evaluate the material para- meters for the skin and this will help to improve the soybean model. Determination of the ultimate tensile strength of the seed coat is another subject for further research. Using elastic thin ring theory Hoki (1973) calculated the ultimate tensile strength of the navy bean's seed coat. The same method is applicable to soybeans, if the viscoelastic behavior of the skin is considered. (6) 108 The effect of initial moisture content, air flow rate and the relative humidity of the drying air on the drying process and the formation of cracks must be studied and incorporated in the model. Two recent studies by Ting et a1. (1978) and White et a1. (1978) indicated that drying damage is highly dependent on initial moisture content, relative humidity and air flow rate. LIST OF REFERENCES LIST OF REFERENCES Alam, A., and G.C. Shove 1972 Hygrosc0picity and thermal properties of soybeans. ASAE paper 72-320, St. Joseph, Michigan. Alam, A., and G.C. Shove 1972 Conditioning soybeans — A review. ASAE paper 72-321, St. Joseph, Michigan. Arora, V.K., S.M. Henderson and T.H. Burkhardt 1973 Rice drying and cracking versus thermal and' mechanical properties. Tans. of the ASAE, l6(2):320-323. Bacchus, N.N. 1971 Mechanical properties of soybeans. Unpub- lished M.S. thesis. Univerisity of Guelph, Ontario, Canada. Becker, H.A. and H.R. Sallans 1955 A study of internal moisture movement in the drying of the wheat kernel. Cereal chemistry, 32(3):212-226. Bellinger, P.L. ed. 1964 Thermal prOperties of grains. Agricultural Engineering Yearbook, American Society of Agricultural Engineers, St. Joseph, Michigan. Bland, D.R. 1960 The theory of linear Viscoelasticity. Pergamon Press, Oxford, England. Brocci, R.A. 1969 Analysis of axisymmetric linear heat conduction problems by finite element method. An ASME paper presented for the winter meetingof ASME, November 20, 1969, L.A., California. Byg, D.M. and W.H. Johnson 1970 Reducing soybean harvest losses. Ohio Report 55(1):17—18. Chen, C.S. and W.H. Johnson 1969 Kinetics of moisture movement in hygroscopic material (1. Theoretical considerations of drying phenomena). Transactions of the ASAE, 12(1):109-113. 109 110 Chirife, J. 1971 Diffusional process in the drying of tapioca root. Journal of Food science, 36(2):327—330. Chittenden, D.H. and A. Hustrulid 1966 Determining drying constants for shelled corn. Transactions of the ASAE, 9(1):52-55. Chizhikov, A. The effect of drying conditions on the strength of maize kernels. Elevat. Prom., 26(9):15—16. Christensen, R.M. 1971 Theory of Viscoelasticity. An introduction. Academic Press, New York and London. Crank, J. 1964 The mathematics of diffusion. Oxford University Press, New York, N.Y. Dal Fabbro, I.M. 1979 Strain failure of apple material. Unpublished PhJ).thesis, Michigan State University, East Lansing, Michigan. De Baerdemaeker, J.G. 1975 Experimental and numerical techniques related to the stress analysis of apples under static loads. Unpublished Ph.D. thesis, Michigan State University, East Lansing, Michigan. Ekstrom, G.A. 65 Coefficient of cubical thermal expansion and some tensile properties of corn kernels and their relationship to stress crack formation. Ph.D. thesis, Purdue University Library. Ekstrom, G.A., J.B. Liljedahl and R.M. Peart 1966 Thermal expansion and tensile properties of corn kernels and their relationship to cracking during drying. Transactions of the ASAE, 9(4):556—561. Esau, K. - Anatomy of seed plants. John Wiley, New York, New York. Ferry, J.D. 1961 Viscoelastic properties of polymers. John Wiley and Sons, Inc., New York, New York. Flfigge, W. 1975 111 Viscoelasticity. Second editon, Springer - Verlag Publishing Co., Berlin/Heidelberg/ New York. Fridley, R.B., R.A. Bradley, J.A. Rumsey and R.A. Adrian. 1968 Goncharova, 1962 Gradowczyk, 1969 Grosh, G.M. 1959 Gupta, K.K. 1974 Gustafson, R 1977 Hall, c.w. 1961 Hammerle, J. 1968 Some aspects of elastic behavior of peaches and pears. Transactions of the ASAE, 11(46): 46-49. 2. Investigation into the effect of hydrothermal treatment of grain on the changes in its structural and mechanical properties. Mukrom. Elevat. Prom., 28(9):8-10. M.H. and F. Moavenzadeh Characterization of linear viscoelastic materials. Transaction of the Society of Rheology, 13(2):173-191. and M. Milner Water penetration and internal cracking in tempered wheat grains. Cereal Chemistry, 36:1—7. and E. Heer Viscoelastic structures. In structural mechanics computer programs. Edited by W. Pilkey, K. Saczalski and H. Schaeffer. University Press of Virginia, Charlotteville, Virginia. .J., D.R. Thompson and S. Sokhansanj Temperature and stress analysis of corn kernel - Finite element analysis. ASAE paper No. 77-5517 American Society of Agricultural Engineers, St. Joseph, Michigan. Drying farm crops. Chapter 49, Agriculture Engineering Handbook, edited by C.B. Richey, McGraw Hill Book Co., Inc. New York, New York. R. Failure in a thin viscoelastic slab subjected to temperature and moisture gradients. Ph.D. dissertation, Pennsylvania State University. University Microfilms International, Ann Arbor, Michigan. 112 Hammerle, J.R. 1970 The determination of transient temperature and moisture gradients in a slab. ASAE paper No.70-335, American Society of Agricultural Engineers, St. Joseph, Michigan. Hammerle, J.R. 1972 Theoretical analysis of failure in a visco- elastic slab subjected to temperature and moisture gradients. Transactions of the ASAE, 5:960-965. Hammerle, J.R. and N.N. Mohsenin 1970 Heer, E. And 1969 Henderson, S. 1954 Henderson, S 1967 Herum, F.L., 1973 Hilton, H.H. 1954 Hoki, M.O. 1973 Hustrulid, A. 1959 Tensile relaxation modulus of corn horny endosperm as a function of time, temperature and moisture content. Transactions of the ASAE, 13(3):372-375. J.C. Chen Finite element formulation for linear thermo- viscoelastic material. Technical report 32- 1381, Jet PrOpulsion Laboratory, Pasadina, California. M. The causes and characteristics of rice check- ing. Rice Journal, 57(5):16-18. .M. and R.L. Perry Agricultural process engineering. Second edition. The AVI Publishing Company, Inc. Westport, Connecticut. H.J. Barre, K. Majidzadeh and N.C. Saxena Viscoelastic behavior of soybeans due to tem- perature and moisture content. ASAE paper No. 73-319, American Society of Agricultural Engineers, St. Joseph, Michigan. Thermal stresses in thick walled cylinders exhibiting temperature dependent viscoelastic properties of the Kelvin type. Proc. of ASME, 2(3):547-553. Mechanical strength and damage analysis of navy beans. Ph.D. thesis, Michigan State University, East Lansing, Michigan. and A.M. Flikke Theoretical drying curve for shelled corn. Transaction of the ASAE, 2(1):112-114. 113 Hyde, E.O.C. 1954 The function of the hilum in some Papilionaceae in relation to the ripening of the seed and the permeability of the testa. Annals of Botany, 18:241-256. Jasansky, A. and W.K. Bilanski 1973 Thermal conductivity of whole and ground soybeans. Transactions of ASAE, 16(1): 100-103. Kunze, O.R. and C.W. Hall 1965 Relative humidity changes that causes brown rice to crack. Transactions of the ASAE, 8(3):396-399, 405. Landau, H.G., J.H. Weiner and E.E. Zwicky 1960 Thermal stress in a viscoelastic plastic with temperature dependent yield stress. Journal of Applied Mech., Transactions of ASME, 82(3):297-300. Leaderman, H. 1943 Elastic creep properties of filamentous materials and the other high polymers. p. 175. Textile foundation, Washington, DC. Lee, E.H. 1964 Some recent developments in linear visco- elastic stress analysis. Proc. of the 11th International Congress of Applied Mechanics: 396-420. Mannapperuma, J.D. 1975 Analysis of thermal and moisture stresses caused during drying of brown rice. Unpublished M.S. thesis, Louisiana State University. Martin, C.R. and L.E. Stephens 1976 Broken cron and dust generated during repeated handling. ASAE paper 76-3026, American Society of Agricultural Engineers, St. Joseph, Michigan. Miles, J A., and G.E. Rehkugler 1971 The development of a failure criterion for apple flesh. ASAE paper no. 71-330, American Society of Agricultural Engineers, St. Joseph, Michigan. Miller, M.K. 1966 Milner, M. 1950 114 and W.T. Guy, Jr. Numerical inversion of the Laplace transform by use of Jacobi Polynomials. SIAM Journal of Numerical Analysis, 3(4):624-635. Biological processes in stored soybeans. Chapter 13, Soybeans and Soybean Products, Volume I, edited by K.S. Markley. Inter- science Publishers, Inc., New York, New York. Milner, M. and J.A. Shellenberger. 1953 Misra, R.N. 1977 Misra, R.N. 1978 Physical prOperties of weathered wheat in relation to internal fissuring detected radio- graphically. Cereal Chemistry, 30:202-212. Failure and stress analysis in soybeans sub— jected to a temperature and moisture gradient. Unpublished Ph.D. thesis. North Carolina State University, Raleigh, North Carolina. and J.H. Young Finite element analysis of simultaneous mois- ture diffusion and shrinkage of soybeans during drying. ASAE paper 78-3056, American Society of Agricultural Engineers, St. Joseph, Michigan. Mohsenin, N.N. 1968 Morland, L.W. 1960 Morrow, C.T. 1965 Muki, R. and 1961 Physical prOperties of plant and animal mate- rials: parts 1 and 2 of volume 1. The Pennsylvania State University, University Park, Pennsylvania. and E.H. Lee Stress analysis for linear viscoelastic materials with temperature variation. Trans- acions of the Society of Rheology, 4:233-249. Viscoelasticity in a selected agricultural product. M.S. thesis. The Pennsylvania State University, University Park, Pennsylvania. E. Sternberg On transient thermal stresses in viscoelastic materials with temperature dependent pro- perties. Journal of Applied Mechanics, 6:193-207. 115 Overhults, D.G., G.M. White, H.E. Hamilton and I.J. Ross 1973 Drying soybeans with heated air. Transactions of the ASAE, l6(l):llZ-113. ' Paulsen, M.R., and G.H. Brusewitz 1976 Coefficient of cubical thermal expansion for spanish peanut kernels and skins. Trans- action of ASAE, 19(3):592-600. Perry, J.S. 1959 Mechanical damage to peabeans as affected by moisture, temperature and impact loading. Unpublished Ph.D. thesis,Michigan State Univer- sity, East Lansing, Michigan. Rao, V.N.M. 1971 Failure in a viscoelastic cylinder subjected to temperature and moisture graidents. M.S. thesis, North Carolina State University, Raleigh, North Carolina. Rao, V.N.M., D.D. Hamann, and J.R. Hammerle 1975 Stress analysis of a viscoelastic sphere sub- jected to temperature and moisture gradients. Journal of Agricultural Engineering Research, 20:283-293. Rosen, B. 1964 Fracture processes in polymeric solids. John Wiley and Sons, Inc., New York, New York. Sabbah, M.A., G.E. Meyer, H.M. Keaner, and W.L. Roller 1976 Reversed-air drying for fixed bed of soybean seed. ASAE paper No. 76-3023, St. Joseph, Michigan. Sarvacos, G.D. and S.E. Charm 1962 A study of mechanism of fruit and vegetable dehydration. Food Technology, 19:78-84. Saxena, N.C. 1972 Viscoelastic characterizations of soybean grain under quasistatic loading. Ph.D. dissertation, The Ohio State University, University Micro- films International, Ann Arbor, Michigan. Schwarzl, F., and A.J. Staverman 1952 Time-temperature dependence of linear visco- elastic behavior. Journal of Applied Physics, 23:858-843. Segerlind, L.J. 1976 Sharma, M.B. 1964 Silberstein, D 1976 Singh, R., H.J. 1972 116 Applied finite element analysis. John Wiley and Sons., Inc., New York, New York. Viscoelasticity and mechanical properties of polymers. Summer institute of applied mechanics and material science. Department of Engineering Mech., Pennsylvania State University, University Park, Pennsylvania .A. and V.N.M. Rao Mechanical prOperties of peanuts at various temperatures and relative humidiites. Unpublished paper, Department of Food Science, University of Georgia, Athens, Georgia. Barre, and M.Y. Hamdy Drying spherical porous bodies. Trans- actions of the ASAE, 15(2):338-341. Smith, A.K. and S.J. Circle 1976 Stewart, B.R. 1964 Soybeans: chemistry and technology. The AVI Publishing Company, Inc., Westport, Connecticut. Effect of moisture content and specific weight on the internal friction prOperties of sorghum grain. ASAE paper N. 64-804, American Society of Agricultural Engineers, St. Joseph, Michigan. Taylor, R.L., K.S. Pister and G.L. Goudreau 1963 Thompson, R.A. 1963 Timbers, G.E. 1964 Thermomechanical analysis of viscoelastic solids. Int. Journal for Num. Meths. in Engineering, 2:45-59. and G.H. Foster Stress cracks and breakage in artificially dried corn. Market Research Report No. 631, USDA, in C00peration with Agricultural Experiment Station, Purdue University. Some mechanical and rheological properties of Netted Gem. potato. Unpublished M.S. thesis, University of British Columbia. 117 Timoshenko, S.P. and J.N. Goodier 1970 Theory of elasticity. McGraw-Hill book Company, New York, New York. Ting , K.C., G.M. White, I.J. Ross and O.J. Loewer 1978 Seed coat damage in deep-bed drying of soy- beans. ASAE technical paper No. 78-3006, St. Joseph, Michigan. Wang, J.K. 1956 Drying pea beans intermittently at 100°F. Unpublished M.S. thesis, Michigan State University, East Lansing, Michigan. Wang, J.K. and C.W. Hall 1961 Moisture movement in hygroscopic materials. Transactions of the ASAE, 4(1):33-36. Watts, K.C. 1967 Thermal properties of soybeans. M.S. thesis, University of Guelph, Guelph, Ontario, Canada. Watts, K.C. and W.K. Bilanski 1970 Calorimetric determination of the specific heat of soybeans. Canadian Agr. Eng. 12(1):4S-47. Watts, K.C. and W.K. Bilanski 1973 Method for estimating the thermal diffusivity of whole soybeans. Transactions of ASAE, 16(6):ll43-1145. Whitaker, T.B., H.J. Barre and M.Y. Hamdy. 1969 Theoretical and experimental studies of diffusion in spherical bodies with a variable diffusion coefficient. Transactions of the ASAE, 12(5):668-672. Whitaker, T.B. and J.H. Young 1972 Simulation of moisture movement in peanut kernels: Evaluation of diffusion equation. Transactions of the ASAE, 15(1):l63-l66. White, G.M., T.C. Bridges, O.J. Loewer and I.J. Ross 1978 Seed coat damage in thin-layer drying of soybeans as affected by drying conditions. ASAE Technical paper No. 78-3052, St. Joseph, Michigan. White, J.L. 1968 118 Finite element in linear Viscoelasticity. Proceedings of the second conference in Matrix Methods in Structural Mechanics, AFFDL-TR-68-150:489-516. Whitney, J.D. and J.G. Porterfield 1968 Wilson, H.K. 1955 Young, J.H. 1969 Young, J.H. 1971 Zienkiewicz 1971 Zoerb, G.C. 1958 Zoerb, G.C. 1960 Moisture movement in porous, hygroscopic solid. Transactions of the ASAE, 11(5): 716-719. Grain creps. McGraw-Hill Book Company, New York, New York. Simultaneous heat and mass transfer in a porous, hygrosc0pic solid. Transactions of the ASAE 12(5):720-725. and T.B. Whitaker Evaluation of the diffusion equation for describing thin-layer drying of peanuts in the hull. Transactions of the ASAE, 14(2):309-312. The finite element method in engineering science. McGraw-Hill, London. Mechanical and rheological prOperties of grain. Ph.D. thesis, Michigan State University, East Lansing, Michgian. and C.W. Hall Some mechanical and rheological properties of grains. Journal of Agri. Eng. Research 5(1):83-93. APPENDIX A COMPUTER PROGRAMS SIMHMTR AND STRESS 119 N I K .3 FF. 00 9 RV. T ET U .55. . D. H V T UT. U NT 0 n C = {bu 1. ND 6 o.— N :- MH 9pd D.) P SC AL. H ”L . Ta. E .2 Cu: 0.! T P TP )0 9 TK UA 35 .A S 9 Pol (1 N IL N 9 Fl‘ TI 0... 7.7. F." S ”N : 1|: 9 9 I P. 0.53,) 0 RC 6.. 933 M 0K 3833‘s 9’s par-(c1 \I FSQC-N AP" 1.9 a“ N 1..—LI IO)" 2 197R” .9953" 9 H/....UT T 9 91.x 9GHTI UAE1FI‘I) RNR TI PT 93 9O “.3 L TLCV B 919‘) 925. 9/N3UI £3530) 9‘ 64.0305 (JCT Doc... L6CONU RPE)Q.6.0L 9/ LOP. TAOIJZ .lT DOEGCF. HTITTT‘I 9T.H .T. H 955Kv1T TF1“. ~351N8997{i I nXK In: EGAKGF K... SPNNNNNLT 960R 9R AOOOOOTAl/ ESU HTITT.T.T./ :xNNeTT A éSSQVSNNKICqNS QiNNNNNOO IUEI G.L....p-BCCBLHHLATNflvo OPNHHNHHHATI . E." R‘I.*¢IT.0C£ABPL . C TDDDOOPJCROT.N..C‘U 1. n. C CCCCCCC IN°UT 0‘ THE TITLE CARD AND THE CONTROL DAQAHETEFS INF :CRH09OG LRHO 10.“ % BHQNIT'IHTRQDT’H’HHQTINF, H‘ 0‘ TITLE AND CATA HLAOINGS 290 9‘9 I». 0.).“ 3C 212.... 123 92.: 95 97. NTUNCNT CED ITNITIT NNNO‘AITATA III : DMTDMDH HHHNARNAD‘AQ EEEN EOOEC E0 .30..Q.ND.FC~RFO.F 3122 0 2 1. CTOR A( ) ’ZATION OF THE COLUMN V POINTERS AND INITIAL OUTPUT CALCULATION OF THE C C GLOBAL STIFFNESS MATRIX AND GLOBAL FOFCE NATRIX C+H"FJ 20 NFCHT 95 9 o 3 666:.1 013 P IJIIJS 2 n- : = N 2 = = 1:1, : M.H.DL3191. .Ls FCHNLIIQITO 8. GGG:.AO(0K¢ P‘ JJIJJsbr-AFOQ.‘ .._o c C; 3 0.0 So 1 7‘ A0 a a E C 9' g. L E 9 YING CF TH o . CCCCCC NT DA C INDUT AND ECHO PRINT OF ELEM 7. = = 2‘ :11 1.9"? (00 (4:. 9 Z a O Zl‘x 9‘ O 97.2 \I” \I 9‘ .162 1.096 ll“ ‘1 O RRZ REX 9: = 9343 0’) D 9 9 7.22 T.C.“ 9“) ’NI .3ng 93 N ONL 9 in; Ex I. 4.9-",— €.11.FE 9 .>_N.U.06.N)3 N O D 9 37¢ 9). 0 3’2 9 11)\I11 9x = 9115 901 NNTTTTVV.‘ NIEZTTTT 9‘ : .. AEEA 70””!9!” A11RIIQ 0.:(‘091KC P.PP7.FHHF 1 3 2 C C NT MATRICES ELEM- OF THE CALCULATION ) 3 1| 2 . ) 1 9‘ R a a6 ‘ Z . 3 9| Dx . ) a. 9‘ z o R \l. A 2 ¢ ( \l P. \l . J ’ ( 7. a... C 2 n“ ‘9 O z I \l I . .91 o. ’ RC In R K. R A = 1 O B 3.0“ 1 R RK o .. § KO 3 9‘ O \l I z 3 ’J \I ¢ 5 3‘ \l \l ‘I .8 3 3 T— Q. l\ C, 9‘ DE, R gala q n. 1 9 4A l\ 10‘ \I ’E ‘ RIB 2))\l)\l\l3;‘ .1. K O‘ (312231(A A R1 RSI-‘Tl‘n-‘(ZC m \IOR §ZZZRRRU O 0 OK Too-9..)2 A 0..“ 1)),)\I\92‘ ¥ 0 0: (7.31312CI3G’JEO’ R“((l\(nfl°\. 90 Q OEJ (ZZZRaK-K=AH1=1N O 9 2 = = = = : = 23 = \9 : NDI Q)\I’))’ARIIJNI‘ A123123—u- -. ECG/999‘" 8.0.1.9!(‘99‘ E.OFOFFS RBBBPuCCAAOEn-TI.= 1 2 x I R T A N h. <~ 7. S 9 E T N N F E F c H I o E T D .L S 2 E 1 L I F A II o 8 .9 0 T. H L I! 1. G R 1 II o .. 2 o 9 t W. 5 I t L L T b. 9 \l H H R c» \l 9 O 0 ‘ x o J ,L ’ 0' 8 I I at )H .9 N R R s 2 K. K I 0 9 T . LL () ( 2 O N H ‘1 ”H R, R S 1 5 N 0 K O O 9 K O E t R ‘ ‘9’ ,‘ C I T N F Z T) JR 3 T ¢ 0 8 T KK (9 9 R L I 1 N 9 (:9 R) .9 E E T 0 0 2 :- RR 9 J J P 8 c T I 9 N 6.! 9‘ l‘ m o. E E O T 9 I ) 0 3R R 10 v D G c I: In J2 ‘99 19 p. f- : H N 8T. LE 9' = 99 9 .9 9 9 90R 0 is 9: EV JHF R) ‘9’), T 1. “a C 09.. o NN (HNO‘J JKKK N 2“ T9 O 90 R=Ioo‘ 9999 E ch F 0L E J0 1. — HT62R JJJK N E A 0 E L H :9: It)... I)...“ E D’E N o ’0 KKTTFO +0HHHN L 0R m \I \a \I 2?. (DON’T 05555 E 9' CA 39.53 1951 9 91R 0 OIJKZEEE J0: I 09121.... 9X 3(00T‘T1: : = .. F. E, T T 93 ETO5 QTEE. FF/’)\a\o 0 no OJ A LNTALDET.‘ 0T 0 oH-EEHJKJK EN 9 L OT... =IDI‘T ERNN. 3 3‘ 9 9 9 9 N 0 ON! U D‘AToSIEAI CONNGTTGJJKK 0 : INT c IDHOISTNCJSNNLJKL‘ls“ T. LT‘N L (ARI‘IIRJ‘ .. (at = (T : HHNH T EFFC A FECOF : R0 .. FGFFLFFLSSSS R 01E 0 IRFDIJNFKILIIHEEHEEE B. S 5 9 2 o N 1 1. I cc 3 11 I .JJ. 12'! 9 1J8 IQTTTJJK ) 9 BD(A( cuo1).NP.Nau) 9N$L+1 SN A(J n+6).Ac#GHn+&) A(JGF+$ L.Fz.xnn.u. 1.NI .IHT .NNNgLRM PPNW 1.11.9931 . D 9 9 1((BEEEEHNA :IO.‘ J11NN can“ .. UUUUC9R 9 ICF:8(9)CHJ::)NNNNDFTN 7SGTOSJJGGJ)’3IIII G TION OF MOISTURE AND TEMPERATURE MASS AVEFAGE 9R(1)‘Z(2).R(2)‘Z(1)-fi(3)'Z(2)-R(1)'Z(3) ’)‘EF(1)9(RC1)920‘R(2)9R(3))‘EF(Z) 9 §*NS.ID,R(1).Z(1).Rl2),Z(2)oR(3)9Z(3) 1 1 1 \9 ’03, P o 1‘ .993 N T T3R‘o 9 o ZICFZ 1 5 9,9035... 3 \I ‘1’; /R1 D I J 3319,2‘0 " 9 O ((R’AES T \r \r c. RR93—EROZS ‘1 I N P #6.).(86 TNT. P. 4| 0 0 )911NAZOUO NTA)R 32‘. AGSN (ITST E (TR. E / AN ON N E, ZRE ZURT1’ (.TDI 9 N03 )UOOOSleO L111. 1 O 9 I,’2’. a F :OFR 8 993.1 1211921205 9 1TJGOTDD3K 1L=)(((9(MME:O E 9‘“ 0 o Kn- 9 OII‘RRNRUU 0"... 0511TIDDD 03L (: :(UOSSNH‘ N 1..-in, .. : NSDTETSTZ as, .. : NTE ONJ11NJJJJ€12JTTTTLJLN¢A.IL10MJ12I8:DNONIAR:112NST 1::J ==(:::JJ(NNNNL:L9.LOHL ARCMMN MA‘ :(EAMIMMCII 011(0JJF123“K0000A1AP6 U0 EAOECFUUEOUEFOIFRBU‘UUFOR DIJADJJIJJJAABCCCCCJCNuC.RCDRFISSRDSRIDIEARS‘SSIHN 1 L 0A 8 6 1 111 3.4 00” 1 1 7. 6 a a a 6 .Lw _ 121 X 3 V 3 I ’ 9’ V T N X 0 B 5 9. N 9. N 9 b 9. p I N N 7 M- I! .o X K r. 9 .H. 9 Q. 2 o. x F. F \I 3 U , o. , L ) L 9 3 Ah ) C t) ) I V T. ) N nu. . a 9 T L. 9 TE A X Y N ) 9 ”W Du... C p... R M: b. L 8 TT 2 9 A T. X o 1 N 9V 9 6 D H . 0 F 9 )3 x 9 N D. ) 1 9. 9.. a .p o. 7 U Y \I S F = N 99$ / 4. 0 D 41 :- 9 9 P 9 / .L 3 3 9 U .. ... K N) S 9 T. L R 3 T .J v o x pl. N l\ A“ .L U 9 P6 3 3 .L I 2 ) V «K T. . . 91. V.) 9 S. c. P. U A X )T LP3 I S ) J I Y T P. 9 NT. )) ANT. 51 2. ) L a 9 a S E R) 8H TP N 9 9 C U I C H L A I P TI NP IN 41X 5 L T N R C D 0 H IN 9V. NT . 1 r. A I 9 L N m M E PG. Pn. 9H ) vJIT R V H H T 9 T 9 9 N8)? 9 P N 9/ P . P B 9 N 0 a: "I l9 9. .h) N .3915 Y N \I 9' a 3 G :. HR N) 9TP 9 I1. 0 F R H 9 ) 1 1 N A G G 9 Hun/PIN C. .. S 9 . C A 0 F 0.. 9 9 9 D .K A I GZNRT N .. NIL) 7. D R N N I I D. ... .n P. HQ 9TT$N O I AT...) N 9. N F 9 9 T T N 3 V ... C 9 )3... 9N 7. I 2‘... 9.0 0 U R 1. 2 D 9 I A V GI) NIQ)X T 9 TX: 0 I 0 0 ) T .. 1 F T I R AAA (9...- B 9 999 9 C Y. ) H 9 7 T T a .. 9. 5T. .L O 9 x C (a M T! '7. N))I\I h D 1 391;..1. A p: .9 I7— 9 2 T T 9 S 4‘ Q. 0:.) 9 NOT 99.15... T. N 9 €71.19"... a. 0 0 .. ) 7. X4.) 3 C .2 .1. ) .A 2H3 0 SH... F 9 9 9 9 H 0 I L 9 9T. 9 0 G E C .. T 909.. 0 .9 9 T D... 0 T. H "EA T. NFL NPFFP ... C T TXTTX L a. T 9 .. HTT T 9.. ) 9 P T T 4 UTM. 3 0.. LT. TNNNN 0.. v .). .1343 ) ) T. 0 J H. MON 0)4. .. 4. O T. 3 S H P. 99.. HTTTT L X5T/SG 9“ E P. 3 I P. G T c. 66" G... 9 9 C E G H 2 /)6 ) TNT CIY.I2 ) A T a / 939 H I) o T C I TT)X ) 9.. I G H P 9 162 1 /)GXPPF 1 I 0)/)XI5 T 3 E B S ) X F 130T 1IT T T T ) Y X 5:9 9 o :. E1 0 T )16171 9/ )I:.L I t. Ck. .. 0 NE 0 9T... R 0 (N C 1 U 9X 0 N 9LTNNNNN0 0 I 5? 9H 9/XP F 352 o .. R oK)) 9T0) PQPT T B F 9K)B T150155 1 IRTAOOOOOH 9: N 95111/5N 0 92 9) ) P 0 9II N 11L:..9:.N6.: : .. ...V O Q I: 1...T.TU a 8 TT/ 111.1919. 9 3T. QTGTET 9.. CTII1.N CTITT 3 = U 9 9U 9 0)) )UL EIT) T1ATTNN 1 UHNNSSSSSL N o 6TTTTT6T G 6T=T9K3 F onI I TDKMNZK1N1191NS G o:V1 A MUEAINI OIOONNNNN N 94F TA=.AEA 9N N TAIBNTZ 0 NIBH T ONK NTI:N 9 9 9T. N NIB 9 MOPTHTN:O Q9MMEEEEEL N1:0 9MTMTH7U I UDMSIKVO (“PP A :I3LN0TINII0ITL I KOOJ QTNIQNTNTDETNMHHHHHA T..T AQIRIQiO T :AQiTgsT G T1IY L THiLTANgTTTTTNL T TZIT CO:5..C0FN0NU.-.CCI.-.IIIE FCOT ...C..RCP.0.CC T NECOFNFO h FC:0 L 0:.0AF cCCFZFOROA T F0..I C GTHFCING...SNCCH.CanJCD. 7.9 .SU a...» H... 5.? 9... U KP“: BIKIG I . 7:13 C <.Rn.CT.RC—.$.F 7.5ch .2 IQJX 4.. a. 4. D v 1— S .9 1 81 N 7.9 6 7 N 33 8 A 2 9 A a. C 19 J. n. 8 .U .31 T. 1. ... 2 1. C. 7. .. C I 2 1.. P 7.. :1 4.. a. 9. C CC CCC CCC CCC CCC C iflW) 2,1, ’I (195) TXIT - u b OUTDUT 0r THE RESULTS uPIT 21 C 122 E I P 9 A ) ) P 7 9 91. N N 1|! 2 9|! \9 9 o 91 ’16 L 1 1.9 )=1 m 6: 7 SS :1 )7. 9 9 1T. 5 FR... IT J 91 2 9 0 9 6 B. CC 9T 9): IN T, 2 0 RR 9'1 ’11 08) 01 U a 00 JR 1 9 9 HNP G 9 0 o FF 9 9 91) R 9" I O O a \I I, I‘J 9' LP 9 ,(7 O . 9” LT 1 (J ‘I 9 K . N1 KT, \1 6 \I 3 AN 9 I 9 (RI N \l 9: KX? 3 1. O 3 1. DE T. RI 70599 9 )NI 0 9 3 U 9 ON T 9T R 91 N 1 PH 9 £10 0 o 2 5 0 NE I \I)))))\I)II 91R R 2 "x, NTT - 9| . I T L 0)))))))))JJJJJJJJTN )l‘s U ) T. 9 90... O 9 P) P) P E: POJJJJJJJJJ 9 9 9 9 9 9 9 9 99! J 9‘ T ) 9 JNT RTON1 0N1 m 0 H T N) 9 9 9 9 9 9 9 9 9244098761T 9.) E K \I Tl I 9" TOG 9 9 G 9 9 G To T. 120957.3“3994 “333350)91610) R N T. T TKN HS 61 6I T N N 69559““#9h4111111118061TSO) T T IP ABX 7.. )3T )3T N ) FA 9 3119.1111111TTTTTTTT o O 9T. 06 ) AK A NN OTT 9)11I 01I U 0 0 1 1TTTTTTTTTTIIIIIIII)1T1R)2o 1 TN T 9UP )D R7.:XI9=XICT9 L : :IIIIIIIIIIERRRR‘RR711:T717T 9 L9 11)N JE) T CITPQITQCDQ NA JIIIP-KRRRRRRRR: .. .. .. .. .. .. .. 9 9E1 .. 9 91D 0 PF... ) . =1 9 ITS HIE : )E .. )L9E 00 5“ PG .. = = .. = .. .. .. .. .. ))))))))1“2 9.1.43.9 EINT :flvSTT . 1 TL 9E 16 01)“ 97.); :T 0 IO 0))o)\9)\9))))))JJJJJJJJ61TOJ6199T O C 51 O 99..."... 3 LAUZU /TN.U1.TN01TTDN TN TDDOQSUJJJJJJJJJJ 9 9 9 9 9 9 9 9TTT5 9TTTD NUTS 01NKKINTHTNN KENU ENC 9ENSN A DNN5TT9... 9 9 9 9 9 9 9 9 9 90123Q567EEAQIEEAS NNI 5.. T: K09 .. EIR KTN1ITN1ITU:N LG :II1001I12345678911111111TTM1TTTH= NINL5)OT2T0I)LTTU TIT TIT TIOTT UN THH AA TTTTTTTTTTTTTTTTTTIIR IIIRT TH=L1IAT5=1zILINTD FRFOIRFOIFCDF CI DEEOEEOIIIIIIIIIIIIIIIIIIRROORFROO FEKAOTEIOLOJTAROEN IMICFHICOHLSI LC SRRORRORRRRRRRRRRRRRRRRRRHHFDSHMFS IRNCOARBCNOIMCNCRE AA 1. 72 63 CE 0 a .U 00 0 5 o 2 0 70 1.1 R o o 5 00 o 5 o 5 0.. 7 d. A. 4. h. 25 3 4. 1 1 1. 1 1. 1 1.1. 1 CCCC C 66=INPUT9TAFE61=OUTFUT9TAPE19 C p — 123 TINPUT90UTPUT9TAP S H E T) S D/ s I 9 \I E No R R OR A T NE P 5 AT 9 F C BF. P C L .. =H L L 9| 5A. ) 9 R A Q. 9P 9.. ) R .L A 3A P 1 P 8 = L .IP .. P T. In E 9.. 9 P 9 US 1 X. I vI ) NN .. s 3 9 T 9 I O F 95 0 L + .I 8 .‘ E . C \J 9 c P \I T E a: T QED 9' P S N N R .5 C N "L T1. E A T 9 T A ND R CC M R NF N P N N A P N 0 II A A E 9 9 9 9 a 9 9 0 F ST R P "X ) E L N I E .0 NS A 9 E3 1 N C 9 T N S L EA P E L 9 T 9 N P D 9 TG A NL L PE. H P 9 N 9 P NN .1 TE. L A .A::. 8 an H 9 S K .9... X 00 0 C ) P! L 9 9 8 ) T 9 MD A C R S E 99A ) R N 1 N R SEA NS T 9 L 53C 1 P 9 J 9 P RnXu Y :11 9 .1 A I?§b T 91 F .T N 1. EAL LL NV 0 T C A 9.. H 9 N A B 9 TL AN I C 0 S CV54 5 S 9) 9 9 3 59M V0 C: 9 9 S39 X 9 T )N ) H T M R E E R R 995 I H N 13 1 G N AT. .E: 920 H P P Di 9 R Ru 9 9N".9 9 9 RON T S T I I K50 T N L P9 P P L A E NUS 9 o 9 9N1 A 9 C JP J K C FL8 I LRN 0 F 1 N HHP H F N TN T 9 N A" XNCC N K : B 309 N 9 A9 A P 9 .LBU IEEIFT.AA 9 E 9 NNUA L 9 P I) u. N P OON MDT T N L M 9K3 AA ) N )1 ) 9 N RL IN PU 0 A 0 ) G EN9 B 1 9 1+ 1) E 9 TC. TIFOL R "TIC .. ) 9 NU. o 9 \l 9" 9) N I: I: N L. OTC A 955 1 C P 9.. S L F 1 MS NM 9 1 ) OFLTEG C.) C c. O P "N K P‘ P G S 9 SC (:3 ) 9 0 COCOHNRA NO) P BE 9 N95 G P GJ CG 1 P N NT’LANFP..C P. ’19 T. “HJO pt 1;! F..J J 97\.JJ T J F. 5R H T 0 L L PF3 4TN N )IS H T T TA TT A T J HE OPNLEP T T N39 ) PAE T 29 T A A AT AA T P A T TB9NTIAL I I) 9C 1) NTJ ZXE T T T TT 2 K P T A HSA RTAX T TH)3E 91P+ 9 P 25M F D X 9 N0 T FLTBEFOCE )A1I o OdNHLi E 9 9I O .L X SE CL P )A 1 9 9X L CNN 7. T50 E 7 30 9E EPOSF: R 1’79 L A IP 8.... N IR .. 10A F E . I . N H 1 920TL oKPGTIC P 6T .. G I M1DN VT 9 TT C .. .. M1T N.M S .=I T 0T6 A ETJJE o A1T s N U1F 0C LE 1IA51TCI)F 0 EH R 6TTTC1P9==SOC 0 T9 I BDEDNO SS :9:ECNY IEDS IPLB . PPE . F NITA AS ..IP 0 .. NEA3C V. N N IP) N N39 N TNEN IKL C I:DM0HT1TNMDL1) LITHI: L LILILL1LL N1LIA:1PLIL I I AR NLARAR R .. SNL I LHIR 9.) B LNLNLL .. LL 9 .. JLN S JLNL N T CA 99 ECECECFPFPGEAOT AEROXT. H AEAEAA 1AA OITAEOTOTAEA I 0 SP W RFQPRFIIIJJJCOA CRHF3N E CRCRCCJCC OJACRCNOACRC F 1 T. E N 3 1 0 2 S 9 h D I 0 Z A 1 1 2 CCCCCCCCCCCCC C CCC C A JP91)9N'9NBH9NCL9NTS9DTI9IPR) 9 )9 T 1) NP9NBH) ),ATJP91)9ATJ1)9NP9NBN9NCL91PRI) 1),NE.NP.KP,GM,BH.NTS,OTI9ATNP+1)91PR9LF,PAR) ILC TATJGSH91)9NP9N8H.GMT1).8MT1),NE) 8 TOTAL DISPLACEMENT C 124 CS: .SP.CP. 90 a a 1 AJINGS AND DATA HE ~- 5 - OUTPUT 0F TITL C RIAL PROPERTIES‘9/913X9‘TIME‘91;)9‘BUL< HOOUL TITLE /.1ox,2LA«.//) 7 6T6TE: T T HI T..IS EAEA 6 THTH 7. IRIR RCPCSO hFHc U9 1 a .9 2 ac CI = BULK MOCULUS = SHEAR MCDULUS BH GM C C C NODE NUMBER X9‘RADIALLY SYMMETRIC BODY’I/910X9’NEL RTZ)‘./) ) E N , Mu ) H ) B A. 9 O H 9 N F G 9 9 2x N ) 2 2 R2 B 2 R R 9T N 9 9 9 12 9 2 1 1 R 9 P T O P R 9x N F 1 9 9 003 9 ) T = D D 57.9 H 3).) :- I I 793 F 9NL IL 9 9 DI T9138...— A 0 D)0N3 SHUTN 9 C N NQT 9 9 TZA 9’ S 9 9O Lx TIPZ L L0 E2 0.3..)NT T E E1 N9 LL3TS o N NF 3 IT9NN 0E 2))1 UISF/ ON) 9:39 BTTTL n‘ldFErIBZ86x ltiSE E1. 9TLZIE 95 EE T 0:0AA 99¢N11 NLNNE szGCCIT oET... ITOOP 1LKTSST RT U T/II/ DA 0! TP TNN UNSSN NCO 12EAIEAIR OONNO ISODRRTHTTNTU QHESN HT7A=:IR IRNTOBHHNH EFOE...ZRCFRC0ENUOIIO RIDRRR NFINPCRESCOCC 2 30 2° 7 MATERIALS PROPERTY MATRIX S GENERATION OF TH CALL MATPRP (SgGNH9BMM) GLOBAL STIFFNESS MATRIX INITALIZATION 0F TH CCC CCC 7 7. 02’8’)))’) p 00:0222122 NHON92999999 98:11N111111 1N’7T1TTTTTT =9K1A1AAAAAA Ii 9: : = 3 : = = 2 = : I’)’)|I\I))\I 3KT122113323 13P999999999 1F121231332 OOTTTTTTTTTT COSAAAAAAAAA 3 1 C0 700 KK=19NE RATION OF ELEMENT STIFFNESS MATRIX . h b GEN 0 COO MAT STIFFNESS OF ELEMENT PROPERTIES INTO THE GLOBAL INSERTION ccc 125 ) J 9 I ( F \I T H S H L B p: 9 0 I: \I In J .b J 9 9. Q. I ( I ( F NI 9') F P3 T T a S A3 = H‘ 7) S 2 19w 5 9 IhJ NN .1). '95:. .10 : J1uIUU TI J( O)INNN US SJJ‘IIR 0N “NJJPTTU D3: ::(¢NNT09H OJJFTOO.CNUI JJJISCCQESC 7kg 8 In U N .u 0 o N E H T M ) 0 \l 3 R 2 C. F O 3 Z ( 2 F ( S R o F t v D T ) 1 S 2 P. H L 7 ’ :- \IE 3 0 T S 1 1| 1,5 (\I a '“Y 2 Q. 0 On: 0 XI\ ) N1... 03 ACTIXV '13: 3/ "N2! LFH I O 5/(2 5.2?! %. ELE‘ N9 h.‘ ENCH 3F . LI 911 T)IC HM” D‘RRR RR NutNNGA/ru. HGRRRRRRQRR IFOONQ‘N G . = .. .. .. = = .. = : T/IIEi O .OHIIJ,\I\I\I\11\J)N UN<.SLI\TI "H122113323R CONN AT "I8 0 9 Q 9 9 Q Q 9 0U RIM-2:.LOMA B:121231332TDD.HHHAARL .. -h(((l\l\‘(l\(.P-NUOII.tC-CU FR $55....SSSSSRESCDDRRFW 2A 0 C sh 2) 9?; 1‘ = C : 8" 11 0C. puN 1 ET BITIAL STRLIN HLTRIXv ’ J RATION OF' TH: 8 MATRIX AND THE GEN 8),) 4. R7: .2 R232 :0 I 9 Q. 0. o F 9 9t 0:21.67.— 2011 01 OR((( R1:=t.1(BBB : RIJ: .. = = = = : Hl‘ )))),’) T = 00J121212 cake-1.9999999 NA 1112233 EBOO“I\‘I\I\( LPVILBBQZHuaazd )lZo 213 C C C C HATQIX HULTI FLICATICK TO OBTAIN C=(BT)(S) C ) J ’ K ( .5 a ) Y. ’ 23 3K 9 9 9‘" 11 18U .. = _. 98 110K”: .0 U, ZZDSSJ 27.. = 2 = 9 H HT. OOUOUl‘ GCSCSC 52 22 TO OBTAIN ELSTF AND TF NGTH KHe (K J) éu 2.0‘3olh159265‘R3A3/L I ea 5.” 11 «LCJ = = .. § '5 IJAUKHIU R A P 1 E N ' P K 9. On p. I 9 S T \I N 2 9 9 L 1 \1 C’l‘ 1 NZH Q in! o. 2 Pv, 0‘ NB... T t 0 9 A F’z )M G.Ll‘ ruF (CG 2:. N ’ l‘ 9 x 9) E) APZ’LZ MN ’7.Tl\ F‘Z‘I- EF‘DT. GCN/l E c=l NNNNLF IOOOTF TIII/I o U‘NN UQZbShN 77°88FIR ONNNOO 22 = 2 2 TYIU D‘s-.E—b-M” H HSNTDBNHNN" 00UOULO.CNUIIICO UCSCSEC? : .SODCCC 8 2 GOO 000 000 GOO 00C) GOO COO 0000 126 RiAL LENGYH INITIALIZflTION OF THE GLOBAL FORCE HATRIX CO 12 J=1vNP OO 12 JH=1gNCL H 0.3 50‘. 20 NT=NTS S t AHA IONNI‘O‘VvPB. o o gx H. 0.01le "in " 'fl'JOnDDnO-izo # V II ll nDIOOflIO'TI-Vh A )OD ( GSNEPAT NHPPufi H baoflzwb o 2 0 fl IIOODI‘ HHNHDZ 0 ' ‘ ‘ 1am II II II "1'!" N X -4 tufltd'u’"), .J " .‘OuNIUVN N. 0 N\ M...u\w REAR 3.3'R2‘R1“2 n+2 J‘Rl“3. c ‘3.;-3..*R24‘ .C’k; 0R1“3.0 ICATICN TO OBTAIN C=(G)(H) C UVVVRI 0 4H 1 b .4 cammboooomgo-czo c: ohn Nm—anme 1.. an'n‘bH o Oxnn m x 3 a-IIOO H AAAAOH“ Nlo NN N N K9K)‘H(K9J) AND HYDRO NOOAL LOADS (ALPHA‘THETA+EETH‘H) 1":th 4 H 0 TO 30 1 65‘- 66.7/(LENGTH“2.0))‘(C(-91)‘BV(- ) 65‘166o7/(LENGTF“2.J))’(C(291)’BV(;) 7 6 U1 .0 \0 IV NNNNO l0 vmvo 333.3h/(LENGTH“2.0))’(C(1y1) ‘1333.3~/(L5NGTH“2.1))‘(C(2.1) Aliflll. canton. A\OA\D—O C+N¢HZVNVHm 0“ Ct”. AV. ”[1!qu Z