Bending Analysis of Functionally Graded Axisymmetric Circular Plates using the Dual Mesh Finite Domain Method

Abstract

The dual mesh finite domain method (DMFDM) introduced by Reddy employs one mesh for the approximation of the primary variables (primal mesh) and another mesh for the satisfaction of the governing equations (dual mesh). The present study deals with the extension and application of the DMFDM to functionally graded circular plates under axisymmetric conditions. The formulation makes use of the traditional finite element interpolation of the primary variables with a primal mesh and a dual mesh to satisfy the integral form of the governing differential equations, the basic premise of the finite volume method. The method is used to analyze axisymmetric bending of through-thickness functionally graded circular plates using the classical plate theory (CPT) and first-order shear deformation plate theory (FST). The displacement model of the FST and the mixed model of the CPT using the DMFDM are developed along with the displacement and mixed finite element models. Numerical results are presented to illustrate the methodology and a comparison of the generalized displacements and bending moments computed with those of the corresponding finite element models. The influence of the extensional-bending coupling stiffness (due to the through-thickness grading of the material) on the deflections is also brought out.

Keywords:
Dual mesh finite domain method; axisymmetric bending of plates; mixed models; functionally graded materials; numerical results

Graphical Abstract

1. Introduction

1.1 Functionally graded structures

Structures whose material properties are continuously varied in one or more coordinate directions are abundant in nature (e.g., sea shells, bones, etc). With increasing demand for novel materials whose properties can be tailored to suit a particular application, such natural structures provided an inspiration for the design and development of an advanced class of composite materials called functionally graded materials (FGMs) (see, e.g., [¹,²]). This is usually done by providing in-depth graded compositions, microstructure and properties [³]. With advent of advanced manufacturing and material processing techniques, the applications of FGMs have only increased [⁴]. For example, FGMs find applications in space structures as thermal barriers and in turbine rotors, flywheels, and gears as thin coatings, among other applications. While any number of materials can be used to achieve functional gradation, two-constituent functionally graded materials are the most widely used in engineering applications.

The variation in material properties of FGM can be expressed by a thickness function describing the spatial variation of volume fractions of FGM constituents. Various such functions are proposed in the literature (see, e.g., [⁵] for a review). For structural elements like beams and plates, the power-law variation of the modulus of elasticity, while keeping the Poisson’s ratio constant, can be used. If the $z$-coordinate is taken along the thickness of the plate, the variation of Young’s modulus is given by (see Reddy [⁶])

$\begin{array}{r}E\left(z\right)=(E_t-E_b)f\left(z\right)+E_b,f\left(z\right)={\left(\frac{1}{2}+\frac{z}{H}\right)}^n\end{array}$ (1.1)

where $E_t$ and $E_b$ are the material properties of the top and bottom faces of the plate, respectively, $n$ is the power-law index, and $H$ is the plate thickness. Note that when $n=0$, we obtain the single-material plate (with modulus $E_t$).

1.2 Dual mesh finite domain method (DMFDM)

Recently, Reddy [⁷] introduced a numerical approach termed the dual mesh finite domain method (DMFDM) for the solution of second-order differential equations in one and two dimensions with a single unknown. In the DMFDM, the domain of the governing equations is represented with a primal mesh of finite elements to introduce an approximation of the dependent variables, and a dual mesh of finite domains is superimposed on the primal mesh such that the finite element nodes are at the center of the finite domains, except for the nodes on the boundary, where half or quarter finite domains cover the boundary (see [⁷] for the details). Then the governing equations are satisfied in an integral sense over the finite (control) domain, as in the finite volume method (FVM). The second-order terms in the differential equation are integrated by parts (and thereby weakening the differentiability) and expressed as quantities (dual variables) on the interfaces of the dual mesh. Thus, the DMFDM can be viewed as a hybrid method that makes use of two best features of the FEM [⁸], namely, (a) the interpolation of the variables and (b) imposition of physical boundary conditions, and two salient features of the FVM [⁹]: (a) satisfaction of the global balance equations over the finite domain and (b) computation of the secondary variables at the boundaries of the finite domains where they are uniquely defined.

1.3 Present study

In [⁷], the DMFDM was introduced for second-order differential equations involving a single variable, with applications to heat transfer-like problems in one and two dimensions (also see [¹⁰]). Recently, Reddy and Nampally [¹¹] and Reddy, Nampally, and Srinivasa [¹²] extended the method to linear and nonlinear analysis, respectively, of beams with multiple variables (i.e., to systems of coupled second-order differential equations in several dependent unknowns). In this paper, the method is extended to the bending of through-thickness functionally graded circular plates under axisymmetric loads and boundary conditions. Numerical results are presented to illustrate the accuracy of the DMFDM in predicting the deflections and stresses when compared with the finite element solutions. The effect of the extensional-bending stiffness on the deflections is also elucidated.

Following this introduction, a review of the governing equations of functionally graded axisymmetric circulate plates is given in section 2. A reformulation of the governing equations of circular plates (see Reddy [¹³, ¹⁴] for the derivation of the governing equations) as applied to functionally graded plates is presented in section 3. For the mixed model of the CPT, the governing equations are expressed in terms of the generalized displacements and bending moment. The dual mesh finite domain discretization of these equations are presented in section 4. Although the formulations are simple, they are not readily available in the literature. Numerical results are presented in section 5 for hinged and clamped boundary conditions, and the numerical results obtained with the FEM and DMFDM models are compared with the exact solutions to assess the relative accuracy in predicting the generalized displacements and generalized forces. Finally, some remarks about DMFDM and its extensions are outlined in section 6.

2. Governing equations

Consider a through-thickness functionally graded circular plate of thickness $H$ and radius $R$ subjected to axisymmetric distributed load $q$ (i.e., independent of the angular coordinate, $\theta$) on the top face. If further the boundary conditions are also selected to be axisymmetric, then the plate can be considered as functionally graded axisymmetric circular plate. We select the cylindrical coordinate system $(r,\theta ,z)$ to analyze the plate, where $r$ is the radial coordinate outward form the center of the plate $(0\le r\le R)$, $z$ denotes the transverse coordinate $(-H/2\le z\le H/2)$, and $\theta$ is the angular coordinate $(0\le \theta \le 2\pi )$ (see figure 1). The displacement field of such a circular plate based on the first-order shear deformation plate theory (FST) are given by

$u_r\left(r,\theta ,z\right)=u\left(r\right)+z\varphi \left(r\right)$ (2.1a)

$u_{\theta }\left(r,\theta ,z\right)=0$ (2.1b)

$u_z\left(r,\theta ,z\right)=w\left(r\right)$ (2.1c)

The non-zero linear strains based on the above displacements would then be

${\epsilon }_{rr}=\frac{du}{dr}+z\frac{d\varphi }{dr}$ (2.2a)

${\epsilon }_{\theta \theta }=\frac{u}{r}+\frac{z}{r}\varphi$ (2.2b)

${\epsilon }_{rz}=\frac{1}{2}\left(\varphi +\frac{dw}{dr}\right)$ (2.2c)

while the constitutive relations of the functionally graded axisymmetric circular plate are given by

$\left\{\begin{array}{c}{\sigma }_{rr}\\ {\sigma }_{\theta \theta }\\ {\sigma }_{rz}\end{array}\right\}=\frac{E\left(z\right)}{1-{\nu }^2}\left[\begin{array}{ccc}1& \nu & 0\\ \nu & 1& 0\\ 0& 0& 1-\nu \end{array}\right]\left\{\begin{array}{c}{\epsilon }_{rr}\\ {\epsilon }_{\theta \theta }\\ {\epsilon }_{rz}\end{array}\right\}$ (2.3)

Figure 1:
Schematic representation of the axisymmetric circular plate and various stress resultants acting on a differential element of the axisymmetric circular plate.

We further define the stress and moment resultants on the circular plate (see figure 1) as follows:

$N_{rr}={\int }_{-H/2}^{H/2}{\sigma }_{rr}dz=A\left[\frac{du}{dr}+\nu \frac{u}{r}\right]+B\left[\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right]$ (2.4a)

$M_{rr}={\int }_{-H/2}^{H/2}{\sigma }_{rr}zdz=B\left[\frac{du}{dr}+\nu \frac{u}{r}\right]+D\left[\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right]$ (2.4b)

$N_{\theta \theta }={\int }_{-H/2}^{H/2}{\sigma }_{\theta \theta }dz=A\left[\nu \frac{du}{dr}+\frac{u}{r}\right]+B\left[\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right]$ (2.4c)

$M_{\theta \theta }={\int }_{-H/2}^{H/2}{\sigma }_{\theta \theta }zdz=B\left[\nu \frac{du}{dr}+\frac{u}{r}\right]+D\left[\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right]$ (2.4d)

$Q_r={\int }_{-H/2}^{H/2}{\sigma }_{rz}dz=S\left(\varphi +\frac{dw}{dr}\right)$ (2.4e)

where

$\begin{array}{rr}& A={\int }_{-H/2}^{H/2}\frac{E\left(z\right)}{1-{\nu }^2}dz,B={\int }_{-H/2}^{H/2}\frac{E\left(z\right)z}{1-{\nu }^2}dz\\ & D={\int }_{-H/2}^{H/2}\frac{E\left(z\right)z^2}{1-{\nu }^2}dz,S={\int }_{-H/2}^{H/2}\frac{K_{s}E\left(z\right)}{2(1+\nu )}dz\end{array}$ (2.5)

The explicit expressions of the above integrals for the power-law variation of $E\left(z\right)$ (see Eq. (1.1)) are given as

$\begin{array}{rr}& A=\frac{(E_t-E_b)H}{(n+1)(1-{\nu }^2)}+\frac{E_{b}H}{1-{\nu }^2}\\ & B=\frac{(E_t-E_b)H^2}{1-{\nu }^2}\left[\frac{1}{n+2}-\frac{1}{2(n+1)}\right]\\ & D=\frac{(E_t-E_b)H^3}{1-{\nu }^2}\left[\frac{1}{n+3}+\frac{1}{4(n+1)}-\frac{1}{n+2}\right]+\frac{E_b{H}^3}{12(1-{\nu }^2)}\\ & S=\frac{K_s(E_t-E_b)H}{2(n+1)(1+\nu )}+\frac{K_s{E}_{b}H}{2(1+\nu )}\end{array}$ (2.6)

The equations of equilibrium of the functionally graded (FG) axisymmetric circular plate can be obtained using the principle of virtual work. The virtual work statement for FG axisymmetric circular plate is

${\int }_0^{2\pi }{\int }_0^R{\int }_{-H/2}^{H/2}\left({\sigma }_{rr}\delta {\epsilon }_{rr}+{\sigma }_{\theta \theta }\delta {\epsilon }_{\theta \theta }+2{\sigma }_{rz}\delta {\epsilon }_{rz}\right)rdzdrd\theta -{\int }_0^{2\pi }{\int }_0^{R}q\delta wrdrd\theta =0$ (2.7)

Now using the stress and moment resultant definitions of Eqs. (2.4a)-(2.4e) and strain definition of Eqs. (2.2a)-(2.2c) we can rewrite the virtual work statement as

$\begin{array}{r}{\int }_0^R\left[N_{rr}\left(\frac{d\delta u}{dr}\right)+M_{rr}\left(\frac{d\delta \varphi }{dr}\right)+\frac{N_{\theta \theta }}{r}\delta u+\frac{M_{\theta \theta }}{r}\delta \varphi +Q_r\left(\delta \varphi +\frac{d\delta w}{dr}\right)-q\delta w\right]rdr=0\end{array}$ (2.8)

The governing differential equations of the FG axisymmetric circular plate are obtained by taking the Euler-Lagrange equations of the above variational statement. The resulting governing equations read

$-\frac{1}{r}\frac{d}{dr}\left(r{N}_{rr}\right)+\frac{N_{\theta \theta }}{r}=0$ (2.9)

$-\frac{1}{r}\frac{d}{dr}\left(r{Q}_r\right)-q=0$ (2.10)

$-\frac{1}{r}\left[\frac{d}{dr}\left(r{M}_{rr}\right)-M_{\theta \theta }\right]+Q_r=0$ (2.11)

The governing equations derived above are based on the first-order shear deformation theory, where the transverse shear strain is non-zero and assumes a constant value (and hence require the shear correction factor, $K_s$ in Eq. (2.5)). However, as the radius to thickness ratio of the circular plate increases the shear strains are known to tend to zero. In such cases we may also use classical plate theory (i.e., based on Kirchhoff hypothesis) to predict the plate behavior. To accommodate zero transverse shear strains the value of $\varphi$ in FST is taken to be $-\frac{dw}{dr}$ and hence the displacement field of FG axisymmetric circular plate based on classical plate theory (CPT) would become

$u_r\left(r,\theta ,z\right)=u\left(r\right)-z\frac{dw}{dr}$ (2.12a)

$u_{\theta }\left(r,\theta ,z\right)=0$ (2.12b)

$u_z\left(r,\theta ,z\right)=w\left(r\right)$ (2.12c)

The stress and moment resultants on the FG axisymmetric circular plate based on the classical plate theory would then be

$N_{rr}={\int }_{-H/2}^{H/2}{\sigma }_{rr}dz=A\left[\frac{du}{dr}+\nu \frac{u}{r}\right]-B\left[\frac{d^{2}w}{d{r}^2}+\frac{\nu }{r}\frac{dw}{dr}\right]$ (2.13a)

$M_{rr}={\int }_{-H/2}^{H/2}{\sigma }_{rr}zdz=B\left[\frac{du}{dr}+\nu \frac{u}{r}\right]-D\left[\frac{d^{2}w}{d{r}^2}+\frac{\nu }{r}\frac{dw}{dr}\right]$ (2.13b)

$N_{\theta \theta }={\int }_{-H/2}^{H/2}{\sigma }_{\theta \theta }dz=A\left[\nu \frac{du}{dr}+\frac{u}{r}\right]-B\left[\nu \frac{d^{2}w}{d{r}^2}+\frac{1}{r}\frac{dw}{dr}\right]$ (2.13c)

$M_{\theta \theta }={\int }_{-H/2}^{H/2}{\sigma }_{\theta \theta }zdz=B\left[\nu \frac{du}{dr}+\frac{u}{r}\right]-D\left[\nu \frac{d^{2}w}{d{r}^2}+\frac{1}{r}\frac{dw}{dr}\right]$ (2.13d)

In this case the transverse shear strain vanishes; however, the governing equations remain the same as those listed in Eqs. (2.9)-(2.11) with the exception that the shear stress resultant $Q_r$ can only be expressed in terms of the moments $M_{rr}$ and $M_{\theta \theta }$ as given by Eq. (2.11) rather than directly computing in terms of shear strain (which is zero in this case).

The governing differential equations of FG axisymmetric circular plate based on the FST can be expressed in terms of the generalized displacements ($u$, $w$, $\varphi$) by substituting the expressions of stress and moment resultants of Eqs. (2.4a)-(2.4e) into Eqs. (2.9)-(2.11). The resulting equations are listed below.

$-\frac{1}{r}\frac{d}{dr}\left[rA\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rB\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]$ $+\frac{A}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{B}{r}\left(\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right)=0$ (2.14)

$-\frac{1}{r}\frac{d}{dr}\left[rS\left(\varphi +\frac{dw}{dr}\right)\right]-q=0$ (2.15)

$-\frac{1}{r}\frac{d}{dr}\left[rB\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rD\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]$ $+\frac{B}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{D}{r}\left(\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right)+S\left(\varphi +\frac{dw}{dr}\right)=0$ (2.16)

Similarly, the governing differential equations of FG axisymmetric circular plate based on CPT can be expressed in terms of the displacements $u$ and $w$ by first expressing Eq. (2.10) in terms of moments $M_{rr}$ and $M_{\theta \theta }$ using Eq. (2.11) as

$-\frac{1}{r}\left[\frac{d}{dr}\left(r{M}_{rr}\right)-M_{\theta \theta }\right]-q=0$ (2.17)

and then substituting the stress and moment resultants of Eqs. (2.13a)-(2.13d) into the above equation and Eq. (2.9) to give

$-\frac{1}{r}\frac{d}{dr}\left[rA\left(\frac{du}{dr}+\nu \frac{u}{r}\right)-rB\left(\frac{d^{2}w}{d{r}^2}+\frac{\nu }{r}\frac{dw}{dr}\right)\right]$ $+\frac{A}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)-\frac{B}{r}\left(\nu \frac{d^{2}w}{d{r}^2}+\frac{1}{r}\frac{dw}{dr}\right)=0$ (2.18)

$-\frac{1}{r}\frac{d}{dr}\left[\frac{d}{dr}\left[rB\left(\frac{du}{dr}+\nu \frac{u}{r}\right)-rD\left(\frac{d^{2}w}{d{r}^2}+\frac{\nu }{r}\frac{dw}{dr}\right)\right]\right.$ $\left.-B\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+D\left(\nu \frac{d^{2}w}{d{r}^2}+\frac{1}{r}\frac{dw}{dr}\right)\right]-q=0$ (2.19)

3. Equations in terms of displacements and bending moment

3.1 Classical plate theory

The governing equations of FG axisymmetric circular plate based on the classical plate theory when expressed in terms of the displacements, $u$ and $w$, would result in fourth order differential equations. However, since the dual mesh finite domain method is applicable only for first-order or second-order differential equations, we will recast the governing equations (2.9)-(2.11) as second-order differential equations such that the primary variables of the resulting equations are $\{u,w,M_{rr}\}$. First, we will consider the equations (2.13a) and (2.13b). From these two equations we can write an expression of $N_{rr}$ in terms of $u$ and $M_{rr}$. Similarly, $N_{\theta \theta }$ can be expressed in terms of $u$ and $M_{\theta \theta }$ using the equations (2.13c) and (2.13d). In the following we list the resulting equations after some algebraic manipulations.

$N_{rr}=\overline{A}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+\overline{B}{M}_{rr}$ (3.1)

$N_{\theta \theta }=\overline{A}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\overline{B}{M}_{\theta \theta }$ (3.2)

where

$\begin{array}{c}\overline{A}=\frac{D^{\mathrm{*}}}{D}=\frac{AD-B^2}{D},\overline{B}=\frac{B}{D}\end{array}$

To eliminate $M_{\theta \theta }$ from Eq. (3.2) we use Eq. (2.13b) and Eq. (2.13d) to obtain

$M_{\theta \theta }=\nu M_{rr}+B(1-{\nu }^2)\frac{u}{r}-\frac{D(1-{\nu }^2)}{r}\frac{dw}{dr}$ (3.3)

Thus $N_{\theta \theta }$ can be written as

$N_{\theta \theta }=\overline{A}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\overline{B}\nu M_{rr}+\frac{B^2(1-{\nu }^2)}{D}\frac{u}{r}-\frac{B\left(1-{\nu }^2\right)}{r}\frac{dw}{dr}$ (3.4)

Finally, Eq. (2.13b) can be expressed as

$-\left(\frac{d^{2}w}{d{r}^2}+\frac{\nu }{r}\frac{dw}{dr}\right)+\overline{B}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)-\frac{1}{D}{M}_{rr}=0$ (3.5)

The governing equations of FG axisymmetric circular plate based on classical plate theory in terms of $\{u,w,M_{rr}\}$ would then be given by Eq. (2.9), Eq. (2.17) and Eq. (3.5); wherein the expressions for $N_{rr}$, $N_{\theta \theta }$ and $M_{\theta \theta }$ are obtained from Eqs. (3.1), (3.4) and (3.3) respectively. The resulting final equations are as follows:

$-\frac{1}{r}\frac{d}{dr}\left[r\overline{A}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+r\overline{B}{M}_{rr}\right]+\frac{\overline{A}}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)$ $+\frac{\overline{B}\nu }{r}{M}_{rr}-\frac{B\left(1-{\nu }^2\right)}{r^2}\frac{dw}{dr}+\frac{B^2\left(1-{\nu }^2\right)}{D{r}^2}u=0$ (3.6)

$-\frac{1}{r}\frac{d}{dr}\left[r\frac{d{M}_{rr}}{dr}+\left(1-\nu \right)M_{rr}-\frac{B\left(1-{\nu }^2\right)}{r}u+\frac{D\left(1-{\nu }^2\right)}{r}\frac{dw}{dr}\right]-q=0$ (3.7)

$-\frac{1}{r}\frac{d}{dr}\left(r\frac{dw}{dr}\right)+\frac{\left(1-\nu \right)}{r}\frac{dw}{dr}+\overline{B}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)-\frac{1}{D}{M}_{rr}=0$ (3.8)

4. Dual mesh finite domain formulations

4.1 Mixed DMFDM for CPT

As noted earlier, the dual mesh finite domain method is only applicable to equations of second order or less. Thus we develop dual mesh finite domain formulation of the FG axisymmetric circular plate based on classical plate theory using the governing equations (3.6)-(3.8). The resulting formulation is called mixed dual mesh finite domain method (see [¹¹]) since the formulation has both generalized displacements $(u,w)$ and the moment $\left(M_{rr}\right)$ as the primary variables. This is similar to the mixed formulations of finite element method.

Figure 2:
Primal mesh (finite elements) and dual mesh (finite domains) on the computational domain.

To obtain the discretized equations, we consider the computational domain $\Omega =(0,R)$ and divide it into $N$ finite elements (primal mesh) such that each node has its associated finite domain (dual mesh). Except for the finite domains corresponding to the boundary nodes, all the interior finite domains encompass two finite elements (primal mesh) such that one half of each of these two finite elements lie within the finite domain (see figure 2). Thus, for a uniform primal mesh the interior nodes lie at the center of their corresponding finite domains.

The dependent variables are approximated on each finite element using Lagrange interpolation functions (see figure 3). Although different sets of interpolation functions can be used for different primary variables, in the present paper we use same set of linear 1-D Lagrange interpolation functions for all the primary variables. For a typical finite element ${\Omega }^{\left(I\right)}=(r_I,r_{I+1})$, the primary variables can thus be approximated using the set of linear 1-D Lagrange interpolation functions $\{L_1^{\left(I\right)},L_2^{\left(I\right)}\}$ as

$\begin{array}{rr}& u\left(r\right)\approx U_I{L}_1^{\left(I\right)}\left(r\right)+U_{I+1}{L}_2^{\left(I\right)}\left(r\right),\\ & w\left(r\right)\approx W_I{L}_1^{\left(I\right)}\left(r\right)+W_{I+1}{L}_2^{\left(I\right)}\left(r\right)\\ & M_{rr}\left(r\right)\approx M_{r{r}_I}{L}_1^{\left(I\right)}\left(r\right)+M_{r{r}_{I+1}}{L}_2^{\left(I\right)}\left(r\right)\end{array}$ (4.1)

Here $U_I$, $W_I$, $M_{r{r}_I}$ represent the values of $u$, $w$ and $M_{rr}$ respectively at node $I$, while $U_{I+1}$, $W_{I+1}$, $M_{r{r}_{I+1}}$ represent the values of $u$, $w$ and $M_{rr}$ respectively at node $I+1$. It should be noted that when Eq. (2.18) and Eq. (2.19) are used in developing finite elements, a minimum of cubic Hermite interpolation is needed on the transverse deflection $w$ and linear Lagrange interpolation is needed on $u$ [⁸]. However, for the mixed formulations, akin to the one considered here, a minimum of linear Lagrange interpolations on all the primary variables would suffice [¹¹]. The linear Lagrange interpolation functions $L_1^{\left(I\right)}$ and $L_2^{\left(I\right)}$ on a typical finite element ${\Omega }^{\left(I\right)}$ are given by

$\begin{array}{r}{L}_1^{\left(I\right)}=\frac{r_{I+1}-r}{h_I},L_1^{\left(I\right)}=\frac{r-r_I}{h_I}\mathrm{w}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{e}{h}_I=r_{I+1}-r_I\end{array}$

Figure 3:
Lagrange interpolation functions on the finite elements Ω^(I−1) and Ω^(I) for variable u.

Figure 4:
A typical finite domain for the mixed model of CPT with the secondary variables at the faces and the primary variables at the node associated with the finite domain.

Now we proceed to develop the discretized equations corresponding to Eqs. (3.6)-(3.8) by writing their integral statements on a typical interior finite domain $(r_a^{\left(I\right)},r_b^{\left(I\right)})$, associated with node $I$ (see figure 4).

${\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{-\frac{1}{r}\frac{d}{dr}\left[r\overline{A}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+r\overline{B}{M}_{rr}\right]+\frac{\overline{A}}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)\right.$ $+\left.\frac{\overline{B}\nu }{r}{M}_{rr}-\frac{B(1-{\nu }^2)}{r^2}\frac{dw}{dr}+\frac{B^2(1-{\nu }^2)}{D{r}^2}u\right\}rdr=0$ (4.2)

${\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{-\frac{1}{r}\frac{d}{dr}\left[\frac{d}{dr}\left(r{M}_{rr}\right)-\nu M_{rr}-\frac{B(1-{\nu }^2)}{r}u+\frac{D(1-{\nu }^2)}{r}\frac{dw}{dr}\right]-q\right\}rdr=0$ (4.3)

${\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{-\left(\frac{d^{2}w}{d{r}^2}+\frac{\nu }{r}\frac{dw}{dr}\right)+\overline{B}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)-\frac{1}{D}{M}_{rr}\right\}rdr=0$ (4.4)

The first step in developing the discretized equations is to carry out the indicated integration of Eqs. (4.2)-(4.4) such that the resulting boundary terms represent physically meaningful secondary variables which are dual to the primary variables. Eqs. (4.2)-(4.4) after such integration are listed below.

$-N_1^{\left(I\right)}-N_2^{\left(I\right)}+{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{\frac{\overline{A}}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{\overline{B}\nu }{r}{M}_{rr}\right.$ $\left.-\frac{B(1-{\nu }^2)}{r^2}\frac{dw}{dr}+\frac{B^2(1-{\nu }^2)}{D{r}^2}u\right\}rdr=0$ (4.5)

$-V_1^{\left(I\right)}-V_2^{\left(I\right)}-{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}qrdr=0$ (4.6)

$-{\Phi }_1^{\left(I\right)}-{\Phi }_2^{\left(I\right)}+{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{\frac{1-\nu }{r}\frac{dw}{dr}+\overline{B}\left(\frac{du}{dr}+\frac{\nu }{r}u\right)-\frac{1}{D}{M}_{rr}\right\}rdr=0$ (4.7)

In the above equations the secondary variables are defined as follows

$N_1^{\left(I\right)}=-{\left[r\overline{A}\left(\frac{du}{dr}+\frac{\nu }{r}u\right)+r\overline{B}{M}_{rr}\right]}_{r=r_a^{\left(I\right)}}$ (4.8a)

$N_2^{\left(I\right)}={\left[r\overline{A}\left(\frac{du}{dr}+\frac{\nu }{r}u\right)+r\overline{B}{M}_{rr}\right]}_{r=r_b^{\left(I\right)}}$ (4.8b)

$V_1^{\left(I\right)}=-{\left[\frac{d}{dr}\left(r{M}_{rr}\right)-\nu M_{rr}-\frac{B\left(1-{\nu }^2\right)}{r}u+\frac{D\left(1-{\nu }^2\right)}{r}\frac{dw}{dr}\right]}_{r=r_a^{\left(I\right)}}$ (4.9a)

$V_2^{\left(I\right)}={\left[\frac{d}{dr}\left(r{M}_{rr}\right)-\nu M_{rr}-\frac{B(1-{\nu }^2)}{r}u+\frac{D(1-{\nu }^2)}{r}\frac{dw}{dr}\right]}_{r=r_b^{\left(I\right)}}$ (4.9b)

${\Phi }_1^{\left(I\right)}={\left[-r\frac{dw}{dr}\right]}_{r=r_a^{\left(I\right)}}$ (4.10a)

${\Phi }_2^{\left(I\right)}={\left[r\frac{dw}{dr}\right]}_{r=r_b^{\left(I\right)}}$ (4.10b)

Here ($N_1^{\left(I\right)}$, $V_1^{\left(I\right)}$, ${\Phi }_1^{\left(I\right)}$) represent the axial force, shear force and rotation respectively at the left face of the finite domain numbered $I$ (i.e., $r=r_a^{\left(I\right)}$), while ($N_2^{\left(I\right)}$, $V_2^{\left(I\right)}$, ${\Phi }_2^{\left(I\right)}$) represent axial force, shear force and rotation respectively at right face of the finite domain numbered $I$ (i.e., $r=r_b^{\left(I\right)}$) (see figure 4). Since all the primary variables are approximated using linear Lagrange interpolation functions on each finite element, the boundary terms given in Eqs. (4.8a)-(4.10b) can be easily evaluated (see Appendix A), while the integrations of Eqs. (4.5)-(4.7) can be either evaluated numerically (see, e.g., Gauss quadrature rule) or exactly. In the present paper we use Gauss quadrature to evaluate the integral expressions. Thus for an interior node $I$ associated with interior finite domain numbered $I$ the discretized equations would read

$U_{I-1}\left[-\frac{r_a^{\left(I\right)}{\overline{A}}_{I-1}}{h_{I-1}}+\frac{{\overline{A}}_{I-1}{\nu }_{I-1}}{2}+{\int }_{r_a^{\left(I\right)}}^{r_I}\left({\overline{A}}_{I-1}{\nu }_{I-1}\frac{d{L}_1^{(I-1)}}{dr}+\frac{{\overline{A}}_{I-1}}{r}{L}_1^{(I-1)}\right.\right.$ $\left.\left.+\frac{B_{I-1}^2(1-{\nu }_{I-1}^2)}{r{D}_{I-1}}{L}_1^{(I-1)}\right)dr\right]+W_{I-1}\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left(\frac{B_{I-1}({\nu }_{I-1}^2-1)}{r}\frac{d{L}_1^{(I-1)}}{dr}\right)dr\right]$ $+M_{r{r}_{I-1}}\left[\frac{r_a^{\left(I\right)}{\overline{B}}_{I-1}}{2}+{\int }_{r_a^{\left(I\right)}}^{r_I}\left({\overline{B}}_{I-1}{\nu }_{I-1}{L}_1^{(I-1)}\right)dr\right]+U_I\left[\frac{r_a^{\left(I\right)}{\overline{A}}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}\overline{A_I}}{h_I}\right.$ $-\frac{{\overline{A}}_I{\nu }_I-{\overline{A}}_{I-1}{\nu }_{I-1}}{2}+{\int }_{r_a^{\left(I\right)}}^{r_I}\left({\overline{A}}_{I-1}{\nu }_{I-1}\frac{d{L}_2^{(I-1)}}{dr}+\frac{{\overline{A}}_{I-1}}{r}{L}_2^{(I-1)}\right)dr$ $+{\int }_{r_I}^{r_b^{\left(I\right)}}\left({\overline{A}}_I{\nu }_I\frac{d{L}_1^{\left(I\right)}}{dr}+\frac{{\overline{A}}_I}{r}{L}_1^{\left(I\right)}\right)dr+{\int }_{r_a^{\left(I\right)}}^{r_I}\left(\frac{B_{I-1}^2(1-{\nu }_{I-1}^2)}{r{D}_{I-1}}{L}_2^{(I-1)}\right)dr$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}\left(\frac{B_I^2(1-{\nu }_I^2)}{r{D}_I}{L}_1^{\left(I\right)}\right)dr\right]+W_I\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left(\frac{B_{I-1}({\nu }_{I-1}^2-1)}{r}\frac{d{L}_2^{(I-1)}}{dr}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}\left(\frac{B_I({\nu }_I^2-1)}{r}\frac{d{L}_1^{\left(I\right)}}{dr}\right)dr\right]+M_{r{r}_I}\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left({\overline{B}}_{I-1}{\nu }_{I-1}{L}_2^{(I-1)}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}\left({\overline{B}}_I{\nu }_I{L}_1^{\left(I\right)}\right)dr-\frac{r_b^{\left(I\right)}{\overline{B}}_I-r_a^{\left(I\right)}{\overline{B}}_{I-1}}{2}\right]$ $+U_{I+1}\left[-\frac{r_b^{\left(I\right)}{\overline{A}}_I}{h_I}-\frac{{\overline{A}}_I{\nu }_I}{2}+{\int }_{r_I}^{r_b^{\left(I\right)}}\left({\overline{A}}_I{\nu }_I\frac{d{L}_2^{\left(I\right)}}{dr}+\frac{{\overline{A}}_I}{r}{L}_2^{\left(I\right)}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}\left(\frac{B_I^2(1-{\nu }_I^2)}{r{D}_I}{L}_2^{\left(I\right)}\right)dr\right]+W_{I+1}\left[{\int }_{r_I}^{r_b^{\left(I\right)}}\left(\frac{B_I({\nu }_I^2-1)}{r}\frac{d{L}_2^{\left(I\right)}}{dr}\right)dr\right]$ $+M_{r{r}_{I+1}}\left[{\int }_{r_I}^{r_b^{\left(I\right)}}\left({\overline{B}}_I{\nu }_I{L}_2^{\left(I\right)}\right)dr-\frac{r_b^{\left(I\right)}{\overline{B}}_I}{2}\right]=0$ (4.11)

$\begin{array}{rr}& U_{I-1}\left[-\frac{B_{I-1}(1-{\nu }_{I-1}^2)}{2r_a^{\left(I\right)}}\right]+W_{I-1}\left[-\frac{D_{I-1}(1-{\nu }_{I-1}^2)}{r_a^{\left(I\right)}{h}_{I-1}}\right]+M_{r{r}_{I-1}}\left[\frac{1}{2}-\frac{r_a^{\left(I\right)}}{h_{I-1}}-\frac{{\nu }_{I-1}}{2}\right]\\ & +U_I\left[\frac{B_I(1-{\nu }_I^2)}{2r_b^{\left(I\right)}}-\frac{B_{I-1}(1-{\nu }_{I-1}^2)}{2r_a^{\left(I\right)}}\right]+W_I\left[\frac{D_{I-1}(1-{\nu }_{I-1}^2)}{r_a^{\left(I\right)}{h}_{I-1}}+\frac{D_I(1-{\nu }_I^2)}{r_b^{\left(I\right)}{h}_I}\right]\\ & +M_{r{r}_I}\left[\frac{r_a^{\left(I\right)}}{h_{I-1}}+\frac{r_b^{\left(I\right)}}{h_I}+\frac{{\nu }_I-{\nu }_{I-1}}{2}\right]+U_{I+1}\left[\frac{B_I(1-{\nu }_I^2)}{2r_b^{\left(I\right)}}\right]+W_{I+1}\left[-\frac{D_I(1-{\nu }_I^2)}{r_b^{\left(I\right)}{h}_I}\right]\\ & +M_{r{r}_{I+1}}\left[-\frac{1}{2}-\frac{r_b^{\left(I\right)}}{h_I}+\frac{{\nu }_I}{2}\right]-{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}qrdr=0\end{array}$ (4.12)

$U_{I-1}\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left({\overline{B}}_{I-1}r\frac{d{L}_1^{(I-1)}}{dr}+{\overline{B}}_{I-1}{\nu }_{I-1}{L}_1^{(I-1)}\right)dr\right]+W_{I-1}\left[-\frac{r_a^{\left(I\right)}}{h_{I-1}}\right.$ $\left.+{\int }_{r_a^{\left(I\right)}}^{r_I}\left((1-{\nu }_{I-1})\frac{d{L}_1^{(I-1)}}{dr}\right)dr\right]-M_{r{r}_{I-1}}\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left(\frac{r}{D_{I-1}}{L}_1^{(I-1)}\right)dr\right]$ $+U_I\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left({\overline{B}}_{I-1}r\frac{d{L}_2^{(I-1)}}{dr}+{\overline{B}}_{I-1}{\nu }_{I-1}{L}_2^{(I-1)}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}\left({\overline{B}}_{I}r\frac{d{L}_1^{\left(I\right)}}{dr}+{\overline{B}}_I{\nu }_I{L}_1^{\left(I\right)}\right)dr\right]+W_I\left[\frac{r_a^{\left(I\right)}}{h_{I-1}}+{\int }_{r_a^{\left(I\right)}}^{r_I}\left((1-{\nu }_{I-1})\frac{d{L}_2^{(I-1)}}{dr}\right)dr\right.$ $\left.+\frac{r_b^{\left(I\right)}}{h_I}+{\int }_{r_I}^{r_b^{\left(I\right)}}\left((1-{\nu }_I)\frac{d{L}_1^{\left(I\right)}}{dr}\right)dr\right]-M_{r{r}_I}\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\left(\frac{r}{D_{I-1}}{L}_2^{(I-1)}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}\left(\frac{r}{D_I}{L}_1^{\left(I\right)}\right)dr\right]+U_{I+1}\left[{\int }_{r_I}^{r_b^{\left(I\right)}}\left({\overline{B}}_{I}r\frac{d{L}_2^{\left(I\right)}}{dr}+{\overline{B}}_I{\nu }_I{L}_2^{\left(I\right)}\right)dr\right]$ $+W_{I+1}\left[-\frac{r_b^{\left(I\right)}}{h_I}+{\int }_{r_I}^{r_b^{\left(I\right)}}\left((1-{\nu }_I)\frac{d{L}_2^{\left(I\right)}}{dr}\right)dr\right]-M_{r{r}_{I+1}}\left[{\int }_{r_I}^{r_b^{\left(I\right)}}\left(\frac{r}{D_I}{L}_2^{\left(I\right)}\right)dr\right]=0$ (4.13)

In the above equations ${\overline{B}}_{I-1}$ represents the value of $\overline{B}$ in the finite element ${\Omega }^{(I-1)}$ and ${\overline{B}}_I$ represents the value of $\overline{B}$ in the finite element ${\Omega }^{\left(I\right)}$. Similar meaning applies to other coefficients $\overline{A}$, $\nu$, $B$ and $D$.

Figure 5:
Secondary variables on the boundary finite domains.

In a similar fashion, the discretized equations for the finite domain 1 (see figure 5) would be

$-N_1^{\left(1\right)}+U_1\left[\frac{{\overline{A}}_1\left(1-{\nu }_1\right)}{2}+{\int }_0^{0.5h_1}\left({\overline{A}}_1{\nu }_1\frac{d{L}_1^{\left(1\right)}}{dr}+\frac{{\overline{A}}_1}{r}{L}_1^{\left(1\right)}+\frac{B_1^2\left(1-{\nu }_1^2\right)}{r{D}_1}{L}_1^{\left(1\right)}\right)dr\right]$ $+W_1\left[{\int }_0^{0.5h_1}\left(\frac{B_1({\nu }_1^2-1)}{r}\frac{d{L}_1^{\left(1\right)}}{dr}\right)dr\right]+M_{r{r}_1}\left[{\int }_0^{0.5h_1}{\overline{B}}_1{\nu }_1{L}_1^{\left(1\right)}dr-\frac{{\overline{B}}_1{h}_1}{4}\right]$ $+U_2\left[-\frac{{\overline{A}}_1(1+{\nu }_1)}{2}+{\int }_0^{0.5h_1}\left({\overline{A}}_1{\nu }_1\frac{d{L}_2^{\left(1\right)}}{dr}+\frac{{\overline{A}}_1}{r}{L}_2^{\left(1\right)}+\frac{B_1^2(1-{\nu }_1^2)}{r{D}_1}{L}_2^{\left(1\right)}\right)dr\right]$ $+W_2\left[{\int }_0^{0.5h_1}\left(\frac{B_1({\nu }_1^2-1)}{r}\frac{d{L}_2^{\left(1\right)}}{dr}\right)dr\right]+M_{r{r}_2}\left[{\int }_0^{0.5h_1}{\overline{B}}_1{\nu }_1{L}_2^{\left(1\right)}dr-\frac{{\overline{B}}_1{h}_1}{4}\right]=0$ (4.14)

$-V_1^{\left(1\right)}+U_1\left[\frac{B_1\left(1-{\nu }_1^2\right)}{h_1}\right]+W_1\left[\frac{2D_1\left(1-{\nu }_1^2\right)}{h_1^2}\right]+M_{r{r}_1}\left[\frac{{\nu }_1}{2}\right]$ $+U_2\left[\frac{B_1(1-{\nu }_1^2)}{h_1}\right]+W_2\left[-\frac{2D_1(1-{\nu }_1^2)}{h_1^2}\right]+M_{r{r}_2}\left[-1+\frac{{\nu }_1}{2}\right]-{\int }_0^{0.5h_1}qrdr=0$ (4.15)

$-{\Phi }_1^{\left(1\right)}+U_1\left[{\int }_0^{0.5h_1}\left({\overline{B}}_{1}r\frac{d{L}_1^{\left(1\right)}}{dr}+{\overline{B}}_1{\nu }_1{L}_1^{\left(1\right)}\right)dr\right]+W_1\left[\frac{1}{2}+{\int }_0^{0.5h_1}\left((1-{\nu }_1)\frac{d{L}_1^{\left(1\right)}}{dr}\right)dr\right]$ $+U_2\left[{\int }_0^{0.5h_1}\left({\overline{B}}_{1}r\frac{d{L}_2^{\left(1\right)}}{dr}+{\overline{B}}_1{\nu }_1{L}_2^{\left(1\right)}\right)dr\right]+W_2\left[-\frac{1}{2}+{\int }_0^{0.5h_1}\left((1-{\nu }_1)\frac{d{L}_2^{\left(1\right)}}{dr}\right)dr\right]$ $-M_{r{r}_1}\left[{\int }_0^{0.5h_1}\left(\frac{r}{D_1}{L}_1^{\left(1\right)}\right)dr\right]-M_{r{r}_2}\left[{\int }_0^{0.5h_1}\left(\frac{r}{D_1}{L}_2^{\left(1\right)}\right)dr\right]=0$ (4.16)

Finally, the discretized equations on ($N+1$)^th finite domain (see figure 5) would be

$-N_2^{(N+1)}+U_N\left[-\frac{(R-0.5h_N){\overline{A}}_N}{h_N}+\frac{{\overline{A}}_N{\nu }_N}{2}\right.$ $\left.+{\int }_{R-0.5h_N}^R\left({\overline{A}}_N{\nu }_N\frac{d{L}_1^{\left(N\right)}}{dr}+\frac{{\overline{A}}_N}{r}{L}_1^{\left(N\right)}+\frac{B_N^2(1-{\nu }_N^2)}{r{D}_N}{L}_1^{\left(N\right)}\right)dr\right]$ $+W_N\left[{\int }_{R-0.5h_N}^R\left(\frac{B_N({\nu }_N^2-1)}{r}\frac{d{L}_1^{\left(N\right)}}{dr}\right)dr\right]$ $+M_{r{r}_N}\left[{\int }_{R-0.5h_N}^R{\overline{B}}_N{\nu }_N{L}_1^{\left(N\right)}dr+\frac{{\overline{B}}_N(R-0.5h_N)}{2}\right]+U_{N+1}\left[\frac{(R-0.5h_N){\overline{A}}_N}{h_N}\right.$ $\left.+\frac{{\overline{A}}_N{\nu }_N}{2}+{\int }_{R-0.5h_N}^R\left({\overline{A}}_N{\nu }_N\frac{d{L}_2^{\left(N\right)}}{dr}+\frac{{\overline{A}}_N}{r}{L}_2^{\left(N\right)}+\frac{B_N^2(1-{\nu }_N^2)}{r{D}_N}{L}_2^{\left(N\right)}\right)dr\right]$ $+W_{N+1}\left[{\int }_{R-0.5h_N}^R\left(\frac{B_N({\nu }_N^2-1)}{r}\frac{d{L}_2^{\left(N\right)}}{dr}\right)dr\right]$ $+M_{r{r}_{N+1}}\left[{\int }_{R-0.5h_N}^R{\overline{B}}_N{\nu }_N{L}_2^{\left(N\right)}dr+\frac{{\overline{B}}_N(R-0.5h_N)}{2}\right]=0$ (4.17)

$-V_2^{(N+1)}+U_N\left[-\frac{B_N(1-{\nu }_N^2)}{2(R-0.5h_N)}\right]+W_N\left[-\frac{D_N(1-{\nu }_N^2)}{h_N(R-0.5h_N)}\right]$ $+M_{r{r}_N}\left[1-\frac{R}{h_N}-\frac{{\nu }_N}{2}\right]+U_{N+1}\left[-\frac{B_N(1-{\nu }_N^2)}{2(R-0.5h_N)}\right]+W_{N+1}\left[\frac{D_N(1-{\nu }_N^2)}{h_N(R-0.5h_N)}\right]$ $+M_{r{r}_{N+1}}\left[\frac{R}{h_N}-\frac{{\nu }_N}{2}\right]-{\int }_{R-0.5h_N}^{R}qrdr=0$ (4.18)

${-\Phi }_2^{(N+1)}+U_N\left[{\int }_{R-0.5h_N}^R\left({\overline{B}}_{N}r\frac{d{L}_1^{\left(N\right)}}{dr}+{\overline{B}}_N{\nu }_N{L}_1^{\left(N\right)}\right)dr\right]$ $+W_N\left[-\frac{R-0.5h_N}{h_N}+{\int }_{R-0.5h_N}^R\left((1-{\nu }_N)\frac{d{L}_1^{\left(N\right)}}{dr}dr\right)\right]-M_{r{r}_N}\left[{\int }_{R-0.5h_N}^R\frac{r}{D_N}{L}_1^{\left(N\right)}dr\right]$ $+U_{N+1}\left[{\int }_{R-0.5h_N}^R\left({\overline{B}}_{N}r\frac{d{L}_2^{\left(N\right)}}{dr}+{\overline{B}}_N{\nu }_N{L}_2^{\left(N\right)}\right)dr\right]+W_{N+1}\left[\frac{R-0.5h_N}{h_N}\right.$ $\left.+{\int }_{R-0.5h_N}^R\left((1-{\nu }_N)\frac{d{L}_2^{\left(N\right)}}{dr}\right)dr\right]-M_{r{r}_{N+1}}\left[{\int }_{R-0.5h_N}^R\frac{r}{D_N}{L}_2^{\left(N\right)}dr\right]=0$ (4.19)

This completes the mixed dual mesh finite domain formulation of FG axisymmetric circular plate based on CPT.

4.2 Displacement DMFDM for FST

Since the governing differential equations of axisymmetric circular plate based on FST in terms of generalized displacements are second-order differential equations, we can developed dual mesh finite domain method directly from the equations (2.14)-(2.16). The integral statements of the governing equations on a typical interior finite domain $(r_a^{\left(I\right)},r_b^{\left(I\right)})$, associated with node $I$, would be

$0={\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{-\frac{1}{r}\frac{d}{dr}\left[Ar\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+Br\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]\right.$ $+\left.\frac{A}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{B}{r}\left(\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right)\right\}rdr$ (4.20)

$0={\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{-\frac{1}{r}\frac{d}{dr}\left[rS\left(\varphi +\frac{dw}{dr}\right)\right]-q\right\}rdr$ (4.21)

$0={\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{-\frac{1}{r}\frac{d}{dr}\left[rB\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rD\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]\right.$ $\left.+\frac{B}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{D}{r}\left(\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right)+S\left(\varphi +\frac{dw}{dr}\right)\right\}rdr$ (4.22)

The indicated integration in the above equations is carried out such that the resulting boundary terms constitute the secondary variable, dual to the primary variable of the equation considered.

$-N_1^{\left(I\right)}-N_2^{\left(I\right)}+{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{\frac{A}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{B}{r}\left(\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right)\right\}rdr=0$ (4.23)

where

$-V_1^{\left(I\right)}-V_2^{\left(I\right)}-{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}qrdr=0$ (4.24)

$-M_1^{\left(I\right)}-M_2^{\left(I\right)}+{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}\left\{\frac{B}{r}\left(\nu \frac{du}{dr}+\frac{u}{r}\right)+\frac{D}{r}\left(\nu \frac{d\varphi }{dr}+\frac{1}{r}\varphi \right)+S\left(\varphi +\frac{dw}{dr}\right)\right\}rdr=0$ (4.25)

$N_1^{\left(I\right)}=-{\left[rA\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rB\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]}_{r=r_a^{\left(I\right)}}$ (4.26a)

$N_2^{\left(I\right)}={\left[rA\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rB\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]}_{r=r_b^{\left(I\right)}}$ (4.26b)

$V_1^{\left(I\right)}=-{\left[rS\left(\varphi +\frac{dw}{dr}\right)\right]}_{r=r_a^{\left(I\right)}}$ (4.27a)

$V_2^{\left(I\right)}={\left[rS\left(\varphi +\frac{dw}{dr}\right)\right]}_{r=r_b^{\left(I\right)}}$ (4.27b)

$M_1^{\left(I\right)}=-{\left[rB\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rD\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]}_{r=r_a^{\left(I\right)}}$ (4.28a)

$M_2^{\left(I\right)}={\left[rB\left(\frac{du}{dr}+\nu \frac{u}{r}\right)+rD\left(\frac{d\varphi }{dr}+\frac{\nu }{r}\varphi \right)\right]}_{r=r_b^{\left(I\right)}}$ (4.28b)

Figure 6:
A typical finite domain for the displacement model of FST with secondary variables at the faces and the primary variables at the node associated with the finite domain.

Here ($N_1^{\left(I\right)}$, $V_1^{\left(I\right)}$, $M_1^{\left(I\right)}$) represent the axial force, shear force and moment respectively at the left face of the finite domain numbered $I$ (i.e., $r=r_a^{\left(I\right)}$), while ($N_2^{\left(I\right)}$, $V_2^{\left(I\right)}$, $M_2^{\left(I\right)}$) represent axial force, shear force and moment respectively at right face of the finite domain numbered $I$ (i.e., $r=r_b^{\left(I\right)}$) (see figure 6). Assuming linear Lagrange interpolations of the primary variables $\{u,w,\varphi \}$, Eqs. (4.23)-(4.25) on an interior finite domain can be evaluated (see Appendix A for the formulas used) to

$U_{I-1}\left[\frac{A_{I-1}{\nu }_{I-1}}{2}-\frac{A_{I-1}{r}_a^{\left(I\right)}}{h_{I-1}}+{\int }_{r_a^{\left(I\right)}}^{r_I}{A}_{I-1}\left({\nu }_{I-1}\frac{d{L}_1^{(I-1)}}{dr}+\frac{L_1^{(I-1)}}{r}\right)dr\right]$ $+{\varphi }_{I-1}\left[\frac{B_{I-1}{\nu }_{I-1}}{2}-\frac{B_{I-1}{r}_a^{\left(I\right)}}{h_{I-1}}+{\int }_{r_a^{\left(I\right)}}^{r_I}{B}_{I-1}\left({\nu }_{I-1}\frac{d{L}_1^{(I-1)}}{dr}+\frac{L_1^{(I-1)}}{r}\right)dr\right]$ $+U_I\left[\frac{r_a^{\left(I\right)}{A}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}{A}_I}{h_I}+\frac{A_{I-1}{\nu }_{I-1}-A_I{\nu }_I}{2}+{\int }_{r_a^{\left(I\right)}}^{r_I}{A}_{I-1}\left({\nu }_{I-1}\frac{d{L}_2^{(I-1)}}{dr}+\frac{L_2^{(I-1)}}{r}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}{A}_I\left({\nu }_I\frac{d{L}_1^{\left(I\right)}}{dr}+\frac{L_1^{\left(I\right)}}{r}\right)dr\right]+{\varphi }_I\left[\frac{r_a^{\left(I\right)}{B}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}{B}_I}{h_I}+\frac{B_{I-1}{\nu }_{I-1}-B_I{\nu }_I}{2}\right.$ $\left.+{\int }_{r_a^{\left(I\right)}}^{r_I}{B}_{I-1}\left({\nu }_{I-1}\frac{d{L}_2^{(I-1)}}{dr}+\frac{L_2^{(I-1)}}{r}\right)dr+{\int }_{r_I}^{r_b^{\left(I\right)}}{B}_I\left({\nu }_I\frac{d{L}_1^{\left(I\right)}}{dr}+\frac{L_1^{\left(I\right)}}{r}\right)dr\right]$ $+U_{I+1}\left[-\frac{A_I{\nu }_I}{2}-\frac{A_I{r}_b^{\left(I\right)}}{h_I}+{\int }_{r_I}^{r_b^{\left(I\right)}}{A}_I\left({\nu }_I\frac{d{L}_2^{\left(I\right)}}{dr}+\frac{L_2^{\left(I\right)}}{r}\right)dr\right]$ $+{\varphi }_{I+1}\left[-\frac{B_I{\nu }_I}{2}-\frac{B_I{r}_b^{\left(I\right)}}{h_I}+{\int }_{r_I}^{r_b^{\left(I\right)}}{B}_I\left({\nu }_I\frac{d{L}_2^{\left(I\right)}}{dr}+\frac{L_2^{\left(I\right)}}{r}\right)dr\right]=0$ (4.29)

$W_{I-1}\left[-\frac{r_a^{\left(I\right)}{S}_{I-1}}{h_{I-1}}\right]+{\varphi }_{I-1}\left[\frac{r_a^{\left(I\right)}{S}_{I-1}}{2}\right]+W_I\left[\frac{r_a^{\left(I\right)}{S}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}{S}_I}{h_I}\right]$ $+{\varphi }_I\left[\frac{r_a^{\left(I\right)}{S}_{I-1}-r_b^{\left(I\right)}{S}_I}{2}\right]+W_{I+1}\left[-\frac{r_b^{\left(I\right)}{S}_I}{h_I}\right]+{\varphi }_{I+1}\left[-\frac{r_b^{\left(I\right)}{S}_I}{2}\right]-{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}qrdr=0$ (4.30)

$U_{I-1}\left[\frac{B_{I-1}{\nu }_{I-1}}{2}-\frac{B_{I-1}{r}_a^{\left(I\right)}}{h_{I-1}}+{\int }_{r_a^{\left(I\right)}}^{r_I}{B}_{I-1}\left({\nu }_{I-1}\frac{d{L}_1^{(I-1)}}{dr}+\frac{L_1^{(I-1)}}{r}\right)dr\right]$ $+{\varphi }_{I-1}\left[\frac{D_{I-1}{\nu }_{I-1}}{2}-\frac{D_{I-1}{r}_a^{\left(I\right)}}{h_{I-1}}+{\int }_{r_a^{\left(I\right)}}^{r_I}{D}_{I-1}\left({\nu }_{I-1}\frac{d{L}_1^{(I-1)}}{dr}+\frac{L_1^{(I-1)}}{r}\right)dr\right]$ $+U_I\left[\frac{r_a^{\left(I\right)}{B}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}{B}_I}{h_I}+\frac{B_{I-1}{\nu }_{I-1}-B_I{\nu }_I}{2}+{\int }_{r_a^{\left(I\right)}}^{r_I}{B}_{I-1}\left({\nu }_{I-1}\frac{d{L}_2^{(I-1)}}{dr}+\frac{L_2^{(I-1)}}{r}\right)dr\right.$ $\left.+{\int }_{r_I}^{r_b^{\left(I\right)}}{B}_I\left({\nu }_I\frac{d{L}_1^{\left(I\right)}}{dr}+\frac{L_1^{\left(I\right)}}{r}\right)dr\right]+{\varphi }_I\left[\frac{r_a^{\left(I\right)}{D}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}{D}_I}{h_I}+\frac{D_{I-1}{\nu }_{I-1}-D_I{\nu }_I}{2}\right.$ $\left.+{\int }_{r_a^{\left(I\right)}}^{r_I}{D}_{I-1}\left({\nu }_{I-1}\frac{d{L}_2^{(I-1)}}{dr}+\frac{L_2^{(I-1)}}{r}\right)dr+{\int }_{r_I}^{r_b^{\left(I\right)}}{D}_I\left({\nu }_I\frac{d{L}_1^{\left(I\right)}}{dr}+\frac{L_1^{\left(I\right)}}{r}\right)dr\right]$ $+U_{I+1}\left[-\frac{B_I{\nu }_I}{2}-\frac{B_I{r}_b^{\left(I\right)}}{h_I}+{\int }_{r_I}^{r_b^{\left(I\right)}}{B}_I\left({\nu }_I\frac{d{L}_2^{\left(I\right)}}{dr}+\frac{L_2^{\left(I\right)}}{r}\right)dr\right]$ $+{\varphi }_{I+1}\left[-\frac{D_I{\nu }_I}{2}-\frac{D_I{r}_b^{\left(I\right)}}{h_I}+{\int }_{r_I}^{r_b^{\left(I\right)}}{D}_I\left({\nu }_I\frac{d{L}_2^{\left(I\right)}}{dr}+\frac{L_2^{\left(I\right)}}{r}\right)dr\right]$ $+\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\frac{S_{I-1}r}{2}dr\right]{\varphi }_{I-1}+\left[{\int }_{r_a^{\left(I\right)}}^{r_I}\frac{S_{I-1}r}{2}dr\right]{\varphi }_I+\left[{\int }_{r_I}^{r_b^{\left(I\right)}}\frac{S_{I}r}{2}dr\right]{\varphi }_I$ $+\left[{\int }_{r_I}^{r_b^{\left(I\right)}}\frac{S_{I}r}{2}dr\right]{\varphi }_{I+1}+\left[{\int }_{r_a^{\left(I\right)}}^{r_I}{S}_{I-1}r\frac{d{L}_1^{(I-1)}}{dr}dr\right]W_{I-1}+\left[{\int }_{r_a^{\left(I\right)}}^{r_I}{S}_{I-1}r\frac{d{L}_2^{(I-1)}}{dr}dr\right]W_I$ $+\left[{\int }_{r_I}^{r_b^{\left(I\right)}}{S}_{I}r\frac{d{L}_1^{\left(I\right)}}{dr}dr\right]W_I+\left[{\int }_{r_I}^{r_b^{\left(I\right)}}{S}_{I}r\frac{d{L}_2^{\left(I\right)}}{dr}dr\right]W_{I+1}=0$ (4.31)

In evaluating Eq. (4.31), the coefficients of $\varphi$ corresponding to the integral

$\begin{array}{c}{\int }_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}rS\left(\varphi +\frac{dw}{dr}\right)dr\end{array}$

are evaluated by considering $\varphi$ to be constant within each element. Thus in element ${\Omega }^{(I-1)}$, $\varphi =\frac{{\varphi }_I+{\varphi }_{I-1}}{2}$ while in the element ${\Omega }^{\left(I\right)}$, $\varphi =\frac{{\varphi }_{I+1}+{\varphi }_I}{2}$. This is done to remedy the shear locking (see [⁸]).

Similarly, the discretized equations on 1^st finite domain would be

$-N_1^{\left(1\right)}+U_1\left[\frac{A_1(1-{\nu }_1)}{2}+{\int }_0^{0.5h_1}{A}_1\left({\nu }_1\frac{d{L}_1^{\left(1\right)}}{dr}+\frac{L_1^{\left(1\right)}}{r}\right)dr\right]$ $+{\varphi }_1\left[\frac{B_1(1-{\nu }_1)}{2}+{\int }_0^{0.5h_1}{B}_1\left({\nu }_1\frac{d{L}_1^{\left(1\right)}}{dr}+\frac{L_1^{\left(1\right)}}{r}\right)dr\right]$ $+U_2\left[-\frac{A_1(1+{\nu }_1)}{2}+{\int }_0^{0.5h_1}{A}_1\left({\nu }_1\frac{d{L}_2^{\left(1\right)}}{dr}+\frac{L_2^{\left(1\right)}}{r}\right)dr\right]$ $+{\varphi }_2\left[-\frac{B_1(1+{\nu }_1)}{2}+{\int }_0^{0.5h_1}{B}_1\left({\nu }_1\frac{d{L}_2^{\left(1\right)}}{dr}+\frac{L_2^{\left(1\right)}}{r}\right)dr\right]=0$ (4.32)

$-V_1^{\left(1\right)}+W_1\left[\frac{S_1}{2}\right]+{\varphi }_1\left[-\frac{h_1{S}_1}{4}\right]+W_2\left[\frac{-S_1}{2}\right]+{\varphi }_2\left[-\frac{h_1{S}_1}{4}\right]-{\int }_0^{0.5h_1}qrdr=0$ (4.33)

$-M_1^{\left(1\right)}+U_1\left[\frac{B_1(1-{\nu }_1)}{2}+{\int }_0^{0.5h_1}{B}_1\left({\nu }_1\frac{d{L}_1^{\left(1\right)}}{dr}+\frac{L_1^{\left(1\right)}}{r}\right)dr\right]+W_1\left[{\int }_0^{0.5h_1}{S}_{1}r\frac{d{L}_1^{\left(1\right)}}{dr}dr\right]$ $+{\varphi }_1\left[\frac{D_1(1-{\nu }_1)}{2}+{\int }_0^{0.5h_1}{D}_1\left({\nu }_1\frac{d{L}_1^{\left(1\right)}}{dr}+\frac{L_1^{\left(1\right)}}{r}\right)dr+{\int }_0^{0.5h_1}\frac{S_{1}r}{2}dr\right]$ $+U_2\left[-\frac{B_1(1+{\nu }_1)}{2}+{\int }_0^{0.5h_1}{B}_1\left({\nu }_1\frac{d{L}_2^{\left(1\right)}}{dr}+\frac{L_2^{\left(1\right)}}{r}\right)dr\right]+W_2\left[{\int }_0^{0.5h_1}{S}_{1}r\frac{d{L}_2^{\left(1\right)}}{dr}dr\right]$ $+{\varphi }_2\left[-\frac{D_1(1+{\nu }_1)}{2}+{\int }_0^{0.5h_1}{D}_1\left({\nu }_1\frac{d{L}_2^{\left(1\right)}}{dr}+\frac{L_2^{\left(1\right)}}{r}\right)dr+{\int }_0^{0.5h_1}\frac{S_{1}r}{2}dr\right]=0$ (4.34)

while the discretized equations on ($N+1$)^th finite domain would be

${-N}_2^{(N+1)}+U_N\left[\frac{A_N{\nu }_N}{2}-\frac{A_N(R-0.5h_N)}{h_N}+{\int }_{R-0.5h_N}^R{A}_N\left({\nu }_N\frac{d{L}_1^{\left(N\right)}}{dr}+\frac{L_1^{\left(N\right)}}{r}\right)dr\right]$ $+{\varphi }_N\left[\frac{B_N{\nu }_N}{2}-\frac{B_N(R-0.5h_N)}{h_N}+{\int }_{R-0.5h_N}^R{B}_N\left({\nu }_N\frac{d{L}_1^{\left(N\right)}}{dr}+\frac{L_1^{\left(N\right)}}{r}\right)dr\right]$ $+U_{N+1}\left[\frac{A_N(R-0.5h_N)}{h_N}+\frac{A_N{\nu }_N}{2}+{\int }_{R-0.5h_N}^R{A}_N\left({\nu }_N\frac{d{L}_2^{\left(N\right)}}{dr}+\frac{L_2^{\left(N\right)}}{r}\right)dr\right]$ $+{\varphi }_{N+1}\left[\frac{B_N(R-0.5h_N)}{h_N}+\frac{B_N{\nu }_N}{2}+{\int }_{R-0.5h_N}^R{B}_N\left({\nu }_N\frac{d{L}_2^{\left(N\right)}}{dr}+\frac{L_2^{\left(N\right)}}{r}\right)dr\right]=0$ (4.35)

$-V_2^{(N+1)}+W_N\left[-\frac{(R-0.5h_N)S_N}{h_N}\right]+{\varphi }_N\left[\frac{(R-0.5h_N)S_N}{2}\right]$ $+W_{N+1}\left[\frac{(R-0.5h_N)S_N}{h_N}\right]+{\varphi }_{N+1}\left[\frac{(R-0.5h_N)S_N}{2}\right]-{\int }_{R-0.5h_N}^{R}qrdr=0$ (4.36)

$-M_2^{(N+1)}+U_N\left[\frac{B_N{\nu }_N}{2}-\frac{B_N(R-0.5h_N)}{h_N}+{\int }_{R-0.5h_N}^R{B}_N\left({\nu }_N\frac{d{L}_1^{\left(N\right)}}{dr}+\frac{L_1^{\left(N\right)}}{r}\right)dr\right]$ $+{\varphi }_N\left[\frac{D_N{\nu }_N}{2}-\frac{D_N(R-0.5h_N)}{h_N}+{\int }_{R-0.5h_N}^R{D}_N\left({\nu }_N\frac{d{L}_1^{\left(N\right)}}{dr}+\frac{L_1^{\left(N\right)}}{r}\right)dr\right]$ $+U_{N+1}\left[\frac{B_N(R-0.5h_N)}{h_N}+\frac{B_N{\nu }_N}{2}+{\int }_{R-0.5h_N}^R{B}_N\left({\nu }_N\frac{d{L}_2^{\left(N\right)}}{dr}+\frac{L_2^{\left(N\right)}}{r}\right)dr\right]$ $+{\varphi }_{N+1}\left[\frac{D_N(R-0.5h_N)}{h_N}+\frac{D_N{\nu }_N}{2}+{\int }_{R-0.5h_N}^R{D}_N\left({\nu }_N\frac{d{L}_2^{\left(N\right)}}{dr}+\frac{L_2^{\left(N\right)}}{r}\right)dr\right]$ $+W_N\left[{\int }_{R-0.5h_N}^R{S}_{N}r\frac{d{L}_1^{\left(N\right)}}{dr}dr\right]+{\varphi }_N\left[{\int }_{R-0.5h_N}^R\frac{S_{N}r}{2}dr\right]$ $+W_{N+1}\left[{\int }_{R-0.5h_N}^R{S}_{N}r\frac{d{L}_2^{\left(N\right)}}{dr}dr\right]+{\varphi }_{N+1}\left[{\int }_{R-0.5h_N}^R\frac{S_{N}r}{2}dr\right]=0$ (4.37)

5. Numerical Results

To illustrate the workings of the dual mesh finite domain models presented in the previous sections, we consider two cases: (a) hinged (or simply supported) and (b) clamped functionally graded circular plates, subjected to uniformly distributed load of intensity $q$. A comparison between the numerical results obtained with FEM (see Appendix B for various models) and DMFDM models for each case considered is presented. There are three models of FEM and two models of DMFDM, as summarized below:

• FE-CP(D) - Displacement finite element model for CPT

• FE-CP(M) - Mixed finite element model for CPT

• FE-FP(D) - Displacement finite element model for FST

• DM-CP(M) - Mixed dual mesh finite domain model for CPT

• DM-FP(D) - Displacement dual mesh finite domain model for FST

We investigate the effect of mesh and power-law index $n$ on the deflections and stresses. We consider a functionally graded circular plate of radius $R=10$ in, thickness $H=0.1$ in and subjected to axisymmetric uniformly distributed load intensity $q=0.5$ lb/in². The functionally graded material properties are taken to be $E_t=3\times {10}^7$ psi, $E_b=3\times {10}^6$ psi, and $\nu =0.3$ (units are stated here only to give a physical meaning; one may replace them with appropriate metric units).

The stresses and stress resultants can be post-computed in the numerical models described. For the displacement models, once the generalized displacements are obtained, their derivatives can be easily determined. Then the stresses can be directly computed from the constitutive relations of Eq. (2.3), while the stress resultants and moments can be computed using Eqs. (2.4a)-(2.4) or Eqs. (2.13a)-(2.13d) depending on whether FST is used or CPT is used. For the mixed models, the second derivative of the transverse deflection needs to be calculated from the moment $M_{rr}$, which is determined as a nodal variable in mixed models. The second derivative of transverse deflection $w$ can be obtained from the following equation:

$\frac{d^{2}w}{d{r}^2}=-\frac{\nu }{r}\frac{dw}{dr}+\overline{B}\left(\frac{du}{dr}+\nu \frac{u}{r}\right)-\frac{1}{D}{M}_{rr}$ (5.1)

Then the stresses can be computed using Eqs. (2.3) while the stress resultants can be computed using Eqs. (3.1) and (3.4).

5.1 Hinged plates

The boundary conditions for the primary variables of hinged axisymmetric circular plate in various models are as follows:

$\mathrm{D}\mathrm{i}\mathrm{s}\mathrm{p}\mathrm{l}\mathrm{a}\mathrm{c}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{n}\mathrm{t}\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}\mathrm{s}:u\left(0\right)=0,\frac{dw}{dr}\left(0\right)\mathrm{o}\mathrm{r}\varphi \left(0\right)=0,u\left(R\right)=w\left(R\right)=0$ (5.2)

$\mathrm{M}\mathrm{i}\mathrm{x}\mathrm{e}\mathrm{d}\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}\mathrm{s}:u\left(0\right)=0,u\left(R\right)=w\left(R\right)=M\left(R\right)=0$ (5.3)

A comparison between the transverse deflections at the center of a hinged homogeneous circular plate obtained from various models is presented in Table 1 When 32 elements are used, all models give the same results as the exact solution (see [¹³-¹⁵]) up to the fourth decimal point. The mixed finite element model has a slow mesh convergence rate at the origin for transverse deflection, $w$, and hence a non-uniform finite element mesh with a finer mesh at the origin is used in reporting the values of $w\left(0\right)$. The non-uniform mesh is chosen such that for two element mesh, the lengths of the elements are in ratio 1:9, with the shortest element being at origin. All the subsequent refinements are made by breaking the finite elements into half. Thus for four element mesh the element lengths would be $\{0.5R,0.5R,4.5R,4.5R\}$ and further refinement is made by breaking these elements into their halves and so on. Since the radius to thickness ratio $R/H=100$, both classical plate theory and first order shear deformation theory predict the same results.

Figure 7(a) contains a comparison of the moment, $M_{rr}$, obtained from displacement and mixed finite element models and mixed dual mesh finite domain model of the CPT when a uniform mesh of 32 elements is used, while Fig. 7(b) includes a comparison of the moment, $M_{rr}$, obtained from displacement finite element model and displacement dual mesh finite domain model of the FST when a uniform mesh of 32 elements is used. While the moments obtained from all the models are very close, the mixed finite element model and mixed dual mesh finite domain model deviate from the displacement finite element model at the origin ($r=0$). Both mixed models (i.e., mixed FEM and mixed DMFDM models) are prone to give erroneous moments at the origin. It should be noted that this anomalous behavior at the origin is not a characteristic of the dual mesh finite domain method; both the mixed finite element model and mixed dual mesh finite domain model exhibit such behavior. However, the displacements models do not exhibit such a behavior.

Table 1:
Center transverse deflection, $w\left(0\right)$, of hinged homogeneous axisymmetric circular plate as predicted by various models. The results are given in inches.

Mesh	FE-CP(D)	FE-CP(M)	FE-FP(D)	DM-CP(M)	DM-FP(D)
2	0.1159	0.1120	0.1076	0.1149	0.1052
4	0.1159	0.1136	0.1143	0.1156	0.1134
8	0.1159	0.1155	0.1156	0.1158	0.1153
16	0.1159	0.1159	0.1159	0.1159	0.1158
32	0.1159	0.1159	0.1159	0.1159	0.1159
64	0.1159	0.1159	0.1159	0.1159	0.1159
Exact	0.1159	0.1159	0.1159	0.1159	0.1159

Figure 7:
Moment of hinged homogeneous circular plate as predicted by various models. (a) Numerical models based on the CPT. (b) Numerical models based on the FST.

5.2 Clamped plates

The boundary conditions on the primary variables of clamped axisymmetric circular plate in various models are as follows:

$\mathrm{D}\mathrm{i}\mathrm{s}\mathrm{p}\mathrm{l}\mathrm{a}\mathrm{c}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{n}\mathrm{t}\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}\mathrm{s}:u\left(0\right)=0,\frac{dw}{dr}\left(0\right)\mathrm{o}\mathrm{r}\varphi \left(0\right)=0,$ $u\left(R\right)=w\left(R\right)=0,\frac{dw}{dr}\left(R\right)\mathrm{o}\mathrm{r}\varphi \left(R\right)=0$ (5.4)

$\mathrm{M}\mathrm{i}\mathrm{x}\mathrm{e}\mathrm{d}\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}\mathrm{s}:u\left(0\right)=0,u\left(R\right)=w\left(R\right)=0$ (5.5)

A comparison of the stress values, ${\sigma }_{rr}$, on the top face of the clamped homogeneous circular plate as predicted by various numerical models based on CPT with the analytical solution [¹⁴] are presented in Fig. 8. Although all the models are in good agreement with the analytical solution for most part of the computational domain, the mixed models (i.e., mixed FEM model and mixed DMFDM model) deviate from the analytical solution at and very close to the origin. This is due to the erroneous prediction of moment $M_{rr}$ near the origin by the mixed models. Further, Fig. 9 contains a comparison of stress ${\sigma }_{rr}$ on the top face of clamped homogeneous circular plate as predicted by various numerical models based on the FST with that of the analytical solution [¹³, ¹⁴]. Since both the finite element model and dual mesh finite domain model of the FST are displacement based, the stress values are in good agreement with the analytical solution through out the computational domain.

Finally, to illustrate the working of the dual mesh finite domain method in analyzing the functionally graded axisymmetric circular plates, we compare the transverse deflection at the center, $w\left(0\right)$, from the analytical solution of the CPT [¹³-¹⁵] with the transverse deflections obtained from various models described earlier. Various values of the power law index, $n$, are considered to bring out the effect of $n$ on the bending deflections of the plate. Table 2 contains a comparison between the various models with that of the analytical solution (of the FST). For this thin plate ($R/H=100$), both the CPT and the FST predict solutions close to each other [¹³-¹⁵]; 32 elements are used for the results reported in Table 2. Uniform mesh is used in all models except for the mixed FEM model of the CPT. A non-uniform mesh described earlier is used for FE-CP(M) to achieve faster mesh convergence. It can be seen that as the power-law index, $n$, increases the plate stiffness decreases (since plate material tends toward $E_b<E_t$).

Figure 8:
Values of ${\sigma }_{rr}(r,H/2)$ as predicted by various models based on the CPT.

Figure 9:
Values of ${\sigma }_{rr}(r,H/2)$ as predicted by various models based on the FST.

Table 2:
Center transverse deflection, w(0), of a clamped FGM axisymmetric circular plate for various values of n, as predicted by various models; 32 elements are used in all the models. The results are given in inches.

n	FE-CP(D)	FE-CP(M)	DM-CP(M)	FE-FP(D)	DM-FP(D)	FST Exact
0.0	0.0284	0.0285	0.0284	0.0285	0.0284	0.0284
1.0	0.0665	0.0666	0.0665	0.0666	0.0665	0.0666
2.0	0.0976	0.0976	0.0975	0.0977	0.0976	0.0976
3.0	0.1152	0.1152	0.1151	0.1153	0.1152	0.1152
4.0	0.1244	0.1244	0.1242	0.1245	0.1244	0.1244
5.0	0.1296	0.1296	0.1294	0.1297	0.1296	0.1296
7.5	0.1370	0.1370	0.1369	0.1371	0.1370	0.1370
10.0	0.1425	0.1425	0.1424	0.1426	0.1425	0.1425
12.0	0.1467	0.1467	0.1466	0.1468	0.1469	0.1467
15.0	0.1527	0.1528	0.1527	0.1529	0.1527	0.1523
20.0	0.1622	0.1622	0.1621	0.1623	0.1622	0.1622

6. Concluding remarks

The dual mesh finite domain method (DMFDM) is presented for the linear, axisymmetric bending analysis of two-constituent through thickness functionally graded circular plates. Both classical (CPT) and shear deformation (FST) theories are used. The mixed model of the CPT and displacement model of the FST are developed using the DMFDM. The mixed model includes the bending moment as a dependent variable in addition to the radial and transverse displacement. Numerical results indicate that the mixed DMFDM model has a better mesh convergence characteristic than the mixed FEM model. In general, the DMFDM models give as accurate results as the FEM, with the DMFDM giving the bending moments directly at the nodes. While both FEM and DMFDM have comparable accuracy, the DMFDM has less formulative steps and computational expense due to the fact that there is no element concept and hence no assembly of element equations; these features will have significant savings in multidimensional problems. Extensions of the DMFDM to higher-order and non-classical theories of beams, plates, and shells, accounting for geometric and material nonlinearities, and the solution of the Navier-Stokes equations of fluid dynamics [¹⁶] for non-Newtonian fluids are awaiting attention.

Acknowledgments

The research reported herein was supported by the O’Donnell Foundation Chair IV held by the second author.

References

• ^[1]
  D. Hasselman, G. Youngblood, Enhanced thermal stress resistance of structural ceramics with thermal conductivity gradient, Journal of the American Ceramic Society 61 (1,2) (1978) 49-53.
• ^[2]
  M. Koizumi, The concept of FGM, Ceramic Transactions, Functionally Gradient Materials 34 (1993) 3-10.
• ^[3]
  A. Kawasaki, R. Watanabe, Concept and P/M fabrication of functionally graded materials, Ceram. Int. 23 (1) (1997) 73-83.
• ^[4]
  M. Naebe, K. Shirvanimoghaddam, Functionally graded materials: A review of fabrication and properties, Appl. Mater. Today 5 (2016) 223-245.
• ^[5]
  S. Nikbakht, S. Kamarian, M. Shakeri, A review on optimization of composite structures Part II: Functionally graded materials, Compos. Struct. 214 (2019) 83-102.
• ^[6]
  J. N. Reddy, Mechanics of Laminated Composite Plates and Shells: Theory and Analysis, 2nd Edition, CRC Press, 2004.
• ^[7]
  J. N. Reddy, A dual mesh finite domain method for the numerical solution of differential equations, Int. J. Comput. Methods Eng. Sci. Mech. 20 (3) (2019) 212-228.
• ^[8]
  J. N. Reddy, An Introduction to the Finite Element Method, 4th Ed, McGraw-Hill, New York, NY, 2019.
• ^[9]
  S. Mazumder, Numerical Methods for Partial Differential Equations. Finite Difference and Finite Volume Methods, Elsevier, New York, NY, 2016.
• ^[10]
  J. N. Reddy, M. Martinez, A dual mesh finite domain method for steady-state convection-diffusion problems, Comput. Fluids, accepted for publication.
• ^[11]
  J. N. Reddy, P. Nampally, A dual mesh finite domain method for the analysis of functionally graded beams, Composite Structures 251 (2020) 112648.
• ^[12]
  J. N. Reddy, P. Nampally, A. R. Srinivasa, Nonlinear analysis of functionally graded beams using the dual mesh finite domain method and the finite element method, International Journal of Non-Linear Mechanics, 127 (2020) 103575.
• ^[13]
  J. N. Reddy, Energy Principles and Variational Methods in Applied Mechanics, 3rd Edition, John Wiley & Sons, New York, NY, 2017.
• ^[14]
  J. N. Reddy, Theory and Analysis of Elastic Plates and Shells, 2nd Edition, CRC Press, Boca Raton, FL, 2007.
• ^[15]
  C.M. Wang, J. N. Reddy, and K.H. Lee, Shear Deformation Theories of Beams and Plates. Relationships with Classical Solution, Elsevier, U.K., 2000.
• ^[16]
  J. N. Reddy, N. Kim, and M. Martinez, A dual mesh control domain method for the solution of nonlinear Poisson's equation and the Navier-Stokes equations for incompressible fluids, Physics of Fluids, 32 (2020) https://doi.org/10.1063/5.0026274
  » https://doi.org/10.1063/5.0026274

Appendix A. Discrete forms of various functions

For a function $a\left(r\right)$, which is assumed to be element wise constant, and an independent variable $u\left(r\right)$, the following relations are useful in deriving the discretized equations of the DMFDM:

${\left[ar\frac{du}{dr}\right]}_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}=\left(\frac{r_a^{\left(I\right)}{a}_{I-1}}{h_{I-1}}\right)U_{I-1}-\left(\frac{r_a^{\left(I\right)}{a}_{I-1}}{h_{I-1}}+\frac{r_b^{\left(I\right)}{a}_I}{h_I}\right)U_I+\left(\frac{r_b^{\left(I\right)}{a}_I}{h_I}\right)U_{I+1}$ (A.1)

${\left[ar\frac{du}{dr}\right]}_{r=0.5h_1}=\frac{a_1}{2}\left(U_2-U_1\right),{\left[ar\frac{du}{dr}\right]}_{r=R-0.5h_N}=\frac{a_N\left(R-0.5h_N\right)}{h_N}\left(U_{N+1}-U_N\right)$ (A.2)

$[rau{]}_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}=-\left(\frac{r_a^{\left(I\right)}{a}_{I-1}}{2}\right)U_{I-1}+\left(\frac{r_b^{\left(I\right)}{a}_I-r_a^{\left(I\right)}{a}_{I-1}}{2}\right)U_I+\left(\frac{r_b^{\left(I\right)}{a}_I}{2}\right)U_{I+1}$ (A.3)

$[rau{]}_{r=0.5h_1}=\frac{a_1{h}_1}{4}(U_1+U_2),[rau{]}_{r=R-0.5h_N}=\frac{a_N\left(R-0.5h_N\right)}{2}\left(U_{N+1}+U_N\right)$ (A.4)

$[au{]}_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}=-\left(\frac{a_{I-1}}{2}\right)U_{I-1}+\left(\frac{a_I-a_{I-1}}{2}\right)U_I+\left(\frac{a_I}{2}\right)U_{I+1}$ (A.5)

$[au{]}_{r=0.5h_1}=\frac{a_1}{2}(U_1+U_2),[au{]}_{r=R-0.5h_N}=\frac{a_N}{2}\left(U_N+U_{N+1}\right)$ (A.6)

$[r^{2}au{]}_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}=-\frac{r_a^{(I{)}^2}{a}_{I-1}}{2}{U}_{I-1}+\left(\frac{r_b^{(I{)}^2}{a}_I-r_a^{(I{)}^2}{a}_{I-1}}{2}\right)U_I+\frac{r_b^{(I{)}^2}{a}_I}{2}{U}_{I+1}$ (A.7)

$[r^{2}au{]}_{r=0.5h_1}=\frac{a_1{h}_1^2}{8}(U_1+U_2),[r^{2}au{]}_{r=R-0.5h_N}=\frac{a_N(R-0.5h_N{)}^2}{2}\left(U_N+U_{N+1}\right)$ (A.8)

${\left[r^{2}a\frac{du}{dr}\right]}_{r_a^{\left(I\right)}}^{r_b^{\left(I\right)}}=\frac{r_a^{(I{)}^2}{a}_{I-1}}{h_{I-1}}{U}_{I-1}-\left(\frac{r_b^{(I{)}^2}{a}_I}{h_I}+\frac{r_a^{(I{)}^2}{a}_{I-1}}{h_{I-1}}\right)U_I+\frac{r_b^{(I{)}^2}{a}_I}{h_I}{U}_{I+1}$ (A.9)

${\left[r^{2}a\frac{du}{dr}\right]}_{r=0.5h_1}=\frac{a_1{h}_1}{4}\left(U_2-U_1\right),{\left[r^{2}a\frac{du}{dr}\right]}_{r=R-0.5h_N}=\frac{a_N(R-0.5h_N{)}^2}{h_N}\left(U_{N+1}-U_N\right)$ (A.10)

Appendix B. Finite element models of FGM axisymmetric circular plates

B.1 Mixed finite element model of the CPT

Equations (3.6)-(3.8) are used in developing the mixed finite element model of an axisymmetric circular plate based on the CPT. The mixed FE model of a typical element can be written as

${\left[\begin{array}{ccc}{\mathbf{K}}^{11}& {\mathbf{K}}^{12}& {\mathbf{K}}^{13}\\ {\mathbf{K}}^{21}& {\mathbf{K}}^{22}& {\mathbf{K}}^{23}\\ {\mathbf{K}}^{31}& {\mathbf{K}}^{32}& {\mathbf{K}}^{33}\end{array}\right]}^{\left(e\right)}{\left\{\begin{array}{c}\mathbf{U}\\ \mathbf{W}\\ {\mathbf{M}}_{\mathbf{r}\mathbf{r}}\end{array}\right\}}^{\left(e\right)}={\left\{\begin{array}{c}{\mathbf{F}}^1\\ {\mathbf{F}}^2\\ {\mathbf{F}}^3\end{array}\right\}}^{\left(e\right)}$ (B.1)

where $\mathbf{U}$ denotes the vector of nodal displacements $u$, $\mathbf{W}$ denotes the vector of nodal displacements $w$, and ${\mathbf{M}}_{\mathbf{r}\mathbf{r}}$ denotes the vector of moments $M_{rr}$ at the nodes of the element. The non-zero coefficients of the above equations are given by

$K_{ij}^{11}={\int }_{r_a}^{r_b}\left\{\overline{A}r\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+\overline{A}\nu \left(\frac{d{\psi }_i}{dr}{\psi }_j+{\psi }_i\frac{d{\psi }_j}{dr}\right)+\frac{\overline{A}}{r}{\psi }_i{\psi }_j+\frac{B^2(1-{\nu }^2)}{Dr}{\psi }_i{\psi }_j\right\}dr$

$K_{ij}^{12}={\int }_{r_a}^{r_b}\frac{B({\nu }^2-1)}{r}{\psi }_i\frac{d{\psi }_j}{dr}dr,K_{ij}^{13}={\int }_{r_a}^{r_b}\left\{\overline{B}\left(r\frac{d{\psi }_i}{dr}{\psi }_j+\nu {\psi }_i{\psi }_j\right)\right\}dr$

$K_{ij}^{21}={\int }_{r_a}^{r_b}\frac{B({\nu }^2-1)}{r}\frac{d{\psi }_i}{dr}{\psi }_{j}dr,K_{ij}^{22}={\int }_{r_a}^{r_b}\frac{D(1-{\nu }^2)}{r}\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}dr$

$K_{ij}^{23}={\int }_{r_a}^{r_b}\left\{r\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+(1-\nu )\frac{d{\psi }_i}{dr}{\psi }_j\right\}dr$

$K_{ij}^{31}={\int }_{r_a}^{r_b}\left\{\overline{B}\left(r{\psi }_i\frac{d{\psi }_j}{dr}+\nu {\psi }_i{\psi }_j\right)\right\}dr$

$K_{ij}^{32}={\int }_{r_a}^{r_b}\left\{r\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+(1-\nu ){\psi }_i\frac{d{\psi }_j}{dr}\right\}dr,K_{ij}^{33}=-{\int }_{r_a}^{r_b}\frac{r}{D}{\psi }_i{\psi }_{j}dr$

$F_i^1={\psi }_i\left(r_a\right)Q_1+{\psi }_i\left(r_b\right)Q_4,F_i^2={\int }_{r_a}^{r_b}qrdr+{\psi }_i\left(r_a\right)Q_2+{\psi }_i\left(r_b\right)Q_5$

$F_i^3={\psi }_i\left(r_a\right)Q_3+{\psi }_i\left(r_b\right)Q_6$

Here $i,j=\{\mathrm{1,2},\dots ,NPE\}$, with $NPE$ being number of nodes in an element. Further, ${\psi }_i$ are the Lagrange interpolation functions used for $(u,w,M_{rr})$ and

$\begin{array}{rr}& Q_1=-[r{N}_{rr}{]}_{r=r_a},Q_4=[r{N}_{rr}{]}_{r=r_b}\\ & Q_2=-{\left[\frac{d}{dr}\left(r{M}_{rr}\right)-M_{\theta \theta }\right]}_{r=r_a},Q_5={\left[\frac{d}{dr}\left(r{M}_{rr}\right)-M_{\theta \theta }\right]}_{r=r_b}\\ & Q_3={\left[r\frac{dw}{dr}\right]}_{r=r_a},Q_6={\left[-r\frac{dw}{dr}\right]}_{r=r_b}\end{array}$

B.2 Displacement finite element model of the CPT

Equations (2.18)-(2.19) are used to develop the displacement finite element model based on the CPT. The displacement finite element model of a typical element is of the form

${\left[\begin{array}{cc}{\mathbf{K}}^{11}& {\mathbf{K}}^{12}\\ {\mathbf{K}}^{21}& {\mathbf{K}}^{22}\end{array}\right]}^{\left(e\right)}{\left\{\begin{array}{c}\mathbf{U}\\ \mathbf{\Delta }\end{array}\right\}}^{\left(e\right)}={\left\{\begin{array}{c}{\mathbf{F}}^1\\ {\mathbf{F}}^2\end{array}\right\}}^{\left(e\right)}$ (B.2)

where $\mathbf{U}$ denotes the vector of nodal displacements $u$, $\mathbf{\Delta }$ denotes the vector of nodal transverse deflections $w$, and rotations $-\frac{dw}{dr}$ within the element. The non-zero coefficients of the above equations are

$K_{ij}^{11}={\int }_{r_a}^{r_b}\left\{Ar\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+A\nu \left(\frac{d{\psi }_i}{dr}{\psi }_j+{\psi }_i\frac{d{\psi }_j}{dr}\right)+\frac{A}{r}{\psi }_i{\psi }_j\right\}dr$

$K_{iJ}^{12}={\int }_{r_a}^{r_b}\left\{-Br\frac{d{\psi }_i}{dr}\frac{d^2{\phi }_J}{d{r}^2}-B\nu \frac{d{\psi }_i}{dr}\frac{d{\phi }_J}{dr}-B\nu {\psi }_i\frac{d^2{\phi }_J}{d{r}^2}-\frac{B}{r}{\psi }_i\frac{d{\phi }_J}{dr}\right\}dr$

$K_{Ij}^{21}={\int }_{r_a}^{r_b}\left\{-Br\frac{d^2{\phi }_I}{d{r}^2}\frac{d{\psi }_j}{dr}-B\nu \frac{d{\phi }_I}{dr}\frac{d{\psi }_j}{dr}-B\nu \frac{d^2{\phi }_I}{d{r}^2}{\psi }_j-\frac{B}{r}\frac{d{\phi }_I}{dr}{\psi }_j\right\}dr$

$K_{IJ}^{22}={\int }_{r_a}^{r_b}\left\{Dr\frac{d^2{\phi }_I}{d{r}^2}\frac{d^2{\phi }_J}{d{r}^2}+D\nu \left(\frac{d^2{\phi }_I}{d{r}^2}\frac{d{\phi }_J}{dr}+\frac{d{\phi }_I}{dr}\frac{d^2{\phi }_J}{d{r}^2}\right)+\frac{D}{r}\frac{d{\phi }_I}{dr}\frac{d{\phi }_j}{dr}\right\}dr$

$F_i^1={\psi }_i\left(r_a\right)Q_1+{\psi }_i\left(r_b\right)Q_4$

$F_I^2={\int }_{r_a}^{r_b}qrdr+{\phi }_I\left(r_a\right)Q_2+{\left.\frac{d{\phi }_I}{dr}\right|}_{r_a}{Q}_3+{\phi }_I\left(r_b\right)Q_5+{\left.\frac{d{\phi }_I}{dr}\right|}_{r_b}{Q}_6$

Here $i,j=\{\mathrm{1,2},\dots ,NPE\}$ and $I,J=\left\{\mathrm{1,2},\dots ,2NPE\right\}$ with NPE being the number of nodes in an element; ${\psi }_i$ are Lagrange interpolation functions used for $u$ and ${\phi }_I$ are the Hermite interpolation functions used for $w$, and

$\begin{array}{rr}& Q_1=-[r{N}_{rr}{]}_{r=r_a},Q_4=[r{N}_{rr}{]}_{r=r_b}\\ & Q_2=-{\left[\frac{d}{dr}\left(r{M}_{rr}\right)-M_{\theta \theta }\right]}_{r=r_a},Q_5={\left[\frac{d}{dr}\left(r{M}_{rr}\right)-M_{\theta \theta }\right]}_{r=r_b}\\ & Q_3=-{\left[r{M}_{rr}\right]}_{r=r_a},Q_6={\left[r{M}_{rr}\right]}_{r=r_b}\end{array}$

B.3 Displacement finite element model of the FST

Equations (2.14)-(2.16) are used to develop the displacement finite element model of a FG axisymmetric circular plate based on the FST. The displacement finite element equations on a typical element are of the form

${\left[\begin{array}{ccc}{\mathbf{K}}^{11}& {\mathbf{K}}^{12}& {\mathbf{K}}^{13}\\ {\mathbf{K}}^{21}& {\mathbf{K}}^{22}& {\mathbf{K}}^{23}\\ {\mathbf{K}}^{31}& {\mathbf{K}}^{32}& {\mathbf{K}}^{33}\end{array}\right]}^{\left(e\right)}{\left\{\begin{array}{c}\mathbf{U}\\ \mathbf{W}\\ \mathbf{\Phi }\end{array}\right\}}^{\left(e\right)}={\left\{\begin{array}{c}{\mathbf{F}}^1\\ {\mathbf{F}}^2\\ {\mathbf{F}}^3\end{array}\right\}}^{\left(e\right)}$ (B.3)

where $\mathbf{U}$ denotes the nodal displacements $u$, $\mathbf{W}$ denotes the vector of nodal transverse deflections $w$, and $\mathbf{\Phi }$ denotes the vector of nodal rotations $\varphi$ within the element. The non-zero coefficients of the above equation are

$K_{ij}^{11}={\int }_{r_a}^{r_b}\left\{Ar\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+A\nu \left(\frac{d{\psi }_i}{dr}{\psi }_j+{\psi }_i\frac{d{\psi }_j}{dr}\right)+\frac{A}{r}{\psi }_i{\psi }_j\right\}dr$

$K_{ij}^{13}={\int }_{r_a}^{r_b}\left\{Br\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+B\nu \left(\frac{d{\psi }_i}{dr}{\psi }_j+{\psi }_i\frac{d{\psi }_j}{dr}\right)+\frac{B}{r}{\psi }_i{\psi }_j\right\}dr$

$K_{ij}^{22}={\int }_{r_a}^{r_b}Sr\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}dr,K_{ij}^{23}={\int }_{r_a}^{r_b}Sr\frac{d{\psi }_i}{dr}{\psi }_{j}dr$

$K_{ij}^{31}={\int }_{r_a}^{r_b}\left\{Br\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+B\nu \left(\frac{d{\psi }_i}{dr}{\psi }_j+{\psi }_i\frac{d{\psi }_j}{dr}\right)+\frac{B}{r}{\psi }_i{\psi }_j\right\}dr$

$K_{ij}^{32}={\int }_{r_a}^{r_b}Sr{\psi }_i\frac{d{\psi }_j}{dr}dr$

$K_{ij}^{33}={\int }_{r_a}^{r_b}\left\{Dr\frac{d{\psi }_i}{dr}\frac{d{\psi }_j}{dr}+D\nu \left(\frac{d{\psi }_i}{dr}{\psi }_j+{\psi }_i\frac{d{\psi }_j}{dr}\right)+\frac{D}{r}{\psi }_i{\psi }_j+Sr{\psi }_i{\psi }_j\right\}dr$

$\begin{array}{rr}& \\ & \\ F_i^1& ={\psi }_i\left(r_a\right)Q_1+{\psi }_i\left(r_b\right)Q_4,F_i^2={\psi }_i\left(r_a\right)Q_2+{\psi }_i\left(r_b\right)Q_5,F_i^3={\psi }_i\left(r_a\right)Q_3+{\psi }_i\left(r_b\right)Q_6\end{array}$

Here $i,j=\{\mathrm{1,2},\dots ,NPE\}$, with NPE being the number of nodes in an element; ${\psi }_i$ are the Lagrange interpolation functions used for ($u,w,\varphi$) and $\begin{array}{l}{Q}_1=-{\left[r{N}_{rr}\right]}_{r=r_a},Q_4={\left[r{N}_{rr}\right]}_{r=r_b}\\ Q_2=-{\left[r{Q}_r\right]}_{r=r_a},Q_5={\left[r{Q}_r\right]}_{r=r_b}\\ Q_3=-{\left[r{M}_{rr}\right]}_{r=r_a},Q_6={\left[r{M}_{rr}\right]}_{r=r_b}\end{array}$

Edited by

Editor:

Marcílio Alves.

Publication Dates

• Publication in this collection
  09 Nov 2020
• Date of issue
  2020

History

• Received
  30 July 2020
• Reviewed
  03 Aug 2020
• Accepted
  19 Aug 2020
• Published
  08 Sept 2020