Development of an Analytical Model to Predict Stress–Strain Curves of Short Fiber-Reinforced Polymers with Six Independent Parameters

Abstract

Mechanical properties of fiber-reinforced polymers are sensitive to environmental influences due to the presence of the polymer matrix but inhomogeneous and anisotropic due to the presence of the fibers. Hence, structural analysis with mechanical properties as a function of loading, environment, design, and material condition produces more precise, reliable, and economic structures. In the present study, an analytical model is developed that can predict engineering values as well as non-linear stress–strain curves as a function of six independent parameters for short fiber-reinforced polymers manufactured by injection molding. These parameters are the strain, temperature, humidity, fiber content, fiber orientation, and thickness of the specimen. A three-point test matrix for each independent parameter is used to obtain experimental data. To insert the effect of in-homogenous and anisotropic distribution of fibers in the analytical model, microCT analysis is done. Similarly, dynamic mechanical thermal analysis (DMTA) is done to insert the viscoelastic effect of the material. The least mean square regression method is used to predict empirical formulas. The standard error of regression for the fitting of the model with experimental stress–strain curves is closely controlled below 2% of the stress range. This study provides user-specific material data for simulations with specific material, loading, and environmental conditions.

1. Introduction

Short fiber-reinforced polymer (SFRP) parts are widely used in industries as they can be easily molded into complex shapes. However, the orientation of the fibers varies from one point to another in composite structures. In complex shapes such as the dome of a pressure vessel or plastic gears, plastic axle, bicycle seats, etc., the fiber angle varies locally due to the geometry of the structure, fabrication process, and type of fiber used in the composite material. This induces a strong heterogeneity throughout the structures, enhancing anisotropic mechanical behavior. For continuous fiber composites, finite element (FE) analysis is well developed to consider variations in fiber angle locally in stress analysis [1,2,3,4,5], whereas for SFRP parts, the micromechanical models in the FE analysis use the same empirical formulas that are used for continuous fiber composites with some modification. Due to the short fiber length, the randomness of the fiber arrangement significantly varies throughout the specimen. Hence, local fiber orientation distribution affects material characterization. Changes in the microstructure of short fibers should be considered in the calculation of mechanical properties in specimen level analysis and then in finite element analysis of the component.

For example, injection-molded short fiber-reinforced plates show fibers aligned towards the molding direction in outer peripheral layers and transversely deviated in the core layer. Therefore, heterogeneity and anisotropy in material properties should be considered in structural analysis. Mechanical properties of an injection-molded specimen are the combination of the mechanical properties of each layer. Such properties can be estimated by developing a model as a function of the degree of anisotropy of each layer. Engineering values created by such a model can insert the variation of local fiber orientation in the FE analysis of an injection-molded component [6]. In addition to fiber contribution, the polymer also contributes to sensitivity in mechanical properties due to environmental conditions such as temperature and humidity. Hence, in designing SFRP parts, the sensitivity of mechanical properties due to both the fiber and polymer matrix should be considered in material data. Such material data inserted in FE analysis will provide an economical, effective, and precise design.

One of the material data components is the stress–strain curve. The aim of this project is to develop a model that can provide the stress–strain curve considering inhomogeneity and anisotropy of the material due to fiber orientation distribution and sensitivity towards environmental conditions in addition to material selection. In stress analysis, the material is pulled with a certain pulling speed, which develops strain in it. Due to the strain, the resistance in the material increases, called stress. Hence, stress σ is dependent on strain ε. Therefore, the mathematical formula for this model should be as follows:

$$\sigma =f\left(\epsilon \right)$$

(1)

Stress varies linearly until a certain point, which indicates elastic deformation. It follows Hook’s law of elasticity. After the yield point, the material no longer follows Hook’s law due to strain hardening or plastic deformation. Ramberg and Osgood used a three-parameter equation to predict stress–strain relationship beyond the yield point. The equation is as follows [7,8]:

$$\epsilon =\frac{\sigma }{E_o}+k{\left(\frac{\sigma }{{\sigma }_k}\right)}^n$$

(2)

E_o is Young’s modulus or the elastic modulus and σ_k is proof stress corresponding to the plastic strain k. Parameter n describes the bend of the stress–strain curve. The elastic part of the stress–strain curve, which is the first part of the equation, follows Hook’s law of elasticity whereas the plastic part of the stress–strain curve, which is the second part of the equation, follows a power law of the non-dimensional stress ratio. This equation was designed initially for metals such as aluminum where k was generally accepted as 0.2% of the plastic strain and n is a material constant, which is calculated based on 0.01% to 0.2% proof stress [9,10]. This value gives the measure of work hardening or plastic deformation. It varies from 0–0.5 [11]. Equation (2) predicts the stress–strain curve until 0.2% of plastic strain for metals, but after that it cannot follow the curvature of the stress–strain curve accurately [10]. To overcome this limitation, the stress–strain curve is divided into elastic and plastic regions. Tayler series expansion of the Ramberg–Osgood (RO) equation is used for fitting the stress–strain curve in the plastic region [12]. However, this equation uses values of the plastic strain at ultimate and yield limit to calculate the value-modified material constant n [13]. Kamaya et al. [9] used the yield and ultimate strength value to develop a modified version of the RO equation with the help of the J integral. The accuracy of the curves varied from 2–10%. The RO equation requires values of k, n, and plastic strain.

To predict the stress–strain curve for a composite material, the fiber orientation distribution should be considered. Several numerical approaches were developed in LS-DYNA by overlapping fiber orientation distribution from an injection molding simulation model to the finite shell element of the anisotropic structural simulation. [14,15,16,17]. Dean et al. [18] tested a macro-mechanical model in which the average of the layerwise fiber orientation tensor in each direction (flow direction of injection molding and transverse direction) was inserted in a macro-mechanical invariant based on the anisotropic constitutive model mentioned in [18,19]. This FE simulation requires a huge computational time and cost.

An analytical model can be a more economical solution in the design of SFRP. Several analytical models were developed to predict the mechanical properties of SFRP using failure criteria usually used for unidirectional (UD) laminates such as the Tsai–Hill criterion [20], theory of linear elasticity for orthotropic material, Halpin-–Tsai–Nielsen criterion [21], etc. [5,22,23] by considering the specimen as a pile of three UD laminates. The fiber orientation distribution within these layers varies according to the thickness of the injection-molded plate [18].

The thickness of the specimen has a significant influence on anisotropy as well as mechanical properties, especially in injection-molded plates. The fiber fraction ratio aligned to the flow direction is high in thin plates (1 mm) as compared to thick plates (>2 mm) [20,24]. Due to the difference in the fiber orientation distribution, there is a significant difference in the normalized modulus as the thickness of the specimen increases. The E modulus of a 0° fiber-oriented specimen decreases whereas at 90°, increases with an increase in the thickness of the specimen [20]. Moreover, temperature increases the difference in the tensile modulus with the thickness of the specimen for fiber angles greater than 30° [20]. This could be due to matrix-dominated behavior at a fiber angle higher than 30° since the matrix (polymer) is sensitive to temperature. Hence, a single model, which can insert synergetic effects of all influential parameters is required to predict the mechanical behavior of SFRP.

Therefore, in this project, an analytical model is developed, which provides material data (stress–strain curve) considering the influence of environmental, loading, material, and design condition.

$$\sigma =f\left(\epsilon , temp, fo,fc,RH, th\right)$$

(3)

The aim of this project is to develop a formula as mentioned in Equation (3), which can predict the stress–strain curve at any arbitrary temperature (temp), fiber orientation (fo) and fiber content (fc), relative humidity (RH), and thickness of the specimen (th). Sensitivity in mechanical properties due to the fibers can be considered by inserting fiber orientation and fiber content in the formula. Temperature and humidity will add more pronounced viscoelastic behavior of the polymer in the formula.

If we compare a stress–strain curve of a metal and composite material, it is evident that in metals, the stress–strain curve transits quickly and steeply from the elastic to plastic region (Figure 1b). However, in SFRP, this transition is gradual and slow (Figure 1a). Therefore, the stress–strain curve for SFRP should be divided into three distinguished parts, which are defined as follows:

• Linear part: the part of stress–strain curve before the elastic limit (blue part in Figure 1a).
• Onset of the bend: This is the part of the stress–strain curve where the graph starts following a curve. The end-limit of this part will be the start of a line. This limit can be calculated through trial and error and can vary with material composition. This can also be referred to as the elasto–plastic region in the case of metals (yellow part in Figure 1a).
• Offset of the bend: This part of the stress–strain curve is approximately linear after the bend. This can also be referred to as the plastic region. This curve starts from the offset of the bend until fracture of the specimen (red part in Figure 1a).

These three parts can be predicted separately by using some mathematical functions.

The Ramberg–Osgood (RO) equation cannot be used here because it requires the value of the plastic strain. To calculate the plastic strain, values of stress and strain are required. However, in this project, both of these values for an arbitrary material condition are not available.

Therefore, a new empirical formula should be developed. The first step in designing the empirical formula is to predict the elastic modulus. Several models were developed to predict the elastic modulus, i.e., rule of mixture (ROM), inverse rule of mixture (IROM), and Halpin Tsai and Bowyer–Bader model [5,23,25]. All these models are designed for fiber laminates where fibers are continuous and compactly packed. All these models use the fiber volume fraction and other geometric parameters of the fiber. They require the elastic modulus of the fiber and matrix separately. This requirement of the measuring volume fraction ratio and fiber length makes the model user unfriendly and complicated. Moreover, equations involved in these models do not consider the change in temperature and moisture, which modifies not only the elastic modulus of the fiber and polymer but also the length of the fiber [20]. Hence, a model using the mass fraction ratio is more friendly and practical. There is no need for extra effort in converting the mass fraction to the volume fraction. A model described in [26] uses the mass fraction of the fiber. However, this model cannot accommodate variation in temperature.

The modulus of elasticity varies with temperature if the fibers are in the transverse direction of the load. The longitudinal modulus of elasticity of the UD laminate GF PP is independent of temperature, whereas the transverse modulus of elasticity has a decaying tendency with increases in temperature [27]. An empirical formula mentioned in [27] to predict the elastic modulus has inserted a factor of the normalized temperature ratio with the melting point of the polymer. It predicts the E modulus quite accurately but cannot accommodate variations in the fiber angle. Zhai et al. tested the famous Mori–Tanka micromechanical model [28] to predict the elastic modulus with different temperatures and fiber angles [29]. The model requires the orientation tensor [30], elastic moduli of fiber, and matrix and fiber length ratio. This model is dependent on the shape factor of the fibers. Neither fiber angle nor temperature is an independent variable in this model, which is the requirement for this study. To add the variation of the elastic modulus of composite material due to temperature, viscous behavior of the polymer must be inserted in the form of mathematical functions.

In this project, a model is designed in which the relationship between the microstructural properties and mechanical properties of composite material is developed. This is done by studying the arrangement of fibers and matrix in the material at its micro-scale level and by studying temperature-dependent behavior of the polymer. Micro-computer tomography (μCT) analysis helps us to understand the fiber arrangement and its influence on mechanical properties. The stiffness of the polymer is dependent on the temperature. DMTA (Dynamic mechanical thermal analysis) describes the stiffness of the material over the whole temperature range. Parameters of these analyses should be included in the formula for the stress–strain curve in the form of some empirical equations, which will be described in Section 3. In Section 2, experimental methodology, observation of tensile tests, and the development of an analytical model will be discussed.

2. Methodology and Experiment

Short glass fiber-reinforced polyamide (PA46 GF) with different fiber contents (15, 30, and 60 percent by weight) was used to develop the stress–strain model. To add sensitivity of the material towards temperature and humidity in the model, DMTA was conducted. All specimens were pre-conditioned for 50% RH following the methodology described in standard DIN EN ISO 1110 for pre-conditioning of polyamide. Small dog bone specimens of type 1BA from standard DIN EN ISO 527-2 were milled from injection-molded plates with a dimension of 80 × 80 mm². The fiber angle of the specimen was varied by rotating the specimen with respect to the molding direction of the injection-molded plate. The number of specimens per plate was determined on the basis of prior fiber orientation distribution analysis.

2.1. Fiber Orientation Distribution Analysis

The fiber orientation in the injection-molded plates is not uniform throughout the thickness, width, and length of the plate [15,18,20]. Fibers at the outer periphery of the plates are aligned towards the molding direction, whereas at the core they are transversely deviated. Hence, fiber orientation distribution analysis is necessary to select the position for the specimens in the plate where maximum fibers are oriented towards the assigned direction (molding direction). Inhomogeneity of fiber distribution varies with the fiber content, thickness, and position of the specimen with respect to the plate.

Therefore, three specimens from each plate were milled as shown in Figure 2. A small section at the center of ex-centric and centric specimens was scanned by µCT. Fiber distribution and orientation tensor [30] were studied with the help of the software VG Studio from Volume Graphics^®. Fiber angles with respect to the molding direction were compared in each layer of the specimen. Orientation tensor, index of anisotropy, and the variation of the core-shell layer thickness were calculated [30,31]. The index of anisotropy is a quantitative way to characterize a sample on a scale of 0 to 1 where 0 is isotropic and 1 is perfectly anisotropic [32]. The optimum number of specimens per plate for the tensile tests was determined by comparing all these factors in ex-centric and centric specimens for each test matrix point. Tensile tests for each ex-centric and centric specimen were also done to quantify the deviation in the mechanical properties within the plate. The specimens with the highest fraction of fibers orientated towards the assigned orientation (i.e., with respect to the molding direction) and highest index of anisotropy were chosen for further tensile tests.

2.2. Test Matrix

Three test points per independent variable were taken as shown in

Table 1. Test matrix to produce experimental data for the stress–strain model.

Fiber Content	Fiber Orientation	Thickness	Temperature	Humidity
15%, 30%, 60%	0°, 30°, 90°	2 mm, 3 mm	−20 °C, RT, 80 °C	50% RH

. With three experimental points, a linear function (or line) can be predicted. However, there are only two testing points for the thickness of the specimen. Hence, the thickness is a nominal variable rather than numerical.

Fiber orientation in Table 1 is the assigned orientation of the specimen with respect to the molding direction. The number of specimens per plate was determined based on fiber orientation analysis mentioned in Section 2.1. In each specimen, most of the fibers are orientated towards the assigned orientation. The local fiber orientation distribution is not considered in this analytical model. Humidity is kept constant at 50% RH to reduce experimental efforts. However, the analytical model has a capability to insert humidity as numerical variable [33] Five tensile tests were conducted for each test matrix point to reduce data scattering. The same procedure was used for DMTA tests [33].

2.3. Experimental Observation

Tensile tests were performed on a uniaxial tensile testing machine with a load cell of 10 kN or 250 kN depending upon the fiber content of the material. Strain was calculated by digital image correlation from the video captured by a high-resolution 2D camera system. A climatic chamber was attached to perform tensile tests at high and low temperature. This chamber cannot control humidity. However, the effect of high temperature on the moisture content in the specimen was considered negligible if the specimens with a fiber content of 30% remained in the chamber for 8 min, those with a fiber content of 60% remained for 5 min, and those with a content of 15% remained for 17 min or less. None of the tests took more time than these. Hence, the relative humidity of the specimens was assumed to be constant throughout the tests. High and low temperature tests were closely controlled for the isothermal condition with a variation of ±2 °C, but the room-temperature test was not that closely controlled. It varied from 17–27 °C. To ensure uniform distribution of the heat, a thermocouple was attached in close proximity at the center of the specimen. All experimental stress–strain curves are plotted in Figure 3. From the graphs, it can be easily concluded that the fracture strain is lower at a lower temperature and vice versa due to the inverse relationship of stiffness and ductility.

As shown in Figure 3, the stiffness of the curves at −20 °C are quite high as compared to 80 °C of the same material. The water content in the specimens is frozen at low temperature such as −20 °C. This increases the stiffness of the material. Moreover, due to the sensitivity of the polymer to temperature, stiffness varies due to the glass transition temperature of the polymer. At low temperature, in this case, −20 °C, PA46 is in the glassy state. Therefore, the curve shows higher stiffness as compared to 80 °C where PA46 stays in a rubbery state.

Curves for the higher fiber content show higher stiffness. However, with lower orientation angles, the stress–strain curves show higher stiffness. Hence, it is concluded from this observation that the mathematical formula must have special parameters or procedures to accommodate the variation in stiffness with respect to the fiber content, fiber orientation as well as humidity and temperature. Since the same observations were noticed in the curves of storage modulus, a similar approach is used to develop an analytical model for the storage modulus [33]. Differences in Young’s modulus between fiber orientations are due to the fiber content. The higher the fiber content of the specimen, the higher the difference in Young’s modulus for the same fiber orientation. Similarly, the strength of the material increases linearly with the fiber content [34]. The variation of Young’s modulus with respect to the fiber orientation within the same fiber content is not linear. This is well explained in the classical laminate theory. Young´s modulus of fiber-reinforced composites shows hyperbolic behavior if the fiber angle changes in regard to the loading direction [5,20]. Similar behavior is seen in the experimental curves of the storage modulus for the same material and specimens through DMTA [33]. Thus, Young´s modulus could be predicted by the storage modulus.

From the mathematical point of view, the stress–strain curves in Figure 3 can be divided into three parts as described earlier. The linear part can be described by the Young´s modulus. Then, it follows a bend. In the case of SFRP, this onset of bending is longer as compared to metals. After this bend, the material shows purely plastic behavior that could be described as a straight line with lower/minimal slope. This approach will be used in designing empirical formulas.

2.4. Analytical Approach

An analytical approach for prediction of the stress–strain curve can be divided into two steps. First, the stiffness of the material should be predicted. This will help to predict the linear part of the curve. Later, the non-linear part of the curve can be assumed based on some mathematical formulas. Young´s modulus from tensile tests and the storage modulus from DMTA tests is similar based on the basic principle of testing [35,36]. In tensile tests, specimens are strained in a quasi-static condition. In DMTA tests, the specimen is subjected to oscillatory strain under varying temperature. Stiffness has been calculated at each temperature, which consists of the storage modulus and loss modulus. Therefore, the stiffness of the composite material should follow the same behavior as the storage modulus curve with respect to the temperature [37,38].

To analyze and validate the hypotheses, all values of storage modulus calculated from the DMTA tests were plotted against Young´s modulus determined by tensile tests for several fiber orientations, fiber contents, and temperatures as shown in Figure 4. A linear relationship can be observed between Young´s modulus and storage modulus by using the least mean square method of regression. Hence, Young´s modulus can be predicted with the help of the storage modulus model.

Experiments performed at room temperature were discarded because room temperature lies in the glass transition temperature range of PA46 (15–40 °C). In this region, the curve of the storage modulus is very sensitive to small temperature deviations. To predict the stiffness of the material, the storage modulus for arbitrary conditions is required. A separate analytical model for the storage modulus is designed that can predict the storage modulus with four independent parameters (i.e., fiber content, fiber orientation, temperature, and humidity). The development of this model called the storage modulus model is described in a conference paper [33]. The next step is to design the formula for the non-linear part of the stress–strain curve. The stiffness of the composite is the combination of fiber and matrix behavior. The fiber contributes towards stiffness of the composite whereas the polymer influences the toughness. As shown in Figure 3, the sensitivity of the stress–strain curve with respect to temperature and humidity is also due to the presence of the polymer. Therefore, the stress curve is a collection of the behavior contributed by both the fiber and polymer. Similar observations can be noticed in the storage modulus model [33]. In other words, the fiber content increases linearity, whereas the polymer influences the commencement of curvature in the stress–strain curve. The stress values in the elastic range are not exactly linear, which can be seen clearly from Figure 1 and Figure 3. It is a combination of a line and a curve. This is due to the combination of the fiber and matrix. The plastic deformation of the matrix happens much faster as compared to the fibers. This has been investigated by plotting the rate of change of the E modulus throughout the test, which is expected to be constant or near constant at least at the elastic range of the stress–strain curve. This happens in metals, but in short fiber-reinforced polymer (especially PA46), this is not the case. Thus, a single mathematical function cannot predict the stress of the composite material. It can be considered as a combination of linear and non-linear functions in which the linear function influences the stiffness and the non-linear function influences the toughness of the curve. Linearity or slope of the stress–strain curve is governed by the fiber content. Each value of stress in Equation (3) can be assumed as weighted moving average of the fiber and polymer contribution biased by the fiber content by weight. Equation (3) can be described as:

$$\sigma =\frac{Fiber content[\%]}{100}{f}_l\left(\epsilon \right)+\frac{100-Fiber content[\%]}{100}{f}_{nl}\left(\epsilon \right)$$

(4)

Here, $f_l\left(\epsilon \right)$ and $f_{nl}\left(\epsilon \right)$ are the mathematical functions for linearity and non-linearity of the curve, respectively. $f_l\left(\epsilon \right)$ is a linear function with E modulus of the composite. This E modulus is predicted by the storage modulus model by simply inserting all arbitrary conditions. However, the non-linear function $f_{nl}\left(\epsilon \right)$ can be a single or a combination of mathematical functions. The non-linear function will be estimated based on observations of experimental data of the material.

Since the experimental data were evaluated on the basis of standard ISO 527-2 [39], Young´s modulus is calculated using the range from 0.05% to 0.25% for the strain. Similarly, the elastic limit for the predicted curve should be 0.25% strain. As mentioned in Section 1, the predicted curve should also be divided into three parts. The offset of the bend in the stress–strain curve could start from a certain fraction of the fracture strain. For example, one-fourth of the fracture strain could be considered as offset of the bend. This value is calculated through trial and error or by observation. An initial estimate of the offset of the bend can be estimated by fitting the Holloman equation of strain hardening [40] to the available experimental data. Hence, the elaborated version of Equation (4) is as follows:

$$\begin{array}{cc}{\sigma }_e=\frac{Fiber content}{100}{f}_l\left(\epsilon \right)+\frac{100-Fiber content}{100}{f}_{nl1}\left(\epsilon \right)\hfill & \hfill \hspace{1em} [\epsilon <0.25\%]\\ {\sigma }_{elpl}=\frac{Fiber content}{100}{E}_{a_{0.25\%\epsilon }}\epsilon +\frac{100-Fiber content}{100}{f}_{n{l}_1}\left(\epsilon \right)\hfill & \hfill [0.25\%<\epsilon <x\% FS]\\ {\sigma }_{pl}=f_{n{l}_2}\left(\epsilon \right)+{\sigma }_{elp{l}_{x\%FS}}\hfill & \hfill \hspace{1em}\hspace{1em} [\epsilon >x\% FS]\\ \hfill f_{nl1}\left(\epsilon \right) \& f_{nl2}\left(\epsilon \right)=any mathematical function,FS=fracture strain\end{array}\}$$

(5)

Equation (5) shows the formulas used to predict three parts of the stress–strain curve, i.e., linear part, onset of bending and offset of bending, as illustrated by Figure 1. $E_{a_{0.25\%\epsilon }}$ is the value of the apparent elastic modulus at 0.25% strain. The first part of the stress is referred to as ${\sigma }_e$, which indicates the elastic part of the stress–strain curve. The second part is ${\sigma }_{elpl}$, which refers to the onset of the bend of the stress–strain curve, which starts from 0.25%ε until x% of the fracture strain. This x value has been estimated by trial and error. The initial value has been estimated by plotting the Hollomon equation [40] on the graph of true stress and plastic strain. The strain values, from where the Hollomon equation fits linearly to the plastic strain is considered as the first initial estimate for the value of x. The last part of the stress–strain curve is the offset of the bend, which is denoted by ${\sigma }_{pl}$. The non-linear equation/s iterate from the strain values at x% of the fracture strain until fracture strain. ${\sigma }_{elp{l}_{x\%FS}}$ is the stress value at x% of fracture strain (FS), which is the initial value for the iteration. In this analytical model, the non-linear function $f_{nl}$ is described by two separate mathematical functions for the onset and offset of the bend, that is $f_{n{l}_1}$ and $f_{n{l}_2}$, respectively. Equation (5) was used to fit all experimental data, taking strain as an independent variable. For the non-linear function, a set of or single mathematical equations are used that can imitate the curvature nature of the onset and offset of the bend in the stress–strain curve. A curve-fitting method from python coding language was used to determine the best-fitting parameters for the mathematical functions mentioned in Equation (5).

To add temperature as an additional independent variable, the value of the fracture strain at the particular temperature is required. Similarly, the fitting parameters of all mathematical functions must be predicted. Hence, another code in python language was written to predict fracture strain and fitting parameters. The algorithm of this code is described in the next Section 2.5.

2.5. Analytical Model to Predict Fracture Strain and Fitting Parameters

To insert additional variables such as temperature, fiber content, fiber orientation, and thickness in Equation (5), fracture strain must be predicted first. Then, the fitting parameters for function $f_{nl1}$ and $f_{nl2}$ in Equation (5) can be predicted. The value of the fracture strain and the fitting parameters are predicted by using the algorithm shown in Figure 5. An analytical model for the storage modulus was developed based on the same test matrix as mentioned in Table 1 [33].

This algorithm as shown in Figure 5 is followed by the code to predict the fracture strain. There are three different models named Wish FO, Wish FC, and Wish FS. As the names says, these models provide the values of the fracture strain or any fitting parameter for the mathematical functions at an arbitrary temperature, fiber orientation or fiber content. From model Wish FS, the value of the fracture strain is predicted for different temperatures, whereas from Wish FO, the fracture strain is predicted for different fiber orientations. Similarly, Wish FC predicts the value of the fracture strain for different fiber contents. The workflow or procedure to predict the values from these models is described in the next paragraph.

Wish FS creates a normalized curve of the storage modulus with respect to temperature. Experimental fracture strain values for −20 °C and 80 °C are projected on this normalized curve. Hence, the normalized value for fracture strain at any temperature can be calculated by projecting on this normalized curve. The experimental fracture strain of room temperature has been excluded. Room temperature is in the range of 17 to 27 °C. The exact temperature was not recorded. The glass transition temperature for PA46 lies in this range. On the other hand, the storage modulus model is sensitive to temperature [33]. The more accurate the temperature values, the more accurate the prediction of the fracture strain. The experimental fracture strain of any two values of the fiber content or fiber orientation in the case of Wish FO will be projected on the normalized curve. Now the normalized value of fracture strain at any arbitrary fiber orientation or content can be calculated by mapping back the normalized value of the storage modulus at that particular fiber orientation or content.

The same algorithm will be run to calculate the fitting parameters for the mathematical functions. The algorithm in Figure 5 will again execute replacing the fracture strain with the required fitting parameter. For the 4P model with four independent variables, the algorithm will be executed separately for all fitting parameters of the mathematical functions for each iteration. Repetition and a combination of executing these three models (Wish FS, Wish FC, and Wish FO) can predict the fracture strain and fitting parameters at any arbitrary fiber orientation, fiber content, temperature, and thickness.

3. Comparison of Analytical Models with Experimental Data

Equation (5) is used to predict the stress–strain curve. First, strain is taken as an independent variable and the curve is predicted with the help of curve fitting by the least square method. The standard error of regression is calculated. This model is called the 1P model. Later, all other additional independent variables are inserted with the help of the algorithm mentioned in Figure 5. This model is called the 4P model. The value of the fracture strain will instruct the iteration code of the 4P model to stop iterating the values of stress at a certain strain value i.e., fracture strain. Apart from the fracture strain, the fitting parameters of all mathematical functions for any arbitrary material condition are also required. The same algorithm mentioned in Figure 5 is used to predict these values as well.

3.1. 1P Model with One Independent Variable

The 1P model predicts the stress–strain curve directly by fitting to experimental data. The fracture strain value of each model curve is taken from experimental data. Hence, the standard error of regression is quite low as shown in Figure 6. The standard error of regression of the model curves (red) is less than 2% of the stress range. To predict the stress–strain curve for any temperature between −20 °C and 80 °C, the 4P model is used as described in the next section.

3.2. 4P Model with Four Independent Variables

The fracture strain and other fitting parameters of mathematical functions are predicted based on the procedure and concept described in the algorithm mentioned in Figure 5. Some of the parameters are predicted by the storage modulus model [33], for example, the stiffness of the material. The stiffness of PA46GF varies with temperature. Hence, the storage modulus model creates variation in the stiffness according to the viscous behavior.

Figure 7 shows verification of the 4P model with experimental data. Here, green curves are experimental stress strain curves whereas red curves are the 1P model, and colorful dotted lines are the 4P model. It shows that there is a cluster of curves near −20 °C and 80 °C. This verifies that the stress–strain model is working according to the viscous behavior of the material. For 50% RH conditioned PA46, the glass transition region lies between 15° and 40 °C. The variation in the stiffness below and above the glass transition temperature is less compared to the glass-transition temperature region. Therefore, the stress curves at temperatures ranging from 10 to 40 °C are uniformly spread out. Hence, incorporating the storage modulus model to predict the E modulus is a very important part of this analytical model.

In the glass-transition region, the stress values are highly sensitive to the temperature as Young’s modulus varies strongly with temperature. Room temperature lies in this region. However, the exact temperature of the laboratory was not measured. Hence, the standard error of regression of the 4P model at room temperature, which is assumed to be 25 °C, shows abnormally high values. The error has been accumulated from the 1P model to the 4P model due to four occurrences of the prediction of the fitting parameters and fracture strain. Some non-physical behavior of the curves is noticed at the higher temperature. For example, at the junction of two mathematical functions, there is a sudden drop of the curve, which is mathematically correct but physically incorrect with respect to the stress behavior of the material. Excessive bending of the stress–strain curve could be mathematically true based on the equations used but physically incorrect with respect to the material behavior. To reduce such non-physical behavior of the mathematical functions, some check functions were inserted in the 4P model. However, the other stress–strain curves for −20 °C and 80 °C show standard error of regression less than 2% of the stress ranges. The same model has been verified with experimental data for 2-mm-thick specimens and with the experimental data for dry specimens (room temperature condition).

4. Discussion and Conclusions

The paper proposed an analytical model that predicts the stress–strain curve for injection-molded SFRP considering five extra independent variables, namely, the fiber content, fiber orientation, temperature, humidity, and thickness of the specimen. The empirical formulas accommodate synergetic effects of each parameter by incorporating fiber properties with the help of µCT analysis and polymer properties with the help of DMTA analysis. This is the reason behind achieving an accuracy of 2% of the stress range in the prediction of mechanical properties. Mechanical characterization shows that the mechanical properties of the SFRP are strongly dependent on viscous behavior of the polymer. At temperatures higher than glass temperature, the polymer is in a rubbery phase. Hence, SFRP shows ductile behavior in the stress–strain curve whereas at temperatures lower than the glass transition temperature, it shows brittle behavior. In the glass transition region, the storage modulus varies almost linearly, and the stress–strain curves are equally spread out.

This model is developed for injection-molded short glass fiber-reinforced polyamide PA46. The following assumptions are considered:

• The stiffness of the material is predicted through the storage modulus. The frequency used for DMTA analysis is 10 Hz. It is assumed that the DMTA test at this frequency can be similar to the quasi-static tensile test.
• The variation of the stiffness of the material with respect to fiber orientation is assumed to be linear due to a lack of experimental points.
• It is assumed that the fracture strain and fitting parameters of all mathematical functions follow the same storage modulus curve from DMTA [33].

Further research can be done to determine the relationship of the frequency of the DMTA test to the tensile test. In this way, another model can also be designed to predict the fatigue behavior of SFRP. It would be interesting to develop a similar model for continuous fiber-reinforced composites and other matrices. The approach of the model will remain the same if the matrix material changes except for some modifications in the mathematical functions. Similarly, an analytical model for other loading conditions can also be developed, for example, compression or shear loading. Since this model has thickness as an independent variable, layer-wise mechanical properties for injection-molded specimens can also be extracted [6].

In engineering design, the stiffness, strength, yield strength, and ultimate strain are more important than the fracture strain. In most of the cases, a component is designed until the yield strain of the material with a factor of safety. Hence, this model provides these values with high accuracy, which is suitable for component development. This model can be inserted in FEA. In-situ material data provided from this analytical model according to the change in the testing environment can be mapped in each finite element. This will help the simulation engineer to design economical and reliable SFRP parts.