# Multibody Modeling Method for UHV Porcelain Arresters Equipped with Lead Alloy Isolation Device

Authors: Xiaochao Su; Lei Hou; Zhubing Zhu; Yushu Chen
Journal: Shock and Vibration
 (2021)

---

## Abstract

This paper presents a multibody modeling method for seismic analysis of UHV porcelain surge arresters equipped with a kind of seismic isolation device. An UHV arrester is modeled as a planar multibody system whose number of DOF is equal to the number of the arrester units. Joint coordinate method is adopted to construct the governing equations of motion. The seismic isolation device utilizing a number of lead alloy dampers as its core energy dissipation components is also investigated. An analytical model of this device is given by modeling each lead alloy damper as a hysteretic spring and reducing the entire device to a planar system consisting of a range of hysteretic springs. Its mechanical characteristic is derived theoretically, and the obtained moment-angle relationship is expressed as a system of differential algebraic equations. The initial rotational stiffness of the device is formulated in terms of the structural and mechanical parameters of the device. This analytic expression is used in estimating the fundamental frequency of the isolated equipment. By this modeling method, it is easy to construct the governing equations of motion for the isolated system. An UHV arrester specimen is analyzed by this proposed method. The effectiveness of the isolation device in terms of reducing the internal base moment is significant and the influence of system parameters on the effectiveness is also discussed. The proposed method shows its potential usefulness in optimal design of the isolation device.


---

## Body

## 1. Introduction

Surge arrester is one of the core components of an electrical substation. Nowadays, porcelain-housed metal oxide varistor (MOV) surge arresters are commonly equipped in substations. A typical ultra-high voltage (UHV) porcelain-housed MOV surge arrester is usually composed of several units and one grading ring, which are connected to each other in series with flanges. The structural schematic diagram of an UHV porcelain-housed surge arrester consisting of four units is depicted in Figure1. Due to their structural characteristics and material property, porcelain surge arresters are easy to be damaged during earthquakes. Thus, improving the seismic performance of the surge arresters is one of the keys to increase the seismic resistance of the entire substation. This issue is faced especially by those UHV substations, because UHV porcelain arresters are almost 12 meters in height and 10 tons in weight and usually have fundamental frequencies ranging from 1 to 10 Hz that match the predominant frequencies of the most of earthquake ground motion [1].Figure 1
Structural schematic of an UHV porcelain arrester.In addition to porcelain surge arresters, many other porcelain post electrical equipment, such as post insulators and disconnect switches, also have high seismic vulnerability. Adopting seismic protective devices is one of the ways to improve their seismic performance. In recent research studies, the use of wire rope in base isolation for circuit breakers has been investigated in [2, 3]. The wire rope isolator was modeled by a modified Bouc–Wen model, and the isolated equipment was analyzed by finite element (FE) method. An innovative multiple tuned mass damper (MTMD) seismic mitigation strategy for porcelain post electrical equipment was investigated in [4, 5]. In [6] Xie et al. studied a new kind of assembled isolation device combining wire rope isolators and linear viscous dampers for bypass switch. The isolation device was modeled by a simplified system of 2-DOF, and its effectiveness was validated numerically by using FE analysis and experimentally through shaking table test. This kind of isolation device is also studied by Puff et al. in [7] for circuit breaker. A restraint system consisting of rods, springs, and viscous dampers was investigated in [8] for seismic protection of an 800 kV thyristor valve. A 4-DOF simplified model was developed to investigate the optimal parameters for the restraint system. In [9], friction springs damper used to retrofit tall high voltage equipment for seismic performance improving was discussed. The use of this kind of damping device on seismic protection for electrical equipment can also be found in [10, 11]. A novel seismic isolation device utilizing lead alloy dampers (LADs) as its core energy dissipation components has been studied by Cheng et al. in [12]. The isolation device is equipped between the base unit and the support structure. It is composed primarily of two steel connection plates, a cushion block (CB), and a circle of uniformly distributed LADs. The upside steel plate (USP) is connected to the base unit, and the downside steel plate (DSP) is mounted to the support. The schematic diagram of the assembly is depicted in Figure 2. When an isolated arrester undergoes rocking motion during seismic events, parts of the LADs are subjected to tension, and the others are subjected to compression. Seismic energy is dissipated through the LADs undergoing cyclic yielding deformation, so that the energy transmitted to the above structure is reduced. In [12], the LAD was modeled by the classical Bouc–Wen model, and the isolated electrical equipment was analyzed by using FE method.Figure 2
Schematic diagram of the seismic isolation device: (a) lateral view, (b) plan view.
(a)(b)It is common to adopt FE modeling method in seismic analysis for porcelain post electrical equipment, such as in [13–15], while simplified modeling methods are also adopted in research on electrical equipment in substations. In [16, 17], generalized single degree of freedom system was used to model a post type of electrical equipment. A 4-DOF system was introduced in [18] to model support-equipment structures in substations and the influence of structure parameters on the dynamic properties of the system was analyzed by using the proposed simplified model. In optimal design of an isolation device, an equivalent model having small number of DOF was proposed in [19] for isolated UHV bypass switches, because FE method is computational time consuming.For porcelain postelectrical equipment, the current Chinese Code [20] in China recommends a simplified mass-spring modeling method with special care on the stiffness property of the flange joint. An arrester unit in Figure 1 mainly consists of one porcelain housing and core MOV, and each end of the porcelain housing is cemented to joint with metal flange, which is schematically shown in Figure 3. Due to the fact that the cement is generally mechanically weaker than the porcelain housing, an arrester has unevenly distributed stiffness in the axial direction. The third of the authors Zhu has verified in [21] that the simple cantilever beam model is unsuitable in modeling UHV porcelain arresters. Thus, more attention should be paid to the cement joints. An empirical formula to calculate the rotational stiffness at the cement joint was introduced in the Japanese Standard [22] and also the Chinese Code [20], which is expressed as follows:(1)Kc=ξdchc2te,where Kc is the rotational stiffness of the porcelain arrester unit at the cement joint, dc represents the outer diameter of the porcelain housing at the cement joint, hc and te denote the height and the thickness of cement, respectively, and ξ is considered as the coefficient to compute the rotational stiffness. The value of ξ is given by the Chinese Code [20] as 6.54×107 with dc, hc, and te having the unit of meter (m) and Kc having the unit of N⋅m⋅rad−1. In recent researches [23, 24] on the coefficient ξ for UHV porcelain electrical equipment, a new formulation is given as follows:(2)ξ=107×6.54,dc≤0.275,−15.4dc+10.775,0.275<dc<0.375,5.00,dc≥0.375.Figure 3
Structural schematic of one porcelain arrester unit.Taking account of the rotation stiffness at cement joint, the Chinese Code [20] recommends one can model the flange joint as an equivalent beam element in FE analysis or as a rotational spring in mass-spring modeling. When following the recommendations, the UHV porcelain arrester shown in Figure 1 would be modeled as a nonuniform beam or a mass-spring system as depicted in Figures 4(a) and 4(b), respectively. An analytical model given in [25] modeling each unit as one beam element connected to each other by two rotational springs is shown in Figure 4(c).Figure 4
Three kinds of models of an UHV arrester: (a) nonuniform beam model, (b) mass-spring model, (c) analytical model in [25].
(a)(b)(c)The aim of this paper is to present a new modeling method for UHV porcelain surge arresters equipped with the novel seismic isolation device investigated in [12]. The contributions of this study are summarized as follows. Firstly, an UHV porcelain arrester is modeled as a planar multibody system, whose number of DOF is equal to the number of the porcelain units which is usually less than 6. The governing equations of motion are constructed through joint coordinate method. Secondly, an analytical model of the seismic isolation device is given and its mechanical characteristics are derived theoretically. On characterizing the moment-angle relationship of the device, the work in this paper goes a little further than that in [12]. A new parameter indicating the relation between ground motion and LADs distribution is introduced, which was not considered in the previous research work. An analytic expression for the initial rotational stiffness of the isolation device is obtained, which is formulated in terms of the structural and mechanical parameters of the device. Finally, a numerical analysis on a specimen of arrester shows the application of the proposed method in analyzing the effectiveness of the isolation device and the influence of the system parameters on it.


## 2. Multibody Model of an UHV Porcelain Arrester

Because the distributed stiffness decreases at cement joints and the structure of an arrester is planar symmetric, we model each arrester unit as a planar rigid rod and characterize its resistance to the deflection by two torsional springs separately attached to the ends of the rod. Practically, it is logical to assume uniform mass distribution for the rod, so that the rod’s mass center is located at its geometric center. Figure5(a) shows the rigid body model of an individual arrester unit. The flange joint between arrester units is modeled as a pin joint. The torsional springs at the linked ends are considered to be connected in series and can be represented by a single equivalent torsional spring located at the pin joint. The grading ring is modeled as a lumped mass fixed to the top flange of the uppermost unit. Figure 5(b) shows the multibody model of the UHV porcelain arrester depicted in Figure 1. Damping elements can be added to the pin joints when damping effect is considered. So in this paper, an UHV porcelain arrester is modeled as a planar multibody system, more specifically, a planar open-chain system containing only a single one-branch tree.Figure 5
Diagram of the multibody model of an UHV porcelain arrester: (a) rigid body model of one unit; (b) multibody model of the entire arrester.
(a)(b)


### 2.1. Constructing the Governing Equations of Motion

With regard to an UHV porcelain arrester consisting ofn units and one grading ring, a (n+1)-body system can be constructed by using the above modeling method. For purpose of simplifying analysis, assumptions that linear torsional springs and linear viscous dampers are used at the pin joints are adopted in the following sections. The input ground motion is also reduced to a horizontal ground motion in the following analysis, but it is easy to extend to include vertical ground motion.To construct the equations of motion for the (n+1)-body system, the joint coordinate method in [26] is mainly used. The main processes are summarized as follows:(1)
Index the bodies of the system(2)
Construct a global nonmoving Cartesian frame of reference(3)
Define a set of joint coordinatesθ equaling in number to the number of system’s DOF(4)
Express the absolute body coordinatesc in terms of the joint coordinates θ(5)
Derive the velocity transformation matrixB, i.e., the Jacobian of c with respect to θ(6)
Construct the mass matrixM˜ and applied force array h of the multibody systemOnce the matricesM˜ and B and vectors h and c have been constructed, the equations of motion for the system can be derived by substituting them into the following equation:(3)BTM˜c¨−h=0.Detailed processes for constructing the equations are presented as follows. Assign numbers from 1 ton to the rods from the base to the top one, and assign index “top” to the lumped mass. Denote by joint (i) the one linking body (i−1) and body (i). Choose a global nonmoving Cartesian frame of reference at initial rest state, whose origin O is located at the root joint, x-axis is along the direction of horizontal ground motion, and y-axis is along the axial direction of the arrester. Denote by uGt the horizontal ground displacement history. Denote by θi the absolute rotation angle of body (i). At some instant of time, the configuration of the multibody system and the reference frame is depicted in Figure 6.Figure 6
Reference frame and configuration of the multibody system at some instant of time.Obviously, the system hasn-DOF. Denote by θ=θ1,…,θnT the joint coordinate of the system. The global coordinate of body (i), denoted by ci=xi,yi,ϕiT, can be easily expressed in terms of θ as follows:(4)xi=uGt−l1sinθ1−⋯−li−1sinθi−1−li2sinθi,yi=l1cosθ1+⋯+li−1cosθi−1+li2cosθi,ϕi=θi.where i=1,…,n and li is the length of the rod (i), while for the top body(5)xtop=uGt−l1sinθ1−⋯−lnsinθn,ytop=l1cosθ1+⋯+lncosθn,ϕtop=θn.Then, the configurationc=c1T,…,cnT,ctopTT of the system can be expressed in terms of the joint coordinate θ, and the velocity transform matrix B and acceleration vector c¨ can be derived accordingly.To construct the mass matrixM˜ of the system, denote by mi the mass of rod (i) and mtop the mass of the top lumped mass. Then the moment of inertia for rod (i) about the mass center in the plane can be represented as Ji=mili2/12. For the lumped mass, moment of inertia is ignored, represented by Jtop=0. Then, the mass matrix M˜ is given as a block diagonal matrix of the matrices M˜1,…,M˜nandM˜top in which M˜i=mi000mi000Ji is a 3×3 diagonal matrix for i=1,…,n,top.To construct the array of applied forces for the system, denote byhi=fix,fiy,ηiT the applied force array on body (i), in which fix and fiy represent, respectively, the horizontal and vertical resultant force and ηi represents the resultant moment. Benefitting from adopting joint coordinate, there is no need to consider the reaction forces between linked bodies in this dynamic equations construction. The forces and moments applied on body (i) resulting from gravity and springs are depicted in Figure 7. In horizontal direction, there is no applied force acting on mass center so that fix=0. In vertical direction, only gravitational force oriented in negative y-axis is applied through mass center, which can be expressed as fiy=−mig. Both springs located at joint (i) and (i+1) contribute moments to body (i). Define Δθi=θi−θi−1 as the relative rotational angle of body (i) to body (i−1), and denote by si the stiffness of the linear torsional spring located at joint (i). Define Mis=siΔθi, and then the moment associated with springs acting on the body (i) as shown in Figure 7 can be expressed as ΔMis=Mi+1s−Mis. In a similar way, when considering the moments due to the linear dampers, denote by di the damping coefficient of the damper located at joint (i), and define Δθ˙i=θ˙i−θ˙i−1 as the relative rotational velocity of body (i) to body (i−1), and Mid=diΔθ˙i. Then, the moment associated with rotational dampers acting on body (i) can be represented by ΔMid=Mi+1d−Mid.Figure 7
Forces applied on body (i).To sum up, for body (i) the force array is hi=0,−mig,ΔMis+ΔMidT, and for the top lumped mass the force array is htop=0,−mtopg,0T. Then, the applied force array for the system is derived as h=h1T,…,hnT,htopTT.The equations of motion for the system are derived by substituting the obtainedc, B, M˜, and h into equation (3). The result can be written in a compact form by using matrix and vector notation as follows:(6)Mθθ¨+Dθ˙+Sθ+Fθ,θ˙=u¨Gtmιθ.It is worth noting that the above ordinary differential equations are nonlinear. The matrixMθ is not constant. It includes terms in form of cosθi−θj that are nonlinear functions of the time dependent variable θi. The vector Fθ,θ˙ represents an array of nonlinear forces, and it contains terms in form of θ˙i2 and cosθi−θj. It also contains nonlinear terms of form sinθi because the effect of gravity has been considered in the construction of the equations of motion. On the right-hand side, the vector mιθ includes cosθi−θj form of terms.From a practical point of view, the angle of rotationθi is so small that small angle approximation can be used to linearize equation (6). The linearized equations can be expressed as follows by adopting the notation of m^i=mi+1+⋯mn+mtop(7)Mθ¨+Dθ˙+S−Gθ=u¨Gtmι,where M is the general mass matrix having the following form:(8)M=m13+m^1l12m22+m^2l1l2…mn2+m^nl1lnm22+m^2l2l1m23+m^2l22…mn2+m^nl2ln⋮⋮⋱⋮mn2+m^nlnl1mn2+m^nlnl2…mn3+m^nln2, representing the general damping matrix having form as follows:(9)D=d1+d2−d2…00−d2d2+d3…00⋮⋮⋱⋮⋮00…dn−1+dn−dn00…−dndn,and S−G can be considered as the general stiffness matrix with S and G having the following forms, respectively:(10)S=s1+s2−s2…00−s2s2+s3…00⋮⋮⋱⋮⋮00…sn−1+sn−sn00…−snsn,G=gm12+m^1l10…00gm22+m^2l2…0⋮⋮⋱⋮00…gmn2+m^nln.The vectormι on the right-hand side of equation (7) has form as follows:(11)mι=m12+m^1l1m22+m^2l2⋮mn2+m^nln.From the perspective of vibration mechanics, equation (7) represents a linear vibrating system that has a fully coupling mass matrix and excitation force induced by ground acceleration. At the present stage, an UHV porcelain arrester subjected to earthquake induced ground motion is modeled as a multibody system, and its governing equations of motion have been derived through joint coordinate method.


### 2.2. Analysis with the Equations of Motion

After the linear equation (7) has been established, the natural frequencies of the system can be derived by solving the generalized eigenvalue problem of matrix S−G and M.Seismic response history analysis can be taken by numerically integrating equation (7) in time-stepping methods. The methods commonly used include central difference method and Newmark’s method [27]. Of greatest interests in the investigation of arrester’s seismic performance are the responses of the internal force and the relative displacement and absolute acceleration at the top end. As soon as the response history of θ, θ˙, and θ¨ is evaluated, the top horizontal displacement relative to the base can be expressed by applying small angle approximation as(12)x˜top=−∑i=1nliθi.Similarly, the total absolute acceleration response at the top end can be formulated as the following linear approximation:(13)x¨top=u¨G−∑i=1nliθ¨i.The response of the internal moment at joint (i) can be expressed as follows:(14)Mi=−Mis=−siθi−θi−1.Static analysis can also be carried out to the system. In order to construct the equilibrium equations for the system subjected to external static load, the applied forces arrayh should be modified by removing the damping terms and adding external forces at associated position. Then substitute it into equation (3) and drop the inertia and damping terms; the static equilibrium equations can be derived as the following expression:(15)BTh=0.The above equation is a system of nonlinear algebraic equations with respect to the joint coordinateθ. Newton–Raphson method can be used to compute the value of θ. Then the other concerned response quantities can be evaluated. For example, the relationship between the lateral force and displacement at the top end can be investigated analytically.


### 2.3. Parameters Specification

In a practical analysis, specific values of the parameters in the system equation (7) should be given. The mass and length of the arrester units can be obtained from the manufacturer. However, the stiffness at cement joint is usually unknown. As described at the start of this section, the torsional spring located at pin joint is considered as two associated torsional springs connected in series. Denote the rotational stiffness of the spring at the bottom end of rod (i) as si− and that at the top end of rod (i−1) as si−1+, and then the equivalent rotational stiffness si can be evaluated as follows:(16)si=si−si−1+si−+si−1+.The values ofsi− and si−1+ can be calculated through equations (1) and (2), in which the cementing parameters can be obtained from manufacturer. Thus, the necessary stiffness parameters of the model can be obtained through the structural parameters of the equipment.When performing numerical calculation in response history analysis, the damping matrixD in the system equation (7) needs to be specified. For lightly damped system that has similar damping mechanisms distributed throughout the structure, it is appropriate to assume that the system has classical damping [27]. Rayleigh damping assumption is usually adopted to construct the proportional damping matrix. Accordingly, the damping matrix D can be constructed from the following free vibration system induced by equation (7) without gravity effect:(17)Mθ¨+Dθ˙+Sθ=0.In this paper, stiffness proportional dampingD=bS is used to specify the damping matrix. The coefficient b is determined by the first natural circular frequency ω and the first-order modal damping ratio ζ of the system (equation (17)).


## 2.1. Constructing the Governing Equations of Motion

With regard to an UHV porcelain arrester consisting ofn units and one grading ring, a (n+1)-body system can be constructed by using the above modeling method. For purpose of simplifying analysis, assumptions that linear torsional springs and linear viscous dampers are used at the pin joints are adopted in the following sections. The input ground motion is also reduced to a horizontal ground motion in the following analysis, but it is easy to extend to include vertical ground motion.To construct the equations of motion for the (n+1)-body system, the joint coordinate method in [26] is mainly used. The main processes are summarized as follows:(1)
Index the bodies of the system(2)
Construct a global nonmoving Cartesian frame of reference(3)
Define a set of joint coordinatesθ equaling in number to the number of system’s DOF(4)
Express the absolute body coordinatesc in terms of the joint coordinates θ(5)
Derive the velocity transformation matrixB, i.e., the Jacobian of c with respect to θ(6)
Construct the mass matrixM˜ and applied force array h of the multibody systemOnce the matricesM˜ and B and vectors h and c have been constructed, the equations of motion for the system can be derived by substituting them into the following equation:(3)BTM˜c¨−h=0.Detailed processes for constructing the equations are presented as follows. Assign numbers from 1 ton to the rods from the base to the top one, and assign index “top” to the lumped mass. Denote by joint (i) the one linking body (i−1) and body (i). Choose a global nonmoving Cartesian frame of reference at initial rest state, whose origin O is located at the root joint, x-axis is along the direction of horizontal ground motion, and y-axis is along the axial direction of the arrester. Denote by uGt the horizontal ground displacement history. Denote by θi the absolute rotation angle of body (i). At some instant of time, the configuration of the multibody system and the reference frame is depicted in Figure 6.Figure 6
Reference frame and configuration of the multibody system at some instant of time.Obviously, the system hasn-DOF. Denote by θ=θ1,…,θnT the joint coordinate of the system. The global coordinate of body (i), denoted by ci=xi,yi,ϕiT, can be easily expressed in terms of θ as follows:(4)xi=uGt−l1sinθ1−⋯−li−1sinθi−1−li2sinθi,yi=l1cosθ1+⋯+li−1cosθi−1+li2cosθi,ϕi=θi.where i=1,…,n and li is the length of the rod (i), while for the top body(5)xtop=uGt−l1sinθ1−⋯−lnsinθn,ytop=l1cosθ1+⋯+lncosθn,ϕtop=θn.Then, the configurationc=c1T,…,cnT,ctopTT of the system can be expressed in terms of the joint coordinate θ, and the velocity transform matrix B and acceleration vector c¨ can be derived accordingly.To construct the mass matrixM˜ of the system, denote by mi the mass of rod (i) and mtop the mass of the top lumped mass. Then the moment of inertia for rod (i) about the mass center in the plane can be represented as Ji=mili2/12. For the lumped mass, moment of inertia is ignored, represented by Jtop=0. Then, the mass matrix M˜ is given as a block diagonal matrix of the matrices M˜1,…,M˜nandM˜top in which M˜i=mi000mi000Ji is a 3×3 diagonal matrix for i=1,…,n,top.To construct the array of applied forces for the system, denote byhi=fix,fiy,ηiT the applied force array on body (i), in which fix and fiy represent, respectively, the horizontal and vertical resultant force and ηi represents the resultant moment. Benefitting from adopting joint coordinate, there is no need to consider the reaction forces between linked bodies in this dynamic equations construction. The forces and moments applied on body (i) resulting from gravity and springs are depicted in Figure 7. In horizontal direction, there is no applied force acting on mass center so that fix=0. In vertical direction, only gravitational force oriented in negative y-axis is applied through mass center, which can be expressed as fiy=−mig. Both springs located at joint (i) and (i+1) contribute moments to body (i). Define Δθi=θi−θi−1 as the relative rotational angle of body (i) to body (i−1), and denote by si the stiffness of the linear torsional spring located at joint (i). Define Mis=siΔθi, and then the moment associated with springs acting on the body (i) as shown in Figure 7 can be expressed as ΔMis=Mi+1s−Mis. In a similar way, when considering the moments due to the linear dampers, denote by di the damping coefficient of the damper located at joint (i), and define Δθ˙i=θ˙i−θ˙i−1 as the relative rotational velocity of body (i) to body (i−1), and Mid=diΔθ˙i. Then, the moment associated with rotational dampers acting on body (i) can be represented by ΔMid=Mi+1d−Mid.Figure 7
Forces applied on body (i).To sum up, for body (i) the force array is hi=0,−mig,ΔMis+ΔMidT, and for the top lumped mass the force array is htop=0,−mtopg,0T. Then, the applied force array for the system is derived as h=h1T,…,hnT,htopTT.The equations of motion for the system are derived by substituting the obtainedc, B, M˜, and h into equation (3). The result can be written in a compact form by using matrix and vector notation as follows:(6)Mθθ¨+Dθ˙+Sθ+Fθ,θ˙=u¨Gtmιθ.It is worth noting that the above ordinary differential equations are nonlinear. The matrixMθ is not constant. It includes terms in form of cosθi−θj that are nonlinear functions of the time dependent variable θi. The vector Fθ,θ˙ represents an array of nonlinear forces, and it contains terms in form of θ˙i2 and cosθi−θj. It also contains nonlinear terms of form sinθi because the effect of gravity has been considered in the construction of the equations of motion. On the right-hand side, the vector mιθ includes cosθi−θj form of terms.From a practical point of view, the angle of rotationθi is so small that small angle approximation can be used to linearize equation (6). The linearized equations can be expressed as follows by adopting the notation of m^i=mi+1+⋯mn+mtop(7)Mθ¨+Dθ˙+S−Gθ=u¨Gtmι,where M is the general mass matrix having the following form:(8)M=m13+m^1l12m22+m^2l1l2…mn2+m^nl1lnm22+m^2l2l1m23+m^2l22…mn2+m^nl2ln⋮⋮⋱⋮mn2+m^nlnl1mn2+m^nlnl2…mn3+m^nln2, representing the general damping matrix having form as follows:(9)D=d1+d2−d2…00−d2d2+d3…00⋮⋮⋱⋮⋮00…dn−1+dn−dn00…−dndn,and S−G can be considered as the general stiffness matrix with S and G having the following forms, respectively:(10)S=s1+s2−s2…00−s2s2+s3…00⋮⋮⋱⋮⋮00…sn−1+sn−sn00…−snsn,G=gm12+m^1l10…00gm22+m^2l2…0⋮⋮⋱⋮00…gmn2+m^nln.The vectormι on the right-hand side of equation (7) has form as follows:(11)mι=m12+m^1l1m22+m^2l2⋮mn2+m^nln.From the perspective of vibration mechanics, equation (7) represents a linear vibrating system that has a fully coupling mass matrix and excitation force induced by ground acceleration. At the present stage, an UHV porcelain arrester subjected to earthquake induced ground motion is modeled as a multibody system, and its governing equations of motion have been derived through joint coordinate method.


## 2.2. Analysis with the Equations of Motion

After the linear equation (7) has been established, the natural frequencies of the system can be derived by solving the generalized eigenvalue problem of matrix S−G and M.Seismic response history analysis can be taken by numerically integrating equation (7) in time-stepping methods. The methods commonly used include central difference method and Newmark’s method [27]. Of greatest interests in the investigation of arrester’s seismic performance are the responses of the internal force and the relative displacement and absolute acceleration at the top end. As soon as the response history of θ, θ˙, and θ¨ is evaluated, the top horizontal displacement relative to the base can be expressed by applying small angle approximation as(12)x˜top=−∑i=1nliθi.Similarly, the total absolute acceleration response at the top end can be formulated as the following linear approximation:(13)x¨top=u¨G−∑i=1nliθ¨i.The response of the internal moment at joint (i) can be expressed as follows:(14)Mi=−Mis=−siθi−θi−1.Static analysis can also be carried out to the system. In order to construct the equilibrium equations for the system subjected to external static load, the applied forces arrayh should be modified by removing the damping terms and adding external forces at associated position. Then substitute it into equation (3) and drop the inertia and damping terms; the static equilibrium equations can be derived as the following expression:(15)BTh=0.The above equation is a system of nonlinear algebraic equations with respect to the joint coordinateθ. Newton–Raphson method can be used to compute the value of θ. Then the other concerned response quantities can be evaluated. For example, the relationship between the lateral force and displacement at the top end can be investigated analytically.


## 2.3. Parameters Specification

In a practical analysis, specific values of the parameters in the system equation (7) should be given. The mass and length of the arrester units can be obtained from the manufacturer. However, the stiffness at cement joint is usually unknown. As described at the start of this section, the torsional spring located at pin joint is considered as two associated torsional springs connected in series. Denote the rotational stiffness of the spring at the bottom end of rod (i) as si− and that at the top end of rod (i−1) as si−1+, and then the equivalent rotational stiffness si can be evaluated as follows:(16)si=si−si−1+si−+si−1+.The values ofsi− and si−1+ can be calculated through equations (1) and (2), in which the cementing parameters can be obtained from manufacturer. Thus, the necessary stiffness parameters of the model can be obtained through the structural parameters of the equipment.When performing numerical calculation in response history analysis, the damping matrixD in the system equation (7) needs to be specified. For lightly damped system that has similar damping mechanisms distributed throughout the structure, it is appropriate to assume that the system has classical damping [27]. Rayleigh damping assumption is usually adopted to construct the proportional damping matrix. Accordingly, the damping matrix D can be constructed from the following free vibration system induced by equation (7) without gravity effect:(17)Mθ¨+Dθ˙+Sθ=0.In this paper, stiffness proportional dampingD=bS is used to specify the damping matrix. The coefficient b is determined by the first natural circular frequency ω and the first-order modal damping ratio ζ of the system (equation (17)).


## 3. Model of an UHV Arrester Equipped with Seismic Isolation Device

### 3.1. Modeling the Isolation Device

In this section, a model of the seismic isolation device introduced in [12] which is depicted in Figure 2 is presented. Denote by nd the number of the LADs and by r the radius of the LAD circle. Construct two axes as shown in Figure 8(a), such that the x-axis coincides with the direction of the horizontal ground motion. Starting from the positive x-axis, number the LADs successively in a counterclockwise direction. Define the initial angle α as the angle between x-axis and the first LAD. Then, the location of LADi can be determined by the following angular expression:(18)αi=α+i−12πnd.Figure 8
Model of the isolation device: (a) effective arm of force, (b) planar isolation system.
(a)(b)Under the assumption of connection steel plates being rigid, the effective arm of force due to LADi about x′-axis can be expressed as ai=rcosαi. Reduce the whole isolation device to a planar system by means of regarding the connection plates as rigid beams, treating the cushion block as a hinged support fixed to the downside beam and pinned to the midpoint of the upside beam, and modeling each LAD as a hysteretic spring placed at its projection on x-axis that connects the two beams. The schematic diagram of the planar isolation system is depicted in Figure 8(b).Denote byM0 the external bending moment applied to the upside beam and θ0 the associated rotational angle responding to the moment, which is shown in Figure 9. The moment-angle relationship of the planar isolation system is derived theoretically in the following.Figure 9
Moment and associated angle of the planar isolation system.Denote the deformation of LADi caused by upside beam rotation as xi. Then it can be expressed as a multiplication of the rotational angle and the effective arm of force associated with LADi as follows:(19)xi=θ0rcosαi.The sign ofxi indicates the deformation state of LADi in a manner that positive means it is in tension, while negative means it is in compression. Denote by fi the applied force that induces the deformation xi of LADi. Then the moment distributed on LADi is fircosαi, and the external moment M0 is the summation of those distributions as follows:(20)M0=∑i=1ndfircosαi.Once the relationship betweenfi and xi for an individual LAD is given, the relationship between M0 and θ0 can be determined through equations (18)–(20).For LADi assembled in the isolation device, the classical Bouc–Wen model from [28] is used to describe its hysteretic characteristic. The relationship between its hysteretic restoring force history fit and deformation displacement history xit can be described by an operator ΦBW:xit⟶fit, which is governed by the following equations:(21)fi=λFyδxi+1−λFyzi,z˙i=x˙iδA−β+γsgnx˙izizp,where Fy is denoted as yield force, δ represents yield displacement, λ is the parameter of postyielding stiffness ratio, zi is a dimensionless quantity served as an auxiliary internal state variable that obeys a single nonlinear differential equation, the model parameters A, β, γ, and p are dimensionless quantities controlling the shape of the hysteresis loop, and sgn (⋅) is the signum function. For mathematical and physical consistency discussed in [29], two constraints, A=1 and β+γ=1, are imposed to reduce the Bouc–Wen model to a strain-softening formulation with well-defined mechanical properties.Combining equations (19)–(21) and applying the following equality for positive integer n as n>2(22)∑i=1ncos2α+i−12πn=n2,it can be derived that the relationship between M0 and θ0 of the isolation device with LADs in Bouc–Wen type is governed by the following equations:(23)M0=nd2λFyδr2θ0+1−λFyr∑i=1ndzicosαi,z˙1=rδθ˙0cosα11−β+γsgnθ˙0z1cosα1z1p,⋮z˙nd=rδθ˙0cosαnd1−β+γsgnθ˙0zndcosαndzndp.From a mechanical point of view, the above equations can be considered as the constitutive equations of a hysteretic torsional spring with a group ofnd internal state variables. The initial rotational stiffness of the isolation device, denoted by s0, is defined as(24)s0=dM0dθ0|t=0=M˙0θ˙0|t=0.Based on equations (22), (23) and the initial conditions zi0=0 for i=1,…,nd, we can derive that(25)s0=nd2Fyδr2.The above equation shows that the initial stiffness of the isolation device can be determined by the structural parameters of the device and the mechanical characteristics of the LAD. This analytic expression can be used in estimating the fundamental frequency of the isolated equipment. It is worth noting that the above derivation can be easily extended to the isolation device with LADs modeled by other types of, for example, bilinear or multilinear, hysteresis model.


### 3.2. Modeling the Isolated UHV Porcelain Arrester

The model of an UHV porcelain arrester equipped with the isolation device is depicted in Figure10. The downside rigid beam of the isolation system is fixed to the ground. The upside rigid beam is indexed as body (0) and pinned to the body (1). The spring and damper attached to the bottom end of body (1) are connected to the body (0). So the isolated arrester consisting of n units and one grading ring is modeled as a (n + 2)-body system.Figure 10
Reference frame and configuration of the isolated multibody system at some instant of time.Choose the absolute rotational anglesθ0,θ1,…,θn as the joint coordinates of the isolated system. Then the coordinate of body (0) can be formulated as c0=uGt,0,θ0T. Compared with those arrester units, the mass of the upside plate in the isolation device is so small that it can be neglected. Accordingly, the mass matrix M˜0 is a 3×3 zero matrix, and the applied force array on body (0) can be expressed as h0=0,0,M1s+M1d−M0T, where M1s=s1θ1−θ0, M1d=d1θ˙1−θ˙0 and M0 is the moment due to the isolation system having expression of equation (23). Then, the body configuration, the mass matrix, the applied forces array, and the velocity transformation matrix of the isolated system can be constructed easily.The governing equations of motion for the isolated system are derived, and by applying small angle approximation, they can be expressed in matrix and vector notation as follows:(26)000Mθ¨0θ¨+d1−d1e1T−d1e1Dθ˙0θ˙+s1−s1e1T−s1e1S−G,θ0θ+M00=u¨Gt0mιwhere e1 is a column vector of n elements with the first one being 1 and the others being 0. Equation (26) has a singular mass matrix. From a physical point of view, it represents that a hysteretic spring is directly connected in series to the base spring of the nonisolated system. Thus, the estimation of the fundamental frequency of the isolated system can be obtained through replacing s1 in S by s0s1−/s0+s1− and solving the generalized eigenvalue problems of S−G and M. To solve equations (26), (23) has to be combined with; that results in a system of n second-order equations and nd+1 first-order equations. Convert the equations into first-order form and choose θ,θ˙,θ0,z1,…,znd as the state variables. Because of equation (23), the first-order system can be solved by a solver for fully implicit ordinary differential equations.


## 3.1. Modeling the Isolation Device

In this section, a model of the seismic isolation device introduced in [12] which is depicted in Figure 2 is presented. Denote by nd the number of the LADs and by r the radius of the LAD circle. Construct two axes as shown in Figure 8(a), such that the x-axis coincides with the direction of the horizontal ground motion. Starting from the positive x-axis, number the LADs successively in a counterclockwise direction. Define the initial angle α as the angle between x-axis and the first LAD. Then, the location of LADi can be determined by the following angular expression:(18)αi=α+i−12πnd.Figure 8
Model of the isolation device: (a) effective arm of force, (b) planar isolation system.
(a)(b)Under the assumption of connection steel plates being rigid, the effective arm of force due to LADi about x′-axis can be expressed as ai=rcosαi. Reduce the whole isolation device to a planar system by means of regarding the connection plates as rigid beams, treating the cushion block as a hinged support fixed to the downside beam and pinned to the midpoint of the upside beam, and modeling each LAD as a hysteretic spring placed at its projection on x-axis that connects the two beams. The schematic diagram of the planar isolation system is depicted in Figure 8(b).Denote byM0 the external bending moment applied to the upside beam and θ0 the associated rotational angle responding to the moment, which is shown in Figure 9. The moment-angle relationship of the planar isolation system is derived theoretically in the following.Figure 9
Moment and associated angle of the planar isolation system.Denote the deformation of LADi caused by upside beam rotation as xi. Then it can be expressed as a multiplication of the rotational angle and the effective arm of force associated with LADi as follows:(19)xi=θ0rcosαi.The sign ofxi indicates the deformation state of LADi in a manner that positive means it is in tension, while negative means it is in compression. Denote by fi the applied force that induces the deformation xi of LADi. Then the moment distributed on LADi is fircosαi, and the external moment M0 is the summation of those distributions as follows:(20)M0=∑i=1ndfircosαi.Once the relationship betweenfi and xi for an individual LAD is given, the relationship between M0 and θ0 can be determined through equations (18)–(20).For LADi assembled in the isolation device, the classical Bouc–Wen model from [28] is used to describe its hysteretic characteristic. The relationship between its hysteretic restoring force history fit and deformation displacement history xit can be described by an operator ΦBW:xit⟶fit, which is governed by the following equations:(21)fi=λFyδxi+1−λFyzi,z˙i=x˙iδA−β+γsgnx˙izizp,where Fy is denoted as yield force, δ represents yield displacement, λ is the parameter of postyielding stiffness ratio, zi is a dimensionless quantity served as an auxiliary internal state variable that obeys a single nonlinear differential equation, the model parameters A, β, γ, and p are dimensionless quantities controlling the shape of the hysteresis loop, and sgn (⋅) is the signum function. For mathematical and physical consistency discussed in [29], two constraints, A=1 and β+γ=1, are imposed to reduce the Bouc–Wen model to a strain-softening formulation with well-defined mechanical properties.Combining equations (19)–(21) and applying the following equality for positive integer n as n>2(22)∑i=1ncos2α+i−12πn=n2,it can be derived that the relationship between M0 and θ0 of the isolation device with LADs in Bouc–Wen type is governed by the following equations:(23)M0=nd2λFyδr2θ0+1−λFyr∑i=1ndzicosαi,z˙1=rδθ˙0cosα11−β+γsgnθ˙0z1cosα1z1p,⋮z˙nd=rδθ˙0cosαnd1−β+γsgnθ˙0zndcosαndzndp.From a mechanical point of view, the above equations can be considered as the constitutive equations of a hysteretic torsional spring with a group ofnd internal state variables. The initial rotational stiffness of the isolation device, denoted by s0, is defined as(24)s0=dM0dθ0|t=0=M˙0θ˙0|t=0.Based on equations (22), (23) and the initial conditions zi0=0 for i=1,…,nd, we can derive that(25)s0=nd2Fyδr2.The above equation shows that the initial stiffness of the isolation device can be determined by the structural parameters of the device and the mechanical characteristics of the LAD. This analytic expression can be used in estimating the fundamental frequency of the isolated equipment. It is worth noting that the above derivation can be easily extended to the isolation device with LADs modeled by other types of, for example, bilinear or multilinear, hysteresis model.


## 3.2. Modeling the Isolated UHV Porcelain Arrester

The model of an UHV porcelain arrester equipped with the isolation device is depicted in Figure10. The downside rigid beam of the isolation system is fixed to the ground. The upside rigid beam is indexed as body (0) and pinned to the body (1). The spring and damper attached to the bottom end of body (1) are connected to the body (0). So the isolated arrester consisting of n units and one grading ring is modeled as a (n + 2)-body system.Figure 10
Reference frame and configuration of the isolated multibody system at some instant of time.Choose the absolute rotational anglesθ0,θ1,…,θn as the joint coordinates of the isolated system. Then the coordinate of body (0) can be formulated as c0=uGt,0,θ0T. Compared with those arrester units, the mass of the upside plate in the isolation device is so small that it can be neglected. Accordingly, the mass matrix M˜0 is a 3×3 zero matrix, and the applied force array on body (0) can be expressed as h0=0,0,M1s+M1d−M0T, where M1s=s1θ1−θ0, M1d=d1θ˙1−θ˙0 and M0 is the moment due to the isolation system having expression of equation (23). Then, the body configuration, the mass matrix, the applied forces array, and the velocity transformation matrix of the isolated system can be constructed easily.The governing equations of motion for the isolated system are derived, and by applying small angle approximation, they can be expressed in matrix and vector notation as follows:(26)000Mθ¨0θ¨+d1−d1e1T−d1e1Dθ˙0θ˙+s1−s1e1T−s1e1S−G,θ0θ+M00=u¨Gt0mιwhere e1 is a column vector of n elements with the first one being 1 and the others being 0. Equation (26) has a singular mass matrix. From a physical point of view, it represents that a hysteretic spring is directly connected in series to the base spring of the nonisolated system. Thus, the estimation of the fundamental frequency of the isolated system can be obtained through replacing s1 in S by s0s1−/s0+s1− and solving the generalized eigenvalue problems of S−G and M. To solve equations (26), (23) has to be combined with; that results in a system of n second-order equations and nd+1 first-order equations. Convert the equations into first-order form and choose θ,θ˙,θ0,z1,…,znd as the state variables. Because of equation (23), the first-order system can be solved by a solver for fully implicit ordinary differential equations.


## 4. Analysis of an UHV Porcelain Arrester Specimen

In this section, a 1000 kV porcelain arrester specimen is analyzed by using the proposed method. The specimen is composed of five units of the same type and one different base unit. The cementing parameters of these units are listed in Table1. The rotational stiffness parameters are calculated and listed in Table 2 together with the mass and length parameters.Table 1
Cementing parameters of the porcelain arrester specimen.
UnitEndhcmtemdcmBase unit±0.2000.0100.495Units 1∼5±0.2000.0100.510Table 2
Multibody model parameters of the porcelain arrester specimen.
BodyIndexmikglimsiN⋅m⋅rad−1Base unit13000.6239.9×107Unit 129122.1155.0×107Units 2∼53∼69122.1155.1×107Grading ringTop380——The specimen of arrester was experimentally investigated in [21] for its flexural performance. Two levels of static testing were carried out to investigate the relationship between the force and horizontal displacement at the top free end. In the following computation, the gravitational acceleration is fixed to be g=9.8m⋅s−2. The analytical results compared with the testing results are depicted in Figure 11. It shows very good agreement between them when loading does not exceed 35 kN. The analytical displacement to the first maximum 26.99 kN force is 167.33 mm, which has relative error 1% to the testing result of 169.02 mm, while under second maximum 43.10 kN lateral force, the analytical displacement is 267.15 mm, which has a relative error of 29.02% to the testing result of 376.39 mm. The reason why relative error increases can be summarized as the following statement. The force of 43.10 kN is nearly the lateral static load which could cause collapse of this specific arrester. The stiffness would decrease when cement joint experiences large deformations exceeding the range of linear stiffness. However, the computation has been performed without considering this nonlinear phenomenon of stiffness reduction.Figure 11
Force-displacement relationship at the top free end of the arrester specimen.The following presents a numerical analysis on the effectiveness of the seismic isolation device. The number of the LADs in the adopted device isnd=8 and the radius of the circle is r=0.305m. The LAD has the following parameters: Fy=22.0kN, δ=0.2mm, λ=0.02, β=γ=0.5, and p=1.The fundamental frequency of the nonisolated arrester is calculated to be 1.57 Hz. After equipment with the isolation device, the fundamental frequency decreases to 1.24 Hz. For response history analysis, the input ground motion chosen here is an artificial acceleration history generated according to the required response spectrum (RRS) in the Chinese Code [20]. The artificial accelerogram is shown in Figure 12, and its response spectrum and the RRS are shown in Figure 13. It is normalized to have 1g peak ground acceleration (PGA). This motion is scaled with desired PGA to fit the required seismic qualification level and scaled up by a factor of 1.4 to represent the dynamic amplification of the supporting structure.Figure 12
Artificial earthquake accelerogram.Figure 13
Response spectra of the artificial earthquake and the RRS.At seismic qualification level of0.2g, with assumed damping ratio ζ=2%, under situation of α=0 (see Figure 8), numerical computations are carried out. The response histories of the internal moment at each joint are then calculated. Comparison of the internal base moment between nonisolated and isolated system is shown in Figure 14. It is evident that there is a noticeable decrease of the base moment over the whole duration of the ground motion. The peak value of the base moment shifted from 2.52×105N⋅m to 0.42×105N⋅m with a reduction of 83.2%. The response histories of the absolute acceleration at the top free end of the nonisolated and isolated system are shown in Figure 15. Similar to the base moment, the response of the top absolute acceleration has a noticeable reduction in the whole duration. The peak value changes from −1.26g to −0.45g with a 64.3% decrease in magnitude. In contrast to the base moment and the top absolute acceleration, the response of the relative displacement at the top free end increases in some duration of the motion, which is shown in Figure 16. The peak value changes from 126 mm to 166 mm with an increase of 31.4%. It is worth noting in Figure 16 that there is residual displacement at the end of motion, which is due to the hysteretic property of the isolation device.Figure 14
Response histories and corresponding Fourier spectra of the base moment of the nonisolated and isolated system with 2% damping at 0.2 (g) seismic qualification level underα = 0.Figure 15
Response histories and corresponding Fourier spectra of the absolute top acceleration of the nonisolated and isolated system with 2% damping at 0.2 (g) seismic qualification level underα = 0.Figure 16
Response histories and corresponding Fourier spectra of the relative top displacement of the nonisolated and isolated system with 2% damping at 0.2 (g) seismic qualification level underα = 0.In order to show the influence of the parameters of the device to the response, the maximum magnitudes of three selected quantities, the base moment (BM), top absolute acceleration (TAA), and top relative displacement (TRD), are calculated for the nonisolated and isolated system with damping of 1%, 2%, and 5% and at two seismic qualification levels of0.2g and 0.4g. As the number of LADs in the isolation device is 8, these quantities are also calculated under α=0 and α=π/8. The results and their changes due to the isolation device are illustrated in Figures 17–19. Firstly, these selected quantities of the nonisolated system are proportional to the seismic qualification level because of the linearity of the system. The isolated system does not have this property due to the nonlinearity introduced by the isolation device. It can be seen from Figure 17 that the effectiveness of the isolation system is significant in terms of reducing the maximum BM. The damping ratio has less influence on the maximum BM of the isolated system than that of the nonisolated system. This is because the isolation device contributes the most damping effect in the isolated system. The parameter α has little influence on the maximum BM of the isolated system, which is expected. Figure 18 shows that the maximum TAA varies with the parameters of α, ζ and seismic level in a similar way to the maximum BM illustrated in Figure 17. The reduction of the TAA is also large but less than that of the BM. In contrast to the BM and TAA, Figure 19 shows that the maximum TRD is increased by the isolation system. The increase is notable especially at 0.4g seismic qualification level. The maximum TRD is slightly influenced by ζ and relatively sensitive to α.Figure 17
Maximum base moments (BM) of the nonisolated and isolated system and their change under different parameters ofζ and α at two seismic qualification levels.Figure 18
Maximum top absolute acceleration (TAA) of the nonisolated and isolated system and their change under different parameters ofζ and α at two seismic qualification levels.Figure 19
Maximum top relative displacement (TRD) of the nonisolated and isolated system and their change under different parameters ofζ and α at two seismic qualification levels.


## 5. Conclusions

In this paper, a novel modeling method has been proposed for the seismic analysis of UHV porcelain arresters equipped with an isolation device. An UHV porcelain arrester has been modeled as a planar multibody system based on the fact that the stiffness distribution is uneven along the axial direction of the arrester. A joint coordinate method has been chosen to construct the governing equations of motion. The number of the DOF of the multibody system equals the number of the arrester units, which is usually less than 6. This is beneficial to analysis requiring less numerical computation for determining the optimal design of the isolation device.With regard to the seismic isolation device, classical Bouc–Wen hysteretic model has been used in modeling its core energy dissipation components. The entire device has been reduced to a planar system, and its mechanical characteristic has been derived theoretically. The governing equations of motion for the isolated arrester have been constructed naturally by embedding the isolation system into the nonisolated multibody system. From a mechanical point of view, the isolation device is represented as a hysteretic rotational spring with a number of internal state variables. Its initial rotational stiffness is formulated in terms of the structural and mechanical parameters of the isolation device. This analytic expression can be used to estimate the fundamental frequency of the isolated equipment, which can guide the design of the isolation device.An UHV porcelain arrester specimen has been analyzed by the proposed method. In the static analysis the relationship between lateral loading and displacement at the top free end shows good agreement with the test results before failure occurs. The relative error increases when loading force exceeds a limit which indicates that the nonlinearity of the stiffness at joints should be considered in the situation of near failure. This behavior could be studied in future research work. The presented numerical analysis shows that the effectiveness of the isolation device is significant in terms of reducing the maximum base moment. The reduction of the maximum top absolute acceleration is also noticeable. However, the enlargement of the top relative displacement at high seismic qualification level is a disadvantage of the isolation device.The proposed method primarily focuses on modeling UHV porcelain arresters, while it can be extended to model other porcelain post electrical equipment, such as insulator and capacitor voltage transformer.



---

*Source: 5549045-2021-04-07.xml*