，’，ZHOU Gugming

State Key Laboratory of Mechanics and Control of Mechanical Structures，Nanjing University of Aeronautics and Astronautics，Nanjing 210016，P.R.China

Abstract: The aim of the present study is to develop an efficient weak form quadrature element for free vibration analysis of arbitrarily shaped membranes. The arbitrarily shaped membrane is firstly mapped into a regular domain using blending functions，and the displacement in the element is assumed as the trigonometric functions. Explicit formulations are worked out for nodes of any type and a varying number of nodes. For verifications，results are compared with exact solutions and data obtained by other numerical methods. It is demonstrated that highly accurate frequencies can be obtained with a small number of nodes by present method.

Key words：arbitrarily shaped membrane；free vibration；quadrature element；blending function

0 Introduction

A membrane is characterized by a negligible re?sistance to bending and a dominating tension. Nu?merous membrane structures exist in practice and their applications are still growing in importance［1］.Therefore，free vibration of membranes has attract?ed the attention of many researchers［2-7］. A review on the study of the membrane vibration was given by Jenkins and Korde recently［1］.

Numerical modeling and analysis of structural elements with irregular shapes have continuously been a popular research topic because of the wide?spread applications of irregular shaped elements in various fields［8］. Besides，the enhancement of com?putational efficiency and accuracy has always been an interesting research topic to the computational mechanics community［9］. Several efficient numerical methods， such as the strong form differential quadrature element method（DQEM）［10］，the local radial basis function-based differential quadrature（LRBFDQ）［7］，the discrete singular convolution（DSC）algorithm［11］and the weak form quadrature element method（QEM）［12］，have been developed recently and are still under developing，since an effi?cient numerical method is always valuable to design?ers in certain engineering applications［13］.

Free vibration of membranes is the simplest two-dimensional vibration problem and thus it is of?ten used to test the performance of a numerical tech?nique，e.g.，the collocation method［2］，the p-version finite element method（p-FEM）［3］，the DSC algo?rithm［4］，the DQEM［5-6］，and LRBFDQ［7］. For arbi?trary shaped membranes，the irregular domain is usually mapped into a regular domain by using ei?ther blending functions［3，6］，or shape functions of a Serendipity element［4-5］. Care should be taken that，however，the mapping accuracy affects the solution accuracy［5-6，14］and that inaccurate mapping may even cause non-convergent results if the irregular shape does not have four corners［15-16］，since zero or nega?tive Jacobian determinant will inevitably occur at the corners［16］. Numerical difficulty may be encountered since the derivative with respect toxand/oryat the corner points does not exist. In such cases，the weak form methods may have some advantages over the strong form methods，since the numerical difficulty may be circumvented using Gauss quadra?ture［16］.

Previous research showed that the QEM was highly accurate［9，12，17-19］and possessed the potential to act as a competitive counterpart to other efficient numerical methods and thus was worth being devel?oped further［9］. An undisputable advantage of the QEM over the p-FEM is that the pre- and post-pro?cessing is much more convenient due to the physi?cal meaning of its DOFs［19］. Therefore，the objec?tive of present paper is to develop an efficient quadrature element for the free vibration analysis of arbitrarily shaped membranes. The arbitrarily shape is firstly mapped into a regular domain using blend?ing functions and the discretization is then done on the regular domain. Trigonometric functions are used as the element displacement to formulate a sub-parametric element. Explicit formulations are given for the element with nodes of any type and a varying number of nodes. Numerical examples are given. The obtained results are compared with ex?act solutions and data obtained by other numerical methods for verifications. Finally conclusions are drawn.

1 Weak Form Quadrature Element Formulations

1.1 Expressions of potential energy and kinetic energy

For the investigation of the free vibration be?havior of membranes by a weak form method，the expressions of potential energyUand kinetic energyTof the element are needed and given as

wherewis the transverse displacement and the over dot on it is the first order derivative with respect to timet；Sthe tension per unit length；ρthe mass per unit area；Athe area of the membrane element；and|J| the determinant of Jacobian matrix.ξandηare coordinates in the regular domain.

1.2 Sub?parametric quadrature membrane ele?ment

To develop a sub-parametric quadrature mem?brane element，the arbitrary shaped membrane ele?ment is firstly mapped into a regular domain shown in Fig.1，and then displacement is assumed as the function ofξandη(-1 ≤ξ，η≤1). To increase the mapping accuracy，blending functions are used，namely［3，6，16］

whereXi(s) andYi(s)(i=1，2，3，4；s=ξorη)are the parametric equations of the four edges of the membrane element.

To calculate the first order derivatives with re?spect toxandy，the chain rule of the partial differ?entiation is used

Fig.1 Sketch of a regular domain

whereJ-1is the inverse of Jacobian matrix.

Jacobian matrixJcan be easily computed by Eq.（3）as

The inverse of Jacobian matrixJ-1is given by

Substituting Eq.（6）into Eq.（4）gives

It is worth noting that ?w/?xand ?w/?ydo not exist numerically at a point(ξi，ηi)if|J(ξi，ηi)|=0.

LetNbe the number of nodes in eitherξorηdirection and (ξj，ηi)(i，j=1，2，…，N) be the coor?dinates of element nodes. The element displacement is assumed as

wherelj(ξ)andli(η)are shape functions.

It is reported that p-FEM with trigonometric functions are numerically more stable than orthogo?nal polynomials as the order is increased［3］. Since both the QEM and the p-FEM are weak form meth?ods，therefore，trigonometric functions are used as the shape functions of the QEM，shown as

wherep（≤π/4）is a control variable. Whenpis very small，Eq.（9）is equivalent to the polynomialbased shape functions.

The shape functions are new and used to develop a quadrature element for the first time，al?though Eq.（9） is widely used to compute the weighting coefficients in the strong form harmonic differential quadrature（HDQ）method［20］.

Substituting Eq.（8）into Eq.（1）and perform?ing the numerical integration by Gauss quadrature yield that

where superscriptsξorηmean that the correspond?ing first order derivative is taken with respect toξorη. More precisely，AξjkandAηilare the weighting co?efficients of the first order derivative with respect toξandη.

The weighting coefficientsAξjkcan be explicitly calculated by

can be calculated in a similar way. In this way，kcan be obtained explicitly for anyNand nodes of any type.

It is worth noting that Eq.（12）is reduced to the formulas of the weighting coefficient in the HDQ method［20］if the integration points are also the nodes. Iflj(ξ) andli(η) are Lagrange interpolation functions with polynomials which are commonly used in the QEM，a similar formula to Eq.（12）is also available to compute the weighting coeffi?cientAξjk［17-18］.

Substituting Eq.（8）into Eq.（2）and perform?ing the numerical integration by Gauss quadrature yield that

wheremis the mass matrix andw? the nodal veloci?ty vector.Lis shown as

It is seen that neither the type of nodes nor the number of nodes is fixed a priori in the derivations of the explicit formulas. Therefore，different types of nodes can be used and the solution accuracy can be easily adjusted by changing the number of the nodes in the developed program.

For free vibration analysis，w(ξ，η，t)=W(ξ，η)sinωtand thusw=Wsinωt，whereωis the circular frequency. After dropping the term sinωt，the equation of motion is given by

If more elements are used，the assemblage procedures are similar to the conventional finite ele?ment method（FEM）. After applying the essential boundary condition，the matrix equation is modi?fied as

Solving Eq.（16）by a generalized eigen-value solver yields the frequencies and mode shapes.

2 Numerical Results and Discus?sion

For demonstrations，the widely used Gauss-Lobatto-Legendre（GLL）points in QEM are adopt?ed as the element nodes. An explicit formula to cal?culate the GLL points does not exist and thus the program reported in Ref.［12］is used to computeξk，ηk(k=1，2，…，N). Mapping an irregular do?main without four corners into a regular one is more difficult than the one with four corners， since|J|=0 at the added corner point［16］. Therefore，ex?amples of irregular shaped membranes with a part of an elliptic shape（Fig.2（a）），full elliptic shape（Fig.2（b）），half circle+triangle（Fig.3），and gen?eral quadrilateral shape，are given. In Figs.2，3，symbols ①— ④ represent the four edges of the membrane element which are needed in Eq.（3）. Be?sides，only oneN×N-node element withp=π/8 is used in the analysis，sincep=π/8 is widely used in the HDQ method and different values ofpaffect on?ly the higher mode frequencies which are not studied in this paper. All edges of the membranes are fixed，i.e.，w=0 on all edges.

The non-dimensional frequency parameterΩis defined as

whereais the semi-major axis of the ellipse（Fig.2）or the radius of the circle（Fig.3）.

Fig.2 Sketches of arbitrary shaped membranes

Fig.3 Sketche of a half-circle+triangle membrane

Example 1First，consider the free vibration of a sectorial membrane shown in Fig.2（a）. To map it into a regular domain with four corners，as shown in Fig.1，one additional corner point is needed. Cur?rently the point，i.e.，point 3 shown in Fig.2（a），is located on the curved edge. Other ways to set the additional corner point exist and the effect on the fi?nal results is negligible if the mapping is accurate enough［16］.

The parametric equations of the four edges of the sectorial membrane，as shown in Fig.2（a），are

whereθequals to arctg(a*tan(φ)/b)；aandbare the semi-major and semi-minor axes of the elliptic curve.

The first 10 non-dimensional frequencies are listed in Table 1.The number of nodes in each direc?tion varies from 7 to 12. Available exact solutions and data obtained by conventional and p-version fi?nite element methods，denoted by c-FEM and p-FEM，are included for comparisons.

Table 1 Comparison of non?dimensional frequency parameter Ω for sectorial membrane (φ=π/2)

It is seen that the present element is highly ac?curate. Exact first three non-dimensional frequency parameters can be obtained by using one 9×9-node element. WhenNis 12，all 10 non-dimensional fre?quency parameters are numerically exact. The con?ventional FEM with fine meshes and the p-version FEM using trigonometric sine functions withp=12 yield slightly lower accurate solutions for higher modes.

The first eight non-dimensional frequencies for sectorial membranes with different sector angles and aspect ratios are listed in Table 2. To ensure the so?lution accuracy，one 14×14-node element is used.Obtained results agree well with the data using the conventional finite element with fine meshes. The QEM results listed in Table 2 should be highly accu?rate according to the observation from Table 1 and thus can serve as references.

Table 2 Non?dimensional frequency parameter Ω for various sectorial membranes (N=14)

Table 3 Comparison of non?dimensional frequency parameter Ω for a circular membrane

Example 2Second，consider the free vibra?tion of an elliptic membrane shown in Fig.2（b）. To map it into a regular domain with four corners，as shown in Fig.1，four corner points are needed. Cur?rently four corner points on the curved edge shown in Fig.2（b）are adopted. Other ways to set the four corner points are available and the final results are independent from the way to place the corner points［16］.

The parametric equations of the four edges of the elliptic membrane shown in Fig.2（b）are

whereaandbare the semi-major and semi-minor axes of the elliptic membrane.

The first ten non-dimensional frequencies are listed in Table 3.The number of nodes in each direc?tion varies from 7 to 14. Available exact solutions and data obtained by DQEM are included for com?parisons.

It is seen that the present element is highly accurate. one 14×14-node element（196 DOFs）can yield ten numerically exact frequencies. Nu?merically exact first ten frequencies are also ob?tained by the DQEM［6］employing the same geo?metric mapping with 17×17 grid points （289 DOFs）. Perhaps due to inaccurate mapping，the DQEM［5］can only yield approximate solutions even with 41×41 grid points（1 681 DOFs）. It is expected that the solutions obtained by lower or?der finite element（1 024 DOFs）［2］are not accu?rate enough. If a very fine mesh is used，the ac?curacy of the recalculated data by the lower order finite element can be improved greatly. The small?est number of DOFs demonstrates the effective?ness of the QEM.

The first eight non-dimensional frequencies for elliptic membranes with different aspect ratios are listed in Table 4. To ensure the solution accuracy，one 14×14-node element is used. Obtained results agree well with the DQEM（17×17）data using the same geometric mapping technique. The results listed in Table 4 should be highly accurate and thus can serve as references.

Table 4 Non?dimensional frequency parameter Ω for various elliptic membranes (N=14)

The accuracy of geometric mapping affects the solution accuracy greatly［6，15］. To demonstrate it，shape functions of G-node serendipity elementfi(ξ，η)are used for geometric mapping

The percentage relative error listed in Table 5 is defined as Err=(A-πab)/(πab)*100%，whereAis the mapped area. The general shape functions of the 12-node serendipity element with non-uniformly distributed nodes are given in Ref.［16］.

The fundamental frequencies of various elliptic membranes are listed in Table 5. Geometric map?ping is performed using Eq.（20）（G=8 or 12）.Available data obtained by the DSC withN=15［4］and the DQEM withN=41［5］are included for com?parisons.

It is seen that both the weak form QEM and the strong form DQEM yield the same convergent results if the same mapping technique is used. The mapping error obviously affects the solution accura?cy obtained by various numerical methods. Employ?ing non-uniformly distributed nodes such as the GLL nodes does improve the mapping accuracy and thus improves the accuracy of solutions. Therefore，accurate mapping technique should be always used to obtain highly accurate results.

Example 3Further consider the free vibration of a half-circle+triangle membrane shown in Fig.3.To map it into a regular domain with four corners shown in Fig.1，one more corner point is needed.Two ways of adding the corner point are shown in Fig.3. If the curved edge becomes straight，it reduc?es to a triangular element.

The parametric equations of the four edges ofthe irregular shaped membrane shown in Fig.3（a）and Fig.3（b）are

whereais the radius of the half circle.

The first eight non-dimensional frequencies ob?tained using Eq.（21）are listed in Table 6. The number of nodes in each direction varies from 7 to 14. Existing solutions obtained by other numerical methods are included for comparisons. Results with Eq.（22）are similar to the ones with Eq.（21）and thus only the ones withN=14 are included for com?parisons.

Since exact solutions are not available，the highly accurate DQEM data cited from Ref.［6］（289 DOFs）are used for comparisons. It is seen that accurate non-dimensional frequency parameters can be obtained by using one 13×13-node（169 DOFs）element with Eq.（21）. The accuracy of the QEM with Eq.（22）is slightly lower. Due to the mapping error，data obtained by the DQEM withN=41 （1 681 DOFs）［5］are still not accurate enough. Some frequencies obtained by the LRBFDQ with 1 675 DOFs［7］are not accurate enough. The data obtained by the lower order FEM with 1 089 DOFs［2］are obviously inaccurate. The accuracy of the FEM data with a fine mesh（5 221 DOFs）is much improved. Based on the total num?ber of DOFs，one may conclude that the QEM is ef?ficient and suitable for analysis of arbitrarily shaped membranes. Perhaps due to the fact of that the curved edge is much longer than the other two edg?es，the solution accuracy of the QEM with Eq.（22）is slightly lower than the one with Eq.（21）whenNis small. The two ways to map the irregular shape without four corners into a regular one demonstrate the flexibility of the sub-parametric element formula?tions.

Example 4Last，consider the free vibration of a general quadrilateral membrane. The Cartesian co?ordinates of four corners are（0.0，0.0），（3.0，0.0），（2.5，2.5）and（0.5，2.0）for comparing with exist?ing data.

The irregular domain can be exactly mapped in?to a regular one by Eq.（20）withG=4，where the shape functions arefi(ξ，η)=(1+ξiξ)(1+ηiη)/4（i=1，2，3，4）. For comparisons with the existing data，a slightly different frequency parameterΩ*from Eq.（17）is used，shown as

Table 6 Comparison of non?dimensional frequency parameter Ω for the half?circle+triangle membrane

The first eight frequencies parameters are list?ed in Table 7. The results are exactly the same as the DQEM ones. Since the geometric mapping is ex?act，the data in Ref.［5］should be very accurate. As expected，the conventional FEM can also yield ac?curate results since all edges are straight.

Table 7 Comparison of the frequency parameter Ω* for the general quadrilateral membrane

3 Conclusions

Free vibration of arbitrarily shaped membranes is analyzed by using the weak form quadrature ele?ment method. A novel sub-parametric element is de?veloped. The element domain is accurately mapped into a regular one using blending functions. Trigono?metric functions are used as the elemental displace?ment for the first time. To ease the programming，explicit formulas of a varying number of nodes are given. For verifications，results are compared with exact solutions and data obtained by other numerical methods. The solution accuracy can be improved by increasing the number of nodes easily. It is shown that highly accurate frequencies can be obtained by the proposed element with a small number of nodes.Accurate mapping plays an important role for suc?cessfully using both strong and weak form methods.

AcknowledgementsThis paper is partially supported by the National Natural Science Foundation of China（Nos.52005256，12072154），the Natural Science Foundation of Jiangsu Province（No.BK20190394），the Jiangsu Post-Doc?toral Research Funding Program（No.2020Z437）and the Priority Academic Program Development of Jiangsu Higher Education Institutions. The present investigators are grateful for helpful discussions with Prof.YUAN Zhangxian of Worcester Polytechnic Institute，USA.

AuthorProf. WANG Xinwei graduated in airplane design from Nanjing Aeronautical Institute，China，in 1975，and re?ceived his M.S. degree in solid mechanics from the same in?stitute，in 1981 and Ph.D. degree in solid mechanics from School of Aerospace and Mechanical Engineering，Universi?ty of Oklahoma，USA，in 1989. He was a professor in struc?tural mechanics at the State Key Laboratory of Mechanics and Control of Mechanical Structures，Nanjing University of Aeronautics and Astronautics，China. His research interests include finite element method and its applications，efficient computational methods and their applications，and mechani?cal behavior of composite materials.

Author contributionsProf. WANG Xinwei designed the study and wrote the manuscript. Dr. CAI Deng’an conduct?ed the analysis and interpreted the results. Prof. ZHOU Guangming contributed to the discussion and background of the study. All authors commented on the manuscript draft and approved the submission.

Competing interestsThe authors declare no competing interests.