Thermoelastic stresses alleviation for two-dimensional functionally graded cylinders under asymmetric loading

Abstract

The goal of this research study is to examine, through the finite element method, the efficiency of using two-dimensional (2D) functionally graded materials (FGMs) in lowering the elastic/thermoelastic stresses acting on cylinders. In 2D-FGMs, the properties are assumed to vary in the radial and tangential directions simultaneously, which is rarely investigated in the literature. The cylinder is subjected to asymmetric inner normal traction with/without asymmetric thermal loading. For the considered case studies, results revealed that 2D-FGMs is beneficial compared to the conventional grading. In case of considering mechanical load only, the tangential stress declines by almost 39%. Similarly, accounting for thermomechanical load resulted in radical falls for the tangential and axial stresses by around 63% and 61%, respectively. Accordingly, the von Mises stress declines dramatically with different values allowing for safe load escalation, and enhancing the cylinder’s durability. Finally, no certain values for the used tangential function’s parameters are preferred to have maximum reduction of stresses under all working circumstances, which necessitates performing optimization.

Keywords:

• Asymmetric loading
• finite element method
• functionally graded materials
• stationary cylinder
• thermoelastic analysis
• von Mises stress

1. Introduction

Functionally graded materials (FGMs) are extensively used in many engineering applications due to their sophisticated properties. Many authors investigated the behaviors of different functionally graded (FG) structures, such as plates, cylinders, and discs. For plates, Cong et al. [1], for example, showed that the porosity had significant impacts on the thermomechanical buckling and post-buckling of FG cylinder (FGCs). Kouider et al. [2] examined the effects of the volume fraction, geometrical parameters, and heterogeneity index on the bending and free vibration responses of FG sandwich plates. Others were concerned with the fracture of such structures (e.g., see Refs. [3, 4]).

Regarding FGCs, many parameters were investigated under different load types (mechanical, thermal, electrical, magnetic, etc.) and conditions (symmetric, asymmetric), see for instance Refs. [5–7]. For example, Tokovyy and Ma [8] examined the performance of heterogeneous 1D-FGCs (1D: one-dimensional) subjected to partial thermal/mechanical loading with material properties radially graded. In addition, Batra and Nie [9] used the polar form of the Fourier expansion and the method of Frobenius series to solve the differential equations of eccentric hollow FGCs. An exponential grading model was used by Loghman et al. [10], as well, while scrutinizing the performance of FGCs subjected to nonsymmetrical magneto-thermomechanical loads. Furthermore, the decelerating behaviors of multilayer 1D-FGCs under thermoelastic loads were explored while considering the materials’ properties temperature dependency [11]. Analogous analyses could be found studying the behaviors of 1D-FG discs [12–15].

On a related front, researchers and industry seek lowering the stresses that encounter any mechanical structure to raise its lifetime, durability, safety level, and loading capacity. This also reduces the failure probability. Such goals could be achieved through different means that include but not limited to: optimization of some parameters 16–18, modifying the microstructure to yield a new material with desired properties [19] that includes higher strength [20], or developing two-dimensional (2D: two-dimensional) FGMs [21, 22]. This article is concerned with the later method, where a property becomes dependent on two directions. This method is rarely examined in contrast to the large number of articles discussing 1D-FG structures. For 2D-FGMs, to name just a few, Nemat-Alla [23] showed that they yielded improved performance than 1D-FGMs. Researches on 2D-FGMs are easily found in the Cartesian coordinates, see for instance Refs. [21, 24, 25]. However, in polar coordinates, limited numbers considered it as a combination between both the radial and axial directions [22, 26–30]. On the contrary, extremely rare articles were devoted to considering the variation with respect to the radial and tangential directions, see for example Ref. [31, 32].

In view of the aforementioned comprehensive literature review, the goal of this article is to examine the effectiveness of using 2D-FGMs in alleviating the stresses encountering FGCs, compared to the traditional 1D-FGCs. The properties would vary in both the radial and tangential directions simultaneously. In addition, a finite element (FE) scheme is built to solve the differential equations. Two examples are presented: one with mechanical load only, and the other includes a thermomechanical load. For each example, the idea of applying 2D-FGMs is presented, and its associated results are discussed.

2. Problem formulation

Using polar coordinates ($r,$ $\theta , z$), a stationary cylinder with inner and outer radii $r_{\mathrm{i}}$ and $r_{\mathrm{o}},$ respectively, is considered. The mechanical equilibrium equations in the radial and tangential directions are written as [33]: (1) $${\left(r{\sigma }_r\right)}_{,r}+{\tau }_{r\theta ,\theta }-{\sigma }_{\theta }=0$$(1) (2) $${\left(r{\tau }_{r\theta }\right)}_{,r}+{\sigma }_{\theta ,\theta }+{\tau }_{r\theta }=0$$(2) where a comma denotes partial differentiation. ${\sigma }_r,$ ${\sigma }_{\theta },$ ${\sigma }_z,$ and ${\tau }_{r\theta }$ (${\epsilon }_r,$ ${\epsilon }_{\theta },$ ${\epsilon }_z,$ and ${\epsilon }_{r\theta }$) are the radial, circumferential, axial, and shear stresses (strains), respectively, which are evaluated via the constitutive equations (Hooke’s law) as follows [34]: (3) $$\mathbf{\sigma }=\mathbf{C}[\mathbf{\epsilon }-\mathbf{\alpha }\mathrm{T}]$$(3) with (4) $$\mathbf{\sigma }=\left\{\begin{array}{c}{\sigma }_r\\ {\sigma }_{\theta }\\ {\sigma }_z\\ {\tau }_{r\theta }\end{array}\right\}, \mathbf{\epsilon }=\left\{\begin{array}{c}{\epsilon }_r=u_{,r}\\ {\epsilon }_{\theta }=\left(u+{\vartheta }_{,\theta }\right)/r\\ {\epsilon }_z=w_{,z}\\ {2\epsilon }_{r\theta }={\vartheta }_{,r}+\left(u_{,\theta }-\vartheta \right)/r\end{array}\right\}$$(4) (5) $$\mathbf{C}=\left[\begin{array}{cccc}{C}_{11}& C_{12}& C_{12}& 0\\ & C_{11}& C_{12}& 0\\ & & C_{11}& 0\\ sym& & & C_{66}\end{array}\right], \mathbf{\alpha }=\alpha \left\{\begin{array}{c}1\\ 1\\ 1\\ 0\end{array}\right\}$$(5) where $\alpha$ is the thermal expansion coefficient, $T$ is the temperature field (considering $0°\mathrm{C}$ as a reference temperature), $u,$ $\vartheta ,$ and $w$ are the radial, tangential, and axial displacements, respectively. For consistency, it should be noted that the natural shearing strain (${\epsilon }_{r\theta }$) is the one used here, and it is equal to half of the engineering shear strain. Besides, for an isotropic material, Eq. (6)(6) $$\{\begin{array}{c}{C}_{11}=\frac{E}{1-v^2}, C_{12}=\frac{Ev}{1-v^2}\to \text{Plane stress}\\ C_{11}=\frac{E\left(1-v\right)}{\left(1+v\right)\left(1-2v\right)}, C_{12}=\frac{Ev}{\left(1+v\right)\left(1-2v\right)}\to \text{Plane strain}\end{array}$$(6) lists the materials’ stiffness constants ($C_{11},$ $C_{12}$ and $C_{66}=\left(C_{11}-C_{12}\right)/2$) that are dependent on $E$ (elastic constant) and $v$ (Poisson’s ratio). (6) $$\{\begin{array}{c}{C}_{11}=\frac{E}{1-v^2}, C_{12}=\frac{Ev}{1-v^2}\to \text{Plane stress}\\ C_{11}=\frac{E\left(1-v\right)}{\left(1+v\right)\left(1-2v\right)}, C_{12}=\frac{Ev}{\left(1+v\right)\left(1-2v\right)}\to \text{Plane strain}\end{array}$$(6)

In terms of $T,$ it is dependent on the thermal conductivity ($k$), and is obtained through solving the steady-state heat conduction equation in the 2D polar coordinates [31]: (7) $$\frac{1}{r}{\left(rk{T}_{,r}\right)}_{,r}+\frac{1}{r^2}{\left(k{T}_{,\theta }\right)}_{,\theta }=0$$(7)

3. Finite element formulation

The finite element method (FEM) is widely used in modeling the behaviors of complicated problems in structural mechanics since it is known for its powerfulness and robustness. So, a FE scheme is developed through Matlab software to obtain $u,$ $\vartheta ,$ and $T.$ The domain is discretized using the eight-node ($n_n=8$) isoparametric 2D elements with three degrees of freedom ($u,$ $\vartheta ,$ and $T$) at each node. In FEM, $u,$ $\vartheta ,$ and $T$ are approximately related to the corresponding nodal values through introducing number of shape functions $N$ [35]: (8) $$u\approx \sum _{i=1}^{n_n}{N}_i^e{u}_i, \vartheta \approx \sum _{i=1}^{n_n}{N}_i^e{\vartheta }_i, T\approx \sum _{i=1}^{n_n}{N}_i^e{T}_i$$(8) where $N_i^e$ is the element $e$ shape function at the $i\mathrm{th}$ node.

Afterward, the standard Galerkin’s procedures are followed to obtain the following symbolic FE equation [35]: $\mathbf{Kd}=\mathbf{R},$ where $\mathbf{K}$ represents the global stiffness matrix, $\mathbf{d}=\left\{\begin{array}{cc}\mathbf{U}& \mathbf{T}\end{array}\right\},$ $\mathbf{U}$ includes both of $u$ and $\vartheta ,$ and $\mathbf{R}$ is the external force vector. In more detail, the FE discretized equation for a system composed of total number of elements $n_e$ is: (9) $$\sum _{e=1}^{n_e}\left(\left[\begin{array}{cc}{\mathbf{K}}_{\mathrm{UU}}^e& {\mathbf{K}}_{\mathrm{UT}}^e\\ {\mathbf{K}}_{\mathrm{TU}}^e& {\mathbf{K}}_{\mathrm{TT}}^e\end{array}\right]\left\{\begin{array}{c}{\mathbf{U}}^e\\ {\mathbf{T}}^e\end{array}\right\}=\left\{\begin{array}{c}{\mathbf{R}}_{\mathrm{U}}^e\\ {\mathbf{R}}_{\mathrm{T}}^e\end{array}\right\}\right)$$(9) with ${\mathbf{K}}^{\mathrm{e}}$ and ${\mathbf{R}}^{\mathrm{e}}$ are the $e\mathrm{th}$ element stiffness matrix and force vector, respectively. They are determined as follows: (10) $$\begin{array}{ll}{\mathbf{K}}_{\mathrm{UU}}^e={\int }_{{\mathrm{\Omega }}_e}r\left[{\mathbf{B}}_{\mathbf{U}}\right]\left[\mathbf{C}\right]\left[{\mathbf{B}}_{\mathbf{U}}\right]\mathrm{d}{\mathrm{\Omega }}_{\mathrm{e}}& {\mathbf{K}}_{\mathrm{UT}}^e=-{\int }_{{\mathrm{\Omega }}_e}r\left[{\mathbf{B}}_{\mathbf{U}}\right]\left[\mathbf{\lambda }\right]\left[\mathbf{N}\right]\mathrm{d}{\mathrm{\Omega }}_{\mathrm{e}}\\ {\mathbf{K}}_{\mathrm{TU}}^{\mathrm{e}}=0& {\mathbf{K}}_{\mathrm{TT}}^{\mathrm{e}}={\int }_{{\mathrm{\Omega }}_e}r\left[{\mathbf{B}}_{\mathbf{T}}\right]\left[\mathbf{k}\right]\left[{\mathbf{B}}_{\mathbf{T}}\right]\mathrm{d}{\mathrm{\Omega }}_{\mathrm{e}}\\ {\mathbf{R}}_{\mathrm{U}}^e=\underset{\mathrm{\Gamma }}{\int }r\left[\mathbf{N}\right]{\mathit{\sigma }}_n\mathrm{d}\mathrm{\Gamma }& {\mathbf{R}}_{\mathrm{T}}^e=0\end{array}$$(10) where ${\mathrm{\Omega }}_e$ is the element’s domain such that $\mathrm{d}{\mathrm{\Omega }}_e=\mathrm{d}r\mathrm{d}\theta ,$ $\mathrm{\Gamma }$ is part of the boundary with specified tractions ${\mathit{\sigma }}_n,$ $\mathbf{\lambda }=\mathbf{C}\mathbf{\alpha },$ and $\mathbf{k}$ is the thermal conductivity matrix. Also, ${\mathbf{B}}_{\mathbf{U}}$ is the strain–displacement matrix, and ${\mathbf{B}}_{\mathbf{T}}$ is the gradient matrix [36]. These integrations are executed using the gauss quadrature method with $3\times 3$ integration points within the element.

At this stage, the vector $\mathbf{d}$ can easily obtained ($\mathbf{d}={\mathbf{K}}^{-1}\mathbf{R}$). Then, strains and stresses are calculated directly through the standard steps of FEM at the gauss points within the element after applying the proper boundary conditions, and then extrapolated to the nodes, where they are averaged based on a node’s location [36].

4. Results and discussion

In this study, the results published by Li and Liu [37] are regenerated to benchmark the accuracy of the developed FE Matlab code. In Ref. [37], a stationary cylinder with $r_{\mathrm{o}}=3r_{\mathrm{i}}$ under plane strain conditions (${\epsilon }_z=0$) was investigated. The power-law model was selected to describe the gradation of some properties [37]: (11) $$\beta \left(r\right)={\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }$$(11) where $\beta$ is a generic material property, and $\eta$ is the heterogeneity index that the impacts of its variation were extensively examined in prior studies (e.g., see Refs. [38–41]). Also, the two subscripts $\mathrm{i}$ and $\mathrm{o}$ refer to the property at $r_{\mathrm{i}}$ and $r_{\mathrm{o}},$ respectively. Li and Liu [37] only proposed that $E$ was varying according to Eq. (11)(11) $$\beta \left(r\right)={\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }$$(11) , with $E_{\mathrm{i}}=200 \mathrm{GPa}$ and $E_{\mathrm{o}}=E_{\mathrm{i}}/3$ (Figure 1 shows its distribution), while $\eta =1,$ and $v$ was kept constant at $0.28.$

Figure 1. Contour plot for the distribution of $E$ using Eq. (11)(11) $$\beta \left(r\right)={\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }$$(11) .

For loading, the inner circumference was subjected to a sinusoidal normal traction: ${\sigma }_{\mathrm{r}}\left(r_{\mathrm{i}},\theta \right)=p_0\cos \hspace{0.17em}\left(n_p\theta \right),$ where $p_0=100 \mathrm{MPa}$ is the pressure amplitude, and $n_p=2$ is the normal traction coefficient. Whereas, the cylinder’s outer edge was kept free of stress.

Regarding FEM, different number of elements along the radial and tangential directions was tested. After checking the mesh dependency, a great agreement is found to occur at $n_e=20,000$ ($100$ and $200$ elements in $r$ and $\theta$ directions, respectively), as depicted in Figure 2 where dimensionless (normalized) stresses are plotted. It should be noted that, henceforth, any stress component is divided by $p_0$ to be normalized.

Figure 2. Dimensionless stresses distribution from the current study and Ref. [37]. ${\overline{\sigma }}_r$ and ${\overline{\sigma }}_{\theta }$ are plotted at $\theta =0°,$ and ${\overline{\tau }}_{r\theta }$ is plotted at $\theta =+45°.$

Then, using the same example, the idea of applying the principle of 2D-FGMs is examined to mitigate the stresses through the cylinder. In other words, $\beta$ would be a function in both $r$ and $\theta .$ This idea has been rarely investigated in the literature as shown previously in Section 2. Mathematically, this is to be accomplished by implementing a known function $f(\theta )$ into Eq. (11)(11) $$\beta \left(r\right)={\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }$$(11) . Such function can take any formula, for instance, exponential [31], trigonometric/polynomial [32], etc. Here, the focus would be directed toward the trigonometric function: cosine, such that: (12) $$f\left(\theta \right)=\cos \hspace{0.17em}\left(n_{\theta }\theta +s\right)$$(12) where the angular (frequency) coefficient $n_{\theta }=1$ or $2$ [32], and $s$ is the shifting (phase) angle that is to be investigated, and its selected values $\in [-180, 180].$ Therefore, Eq. (11)(11) $$\beta \left(r\right)={\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }$$(11) becomes as follows after being multiplied by $f\left(\theta \right)$ [Eq. (12)(12) $$f\left(\theta \right)=\cos \hspace{0.17em}\left(n_{\theta }\theta +s\right)$$(12) ]: (13) $$\check{\beta }\left(r,\theta \right)=f\left(\theta \right)\times \left({\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_i\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }\right)$$(13) where $\check{\beta }$ is the modified property in the $r$ and $\theta$ directions. However, the nature of the cosine function would yield negative values of any property and escalates the positive value; hence, comparing results would not be reasonable. Thus, the following modification is proposed to keep $\check{\beta }$ within the limits (${\beta }_{\mathrm{i}}$ and ${\beta }_{\mathrm{o}}$): (14) $$\beta \left(r,\theta \right)={\beta }_{\mathrm{i}}+\left(\frac{\check{\beta }\left(r,\theta \right)-\min \hspace{0.17em}\left(\check{\beta }\right)}{\max \hspace{0.17em}\left(\check{\beta }\right)-\min \hspace{0.17em}\left(\check{\beta }\right)}\right)\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right)$$(14)

To make sure that there is no confusion occurs for the reader, Figure 3 presents, for instance, the variation of $E$ at $n_{\theta }=2$ and 1 with different values of $s,$ where the harmonic pattern of the $\mathrm{cosine}$ function appears clearly.

Figure 3. Contour plots for $E$ using Eq. (14)(14) $$\beta \left(r,\theta \right)={\beta }_{\mathrm{i}}+\left(\frac{\check{\beta }\left(r,\theta \right)-\min \hspace{0.17em}\left(\check{\beta }\right)}{\max \hspace{0.17em}\left(\check{\beta }\right)-\min \hspace{0.17em}\left(\check{\beta }\right)}\right)\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right)$$(14) .

At the beginning, the impacts of $n_{\theta }=2,$ while $s$ taking either the values of $0°,$ $-30°,$ $-90°,$ and $180°,$ are studied. Using Eq. (11)(11) $$\beta \left(r\right)={\beta }_{\mathrm{i}}+\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right){\left(\frac{r-r_{\mathrm{i}}}{r_{\mathrm{o}}-r_{\mathrm{i}}}\right)}^{\eta }$$(11) produces the contours depicted in Figure 4 with ${{\overline{\sigma }}_r}_{\max }=-{{\overline{\sigma }}_r}_{\min }=1.1,$ ${{\overline{\sigma }}_{\theta }}_{\begin{array}{c}max\\ \end{array}}=-{{\overline{\sigma }}_{\theta }}_{\begin{array}{c}min\\ \end{array}}=2.75$ and ${\overline{\tau }}_{r{\theta }_{\max }}=-{\overline{\tau }}_{r{\theta }_{\min }}=0.67.$ It should be mentioned that, henceforth, subscripts max and min are used to refer to a quantity’s upper and lower values, respectively.

Figure 4. Dimensionless stress contours based on FEM solution for 1D-FGC. (a) Dimensionless radial stress, (b) dimensionless tangential stress, and (c) dimensionless shear stress.

Enhanced results are obtained upon applying Eq. (14)(14) $$\beta \left(r,\theta \right)={\beta }_{\mathrm{i}}+\left(\frac{\check{\beta }\left(r,\theta \right)-\min \hspace{0.17em}\left(\check{\beta }\right)}{\max \hspace{0.17em}\left(\check{\beta }\right)-\min \hspace{0.17em}\left(\check{\beta }\right)}\right)\left({\beta }_{\mathrm{o}}-{\beta }_{\mathrm{i}}\right)$$(14) . The stresses’ contours for the selected values of $s$ are presented in Figures 5–8. Evidently, the three stress components experienced different rates of reduction. However, ${\overline{\sigma }}_{\theta }$ is the component that witnessed the significant drop compared to ${\overline{\sigma }}_r$ and ${\overline{\tau }}_{r\theta }.$

Figure 5. Dimensionless stress contours for $f\left(\theta \right)=\cos \hspace{0.17em}(2\theta +0°).$ (a) Dimensionless radial stress, (b) dimensionless tangential stress, and (c) dimensionless shear stress.

For ${\overline{\sigma }}_{r_{\max }},$ nearly $7.7\%$ reduction is applicable at $s=0$ (Figure 5a). This percentage declines if another value of $s$ is chosen. However, at such instant ($s\ne 0$), ${\overline{\sigma }}_{r_{\min }}$ showed an increase in the decline percentage hitting $7.6\%$ at $s=180°$ as depicted in Figure 8a. Regarding ${\overline{\tau }}_{r\theta },$ its upper and lower values almost remained static around a value of $\pm 0.55$ which resembles a moderate decline. For example, in Figure 6c, both of them declined by nearly $17.3\%$ at $s=0.$ This percentage approached $20\%$ at $s=-90°$ for ${\overline{\tau }}_{r{\theta }_{\min }}$ as shown in Figure 7c.

Figure 6. Dimensionless stress contours for $f\left(\theta \right)=\cos \hspace{0.17em}(2\theta -30°).$ (a) Dimensionless radial stress, (b) dimensionless tangential stress, and (c) dimensionless shear stress.

Figure 7. Dimensionless stress contours for $f\left(\theta \right)=\cos \hspace{0.17em}(2\theta -90°).$ (a) Dimensionless radial stress, (b) dimensionless tangential stress, and (c) dimensionless shear stress.

Then, when it comes to ${\overline{\sigma }}_{\theta },$ it can be stated that 2D-FGMs are absolutely advantageous. Numbers reveal that at $s=0°$ (Figure 5b), ${\overline{\sigma }}_{{\theta }_{\max }}$ went down steeply by $\sim 39\%$ compared to the value of 1D-FGMs. Conversely, according to Figure 8b it declined moderately by $\sim 10\%$ at $s=180°.$ Similarly, at $s=-30°,$ it is depicted in Figure 6b that ${\overline{\sigma }}_{{\theta }_{\max }}$ and ${\overline{\sigma }}_{{\theta }_{\min }}$ declined significantly by approximately $38\%$ and $12.5\%,$ respectively.

Figure 8. Dimensionless stress contours for $f\left(\theta \right)=\cos \hspace{0.17em}(2\theta +180°).$ (a) Dimensionless radial stress, (b) dimensionless tangential stress, and (c) dimensionless shear stress.

Furthermore, Figure 9 presents the reduction values occurring in the extreme values of each stress component at $n_{\theta }=2$ and $n_{\theta }=1$ at different values of $s.$ Generally, it is found that there are no differences obtained when the sign of $s$ changes as $f\left(\theta \right)$ has a similar harmonic pattern to the only present load. Also, in this case of loading, the use of $n_{\theta }=2$ (Figure 9a) produced larger reductions than $n_{\theta }=1$ (Figure 9b).

Figure 9. Percentages of stresses reduction at different values of $s$ for $f\left(\theta \right)=\cos (n_{\theta }\theta +s):$ (a) $n_{\theta }=2,$ and (b) $n_{\theta }=1.$

It can be stated that all the stress values witnessed different percentages of reduction. Since ${\sigma }_{\theta }$ has higher value compared to ${\sigma }_r$ and ${\tau }_{r\theta },$ the focus would be directed toward it. If a designer chooses $s=0$ (Figure 9a), a reduction of nearly $39\mathrm{\%}$ and $10\mathrm{\%}$ is applicable for ${{\overline{\sigma }}_{\theta }}_{\max }$ and ${{\overline{\sigma }}_{\theta }}_{\min },$ respectively. Since, $\left|{{\overline{\sigma }}_{\theta }}_{\min }\right|>\left|{{\overline{\sigma }}_{\theta }}_{\max }\right|;$ therefore, ${{\overline{\sigma }}_{\theta }}_{\min }$ is the cylinder’s load limit decisive parameter if considered individually. However, the small reduction in ${{\overline{\sigma }}_{\theta }}_{\max }$ is critical in most composites since metals have limited tensile strength compared to the huge compressive strength of ceramics that necessitates higher reduction in ${{\overline{\sigma }}_{\theta }}_{\max }.$ In contrast, Figure 9a depicts that a compromise ($\sim 26\%$ reduction) is applicable between the two percentages at $s=\pm 90°,$ which is similar to predictions concluded from Figure 9b at $s=\pm 45°$ and $\pm 135°$ where a decline by almost $9\%$ is attainable. Such points can be considered as the optimum points.

Based on the aforementioned results, it can be stated that using 2D-FGMs would significantly mitigate stresses allowing for load escalation. To make sure that this advantage always happens, the second example is presented, where thermomechanical loading is present. The following material properties are considered in addition to the previously mentioned ones with $\eta =1:$ $k_{\mathrm{i}}=209\mathrm{W}/\mathrm{mK},$ $k_{\mathrm{o}}=2\mathrm{W}/\mathrm{mK},$ ${\alpha }_{\mathrm{i}}=23\times {10}^{-6}/°\mathrm{C},$ and ${\alpha }_{\mathrm{o}}={10}^{-5}/°\mathrm{C}$ [42]. In addition, the following boundary conditions are used: (15) $$\{\begin{array}{c}{\sigma }_r\left(r_{\mathrm{i}},\theta \right)=p_0\cos \hspace{0.17em}\left(n_p\theta \right)Pa, {\tau }_{r\theta }\left(r_{\mathrm{i}},\theta \right)=0\\ T\left(r_{\mathrm{i}},\theta \right)=T_0+50\sin \hspace{0.17em}\left(n_T\theta \right)\mathrm{°C}\\ u_r\left(r_{\mathrm{o}},\theta \right)=u_{\theta }\left(r_{\mathrm{o}},\theta \right)=0\\ T\left(r_{\mathrm{o}},\theta \right)=0\mathrm{°C}\end{array}$$(15) where $T_0$ is the reference temperature, and $n_T$ is the temperature’s coefficient.

Figure 10 presents the resulting temperature profile and stresses, for the 1D-FGC at $n_p=n_T=2,$ $T_0=100°\mathrm{C},$ and $p_0=100 \mathrm{MPa}.$ At this stage, another comparing parameter is introduced: the von Mises stress (${\sigma }_{\mathrm{VM}}$) that is calculated according to Eq. (16)(16) $${\sigma }_{\mathrm{VM}}=\frac{1}{\sqrt{2}}{\left({\left({\sigma }_r-{\sigma }_{\theta }\right)}^2+{\left({\sigma }_{\theta }-{\sigma }_z\right)}^2+{\left({\sigma }_z-{\sigma }_r\right)}^2+6{\tau }_{r\theta }^2\right)}^{1/2}$$(16) [43]. In the dimensionless form, it becomes: ${\overline{\sigma }}_{\mathrm{VM}}={\sigma }_{\mathrm{VM}}/p_0.$ (16) $${\sigma }_{\mathrm{VM}}=\frac{1}{\sqrt{2}}{\left({\left({\sigma }_r-{\sigma }_{\theta }\right)}^2+{\left({\sigma }_{\theta }-{\sigma }_z\right)}^2+{\left({\sigma }_z-{\sigma }_r\right)}^2+6{\tau }_{r\theta }^2\right)}^{1/2}$$(16)

Figure 10. Resulting contours of the thermoelastic problem of 1D-FGC ($f\left(\theta \right)=1$). (a) Temperature, (b) dimensionless radial stress, (c) dimensionless tangential stress, (d) dimensionless axial stress, (e) dimensionless shear stress, and (f) dimensionless von Mises stress.

Figure 10 shows that ${\overline{\sigma }}_r$ varied between $-3.06$ and $1,$ $-0.81\le {\overline{\sigma }}_{\theta }\le 10.5,$ $-9.8\le {\overline{\sigma }}_z\le -0.8,$ and ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ stood at approximately $9.98.$ Also, the upper and lower limits of ${\overline{\tau }}_{r\theta }$ leveled off at $\pm 0.25.$ Such results resemble a confirmation of the previous studies that ${\sigma }_{\theta }$ has noteworthy influences on the cylinder’s performance [14, 44].

For 2D-FGMs, only $n_{\theta }=2$ ($s=+90°$ and $-30°$) and $n_{\theta }=1$ ($s=180°$ and $-120°$) are presented. The resulting temperature profiles are shown in Figures 11a–14a. It is seen that there are significant differences in their distributions as they are interconnected with $k.$ Nevertheless, no impacts on the maximum and minimum $T$ are seen since they are controlled by Eq. (15)(15) $$\{\begin{array}{c}{\sigma }_r\left(r_{\mathrm{i}},\theta \right)=p_0\cos \hspace{0.17em}\left(n_p\theta \right)Pa, {\tau }_{r\theta }\left(r_{\mathrm{i}},\theta \right)=0\\ T\left(r_{\mathrm{i}},\theta \right)=T_0+50\sin \hspace{0.17em}\left(n_T\theta \right)\mathrm{°C}\\ u_r\left(r_{\mathrm{o}},\theta \right)=u_{\theta }\left(r_{\mathrm{o}},\theta \right)=0\\ T\left(r_{\mathrm{o}},\theta \right)=0\mathrm{°C}\end{array}$$(15) , and no heat generation is present in Eq. (7)(7) $$\frac{1}{r}{\left(rk{T}_{,r}\right)}_{,r}+\frac{1}{r^2}{\left(k{T}_{,\theta }\right)}_{,\theta }=0$$(7) . But, it can be grasped that $T$ in almost all regions tends to decrease except a small zone around the inner circumference. Consequently, profound effects are expected to encounter the stresses’ readings that capitalize on $T.$

Figure 11. Resulting contours of the thermoelastic problem of 2D-FGC with $f\left(\theta \right)=\cos \hspace{0.17em}\left(2\theta +90°\right).$ (a) Temperature, (b) dimensionless radial stress, (c) dimensionless tangential stress, (d) dimensionless axial stress, (e) dimensionless shear stress, and (f) dimensionless von Mises stress.

These induced stresses are plotted in Figures 11–14. It can be figured out that the usage of 2D-FGMs brought about drastic reductions in ${\overline{\sigma }}_{r_{\min }},$ ${\overline{\sigma }}_{\theta },$ and ${\overline{\sigma }}_z.$ Until this juncture, 2D-FGMs are far better than 1D-FGMs. Two tiny drawbacks are noticed. First, the maximum tensile ${\overline{\sigma }}_r$ tends to slightly grow. Second, at certain values of $s,$ ${\overline{\tau }}_{r\theta }$ shows small rise. Nevertheless, the declines in ${\overline{\sigma }}_{\theta }$ and ${\overline{\sigma }}_z$ substantially outweigh those increases, which have negligible impacts on ${\overline{\sigma }}_{\mathrm{VM}}.$

In numbers, the values of ${\overline{\sigma }}_{r_{\max }}$ witnessed slight growths. Conversely, ${\overline{\sigma }}_{r_{\min }},$ which is critical than ${\overline{\sigma }}_{r_{\max }},$ dropped significantly. For instance, in Figure 11b, it reached $-2.3$ at $n_{\theta }=2$ and $s=90°$ that resembles about $75.2\%$ of the value of the 1D-FGC. Likewise, ${\overline{\sigma }}_{{\theta }_{\max }}$ and ${\overline{\sigma }}_{{\theta }_{\min }}$ experienced substantial plunges. For example, at $n_{\theta }=2$ and $s=-30°$ (Figure 12c), they decreased by around $59\%$ and $51\%,$ respectively. Despite that, ${\overline{\sigma }}_{{\theta }_{\min }}$ is still the serious one since it has larger absolute value compared to ${\overline{\sigma }}_{{\theta }_{\max }}.$ Generally, for the selected values of $n_{\theta }$ and $s,$ the decline percentage for ${\overline{\sigma }}_{{\theta }_{\min }}$ is confined between the previous percentage ($51\%$) and about half of it ($\sim 25\%$) (Figure 13c: $n_{\theta }=1$ and $s=180°$). Similar behaviors are noticed for ${\overline{\sigma }}_z.$ In Figure 14d, for instance, the critical value of ${\overline{\sigma }}_z$ went down by $\sim 27\%$ if compared to the value at $f\left(\theta \right)=1$ shown in Figure 10d. This percentage rose reaching nearly $38\%$ at $n_{\theta }=2$ and $s=-30°$ (Figure 12d).

Figure 12. Resulting contours of the thermoelastic problem of 2D-FGC with $f\left(\theta \right)=\cos \hspace{0.17em}\left(2\theta -30°\right).$ (a) Temperature, (b) dimensionless radial stress, (c) dimensionless tangential stress, (d) dimensionless axial stress, (e) dimensionless shear stress, and (f) dimensionless von Mises stress.

Figure 13. Resulting contours of the thermoelastic problem of 2D-FGC with $f\left(\theta \right)=\cos \hspace{0.17em}\left(\theta +180°\right).$ (a) Temperature, (b) dimensionless radial stress, (c) dimensionless tangential stress, (d) dimensionless axial stress, (e) dimensionless shear stress, and (f) dimensionless von Mises stress.

Figure 14. Resulting contours of the thermoelastic problem of 2D-FGC with $f\left(\theta \right)=\cos \hspace{0.17em}\left(\theta -120°\right).$ (a) Temperature, (b) dimensionless radial stress, (c) dimensionless tangential stress, (d) dimensionless axial stress, (e) dimensionless shear stress, and (f) dimensionless von Mises stress.

Such declines are considerably beneficial for raising the working loads and avoiding the occurrence of plasticity (failure) as appears in the reduction of ${\overline{\sigma }}_{\mathrm{VM}}$ (Figures 11f–14f). It is evident that there are variations in the values of ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ based on $s$ and $n_{\theta }.$ Using $s=90°$ with $n_{\theta }=2,$ yielded a decline in ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ by $\sim 22\%$ (Figure 11f). It grew to almost $54\%$ if $s$ is switched to $-30°$ (Figure 12f). Moreover, nearly two-thirds of that two values are applicable when $n_{\theta }=1$ prevails with $s=180°$ (Figure 13f) and $s=-120°$ (Figure 14f), respectively.

Furthermore, Figure 15a,b presents the variations’ percentages occurring in each stress component at different values of $s$ at $n_{\theta }=2$ and $n_{\theta }=1,$ respectively. First, it is seen that the sign of $s$ yields dissimilar results in contrast to the previous case (mechanical load only). Second, the change in ${\overline{\sigma }}_{r_{\max }}$ is negligible if compared to the change of ${\overline{\sigma }}_{r_{\min }}.$ Third, in Figure 15a, except at $n_{\theta }=2$ ($s=-90°,$ $-120°,$ $-150°$), the upper and lower values of ${\overline{\tau }}_{r\theta }$ experienced gigantic rises. Despite that, these numbers are tricky since they are less than unity. For example, at $n_{\theta }=1$ and $s=-30°$ (Figure 15b), ${\overline{\tau }}_{r{\theta }_{\max }}$ increased slightly from $\sim 0.25$ to $\sim 0.7.$ Fourth, considering $f(\theta )$ brought about considerable declines in ${\overline{\sigma }}_{z_{\min }}$ which is more important than ${\overline{\sigma }}_{z_{\max }}$ since $\left|{\overline{\sigma }}_{z_{\min }}\right|\gg \left|{\overline{\sigma }}_{z_{\max }}\right|.$ Though the latter also witnessed massive drops, it can be neglected as they are fractions of one.

Figure 15. Variation’s percentage for each stress component for the thermoelastic problem of 2D-FGC with $f\left(\theta \right)=\cos \hspace{0.17em}\left(n_{\theta }\theta +s\right).$ (a) $n_{\theta }=2,$ and (b) $n_{\theta }=1.$

Finally, since ${\overline{\sigma }}_{{\theta }_{\min }}$ as well as ${\overline{\sigma }}_{z_{\min }}$ have the largest absolute values compared to other stress values, the change in ${\overline{\sigma }}_{\mathrm{VM}}$ is considerably influenced by their variation. In Figure 16, it can be deduced that ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ declined with different ranges despite the often misguiding growth of $\overline{\tau }$ and parts of ${\overline{\sigma }}_r.$ This decline hits almost $27\%$ at $n_{\theta }=2$ and $s=60°.$ At that instant as shown in Figure 15a, ${\overline{\tau }}_{r{\theta }_{\max }}$ grew by approximately $125\%$ and ${\overline{\sigma }}_{{\theta }_{\min }}$ fell by $\sim 40\%.$ Moreover, the use of $s=-60°$ produced the largest decrease for ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ ($\sim 63\%$) at $n_{\theta }=2.$ Nearly half of that value is the most applicable reduction percentage if $n_{\theta }=1$ is used. That shows that $n_{\theta }=2$ is more advantageous than $n_{\theta }=1.$

Figure 16. Maximum value of the dimensionless von Mises stress (${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$) reduction’s percentage.

However, until now, the analyses in the study are restricted to $n_{\theta }=1$ or $2.$ This forms the foundation of the next part, where different values of $n_{\theta }$ are examined, i.e., $n_{\theta }=\pm 1,$ $\pm 2,$ $\pm 3,$ $\pm 5,$ and $\pm 10.$ Four cases are explored at $s=-60°.$ In this part, along with the variation of $n_{\theta },$ the values of $n_p,$ $n_T,$ $p_0,$ and $T_0$ are altered and listed in Table 1.

Table 1. The pertinent numerical values of some parameters used while investigating different values of ${\mathit{n}}_{\mathit{\theta }}$ of 2D-FGC under thermomechanical loading conditions.

	Case (1)(1) $${\left(r{\sigma }_r\right)}_{,r}+{\tau }_{r\theta ,\theta }-{\sigma }_{\theta }=0$$(1)	Case (2)(2) $${\left(r{\tau }_{r\theta }\right)}_{,r}+{\sigma }_{\theta ,\theta }+{\tau }_{r\theta }=0$$(2)	Case (3)(3) $$\mathbf{\sigma }=\mathbf{C}[\mathbf{\epsilon }-\mathbf{\alpha }\mathrm{T}]$$(3)	Case (4)(4) $$\mathbf{\sigma }=\left\{\begin{array}{c}{\sigma }_r\\ {\sigma }_{\theta }\\ {\sigma }_z\\ {\tau }_{r\theta }\end{array}\right\}, \mathbf{\epsilon }=\left\{\begin{array}{c}{\epsilon }_r=u_{,r}\\ {\epsilon }_{\theta }=\left(u+{\vartheta }_{,\theta }\right)/r\\ {\epsilon }_z=w_{,z}\\ {2\epsilon }_{r\theta }={\vartheta }_{,r}+\left(u_{,\theta }-\vartheta \right)/r\end{array}\right\}$$(4)
${\mathit{n}}_{\mathit{p}}$	$2$	$1$	$2$	$2$
${\mathit{n}}_{\mathit{T}}$	$2$	$2$	$1$	$1$
${\mathit{p}}_0$ ($\mathrm{MPa}$)	$100$	$100$	$100$	$2$
${\mathit{T}}_0$ ($\mathrm{°C}$)	$100$	$100$	$100$	$1000$

To avoid redundancy of results, only the variation of ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ is presented in lieu of all the stress components. Figure 17 concludes that it is not always ensured to have $n_{\theta }=2$ to attain the greatest reduction in ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}.$ In case (1)(1) $${\left(r{\sigma }_r\right)}_{,r}+{\tau }_{r\theta ,\theta }-{\sigma }_{\theta }=0$$(1) , $n_{\theta }=2$ produced the larger reduction in ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ ($\sim 63\%$) compared to the other values. This matches the patterns of the loads, where $n_p=n_{\theta }=2.$ In the second case, also, $n_{\theta }=2$ is preferable despite changing $n_p$ to $1$ and keeping $n_T=2.$ By moving to the third case, it is seen that $n_{\theta }=1$ is better and yielded $48\%$ drop in ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}.$ From these three cases, it can come to a mind that having $n_{\theta }=n_T$ would produce improved ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}.$ But, this is not constantly true, and this is proven in the fourth case, where $n_p=2n_T=2.$ The favorable $n_{\theta }$ is found to be at $\pm 5.$ This stark disagreement in such case with the three previous cases initiates due to the difference in both $p_0$ and $T_0.$ Therefore, it can be understood that they with $n_p$ and $n_{\theta }$ have strategic roles in determining the ideal value of $n_{\theta }$ to be used.

Figure 17. Variation of ${\overline{\sigma }}_{\mathrm{V}{\mathrm{M}}_{\max }}$ with $n_{\theta }$ using different values of $n_p,$ $n_T,$ $p_0,$ and $T_0.$

5. Conclusion

In this study, a FE scheme was developed to check the efficiency of using 2D-FGMs in alleviating the elastic/thermoelastic stresses of FGCs. The power-law model was used to describe materials’ property variation radially, and it was modified to make properties dependent on both the radial and tangential directions concurrently.

Interpreting results showed the beneficial role of considering 2D-FGMs in reducing most of the stress components through the cylinder, which allows for raising the load limits safely. Eminent findings upon applying 2D-FGMs are listed below:

• The tangential stress has the upper hand in determining the behaviors of cylinders. In the second place comes the axial stress that arises in plane strain conditions.

• In case of mechanical loading: the radial, tangential, and shear stress components experienced declines in their values.

• In case of thermomechanical loading: the radial and shear stresses experienced variations in their readings. Conversely, other stresses decreased significantly. However, a decline in the von Mises stress was achieved with substantial values.

• Temperature in many regions of the cylinder tended to decline when the tangential variation of the thermal conductivity was considered, and this was the major cause for the variation of stresses.

• The sign of the shifting angle almost had no impact in the case of mechanical loading. In contrast, it had a great role when both thermal and mechanical loads were present together.

• The angular coefficient had significant impacts on the stresses’ variation.

Eventually, different percentages of stresses’ reductions are to be obtained if the gradation function is changed with its parameters (i.e., heterogeneity (grading) index). Therefore, it is recommended to optimize the angular coefficient and the shifting angle to obtain the global optimum of the stress components. This would be compulsory in case of complex loading and mixed boundary conditions. However, if 2D-FGMs are used without optimizing its associated parameters, inevitable drops of the stresses would occur.

Data availability

The raw/processed data required to reproduce these findings cannot be shared at this time as the data also forms part of an ongoing study.

Disclosure statement

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Additional information

Funding

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.