Mathematical Modeling of Stress in Circuit Cards Represented by Mechanical Oscillatory Systems

The represented paper is aimed at stress calculation in circuit cards with their representation as a type of mechanical oscillatory systems in purpose of their strength assessment especially in resonance conditions. Three types of oscillatory systems are researched: single-mass; multiple mass and oscillatory system with uniformly distributed mass. In all types the cylindrical bending of circuit cards is considered to be a set of beam-strips with rectangular cross-sections so their stress calculation is performed by conventional methods applied in strength of materials and civil engineering. Mathematical model has been developed for maximal dynamic stress and defl ection estimation in circuit card assemblies represented by unique oscillatory system as prismatic beam set on two oscillating supports under inertial resonance excitation generated by constant dynamic force. Comparative analysis of mathematical modeling, MatLab simulation and experimental determination of maximal dynamic stress and defl ection accomplished for three types of oscillatory systems verifi ed proximity of obtained results. Single-mass oscillatory system is proposed as equivalent to multiple mass or uniformly distributed oscillatory systems on condition of their equal mass, geometric, elastic and dissipation characteristics in resonance frequency correspondent to the main mode of oscillation, so mathematical model designed for single-mass oscillatory system is recommended for strength and stiff ness assessment in engineering calculations where possible diff erence in determination of stress in equivalent systems can used as safety factor.


INTRODUCTION
Variety of modern electronic packages and their parts such as circuit cards (CC) or case walls, etc. are very likely to be exposed to mechanical impacts such as vibration and shocks during their operation. Large number of publications discussing problem of dynamic forces analysis in engineering [1][2][3] and in electronic equipment [4][5][6][7] subjected to vibrations, vibration reduction and suppression design [8,9] testify of need for improving strength and providing reliability to electronics. Results of the previous research published in [10,11] emphasized that mechanic impacts, and especially the dynamic forces in CC assemblies are likely to increase manifold so as to damage their bearing parts and electronic components, to which these forces are transmitted especially in resonant oscillations.
The represented research is aimed at stress calculation in CCs with their representation as a type of mechanical oscillatory systems in purpose of their strength assessment. The most attention is drawn to oscillations in resonance conditions. In all representations of oscillatory systems the cylindrical bending of CCs is considered to be a set of beam-strips with rectangular cross-sections so their stress calculation is performed by conventional methods applied in strength of materials and civil engineering [12].

Mathematical Modeling of Stress in Circuit Cards
Represented by Mechanical Oscillatory Systems 304 The advantage of the mathematical modeling and of analytical formulas would be in their ability to predict the dependence of the results on the parameters of the tested systems, and it has been successfully applied for the single-mass oscillatory system. The attempts were also made to model multiple mass or even uniformly distributed oscillatory systems, however they did not result in analytical formulas capable to adequately calculate stress and deflection because of complexity of considered systems. In order to replace and for the case of single-mass oscillatory system to verify analytical formulas the MatLab simulation and experimental tests were applied.
Analytical estimation of stress and deflection conducted in the research was based on experimental verification of actual physical and mechanical parameters, which are likely to vary depending on technology, temperature, shape etc. For this purpose analytical and experimental method of sample parameters was applied. This method is based on identification of parameters by solving the reverse strength problems. In this method strain and displacement, which are traditionally estimated by calculations, are, instead, measured experimentally and considered as given values for calculation of physical and mechanical parameters, which are then used in analytical modeling.

Single-mass oscillatory system
The single-mass oscillatory system is represented by the beam with concentrated mass m (Fig. 1). Weight of the beam is assumed negligibly small in comparison with concentrated mass.
Such representation is applicable for assembly where mass of electronic components exceeds mass of the substrate and therefore produces imbalance in oscillatory system. The beam is supported by the pinned support O providing with one degree of freedom and roller O' support providing two degrees of freedom. These supports transmit oscillation generated by the shaker.
Mass m, receiving kinematic excitation, undergoes acceleration z 1 ′′ in the inertial frame of reference indicated by 0 ( Fig. 1), which represents a support to the shaker generating acceleration z 0 ′′. Then equation of motion for the mass m: where: F k -elastic force; F c -damping force.
In inertial frame of reference absolute acceleration z 1 ′′ may be expressed by relative ∆z′′ and fictitious z 0 ′′ accelerations as: Thus, equation of motion in non-inertial frame of reference, represented by O and O' supports, can be written as: where: m∆z′′ -relative and mz 0 ′′ -fictitious forces of inertia; F k = k∆z; F с = c⋅∆z′; k -stiffness; c -damping coefficient. Equation (3) written in terms of function ∆z(t) is: where: 0 = � = 2 -natural frequency; Shaker generates harmonic oscillations along Z axis described as: where: Z 0 -amplitude; ω -angular frequency; ttime; ϕ -phase of oscillations.
According to identical differential equation published in [10] the obtained equation (4) has identical solution that gives amplitude of forced oscillation: Function (6) expresses dynamic deflection ∆ dyn of the beam with respect to static equilibrium position ∆ st in the oscillatory system ( Fig. 1), which sum expresses the total deflection Δ t max produced by the net force Р: where: P dyn = -F k ; P st = mg; δ -flexibility of the beam in direction of force P.
Flexibility of the beam where a single force is applied is defined by Mohr's method: where: x -linear coordinate of force application; E -Young's modulus; = Formula (8) considers only internal bending moments while internal transverse forces can be neglected.
Maximal total normal stress in the critical cross-section of the beam defined by net force (9) and bending formula [12] is expressed as: where: = Such assumption results in negligible error in calculations that does not exceed 3%.
Mathematical modeling was conducted for parameters given in Table 1. Average values of physical and mechanical parameters were identified by method of sample parameters, in which physical and mechanical parameters were calculated from experimentally measured values of strain and displacement. Method is described in section of experimental verification, the experimental setup is shown in Figure 11.
The strength condition for CC was specified with respect to lowest ultimate strength in the CC assembly represented by soldered joints [11]: where: s = 2.5 -safety factor compensating instability in design and technology.
In estimation model the dynamic force is applied with constant dynamic force amplitude P 0 = const provided by constant acceleration a 0 = Z 0 ω 2 = 10 m/s 2 . Amplitude frequency responses of maximal stress and deflection given in Figure 2 indicate resonance at first critical Moreover, analysis of function (14) revealed that neglecting static action of gravity force P st , which is considerably lower in comparison with dynamic force P dyn , gives the expression of dynamic stress: which testifies that under given conditions maximal stress does not depend on linear coordinate x of the concentrated mass m. Thus, formula (16) is proposed to use in engineering calculations of stress in CCs represented by oscillatory systems with concentrated mass irrespective to location of its mass.
Besides, function (16) as demonstrated by (Fig. 3) can be used to specify geometric parameters of CC which would ensure compliance with strength condition. In particular, the minimal acceptable length l min is proposed to calculate by the following formula: This formula is validated by relative difference ∆σ between (14) and (16) to comply with threshold tolerance of 15%. For this example (Fig. 3) l min = 273 mm at ∆σ = 10.06%. Formula (14) could be used to increase accuracy if needed.
Mathematical model (16) was initially designed for maximal stress calculation in CC assemblies represented by single-mass oscillatory system with one concentrated mass, however real assemblies consist of multiple electronic components and bearing parts and each of them represents an unique concentrated mass. Thus, such mechanical system subjected to vibration represents a multiple mass oscillatory system. Figure 4 illustrates a model of oscillatory system that consists of 3 masses m 1 , m 2 , m 3 and therefore has three dynamic degrees of freedom since three parameters ∆ 1 , ∆ 2 , ∆ 3 indicate position of all concentrated masses at their vertical displacements.

Multiple mass oscillatory system
According to principle of superposition of forces the displacement in specific direction i correspondent to specific degree of freedom equals to the sum of all displacements in this direction (i = 1..3) produced by all forces acting in directions j (j = 1..3) and is defined by Mohr's method: According to equation of motion (3) and expression (18) equation of displacement for the mass m i : Thus a system of differential equations (i = 1..n) can be written and solution in the following form can be obtained: For representation of only free oscillations where in (19) excitation force (fictitious forces of inertia in this case) and damping force are ignored and so are single displacements (flexibilities) δ ij (i ≠ j), then system of differential equation is falling apart for equations of single displacements δ ii with their natural frequencies: Nevertheless determining basic coordinates ∆ i for the system with more than 2 degrees of freedom encounters certain difficulties.
Symmetric systems with masses symmetrically located may produce either directly or inversely symmetric oscillation forms, where acting forces will be correspondently either directly or inversely symmetrical. In this case single displacements will be calculated as for the group of either directly or inversely symmetric unique forces.
Side single displacements which connect directly and inversely symmetric forces become equal zero: When number of masses is higher than 3, then receiving solution may become impossible especially for the case of considering forced and damped oscillations. Therefore solution to the problem of finding natural frequencies, modes of oscillation and, what matters most, maximal total normal stress and deflection in oscillatory system was offered to obtain by numerical method of simulation by MatLab system.
Simulation conducted by visual programming system Simscape Multibody in Matlab represented PCB as a flexible beam performed by multi-body three-dimensional modeling environment of mechanical systems aimed at forming and solving equations of motion. Since stiffness of supports is considered significantly higher than stiffness of the board, supports are considered absolutely rigid. Modeling flexible plate was performed by using method for approximating flexible body by the set of N solid bodies having concentrated parameters and connected with springs and dampers [13]. Every concentrated mass has one degree of freedom -rotation in ZOX plane (Fig. 5). The mass, spring, and damper elements provide the inertial, restorative, and dissipative forces that collectively account for deformation. Linear spring and damper model was considered.
Stiffness and internal viscous friction are functions of material properties and geometry of flexible elements and are equal for all elements.
The value of the stiffness follows from the equality between the spring torque at the joint and the bending moment on a continuous version of the flexible beam unit. Hooke's law gives the spring torque at the joint: where: k r -rotational spring stiffness; γ -deflection angle.
The bending moment on a continuous beam unit: where: R -the bending radius of curvature.
In the limit of very small deflections γ → l/R, where l is undeformed length of the elements, spring stiffness is: In the represented model k r = 12.03 N•m.
Damping (viscous friction) b = 0.0037 for dampers in the joints was specified for experimentally found free oscillations descent rate n = 5.48 1/s (50) of the beam, as described further in the section of experimental verification. Figure 6 demonstrates simulative approximation of 5-mass oscillatory system. The force application and parameters of oscillatory system were identical to single-mass model with the exception of using 5 symmetrically located masses m i = 0.01 kg. Modeling flexible beam (substrate) was performed by using method for approximating flexible body by the set of N = 21 solid bodies having concentrated parameters and connected with springs and dampers [13,14]. Deflections were measured in simulation by Joint Sensor. Maximal stresses were calculated by using bending formula [12], in which maximal bending moments were measured by Body Sensor. 6 pairs of sensors of both types were attached as shown in Figure 6. The simulation of oscillatory system has resulted in amplitude frequency response graphs read by 6 sensors for maximal total normal stress (Fig. 7a) and deflection (Fig. 7b). Noteworthy is that highest levels of measured parameters are correspondent to the first critical frequency ω 01 = 237 rad/s for the main mode of oscillation.
Obviously the more detailed model is the more accurate calculation results will be. Moreover consideration of all masses in real CC assembly representing numerous electronic components and a substrate laminated with conductive tracks, pads and other features and therefore having no-zero mass introduces complex and cumbersome problem. Such multiple mass mechanical system with actually infinite number of degrees of freedom is desirable to study as a system with mass distributed over its volume and in case of a beam representation -over its length and consider this distribution as uniform.

Oscillatory system with uniformly distributed mass
Equation of motion for the beam with uniformly distributed mass is represented for infinitesimal element with length dx (Fig. 8) arbitrarily selected along the length l of the beam. Element dx undergoes vertical displacement described by the function z(x,t) forced by uniformly distributed dynamic inertial force f(x,t). In this case force f(x,t) is represented by fictitious forces of inertia (analogically to single mass oscillatory system described above in noninertial frame of reference) which is spent to overcome uniformly distributed relative force of inertia dm⋅z(t)'', damping force F с , and concentrated elastic force F k in cross-sections of the beam, which isolate element dx by the left and right faces where forces F k (x,t), F k (x+dx,t, and moments of internal resistance M(x,t), M(x+dx,t) appear as shown in Figure 8.
The equation of motion for dx forced by all forces indicated in Figure 8: where: ρ -density; s -cross-sectional area of the beam.
Division by dx, substitution of ( , ) = ( , ) and application of Solution to the equation of free transverse periodic oscillations represented by equation (33), in which the right side of it equals zero, is considered to be the harmonic function of beam deflection: where: A x (x) -amplitude function of beam axis, which represents the main mode of oscillations.
Then С 1 = С 3 = С 4 = 0; С 2 sin (βl) = 0. When С 2 = 0 solution to (38) is trivial, therefore when С 2 ≠ 0 it is considered that: Whence: Considering (37) the natural frequencies of transverse periodic oscillations will be expressed as: Since strength and stiffness assessment is performed in engineering calculations by maximal magnitudes of stress and deflection then consideration of lower critical frequencies or even of the first one ω 01 , when j = 1, representing the main mode of oscillation is required. This statement was approved by MatLab simulation, in which flexible body of the beam (circuit card assembly) was approximated by the set of 21 elements (solid bodies) joined in between them by springs and dampers. The total mass of the beam m = 0.05 kg was uniformly distributed over all elements. 6 sensors attached as shown in Figure 6 measured total normal stresses and deflections.
Simulation has resulted in amplitude frequency response graphs read by 6 sensors for maximal total normal stress (Fig. 9a) and deflection (Fig.  9b). Noteworthy is that highest levels of measured parameters are correspondent to the first critical frequency ω 01 = 255.3 rad/s for the main mode of oscillation.
Comparative analysis of stress and deflection diagrams obtained by the research of oscillatory systems of all mentioned types (Fig. 10) demonstrated that the main mode of oscillatory system with uniformly distributed mass is similar to multiple mass

Analysis of equivalence for strength and stiffness assessment in oscillatory systems
Similarity of main modes of oscillation indicated by proximity of maximal total normal stresses and deflections obtained for three mentioned types of oscillatory systems, on condition of equal dynamic force applied and mass, geometric, elastic and dissipation characteristics, allowed conjecture of representing these systems as equivalent for strength and stiffness assessment.
There is known [15], that an oscillatory system with uniformly distributed mass m can be 312 represented as oscillatory system with concentrated mass converted into reduced mass β⋅m and applied in a specific place on the beam, in this case in its center, to perform strength and deflection assessment by using mathematic model designed for single-mass oscillatory system. Equivalence of two systems can be reached by equality of their kinetic energies therefore the value of reduced mass β⋅m is assumed such that its kinetic energy: is equal to kinetic energy of mass distributed along the length l of the beam: where: z' -displacement velocity of concentrated mass centered to the beam; z x ' -displacement velocity of infinitesimal element, determined by its position on the axis x along the beam. Whence: Hypothesis about similarity of displacement diagrams produced by the static action of weight of distributed mass m and concentrated mass β⋅m gives the relationship: where: ∆ xst -displacement of infinitesimal element dx produced by static action of weight of distributed mass in oscillatory system; ∆ st -displacement of the beam center produced by static action of weight of concentrated mass (Fig. 11).
Displacement of the section x: and displacement of the beam center: Considering (44) and substituting (45) and (46) into (43) results in β = 17/35. Nevertheless using equivalent model of oscillatory system with reduced mass is only valid when mass of CC assembly is uniformly distributed and in case of multiple mass or single-mass oscillatory systems such approach would bring about lowered result in calculation of stress and deflection.
Analysis of simulation results given as graphs in Figure 10 verified that value of maximal total normal stress in oscillatory systems is that less than more distributed their mass is and, on the contrary, stress is highest for single-mass oscillatory system. Thus, calculation of highest values of maximal total normal stress through representation of single-mass oscillatory system as equivalent to any other systems even to oscillatory system with uniformly distributed mass allows strength assessment by using mathematical models designed for single-mass oscillatory system (16) or (14), whereby possible difference in calculation related to using equivalent system instead of original one can be used as safety factor.

Experimental verification
Research objectives were circuit cards representing three types of oscillatory systems: 1) single-mass oscillatory system with concentrated mass centered to the beam; 2) multiple mass oscillatory system consisting of 5 masses placed uniformly along the beam (Fig. 6); 3) Fig. 11. Main mode of oscillation of the beam with uniformly distributed mass oscillatory system with uniformly distributed mass. Circuit cards were selected on condition of equal mass, geometric, elastic and dissipation characteristics and dynamic force application and identical to those in mathematical modeling and MatLab simulation.
Actual elastic and dissipation parameters were determined by the method of sample parameters based on measuring strain and displacement in the experimental setup shown in Figure 12.
Young's module E and descent rate n where found by static middle-point bend test of CC samples under force threshold limit: During the test the strain was measured by strain gauges along with reading force P and deflection ∆. Two strain gauges with 10 mm base were attached both in the longitudinal (direction of maximal normal stress) and in the transverse directions (indicated no strain). For higher accuracy the Young's module E was found in two ways: The mean value of E = 14.5 GPa was calculated over 5 experiments. Tests were conducted by using method of acoustic emission reading from the sensor attached to PCB surface. Detecting no acoustic emission verified of only elastic strain produced during the test and approved Hook's law application. Descent rate n of oscillatory system was found by free oscillation graph read from vibration sensor attached to the center of PCB by the formula: where: A(t) and A(t+T) -amplitudes measured in time t and t + T correspondently, T -oscillation period.
The mean value of n = 5.48 s -1 . The further vibration tests were conducted for experimental verification of mathematical model. Vibration tests were conducted by constant dynamic load with amplitude P 0 provided by constant vibration acceleration a 0 = 10 m/s 2 . CCs were installed in clamps designed as regular fixtures (Fig. 13) mounted on the shaker. Maximal normal stress