Applications of 2D Padé Approximants in Nonlinear Shell Theory: Stability Calculation and Experimental Justification

1. Introduction

The widest class of shells used in the civil and mechanical engineering is the class of shells with developable principal surface. The stress-strain state of shell structures under loads, which corresponds to buckling, is inhomogeneous, significantly bended, and nonlinear. Permanent interest of researchers in the problem of inhomogeneous compression of shells of zero Gaussian curvature has not led so far to a correct solution. Therefore, there is a need for the development and application of new methods that allow considering the problem in a complex setting, the most appropriate to study real behavior of structures.

The approximate analytic integration of nonlinear differential equations of the theory of flexible elastic shells in most practical cases is based on the method of continuation of solution on the artificially introduced parameter. They can be satisfactorily applied only with an effective method of summation. The most natural analytical continuation method is that using Padé approximants (PAs). PAs effectively solves the problem of analytical continuation of power series, and this is a basis of their successful application in the study of applied problems. Currently, the method of PAs is one of the most promising non-linear methods of summation of power series, and the localization of its singular points. Recently, the method of PAs for single-variable functions has been successfully extended to the approximation of two variable functions (2D PAs).

A method that provides polynomial asymptotics of the exact solution of the general form and its meromorphic continuation based on 2D Padé approximants is proposed in this work. Several examples of displacements, stability and vibration calculations for inhomogeneous loaded shells with developable principal surface are presented. The accuracy of 2D PAs theoretical results are confirmed by experiments with stainless steel specimens based on holographic interferometry. It is shown that the application of PAs provides sufficient accuracy in the studied area that confirms the advantage of our proposed approach.

2. Padé approximants

Let us consider Padé approximants (1-D PAs) which allow us to perform somewhat the most natural continuation of the power series. Below, we are going to define the 1-D PAs [1] for a complex variable$z$:

$$F(z)={\displaystyle \sum _{i=0}^{\infty }{f}_i{z}^i},$$E1

$$F_{nm}(z)={\displaystyle \sum _{i=0}^n{p}_i{z}^i}/{\displaystyle \sum _{i=0}^m{q}_i{z}^i},q_0\equiv 1,$$E2

where coefficients $p_i,q_i$ are determined from the following condition: the first $\left(m+n+1\right)$ components of the expansion of rational function $$F_{nm}(z)$$ in the McLaurin series coincide with the same components of $$F(z)$$ series. Then rational function $$F_{nm}$$ is called [n/m] PAs or 1-D [n/m] PAs. The set of $$F_{nm}$$ functions for different m and n forms the so called Padé table.

PAs make a meromorphic continuation of F due to the following theorem of Montessus de Ballore [10]:

Theorem. Let function $$F(z)$$ be meromorphic in a closed circle $\left|z\right|\le r$ with m different poles $z_i$ of multiplicity ${\mu }_i$ in this circle

$0<\left|z_1\right|\le \left|z_2\right|\le \mathrm{...}\le \left|z_m\right|<r$E3

of total multiplicity M,${\displaystyle \sum _{i=1}^m{\mu }_i}=M$.

Then sequence $$F_{NM}(z)$$ converges uniformly to $$F(z)$$ on compact subsets of this circle without poles, and $z_i$ is attracting zeros of the Padé denominator according to its multiplicity:

$\underset{N\to \infty }{\lim }{F}_{NM}\left(z\right)=F\left(z\right),\left|z\right|\le r,z\ne z_i,i=\overline{1,m}.$E4

The most common generalization of PAs are two-dimensional PAs (2-D PAs). For complex variables $z_1,z_2$ let

$$F(z_1,z_2)={\displaystyle \sum _{i=0}^{\infty }{f}_{ij}{z}_1^i{z}_2^j}$$E5

be a holomorphic function near the origin. For any integer sets $n=\left(n_1,n_2\right)$ and$m=\left(m_1,m_2\right)$, i.e. for any$n,m\in Z_{+}^2$, let

$R(n,m)=\left\{r=\frac{p}{q},p={\displaystyle \sum _{i=0}^{n_1}{\displaystyle \sum _{j=0}^{n_2}{p}_{ij}{z}_1^i{z}_2^j,}q={\displaystyle \sum _{i=0}^{m_1}{\displaystyle \sum _{j=0}^{m_2}{q}_{ij}{z}_1^i{z}_2^j,q_{00}\equiv 1}}}\right\},$E6

be the class of rational functions, i.e. the ratio of 2-D polynomials whose degrees do not exceed $n=\left(n_1,n_2\right)$ and $m=\left(m_1,m_2\right)$ for each variable. It may be written briefly as$\mathrm{deg}(p)\le n,\mathrm{deg}(q)\le m$.

Each rational function $r\in R(n,m)$ may be identified with its power series that converges in some neighborhood of the origin. It should be mentioned that $r=p/q\in R$ depends on ${\tau }_{nm}=(n_1+1)(n_2+1)+(m_1+1)(m_2+1)-1$ parameters (the coefficients of p and q).

The set of integer points $I\left(n,m\right)\subset Z_{+}^2$ for fixed $n=\left(n_1,n_2\right)$ and $m=\left(m_1,m_2\right)$ is called the determinative (interpolation) set, if it has the following properties:

1. $\dim I\left(n,m\right)={\tau }_{nm}$E7
2. $\left(n_1+m_1,0\right),\left(0,n_2+m_2\right)\in I\left(n,m\right)$$z_1=0$$z_2=0$
3. $n=\left(n_1,n_2\right)\in I\left(n,m\right)$E8
4. $\left(k_1,k_2\right)\in I\left(n,m\right)$$\left[0,k\right]\subset I\left(n,m\right)$$\left[0,k\right]=\left\{\left(s_1,s_2\right)\in Z_{+}^2:0\le s_j\le k_j,j=1,2\right\}$
5. $\left(n_1+m_1,m_2\right)\in I(n,m)$

$$I_1\left(n,m\right)=\{\left(i,j\right):$$E10

$$\left[0\le i\le n_1,0\le j\le n_2\right]\cup \left[n_1+1\le i\le n_1+m_1,0\le j\le m_2\right]\cup \left[i=0,n_2+1\le j\le n_2+m_2\right]\},$$E11

$$I_2\left(n,m\right)=\{\left(i,j\right):$$E12

$$\left[0\le i\le n_1,0\le j\le n_2\right]\cup \left[0\le i\le m_1,n_2+1\le j\le n_2+m_2\right]\cup \left[n_1+1\le i\le n_1+m_1,j=0\right]\},$$E13

The generalized PAs for given $n=\left(n_1,n_2\right)$ and $m=\left(m_1,m_2\right)$ are defined as the rational function $F_{nm}\in R\left(n,m\right)$ for which $T_{ij}\left(F-F_{nm}\right)=0$ for all$\left(i,j\right)\in I\left(n,m\right)$, where $T_{ij}\left(\varphi \right)$ are Taylor’s coefficients of the power series for function F. The rational function $F_{nm}$ is called the 2-D PAs $\left[\left(n_1,m_1\right)/\left(n_2,m_2\right)\right]$ of $F\left(z_1,z_2\right)$ which corresponds to the determinative set$\left(i,j\right)\in I\left(n,m\right)$.

As in the 1-D case, the existence and uniqueness of PAs (in the sense of the above given definition) for $C^2$ require special type of analysis. It should be mentioned that PAs do not always exist in the sense of the given definition.

Let $m=\left(m_1,m_2\right)\in Z_{+}^2$ be fixed and let the class

$M_m=M_m\left(C^2\right)=\left\{F:F\left(z_1,z_2\right)=\frac{P\left(z_1,z_2\right)}{Q_m\left(z_1,z_2\right)}\right\}$E14

be defined as a class of functions with the properties:

• $P\left(z_1,z_2\right)$
• $\mathrm{deg}{Q}_m=m$

• functions $P\left(z_1,0\right),P\left(0,z_2\right)$ and polynomials $Q_m\left(z_1,0\right),Q_m\left(0,z_2\right)$ are not equal to zero simultaneously.

The most important theorem for using 2-D PAs for meromorphic continuation is the following Montessus de Ballore – type theorem [2]:

Theorem. Let $F\left(z_1,z_2\right)\in M_m$ be given by the power series, $m=\left(m_1,m_2\right)\in Z_{+}^2$be fixed and$n=\left(n_1,n_2\right)\in Z_{+}^2$. Then:

1. For all $n\text{'}=\min \left(n_1,n_2\right)$ that are large enough, there is a unique Padé approximant $$F_{nm}=P_n/q_n$$ for each of the determinative sets$I_j\left(n,m\right),j=1,2$;

2. The sequence $F_{nm}$ for $n\text{'}=\min \left(n_1,n_2\right)\to \infty$ converges uniformly to function $F\left(z_1,z_2\right)$ inside the compact subsets of$G=C^2\backslash \left\{Q_m=0\right\}$. For any compact $E\subset C^2$ the following relationships are true:

$\begin{array}{l}\underset{n\text{'}\to \infty }{\lim }{\Vert Q_m-q_n\Vert }_E^{1/n\text{'}}=0,\\ \underset{n\text{'}\to \infty }{\lim }{\Vert F-F_{nm}\Vert }_E^{1/n\text{'}}=0,\end{array}$E17

where j = 1, 2 and${\Vert *\Vert }_E=\underset{z\in E}{\sup }|*|$.

This is an analog of the classical Montessus de Ballore theorem for the convergence of the rows of Padé tables.

3. Modified method of parameter continuation

Let us introduce a formal definition of the proposed modified method of parameter continuation (MMPC) for systems of ODEs using the terminology of the perturbation method. It is known [3] that any ODE or system of ODEs may be represented by a normal system of ODEs of the first order in respect to unknown functions ${\left\{u_i=u_i\left(\xi \right)\right\}}_{i=1}^n$ in the vicinity of regular point in the interval$\Omega :\xi \in ]0,1[$:

$L{u}_i+R_i\left(\xi ,u_1,\mathrm{...},u_n\right)+N_i\left(\xi ,u_1,\mathrm{...},u_n\right)=g_i\left(\xi \right),L=\frac{d}{d\xi },i=\overline{1,n},$E18

with the BCs on the bounds $\partial \Omega :\xi =0\cup 1$

$${G_j(u_1,\mathrm{...},u_n)|}_{\partial \Omega }=0,j=\overline{1,n}$$E19

Here L and $R_i$ are the linear differential operators, whereas $N_i$ and $G_j$ are the non-linear differential operators. We assume also that point ${\xi }_0=0$ belongs to closure$\Omega$, and$R_i$, $N_i$and $G_j$ are the holomorphic functions for$${\left\{u_i\right\}}_{i=1}^n$$.

Considering ${\left\{u_i=u_i\left(\xi \right)\right\}}_{i=1}^n$ and their derivatives as independent arguments, we introduce operators$R_i$, $N_i$, $F$and $G_j$ as the multidimensional Taylor series

$$R_i+N_i={\displaystyle \sum _{j=1}^n\left(N_{ij}{u}_j+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{ijp}{u}_j{u}_p+\mathrm{...}}\right)},i=\overline{1,n},$$E20

$$F=\left(F_{0}L{u}_1+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{F}_{0p}{u}_{p}L{u}_1+\mathrm{...}}\right)+{\displaystyle \sum _{j=1}^n\left(F_j{u}_j+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{F}_{jp}{u}_j{u}_p+\mathrm{...}}\right)},$$E21

$$G_j={\displaystyle \sum _{q=1}^n\left(G_{jq}\left(u_q-{u_q|}_{\partial \Omega }\right)+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{G}_{jqp}\left(u_q-{u_q|}_{\partial \Omega }\right)\left(u_p-{u_p|}_{\partial \Omega }\right)+\mathrm{...}}\right)},j=\overline{1,n}.$$E22

We also introduce the following power series

$$N_{ij}={\displaystyle \sum _{r=0}^{\infty }{N}_{ij}^r{\xi }^r,}{N}_{ijp}={\displaystyle \sum _{r=0}^{\infty }{N}_{ijp}^r{\xi }^r,\mathrm{...}{F}_j={\displaystyle \sum _{r=0}^{\infty }{F}_j^r{\xi }^r,}{F}_{jp}={\displaystyle \sum _{r=0}^{\infty }{F}_{jp}^r{\xi }^r,\mathrm{...}{g}_i={\displaystyle \sum _{j=0}^{\infty }{g}_{ij}{\xi }^i}}}i,j,p=\overline{1,n}.$$E23

To implement the MMPC, we introduce parameter $$\epsilon$$ as follows

$u_i={\displaystyle \sum _{j=0}^{\infty }{u}_{ij}^M{\epsilon }^j},$E24

$L{u}_i=\epsilon \left(g_i-R_i\left(u_1,\mathrm{...},u_n\right)-N_i\left(u_1,\mathrm{...},u_n\right)\right),i=\overline{1,n},$E25

${G_j({u_1|}_{\partial \Omega },\mathrm{...},{u_n|}_{\partial \Omega })|}_{\partial \Omega }=0,j=\overline{1,n}$E26

where $R_i$ and $N_i$ are always the algebraic operators in this case.

Substituting power series (6) into (7) and splitting it with respect to the powers of$$\epsilon$$, we get

$\begin{array}{l}{\epsilon }^0:L{u}_{i0}^M=0,{G_i(u_{10}^M,\mathrm{...},u_{n0}^M)|}_{\partial \Omega }=0\Rightarrow u_{i0}^M={u_i^{}|}_{\partial \Omega },i=\overline{1,n},\\ {\epsilon }^1:L{u}_{i1}^M-g_i+{\displaystyle \sum _{r=1}^n\left(N_{ir}{u}_{r0}^M+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}{u}_{r0}^M{u}_{p0}^M+\mathrm{...}}\right)}=0,{u_{i1}^H|}_{\partial \Omega }=0\Rightarrow \\ u_{i1}^M={\displaystyle \sum _{j=0}^{\infty }\frac{{\xi }^{j+1}}{\left(j+1\right)}}\left(g_{ij}-{\displaystyle \sum _{r=1}^n\left(N_{ir}^j{u_r^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N_{irp}^j{u}_r^{}|}_{\partial \Omega }{u_p^{}|}_{\partial \Omega }+\mathrm{...}}\right)}\right),\\ \\ u_{i1}^M={\displaystyle \sum _{j=0}^{\infty }\frac{{\xi }^{j+1}}{\left(j+1\right)}}\left(g_{ij}-{\displaystyle \sum _{r=1}^n\left(N_{ir}^j{u_r^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N_{irp}^j{u}_r^{}|}_{\partial \Omega }{u_p^{}|}_{\partial \Omega }+\mathrm{...}}\right)}\right),i=\overline{1,n},\\ {\epsilon }^2:L{u}_{i2}^M+{\displaystyle \sum _{r=1}^n(N_{ir}{u}_{r1}^M}+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}\left(u_{r1}^M{u}_{p0}^M+u_{r0}^M{u}_{p1}^M\right)+\mathrm{...}})=0,{u_{i2}^M|}_{\partial \Omega }=0\Rightarrow \\ u_{i2}^M=-{\displaystyle \sum _{s=0}^{\infty }{\displaystyle \sum _{r=1}^n(N_{ir}^s({\displaystyle \sum _{j=0}^{\infty }\frac{{\xi }^{s+j+2}}{\left(j+1\right)\left(s+j+2\right)}}}}\left(g_{rj}-{\displaystyle \sum _{l=1}^n\left(N_{rl}^j{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{rlq}^j{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }+\mathrm{...}}\right)})\right)+\\ +\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}^s}({u_p^{}|}_{\partial \Omega }({\displaystyle \sum _{j=0}^{\infty }\frac{{\xi }^{s+j+2}}{\left(j+1\right)\left(s+j+2\right)}}\left(g_{rj}-{\displaystyle \sum _{l=1}^n\left(N_{rl}^j{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{rlq}^j{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)})\right)+\\ +{u_r^{}|}_{\partial \Omega }({\displaystyle \sum _{j=0}^{\infty }\frac{{\xi }^{s+j+2}}{\left(j+1\right)\left(s+j+2\right)}}\left(g_{pj}-{\displaystyle \sum _{l=1}^n\left({N_{pl}^j{u}_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{plq}^j{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)}))+\mathrm{...})\right)+\mathrm{...},i=\overline{1,n}.\end{array}$E27

Summation in (8) of the coefficients with the same degrees of $\xi$ for $$\epsilon =1$$ gives

$$\begin{array}{l}{u}_i={\xi }^0\left(\left[{u_i|}_{\partial \Omega }\right]+\left[0\right]+\left[0\right]+\mathrm{...}\right)+\\ +{\xi }^1\left(\left[0\right]+\left[g_{i0}-{\displaystyle \sum _{r=1}^n\left({N_{ir}^0{u}_r^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}^0{u_r^{}|}_{\partial \Omega }{u_p^{}|}_{\partial \Omega }+\mathrm{...}}\right)}\right]+\left[0\right]+\mathrm{...}\right)+\\ +{\xi }^2(\left[0\right]+\left[\frac{g_{i1}}{2}-\frac{1}{2}{\displaystyle \sum _{r=1}^n\left(N_{ir}^1{u_r^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}^1{u_r^{}|}_{\partial \Omega }{u_p^{}|}_{\partial \Omega }+\mathrm{...}}\right)}\right]+\\ +[-{\displaystyle \sum _{r=1}^n(N_{ir}^0(\frac{g_{r0}}{2}}-\frac{1}{2}{\displaystyle \sum _{l=1}^n\left(N_{rl}^0{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{rlq}^0{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }+\mathrm{...}}\right)})+\\ +\frac{1}{2!}{\displaystyle \sum _{p=1}^n}{N}_{irp}^0({u_p^{}|}_{\partial \Omega }(\frac{g_{r0}}{2}-\frac{1}{2}{\displaystyle \sum _{l=1}^n\left(N_{rl}^0{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{rlq}^0{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)})+\\ +{u_r^{}|}_{\partial \Omega }(\frac{g_{p0}}{2}-\frac{1}{2}{\displaystyle \sum _{l=1}^n\left(N_{pl}^0{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{plq}^0{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)}))+\mathrm{...})]+\mathrm{...})+\mathrm{...},i=\overline{1,n}.\end{array}$$E28

Analysis of the obtained approximation suggests that it gives the exact value of coefficients in the power of the independent variable to the extent equal to the order of approximation (taking into account the expansion in power series of expressions in the equation). This guarantees stability of the computation with a limit-order approximation of the independent variable.

One of the possible fields for application of the proposed approach are the nonlinear problems of plates and shells theory. The equations of statics of the geometrically nonlinear thin-walled structures can be reduced to the resolving equations which contain the products and squares of the desired functions and their derivatives [4]. In this case the solution of equation (9) becomes

$$\begin{array}{l}{u}_i={\xi }^0{u_i|}_{\partial \Omega }+{\xi }^1\left(g_{i0}-{\displaystyle \sum _{r=1}^n\left({N_{ir}^0{u}_r^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}^0{u_r^{}|}_{\partial \Omega }{u_p^{}|}_{\partial \Omega }}\right)}\right)+\\ +{\xi }^2(\left(\frac{g_{i1}}{2}-\frac{1}{2}{\displaystyle \sum _{r=1}^n\left(N_{ir}^1{u_r^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{p=1}^n{N}_{irp}^1{u_r^{}|}_{\partial \Omega }{u_p^{}|}_{\partial \Omega }}\right)}\right)-\\ -{\displaystyle \sum _{r=1}^n(N_{ir}^0(\frac{g_{r0}}{2}}-\frac{1}{2}{\displaystyle \sum _{l=1}^n\left(N_{rl}^0{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{rlq}^0{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)})+\\ +\frac{1}{2!}{\displaystyle \sum _{p=1}^n}{N}_{irp}^0({u_p^{}|}_{\partial \Omega }(\frac{g_{r0}}{2}-\frac{1}{2}{\displaystyle \sum _{l=1}^n\left(N_{rl}^0{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{rlq}^0{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)})+\\ {+u_r^{}|}_{\partial \Omega }(\frac{g_{p0}}{2}-\frac{1}{2}{\displaystyle \sum _{l=1}^n\left(N_{pl}^0{u_l^{}|}_{\partial \Omega }+\frac{1}{2!}{\displaystyle \sum _{q=1}^n{N}_{plq}^0{u_l^{}|}_{\partial \Omega }{u_q^{}|}_{\partial \Omega }}\right)}))))+\mathrm{...},i=\overline{1,n}.\end{array}$$E29

The approximation thus obtained is converted to 1-D PAs in respect to $\xi$ or 2-D PAs. 2-D PAs in the form proposed by V. Vavilov [2] is very promising for the use as an analytical continuation. This technique allows us to choose the coefficients of 2-D Taylor series for the construction of unambiguous 2-D PAs with a given structure of the numerator and denominator. It also ensures optimal PAs features in the sense of the theorem of Montessus de Ballore-type. This means homogenous convergence of PAs to the approximated function with an increase of the degree of the numerator and denominator in all points of its meromorphic area. It should be noted that direct application of 2-D PAs does not lead to the anticipated merging of 1-D approximations. This is due to the initial requirements of the 2-D approximation to ensure its transition to 1-D in the case when the second variable is equal to zero [1]. At the same time it is necessary to ensure such a transition, when the parameter is equal to one. This can be achieved by combining this method with 2-D PAs from a converted parameter which maps the unit to zero.

4. Stability investigation

For the analytic continuation of meromorphic solutions to the region and acceleration of the convergence we use 2-D PAs. To do this, the resulting series are reconstructed in rational functions of the form

$$u_i={\displaystyle \sum _{j=0}^{m_1}{\displaystyle \sum _{k=0}^{n_1}{P}_{ijk}\left(P\right){\epsilon }^j{\xi }^k}}/{\displaystyle \sum _{j=0}^{m_2}{\displaystyle \sum _{k=0}^{n_2}{Q}_{ijk}\left(P\right){\epsilon }^j{\xi }^k}},Q_{i00}\equiv 1,$$E30

where $P$ is the parameter of loading.

Since the proposed method is modified, all functions belonging to the boundary value problems can be expanded in powers of the independent variables and parameters. For $$\epsilon =1$$ we obtain

$$u_i={\displaystyle \sum _{j=0}^{m_3}{\displaystyle \sum _{k=0}^{n_1}{\overline{P}}_{ijk}{P}^j{\xi }^k}}/{\displaystyle \sum _{j=0}^{m_4}{\displaystyle \sum _{k=0}^{n_2}{\overline{Q}}_{ijk}{P}^j{\xi }^k}}.$$E31

Rational functions have singular points which are determined by equating of the denominator to zero. According to the theory of bifurcation of solutions of ordinary differential equations, at these points either bifurcations or limit states are achieved. So we can get the estimation of critical point localization by solving equations of the form

$$\underset{i,\xi }{\min }P:{\displaystyle \sum _{j=0}^{m_4}{\displaystyle \sum _{k=0}^{n_2}{\overline{Q}}_{ijk}{P}^j{\xi }^k}}=0$$E32

.

In practice this equation is transformed to its counterpart simpler form with respect to its characteristic point $\xi ={\xi }_0$ (which usually corresponds to the point of maximum transverse displacement for thin shells)

$${\displaystyle \sum _{j=0}^{m_4}{\displaystyle \sum _{k=0}^{n_2}{\overline{Q}}_{ijk}{P}^j{\xi }_0^k}}=0$$E33

.

Example

Let us consider the computational aspects of the proposed approach. We consider three types of PAs with respect to the independent variable, on the specified parameters, and 2-D.

A typical behavior of the approximations for the BVP is governed by the following problem

$$\begin{array}{c}\epsilon z\text{'}+z=1,\\ z(0)=0,0<\epsilon <<1,x\ge 0.\end{array}$$E34

where natural small parameter $$\epsilon$$ is the factor at the highest derivative, as shown in Fig. 1 for ε = 0.1. The exact solution of this BVP follows

$$z=1-\exp \left(-\frac{x}{\epsilon }\right)=\frac{x}{\epsilon }-\frac{x^2}{2\epsilon }+\mathrm{...}+{(-1)}^{n+1}\frac{1}{n!}{\left(\frac{x}{\epsilon }\right)}^n+\mathrm{...}.$$E35

Eq. (12) shows that the exact solution is regular for all real positive $x$ for$\epsilon \ne 0$. But the general term of power series (12) grows rapidly when$x>\epsilon$, and we have to take into account many terms in (2) to obtain an acceptable and reliable approximation. Thus, the accuracy of the used truncated Taylor series is not uniform according to $x$ value.

Figure 1.

The exact solution (solid line) of Eq. (1) for ε = 0.1 and approximate solutions (1 – three terms ADM, 2 –$z^{\left({\epsilon }_1\right)}$

for ADM, 3 – three terms HAM, 4 –$z^{\left(x\right)}$for HAM, 5 – 2-D PAs for MMPC, ADM and HAM).

Let us introduce parameter $${\epsilon }_1$$ as follows

$$z\text{'}=\frac{1}{\epsilon }-{\epsilon }_1\frac{z}{\epsilon },$$E36

and suppose

$$z={\displaystyle \sum _{i=0}^{\infty }{z}_i{\epsilon }_1^i}.$$E37

This way of introducing the parameter leads to a system of successive approximations of Adomian decomposition method (ADM, see [5,6]). Substituting power series (4) into Eq. (3) and splitting it with respect to the powers of$${\epsilon }_1$$, yields

$${\epsilon }_1^0:z_0^{\text{'}}=\frac{1}{\epsilon },z_0(0)=0\Rightarrow z_0^{}=\frac{x}{\epsilon };$$E38

$${\epsilon }_1^1:z_1^{\text{'}}=-\frac{z_0}{\epsilon }=-\frac{x}{{\epsilon }^2},z_1(0)=0\Rightarrow z_1^{}=-\frac{x^2}{2{\epsilon }^2}={(-1)}^1\frac{1}{2!}{\left(\frac{x}{\epsilon }\right)}^2;$$E39

$${\epsilon }_1^2:z_2^{\text{'}}=-\frac{z_1}{\epsilon }=-\frac{x^2}{{\epsilon }^3},z_2(0)=0\Rightarrow z_2^{}=\frac{x^3}{6{\epsilon }^3}={(-1)}^2\frac{1}{3!}{\left(\frac{x}{\epsilon }\right)}^3;$$E40

$${\epsilon }_1^n:z_n^{\text{'}}=-\frac{z_{n-1}}{\epsilon }={(-1)}^n\frac{1}{(n-1)!}{\left(\frac{x}{\epsilon }\right)}^n\left(\frac{1}{\epsilon }\right),z_n(0)=0\Rightarrow z_n^{}={(-1)}^n\frac{1}{n!}{\left(\frac{x}{\epsilon }\right)}^{n+1};$$E41

$$z=\frac{x}{\epsilon }-\frac{x^2}{2{\epsilon }^2}{\epsilon }_1+\frac{x^3}{6{\epsilon }^3}{\epsilon }_1^2+\mathrm{...}+{(-1)}^n\frac{1}{n!}{\left(\frac{x}{\epsilon }\right)}^{n+1}{\epsilon }_1^n+\mathrm{...}.$$E42

For $${\epsilon }_1=1$$ one gets power series (2), which corresponds to the results reported in [5]. To accelerate the convergence we use 1-D and 2-D PAs. One can use 1-D PAs [1/1] ($z^{({\epsilon }_1)},$when $x=const\ne 0$ and $z^{(x)},$ when${\epsilon }_1=const\ne 0$) or 2-D PAs [(1,1)/(1,1)] ($z^{({\epsilon }_1,x)}$). Using PAs one obtains (for$$\epsilon =1$$)

$z^{({\epsilon }_1)}=\frac{x}{\epsilon }\left(1-\frac{3x}{6\epsilon +2x}\right),$E43

$z^{(x)}=z^{({\epsilon }_1,x)}=\frac{2x}{2\epsilon +x}.$E44

R.h.s. of Eq. (15) contains singularity at point $$\epsilon =0$$ in contrast to Eqs. (16). Thus, Eqs. (16) give the approximation with uniform accuracy when $x$ grows.

Let us rewrite Eq. (11) in the following form

$$z\text{'}={\epsilon }_1(1-z+(1-\epsilon )z\text{'}).$$E45

After substitution of the power series (14) in Eq. (17) one obtains a successive approximation of the homotopy analysis method (HAM) in the so called homotopy perturbation form (HPM, see [7])

$${\epsilon }_1^0:z_0^{\text{'}}=0,z_0(0)=0\Rightarrow z_0^{}=0;$$E46

$${\epsilon }_1^1:z_1^{\text{'}}=1-z_0+(1-\epsilon )z_0^{\text{'}}=1,z_1(0)=0\Rightarrow z_1^{}=x;$$E47

$${\epsilon }_1^2:z_2^{\text{'}}=-z_1+(1-\epsilon )z_1^{\text{'}}=-x+(1-\epsilon ),z_2(0)=0\Rightarrow z_2^{}=-\frac{x^2}{2!}+(1-\epsilon )x;$$E48

$$\begin{array}{l}{\epsilon }_1^3:z_3^{\text{'}}=-z_2+(1-\epsilon )z_2^{\text{'}}=\frac{x^2}{2}-(1-\epsilon )x+(1-\epsilon )(-x+(1-\epsilon )),z_3(0)=0\Rightarrow \\ \Rightarrow z_3^{}=\frac{x^3}{3!}-(1-\epsilon )x^2+{(}^{1}x;\end{array}$$E49

$$\begin{array}{l}{\epsilon }_1^4:z_4^{\text{'}}=-z_3+(1-\epsilon )z_3^{\text{'}}=-\frac{x^3}{6}+(1-\epsilon )x^2-{(}^{1}x+(1-\epsilon )\left(\frac{x^2}{2}-2(1-\epsilon )x+{(1-\epsilon )}^2\right),\\ z_4(0)=0\Rightarrow z_4^{}=-\frac{x^4}{4!}+(1-\epsilon )\frac{x^3}{2}-3{(}^1\frac{x^2}{2}+{(}^{1}x;\end{array}$$E50

We obtain the HAM approximation in the following form

$$z=x{\epsilon }_1+\left(-\frac{x^2}{2!}+(1-\epsilon )x\right){\epsilon }_1^2+\left(\frac{x^3}{3!}-(1-\epsilon )x^2+{(1-\epsilon )}^{2}x\right){\epsilon }_1^3+\mathrm{...}.$$E51

For $${\epsilon }_1=1$$ one obtains

$$\begin{array}{l}z=(1+(1-\epsilon )+{(}^1+{(}^1+\mathrm{...})x+\left(-\frac{1}{2!}+(1-\epsilon )-3{(1-\epsilon )}^2\frac{1}{2}+\mathrm{...}\right)x^2+\\ +\left(\frac{1}{3!}+(1-\epsilon )\frac{1}{2}+\mathrm{...}\right)x^3-\frac{x^4}{4!}+\mathrm{...}.\end{array}$$E52

Eq. (18) coincides with Eq. (12) after expanding the coefficients of Eq. (12) in the vicinity of $$\epsilon =1.$$

For the obtained approximations we use PAs as described above and get the following results:

$z^{({\epsilon }_1)}=z^{({\epsilon }_1,x)}=\frac{2x}{2\epsilon +x},$E53

$z^{(x)}=\frac{2{\left(2-\epsilon \right)}^{2}x}{2\left(2-\epsilon \right)+x}.$E54

Let us introduce parameter $${\epsilon }_1$$ in such a way that

$$z\text{'}={\epsilon }_1\frac{1-z}{\epsilon }.$$E55

After substitution of the power series (14) in Eq. (19), one obtains a new system of successive approximations:

$${\epsilon }_1^0:z_0^{\text{'}}=0,z_0(0)=0\Rightarrow z_0^{}=0;$$E56

$${\epsilon }_1^1:z_1^{\text{'}}=\frac{1-z_0}{\epsilon }=\frac{1}{\epsilon },z_1(0)=0\Rightarrow z_1^{}=\frac{x}{\epsilon }={(-1)}^2\frac{1}{1!}{\left(\frac{x}{\epsilon }\right)}^1;$$E57

$${\epsilon }_1^2:z_2^{\text{'}}=-\frac{z_1}{\epsilon }=-\frac{x}{{\epsilon }^2},z_2(0)=0\Rightarrow z_2^{}=-\frac{x^2}{2{\epsilon }^2}={(-1)}^3\frac{1}{2!}{\left(\frac{x}{\epsilon }\right)}^2;$$E58

$${\epsilon }_1^n:z_n^{\text{'}}=-\frac{z_{n-1}}{\epsilon }={(-1)}^{n+1}\frac{1}{(n-1)!}{\left(\frac{x}{\epsilon }\right)}^{n-1}\left(\frac{1}{\epsilon }\right),z_2(0)=0\Rightarrow z_n^{}={(-1)}^{n+1}\frac{1}{n!}{\left(\frac{x}{\epsilon }\right)}^n;$$E59

$$z=\frac{x}{\epsilon }{\epsilon }_1-\frac{1}{2!}{\left(\frac{x}{\epsilon }\right)}^2{\epsilon }_1^2+\frac{1}{3!}{\left(\frac{x}{\epsilon }\right)}^3{\epsilon }_1^3+\mathrm{...}.$$E60

For $${\epsilon }_1=1$$ one obtains PAs as follows:

$$z^{({\epsilon }_1)}=z^{(x)}=z^{({\epsilon }_1,x)}=\frac{2x}{2\epsilon +x}$$E61

.

The coincidence of these approximations demonstrates their proximity to a unique function representing the exact solution in the fractional-rational form.

The ADM approximation describes well the exact solution only for a distance which is comparable with the value of natural small parameter$$\epsilon$$. Despite the fact that the error of solutions of HAM is substantially less than the ADM, HAM does not accurately reflect the nature of solutions, namely the phenomenon of boundary layer in the vicinity of zero. At the same time, PAs for the ADM approximations for independent variable and PAs for the MMPC (1-D and 2-D) give satisfactory qualitative and quantitative results.

Similar results are provided by the analysis of approximations of BVP of the following problem

$$\begin{array}{c}\epsilon z\text{'}+xz=x,\\ z(0)=2,0<\epsilon <<1,x\ge 0,\end{array}$$E62

whose coefficients are given depending on the variable for ε = 0.2. Its graphs are presented in Fig. 2.

Figure 2.

The exact solution (solid line) of Eq. (20) for ε = 0.2 and approximate solutions (1 – three terms ADM, 2 –$z^{\left({\epsilon }_1\right)}$

for ADM, 3 – three terms HAM, 4 –$z^{\left(x\right)}$for HAM, 5 – 2-D PAs for MMPC, ADM and HAM).

Fig. 3 shows the graphs of approximations for strongly non-linear BVP of the form

$\begin{array}{c}\epsilon z\text{'}=z^2+x,\\ z\left(1\right)=1,0<\epsilon <<1,x\ge 0.\end{array}$E63

The graphs show that the solution is well described by the HAM approximation and MHAM-Padé «in average», and badly – in the boundary layer. The ADM approximation and MADM-Padé, on the contrary, is in good agreement with the behavior of solution in the vicinity of zero and in the bad one – on the stationary part. At the same time, 1-D and 2-D PAs, based on approximations of the MMPC, well describe the solution in the whole interval.

5. Calculation of nonlinear deformation and stability of shells

The proposed MMPC method has been applied to calculate the deformation and stability of a long flexible elastic circular cylindrical shell of radius $$R$$ with half the central angle $${\beta }_0$$ in the case of cylindrical bending under uniform external pressure with a simple support of the longitudinal edges. The corresponding system of resolving equations in the normal form is given in [8]. Dependences of "dimensionless intensity of pressure $$P$$ – deflection $w/R$" for the top cross-section of the shell at different angles and dimensionless flexibility $$С={10}^{-4}$$ are shown in Fig. 4. The dependence of the dimensionless intensity of limit load $$P_{*}$$ on the size of half angle $${\beta }_0$$ is shown in Fig. 5.

Figure 3.

Approximate solutions of Eq. (18) for ε = 0.2 (1 – three terms ADM, 2 –$z^{\left({\epsilon }_1\right)}$for ADM, 3 –$z^{\left(x\right)}$

for ADM,, 4 – three terms HAM, 5 –$z^{\left(x\right)}$for HAM, 6 –$z^{\left({\epsilon }_1\right)}$for HAM and MMPC, 7 –$z^{\left(x\right)}$and 2-D Padé for MMPC).

For comparison, Fig. 4b also shows the dependence of the critical loads for inextensible shell obtained by S.P. Timoshenko [8]. We see that dependences are in good agreement, while consideration of deformation of the longitudinal axis substantially affects the value of critical loads of the construction.

Figure 4.

The dependence of the intensity of pressure $$P$$ versus deflection $$w/R$$ for different values of $${\beta }_0$$ (the value of $${\beta }_0$$ is indicated by curves).

The proposed method can be used in a combination with the known asymptotic method. Consider free vibrations of a flexible elastic circular cylindrical shell of radius$$R$$, thickness $h$ and length$L$, backed by a set of uniformly distributed stringers having a simple support at the ends.

The calculation is based on mixed dynamical equations of the theory of shells after splitting them in powers of natural small parameters [9]. The shape of radial deflection $$w$$ satisfies the boundary conditions given in the form

$$w=f_1(t)\sin s_1{x}_1\cos s_2{x}_2+f_2(t){\sin }^2{s}_1{x}_1$$E64

.

Here $$f_1,f_2$$ functions depend on time and are related by the condition of continuity of displacements

$$f_2=0,25R^{-1}{s}_2^2{f}_1^2$$E65

,

where $$s_1=\pi m{l}^{-1},s_2=n$$ are the parameters characterizing the wave generation along the generator and directrix, respectively.

The governing equations can be reduced by the Bubnov – Galerkin method to the Cauchy problem with respect to $$\xi =f_1/R$$ on $$t_1=t\sqrt{B_1/\rho R^2}$$ (all symbols are taken in accordance with [9])

$$\begin{array}{l}\ddot{\xi }+\alpha \xi \left[{\left(\dot{\xi }\right)}^2+\xi \ddot{\xi }\right]+A_1\xi +A_2{\xi }^3+A_3{\xi }^5=0,\\ \begin{array}{ccc}& & \left(\stackrel{·}{}\right)\equiv \frac{d\left(\right)}{d{t}_1},t_1=0:\xi =f,\dot{\xi }=0\end{array}\end{array}$$E66

Figure 5.

The dependence of limit loads $$P_b$$ versus $${\beta }_0$$ (1 - data [8], 2 - calculation).

The application of the proposed method of parameter continuation to the Cauchy problem (22) gives approximation of the second order for the artificial parameter for frequency $\Omega$ of nonlinear oscillations in the form

$$\Omega =\sqrt{\frac{1+f^2\left(A_2/A_1\right)+f^4\left(A_3/A_1\right)}{\left(1+\alpha f\right)}}.$$E67

.

It is seen that the oscillations are not isochronous. This agrees well with previous results reported in reference [9] (Fig. 6). However, our approach allows for a significant reduction of the computation time (in [9] to obtain similar results the approximation of the fourth order is taken).

Figure 6.

Amplitude of the initial disturbance versus oscillation frequency of stringer shell (1 – according the proposed method, 2 – data [9]).

6. Experimental technique

In solving various kinds of problems of modern development and improvement of thin-walled machine elements operating under the conditions of intensive manufacturing process, the use of a holographic interferometry method should be emphasized [10,11,14]. It allows for a more accurate and complete investigation of shell structures under complex stress-strain state. The accuracy of interpretation of holographic interferograms is mainly determined by the number of support points of the design used for the construction regarding displacements and stresses. Improvement of the accuracy requires a large amount of routine preparations for writing the coordinates of points and their corresponding numbers of lines when developing data on a computer, which is particularly important in the case of an experiment. The existing methods of automated data entry and processing of interferograms yield, as a rule, the specific configuration of the optical system and the types of strain state (flat, one-dimensional, etc.), making them difficult to use in this case. In addition, although several authors proposed methods of interpretation [10], they did not fully take into account the statistical nature of input data. For cylindrical shells a method for automated processing of the results of the holographic research has been proposed, which eliminates the above drawbacks [12]. Next, we have extended this technique to study the motion of shell structures of zero Gaussian curvature which is based on modern means of an interactive data processing. The surface of zero Gaussian curvature can be approximated with sufficient accuracy with respect to the system of flat rectangular panels whose sides are segments close to the case which occurred during the analysis of generators. To determine all components (points) of the displacement vector, three holograms of a circuit design interferometer based on a reference beam is used. The interferometer is shown schematically in Fig 7.

Figure 7.

The scheme of the interferometer (1 – laser generator, 2 – mirror, 3 – expanding lens, 4 – studied object, 5 – camera)

After registering the two exposures, i.e. unloaded and loaded state of the object, we get a flat image of the interference pattern corresponding to the observation of points$$M_i\left(x_{нi},y_{нi},z_{нi}\right)$$,$$i=\overline{1,3}$$. Let us enter the order line using a computer in the following manner. The photos of interferograms are scanned and entered into the computer memory in the form of graphic files with the extension, for example, jpg, which is the most popular choice of compression of graphic information on all platforms, or equivalently in other file formats. Next, the file is displayed on the screen in a specially designed box on the toolbar image processing. The information produced is removed by a successive mouse click on the corresponding image points at the request of a specially created database. Algorithms for further processing of the data are widely described in [12]. In the $$X^{\prime }{O}^{\prime }{Y}^{\prime }$$ coordinate system (Fig. 8) associated with the imaging plate, base point $$M_B\left({x^{\prime }}_B,{y^{\prime }}_B\right)$$ and a segment of the $$OY$$ axis of the $$XOY$$coordinate system, whose direction coincides with the vertical axis of the projection, are given. Further calculations are performed in the $$XOY$$ system in which the entered coordinates of the points of lines of equal order are transformed by the formulas

$$x=x^{\prime }\cos \phi +y^{\prime }\sin \phi -x_B,y=y^{\prime }\cos \phi -x^{\prime }\sin \phi -y_{B,}$$E68

where $$\phi$$is the angle of rotation of the $$XOY$$ system with respect to $$X^{\prime }{O}^{\prime }{Y}^{\prime }$$

$$x_B={x^{\prime }}_B\cos \phi +{y^{\prime }}_B\sin \phi ,y_B=-{x^{\prime }}_B\sin \phi +{y^{\prime }}_B\cos \phi .$$E69

Figure 8.

The scheme of approximation of the surface by the system of folds

While computing the physical coordinates, the approximation of the mentioned shell surface by the system of folds is applied (Fig. 8). In this case, the physical coordinates of point $$M(r,z,\theta )$$ in the $$i$$ fold $$(i=\overline{1,n)}$$ are defined by

$$\theta =\frac{{\theta }_i-{\theta }_{i-1}}{x_i-x_{i-1}}\left(x-x_{i-1}\right)+{\theta }_{i-1}$$E70

$$z=z_1+\frac{y-\left[y_{2,i-1}+\left(y_{2,i}-y_{2,i-1}\right)\left(x-x_{2,i-1}\right)/\left(x_{2,i}-x_{2,i-1}\right)\right]}{y_{1,i}-y_{2,i}}\left(z_2-z_1\right),$$E71

where $$r,\theta ,z$$ are the coordinates of a point in the cylindrical coordinate system associated with the axis of the shell, and the shape of the surface is analytically given by equation$$r=r(\theta ,z)$$; $$r_i,{\theta }_j,z_k$$are the physical coordinates of corner points of the considered fold;$$\left(x_{1,i};y_{1,i}\right)$$, $$\left(x_{2,i};y_{2,i}\right)$$are the coordinates of corner points of the projection of the folds into the $$XOY$$ system (Fig. 8).

Therefore, the so far obtained arrays of point coordinates of the lines of equal order corresponding to the three noncoplanar directions of observation allow us to approximate the surface bands. The most appropriate method to do this is the structural analysis of extrapolation (MSEA) [6] using step by step the best choice of the model. Indeed, the formalization of the input source data for inhomogeneous stress-strain state requires a large number of points.

The above method allows us to determine the coordinates of the centers of bands up to 0.1 mm without any additional devices. This procedure is used to significantly increase the number of input points (up to 200-400 for each direction of observation). The use of spline functions for smoothing requires the enumeration of all coordinates of control points for each calculation of the order of the band. This slows down the calculation and requires a significant memory space. In addition, these disadvantages are compounded by the increasing number of control points. The use of MSEA allows each step to obtain unbiased estimates of the effective coefficients of the model to ensure a maximum plausible value of the order of the reference points [11]. This eliminates the problem of choosing a smoothing parameter, with the number needed to calculate the coefficients one order of magnitude smaller than the number of coordinates of reference points. In addition, the incremental method allows us to formalize the process of selecting the optimal order of approximating polynomial based on the assessment of the significance of the model and the adequacy of its source data. Note that in this case the number of points is much larger than the number of estimated parameters, which suggests a considerable power of the statistical tests (like those of Student's, Fisher and Durbin-Watson), and indicates the validity of hypotheses taken in selecting the best model. An increase of the number of points improves a regression model, and a loss of accuracy in the summation can be successfully overcome by standardizing the original data according to the known methods.

Because the shape of the surface is analytically given by equation$$r=r(\theta ,z)$$, it is possible to obtain two-dimensional regression models for the line of the $i$-order for the direction of observation in the form

$$N_i={\displaystyle \sum _{j=0}^{n(i)}{\displaystyle \sum _{k=0}^{n(j)}{b}_{jki}{\theta }^j{z}^k}}.$$E72

Displacements are defined by the equation [4]:

$$MU=\lambda N,$$E73

where:$$M$$– optic matrix;$$U=\left\{u,v,w\right\}$$– vector of displacement;$$\lambda$$– length of the laser wave;$$N={\left\{N_i\right\}}_{i=1}^3$$– vector of lines order.

Further transformation of movements, according to the Cauchy relations and equations of state of the environment, can also yield the stress state at the point. Performing the calculation of the stress-strain state parameters to form and direct the shell with a certain step, it is possible to obtain data for plotting the distribution of displacements and stresses.

7. Experimental justification of MMPC

Loading capacity of cylindrical shells is significantly affected by the unevenness of deformation caused by the ovality ends of the shell [13]. Imperfections in face of shells usually occur as a result of their deformation either under their own weight or during the mechanical handling, storage, as well as installation and assembly of the shells as individual elements. In the case of welded shells, the end face has the form of an oval with $$a$$ and $$b$$ axes and $$a/b$$ compression ratio (or the actual ovality). At the same roundness of the upper end (a/b)_B may be different from the roundness of the bottom one (a/b)_H because of the conditions introduced by a collection with other elements of the design. In all cases the shape of each end should be within the required tolerances, and roundness introduced by the collection process should not be reduced by more than 0.8.

Another form of imperfections arising from the inaccuracy of the assembly is associated with a weak taper angle characterized by forming the membrane to its axle$$\alpha$$. A number of studies [14] consider that a small taper with $$\alpha \le {13}^{\circ }$$ has no significant effect on the magnitude of critical loads of axial compression. However, the results of stability studies of technologically imperfect cylindrical shells based on the multivariate approach [15,16] suggest that an increase of $$\alpha$$ to the value of $$3^{\circ }$$ often leads to a significant change in carrying capacity, and in some cases the interaction with the oval and other factors yields an increase of the critical loads. These studies have shown the need for a more correct approach in establishing the correspondence between the magnitude of these abnormalities and the level of carrying capacity. As a consequence, it is necessary to study the nature of deformation of shells with different ratios of the parameters of roundness and taper.

In order to solve this problem, two-factor second-order experiments on two levels of $$\alpha$$ and a/b, and also on two levels of a/b at lower and upper ends when α° = 1 (when the taper is small the difference between the lower and upper end is missing) have been implemented. Welded specimens with radius R = 71.5 mm and length L = 200 mm, made from plate steel of the mark $$H18N9-n$$ and with thickness δ = 0.25 mm have been tested. The use of a multi-factor approach allows one to solve correctly the problem of the nonlinear joint influence of defects on the loading capacity of the shell.

Tests on the stability of prototypes carried out on a UME-10TM machine showed that the exhaustion of loading capacity of the shell took place at one stage by reaching a limit point.

The loss of stability of a conical shell with the same low ovality ends (Fig. 9a) is in general related to a form close to its own form of stability loss of oval cylindrical shells under the action of uniform axial compression [17], but shifted to a larger shell butt. On one side of the shell there are two or three belt dents located at the larger end. They cover the smaller curvature of the plate and are shifted to the side of panel larger curvature. Local dents have a relatively large size and do not form a regular closed form buckling.

The increase of taper and roundness of the ends leads to a shift of the zone of wave generation into the longitudinal direction to a lower end (Fig. 9b) while maintaining the overall character of buckling.

Figure 9.

Forms of supercritical wave generation of the shell with a small taper and the same low ovality of ends (a); with a large taper, and the same large oval ends (b); with a large taper, and a large oval of the lower extremity (c); with a large taper, and a large oval upper end (d)

At high cone (within a given experiment) increased roundness of the lower end, while maintaining the shape of the upper longitudinal, increases the localization of buckling (Fig. 9c) shifting the dents closer to the lower end, while maintaining the variability in the circumferential direction. Conversely, the prevalence of high cone-roundness of the upper end leads to a significant shift of dents to the end of a large oval (Fig. 9d).

Results of the experiment allow us to derive mathematical models of the form

$${\left(a/b\right)}_{Н}^{\circ }={\left(a/b\right)}_{В}^{\circ }={\left(a/b\right)}_{}^{\circ }:$$

$$K=0,379+0,0029{\alpha }^{\circ }-0,012{\left(a/b\right)}_{}^{\circ }-0,014{\alpha }^{\circ }{\left(a/b\right)}_{}^{\circ };$$E74

$${\alpha }^{\circ }=1:$$

$$K=0,346-0,017{\left(a/b\right)}_{Н}^{\circ }-0,0033{\left(a/b\right)}_{В}^{\circ }-0,013{\left(a/b\right)}_{Н}^{\circ }{\left(a/b\right)}_{В}^{\circ },$$E75

where (...)° is the standardized value; $$K=\frac{{Т}_{cr}}{2\pi E{\delta }^2}$$is the dimensionless ratio of the critical stress; $${Т}_{cr}$$is the critical compressive load; $$E$$is the Young's modulus.

The resulting models are adequate to the experimental data by Fisher criteria at the 5% significance level. The presence of significant second-order terms indicates a significant non-linearity of the relationship between the parameters, and, therefore, incorrect to separate consideration of the parameters and the placement of single-factor experiments.

Let us investigate the derived mathematical models. The corresponding surface of the pair interactions are shown in Fig. 10 and 11. They demonstrate good agreement between calculation results of MMPC and experimental data.

Analysis of the surface in Fig. 10 shows that the increase in single imperfections significantly reduces the carrying capacity of the shell. In addition, in these limits roundness has a greater impact on the setting than the taper. This is consistent with the single-factor experiments reported in [13,14]. But the analysis of Fig. 10 also shows that the simultaneous increase in taper and ovality can lead to an increase in carrying capacity to a level corresponding to the defect of a free shell. This is essentially a nonlinear effect, which could not be found by single-factor experiments.

Further study of the nonlinear interaction of defects (Fig. 11) showed that in the developed cone-of-roundness of the lower shell end has a more significant impact on the setting of critical effort than the roundness of the upper end. The joint increase in roundness of ends leads to an increase in carrying capacity, which is also an essentially nonlinear effect and is in good agreement with the results shown in Fig. 10.

Subcritical deformation has been studied in thin-walled shells with an oval on the lower and upper end being equal to 0.84 and 0.96, respectively, and taper equal to 0^o56' and 2^o16', respectively. The selected values$$\alpha$$, $$\left(a/b\right)$$, $${\left(a/b\right)}_{В}$$and $${\left(a/b\right)}_{Н}$$ correspond to characteristic points of the models [13,14].

A qualitative analysis of the effect of displacement fields on the results of the holographic experiment suggests an important role played by the strain state of shells under non-uniform roundness in the district and in the longitudinal direction (Fig. 12).

With the increase of $$\alpha$$ up to 2^o16' heterogeneity of the radial deflection is shifted to the lower end. A comparison of interferograms obtained at different load levels shows that an increase in the last number of fringes decreases with equal values of the additional load, and this indicates the hardening of structures, possibly caused by high deformability of the shell at the beginning of loading. The deformation of the shell in the experiment depends on the character of ends, and imperfections differ significantly on the panels of varying curvature.

Figure 10.

The surfaces of the pair interactions of imperfections $\alpha$ and $a/b:$ (two-factor model (a), comparison between model, single-factor experiment [13] (rhombus) and MMPC calculations (circles) (b).

Figure 11.

Surfaces of pair interactions with imperfections$${\left(a/b\right)}_{Н}^{\circ },{\left(a/b\right)}_{В}^{\circ }$$.

Figure 12.

Interferogram envelope with the taper and ovality of the small curvature of the panel (a), the joint panel zone (b), and the larger curvature of the panel (c).

Figure 13.

The deformation of the shell with ovality and taper. Solid curves correspond to the middle of panels: black – large curvature, gray – small curvature, dash – forming at the junction of the panels; positive direction goes toward the center of curvature

The change of $${\left(a/b\right)}_{Н}$$ from 0.96 to 0.84 significantly (1.2-1.4 fold) increases compliance of the membranes, while maintaining the overall picture of the distribution of displacements in the circumferential direction and increasing heterogeneity in the longitudinal direction.

The field of displacements was explained semi-automatically on the basis of the above algorithm. The forms of the radial deflection of some shell generatrixes are shown in Fig. 13.

8. Conclusions

A modified method of the parameter continuation (MMPC) is proposed. This method enables simplification of the calculations both at the stage of constructing the model, and also within its continued use due to precise values of the Taylor coefficients for the solution of the degree not exceeding the number of approximation.

The expression to calculate approximations by the MMPC in the general case and with the nonlinearity type of products and squares of the desired functions is presented.

The application of fractional-rational transformation for the polynomial approximation in the form of the 1-D and 2-D PAs used for increasing the degree of convergence and for the analytical continuation of the approximation in the region of its meromorphy was analyzed. It was concluded that such a transformation is justified if it is applied to polynomials which depend on the variable of integration. We used 2-D PAs for the independent variable and for the artificial parameter applying the scheme proposed by V. Vavilov. In this paper it is shown that this transformation provides a satisfactory quality for the approximation behavior and minimizes its error, in spite of the fact that the use of 2-D PAs requires a further theoretical justification.

The estimation of stability using MMPC approximation is also proposed. A study of numerical results was conducted by applying the methods for three model examples which were perturbed with a natural small parameter. It is shown that the application of PAs provides them with sufficient accuracy in the studied area. This paper shows the advantage of approximations which were obtained based on the MMPC.

Calculations of nonlinear deformation and stability of elastic flexible circular cylindrical shell under uniform external pressures and of free oscillations of simply supported stringer shell demonstrated the efficiency and accuracy of the proposed method.

The methodology and results of a holographic experiment with thin low-conical shells having oval ends are presented. They show good agreement with calculation results.

Acknowledgement

J. Awrejcewicz input to this chapter was supported by the Alexander von Humboldt Award.