# Center Impedance Method for Estimating Complex Modulus with Wide Frequency Range and Large Loss Factors

Authors: Ryuzo Horiguchi; Yoshiro Oda; Keito Sato; Hiroto Kozuka; Takao Yamaguchi
Journal: Shock and Vibration
 (2021)

---

## Abstract

A simple method for determining viscoelasticity over a wide frequency range using the frequency response function (FRF) mobility obtained by the center impedance method is presented. As user data comprise the FRF between the velocity of the excitation rod and excitation force, it is challenging to separate the signal and noise. Our proposed method is based on the FRF obtained from the analytical solution of the equation of motion of the viscoelastic beam and relationship between the complex wavenumber (real wavenumber and attenuation constant) of flexural wave and viscoelasticity. Furthermore, a large loss factor can be handled over a wide frequency range without using the half-power bandwidth. In this study, actual FRF mobility data containing noise were processed using preprocessing, inverse calculation, and postprocessing. Preprocessing removed low-coherence data, compensates for the effects of instrument gain, and transformed the FRF into its dimensionless equivalent. Then, inverse calculations were used to solve the mobility equation and determine the complex wavenumber. In postprocessing, the complex wavenumber obtained by the inverse calculation was curve fitted using functions with mechanical significance. Consequently, the storage modulus based on the curve-fitted complex wavenumber was a monotonically increasing frequency function. The loss factor had a smooth frequency dependence such that it has the maximum value at a single frequency. The proposed method can be applied to composite materials, where the application of time-temperature superposition is challenging. We utilized the measured FRF mobility data obtained over a duration of several seconds, and this method can also be applied to materials with large loss factors of 1 or more.


---

## Body

## 1. Introduction

Acoustic radiations from the surface of plates at the floor of a car decrease with its bending vibrations, which is effective for reducing noise from interiors of automobiles, railway trains, and aircrafts. In addition, vibration damping is necessary for passenger comfort. Recently, polymeric composites have been widely used for this purpose, and therefore, vibration damping properties of various plate materials must be evaluated. To this end, tests are performed at the point of beam vibrations, and the storage modulus and loss factor of the specimen are estimated from the responses. The Oberst beam method is used to measure the responses of composite beams with the material attached to an aluminum or steel substrate for damping materials with low storage moduli [1, 2]. According to ASTM E756, the Oberst beam is magnetically vibrated, and the beam responses are measured using a noncontact magnetic transducer. However, the beam center can be vibrated without using a magnetic substrate, and this method is called the center impedance or central excitation method as specified in JIS K 7391 [2]. Although this is increasingly used in Japan than elsewhere, its usage worldwide is expanding because of the relatively convenient installation of the equipment [3]. In the center impedance method, the vibration damping property of the specimen is evaluated using the transfer function between the specimen displacement (velocity or acceleration) and vibration force applied to the specimen center. A piezo element and force gauge are constructed into the impedance head attached to the tip of the vibration rod to measure the vibration force and acceleration.In the extant measurement standards [1, 2], Young’s modulus for a specimen regarded as a purely elastic material is estimated using the resonant frequency. However, that of a viscoelastic material varies because of its large loss factor. These standards indicate that the half-power bandwidth can be used to estimate this factor near the resonant frequency when assuming structural loss factors to be low. When the half-power bandwidth is used, the system comprising viscoelastic beams is approximated by an equivalent lumped constant system consisting of a mass point, spring, and dashpot. Here, it is assumed that the square of the damping ratio is much less than 1 and the structural loss factor is twice the damping ratio. However, in the half-power bandwidth method, a large damping ratio cannot be handled because the second and higher order terms of the damping ratio are ignored. Papagiannopoulos and Hatzigeorgiou [4] increased the upper limit of the damping ratio to ensure that the half-power bandwidth can be applied by incorporating higher order terms; however, only the loss factor near the resonant frequency could be obtained, which is a limitation. Wojtowicki et al. [5] noted that the frequency response function (FRF) of the displacement amplitude expressed as a function of a complex wavenumber contains all requisite information regarding viscoelasticity, and their method does not require the mode-dependent half-power bandwidth. The equivalent viscoelasticity of the Oberst beam was associated with the complex modulus of the damping layer [5]; however, the complex wavenumber and viscoelasticity were not explicitly related and numerical processing using the nonlinear least-squares method was required.For determining the viscoelasticity for a beam specimen conventionally, Young’s modulus and loss factor of the composite beam are estimated at a certain resonant frequency. In the Oberst beam method, the FRF is measured at a constant temperature in an environmental chamber and the viscoelasticity of the damping layer is determined in several incremental steps for each measurement under a constant temperature. Furthermore, the frequency dependence of the viscoelasticity is estimated at an ambient temperature of 20°C by applying a time-temperature superposition (TTS) to the results [5]. However, applying TTS to composite materials having complex configurations is challenging. As an alternative, the FRF obtained from the analytical solution is employed for the equation of motion of the viscoelastic beam over a wide frequency range. This could be applicable to the FRF used by Wojtowicki et al. [5]. However, they applied TTS to the viscoelasticity of polyvinyl chloride (PVC) obtained from the FRF of the Oberst beam near the resonance point for determining the frequency dependence of viscoelasticity at ambient temperature.The material viscoelasticity has been identified previously using the relationship between the complex wavenumber and viscoelasticity in applications other than vibration damping tests. Benatar et al. [6] derived an analytical solution of the Pochhammer frequency equation for complex wavenumbers in a split Hopkinson pressure bar for wavelengths comparable to its diameter of the bar. They identified viscoelasticity by analytically solving the Pochhammer frequency equation for small damping, while Zvietcovich et al. [7] determined it by optically measuring the responses of acoustically vibrating living tissues such as the cornea using Gaussian beam pulses. Owing to the complexity of the governing equations of wave propagation in these studies, the relationship between complex wavenumbers and frequencies was simplified using linearization. Conversely, the solution must be strictly managed to express the analytical solution of the flexural wave in the thin beam specimen in a simple form. We demonstrated that, for thin viscoelastic rods and beams (plates), a simple relationship between the viscoelasticity and complex wavenumber exists [8, 9]. Upon combining these with an analytical solution for the FRF (mobility (V/F)) between the velocity and excitation force, the viscoelasticity of the material can be determined over a wide frequency range from FRF data using complex wavenumbers [9].The remainder of this paper is organized as follows. In Section2, we outline the measurements using the center impedance method. In Section 3, the complex wavenumber for the flexural wave in the viscoelastic beam is expressed as an explicit function of viscoelasticity; the inverse formula for determining the viscoelasticity for a given complex wavenumber (real wavenumber and attenuation constant) is displayed along with the FRF mobility expressed as an analytical function of the complex wavenumber. Consequently, the complex wavenumber is obtained from the FRF by solving this complex function in reverse. The viscoelasticity can be determined at an arbitrary frequency in the measurement range via the complex wavenumber obtained from the measured FRF mobility data and the above analytical treatment [8, 9].The complex wavenumber and viscoelasticity (storage modulus and loss factor) were determined for arbitrary frequencies by inverse calculations using dimensionless mobility data obtained from preprocessing. The obtained viscoelasticity has a one-to-one correspondence with the measured FRF mobility; however, it considerably fluctuates with frequency because of the noise in the measured data. Thus, noise from the FRF was removed by fitting the relationship between the complex wavenumber and frequency with functions of mechanical significance wherein the real wavenumbers were treated as monotonously increasing frequency functions. The ratio of attenuation constant to real wavenumber was expressed as a function of the loss factor and fitted as that of the frequency using a regression curve proportional to the density function of the gamma distribution. A relationship similar to this curve can be observed in the viscoelastic model of Córtes and Elejabarrieta [10]. Finally, the FRF mobility was calculated using the curve-fitted complex wavenumber as compared with the measured results. The application of the inverse complex wavenumber formula on the curve-fitted complex wavenumber transformed the viscoelasticity (storage modulus and loss factor) to a smooth frequency function.


## 2. Experimental System

A schematic of the measurement system configuration is presented in Figure1; the beam center is vibrated with white noise. As per convention, the measurement system was used to evaluate the loss factor using the half-power bandwidth method and estimate the elastic modulus based on the resonant frequency. However, FRF was used differently. Although the FRF mobility (V/F) was calculated based on the measured acceleration and excitation force of the excitation rod, it was expressed as an analytical function of a complex wavenumber using the beam equation of motion, as discussed in Section 3.Figure 1
Experimental setup.


## 3. Identification of Viscoelasticity by Analytical Solutions

This section briefly summarizes the descriptions [8, 9] regarding the governing equations and analytical expression of FRF mobility.


### 3.1. Mechanical Response of a Viscoelastic Beam

#### 3.1.1. Governing Equations

The equation of motion of the free vibration of a viscoelastic beam can be expressed as(1)ρA∂2w∂t2+EI∗∂4w∂x4=0,where w, t, x, ρ, A=bh, and h denote the lateral displacement m, time s, location (m), mass density kg⋅m−3, cross-sectional area m2, and thickness of the beam (m), respectively; furthermore, EI∗=EI1+iη represents the complex flexural rigidity N⋅m2, in which η denotes the loss factor (–) (loss tangent for a uniform beam) and i=−1 represents the imaginary unit.The general solution of equation (1) can be written as(2)w=expiωtC1coshβ∗x+C2sinhβ∗x+C3cosβ∗x+C4sinβ∗x,(3)β∗4=ρAω2EI1+iη,where ω, Ckk=1,2,3,4, and β∗=βR−iα denote the radial frequency rad⋅s−1, constants (m), and complex wavenumber m−1, respectively. Variables βR and α represent the real wavenumber and attenuation constant (factor), respectively, which were derived from equation (3) and expressed as in [8].


#### 3.1.2. Complex Wavenumber Equations

Real wavenumber:(4)βR=ω2ρA8EI1/411+η2+11+η21/21/421+1/1+η21/21/2+11/2.Attenuation constant (factor):(5)α=ω2ρA8EI1/411+η2+11+η21/21/421+1/1+η21/21/2−11/2.


#### 3.1.3. Inverse Equations

The following inverse equations were obtained from equations (4) and (5) [9]:(6a)EI=ω2ρA8P,η=Q22−Q22−1,where(6b)Q=1+α/βR21−α/βR2,P=Q422−Q2Q+12βR4,under the conditions 0≤α/βR<2−1,1≤Q<2.Using the inverse equation (equation (6a)), the storage modulus E and loss factor η were obtained for a given complex wavenumber; equation (6b) denotes conditions for the wave solution (Appendix A).


### 3.2. FRF Mobility at Excitation Point

Using the solution of the equation of beam motion, the dimensionless FRF mobilityW can be expressed as [8, 9](7a)W=u+iv=−iρbhωL2Mob∗=ρbhωL2ImMob∗−iReMob∗,where(7b)Mob∗=−iωyd0Fx=L/2=iβ∗ρbhω1+1/coshβ∗L/2cosβ∗L/2tanβ∗L/2+tanhβ∗L/2.Mob∗ denotes the dimensional FRF mobility (velocity amplitude at the excitation point)/(amplitude of the excitation force) ms−1N−1, Fx=L/2=EI∗∂3yd/∂x3x=L/2 represents the shearing force on the rod (N), L denotes the beam length (m), and yd0 denotes the amplitude of rod displacement (m) [9].W defined using equation (7a) can be expressed as a function of the dimensionless complex wavenumber Z [9] as follows.


#### 3.2.1. Dimensionless Mobility Equation

(8a)W=Z1+1/coshZcosZtanZ+tanhZ.Equation (8a) is the dimensionless mobility equation, where(8b)Z=β∗L2=X+iY,(8c)X=βRL2,(8d)Y=−αL2,(8e)Ymin≤Y≤0,0≤−YX<2−1.Equations (8b)–(8d) are the dimensionless complex wavenumber, real wavenumber, and attenuation constant ×−1, respectively. Equation (8e) represents the constraints, where the first expression is an additional condition for excluding evanescent waves. Figure 2 portrays the behavior of the complex function W with Y as a parameter, where the following properties of W are portrayed as well [9].Figure 2
Dimensionless mobility in complexZ-plane. (a) Absolute value Wvs.X/π. (b) Phase angle argWvs.X/π. When X>1, the existence condition 0≤−Y/X<2−1 holds for Y<2−1.
(a)(b)


#### 3.2.2. Properties of Dimensionless Mobility (W)

(a)
If viscoelasticity is constant, the dimensionless real wavenumberX increases with frequency.(b)
TheX-axis is delimited with resonance points at intervals of approximately π.(c)
The phase angle ofW ranges from 0 to π (180°).(d)
The relationship between the phase angle of dimensional mobilityMob∗ and W is given by argMob∗=argW+π/2.(e)
The phase angle increases from a frequency slightly lower than the antiresonance point and starts to decline from that slightly lower than the resonance point. Although these alterations are stepped at small attenuation, they become gentle as the attenuation increases.When the FRF mobility is measured, the real wavenumberβR and attenuation constant α are obtained for each frequency via the inverse calculation of equation (8a) for Z. The ramp response method developed in a previous study was used for the inverse calculation [9] (Appendix B). In general, properties (b) and (e) are vital for examining the dimensionless mobility data. The real part of Z must be between the dimensionless real wavenumbers for the nR−1th and nRth resonance points. To ensure this, the information regarding the closest resonance point on the high-frequency side of the measurement frequency f designated as the nRth resonance point is used as the table value. Furthermore, viscoelasticity (E, η) was determined for each frequency by applying the inverse formula to the complex wavenumber data.


#### 3.2.3. Resonance

The dimensionless real wavenumber of thenRth resonance point is the nRth root of(9)tanXR+tanhXR=0.


#### 3.2.4. Antiresonance

The dimensionless real wavenumber of thenth antiresonance point represents the nth root of the following equation. In general, almost no vibration is seen near the exciting rod; however, the free end vibrates vigorously at the antiresonant frequency.(10)1+coshXn′cosXn′=0.


#### 3.2.5. Viscoelasticity by Inverse Analysis: For FRF Mobility without Noise

The numerical solution of equations (8a)–(8e) with respect to Z yields the complex wavenumber β∗=βR−iα vs. frequency relationship [9]. Moreover, the viscoelasticity E,η of the beam specimen was determined by applying the complex wavenumber to the inverse formula (equations (6a)–(6b)).We calculated the dimensionless mobilityW of the Oberst beam (Figure 3) based on the viscoelasticity of PVC measured by Wojtowicki et al. [5]. Thereafter, we regarded this as virtual measurement data and performed inverse calculations [9]. The obtained complex wavenumber was applied to equations (6a)–(6b) for obtaining the storage modulus and structural loss factor of the composite beam. As presented in Figure 4, the virtual measurement data (“measured”) for viscoelasticity were reproduced via inverse calculations (“estimated”). In their study [5], a 6 mm-thick and 400 mm-long PVC-based material was glued onto a 1.6 mm-thick aluminum substrate.Figure 3
Dimensionless mobility vs. frequency (Oberst beam). (a) Absolute value. (b) Phase angle.
(a)(b)Figure 4
Viscoelasticity by inverse calculation. (a) Storage modulus. (b) Loss factor.
(a)(b)


## 3.1. Mechanical Response of a Viscoelastic Beam

### 3.1.1. Governing Equations

The equation of motion of the free vibration of a viscoelastic beam can be expressed as(1)ρA∂2w∂t2+EI∗∂4w∂x4=0,where w, t, x, ρ, A=bh, and h denote the lateral displacement m, time s, location (m), mass density kg⋅m−3, cross-sectional area m2, and thickness of the beam (m), respectively; furthermore, EI∗=EI1+iη represents the complex flexural rigidity N⋅m2, in which η denotes the loss factor (–) (loss tangent for a uniform beam) and i=−1 represents the imaginary unit.The general solution of equation (1) can be written as(2)w=expiωtC1coshβ∗x+C2sinhβ∗x+C3cosβ∗x+C4sinβ∗x,(3)β∗4=ρAω2EI1+iη,where ω, Ckk=1,2,3,4, and β∗=βR−iα denote the radial frequency rad⋅s−1, constants (m), and complex wavenumber m−1, respectively. Variables βR and α represent the real wavenumber and attenuation constant (factor), respectively, which were derived from equation (3) and expressed as in [8].


### 3.1.2. Complex Wavenumber Equations

Real wavenumber:(4)βR=ω2ρA8EI1/411+η2+11+η21/21/421+1/1+η21/21/2+11/2.Attenuation constant (factor):(5)α=ω2ρA8EI1/411+η2+11+η21/21/421+1/1+η21/21/2−11/2.


### 3.1.3. Inverse Equations

The following inverse equations were obtained from equations (4) and (5) [9]:(6a)EI=ω2ρA8P,η=Q22−Q22−1,where(6b)Q=1+α/βR21−α/βR2,P=Q422−Q2Q+12βR4,under the conditions 0≤α/βR<2−1,1≤Q<2.Using the inverse equation (equation (6a)), the storage modulus E and loss factor η were obtained for a given complex wavenumber; equation (6b) denotes conditions for the wave solution (Appendix A).


## 3.1.1. Governing Equations

The equation of motion of the free vibration of a viscoelastic beam can be expressed as(1)ρA∂2w∂t2+EI∗∂4w∂x4=0,where w, t, x, ρ, A=bh, and h denote the lateral displacement m, time s, location (m), mass density kg⋅m−3, cross-sectional area m2, and thickness of the beam (m), respectively; furthermore, EI∗=EI1+iη represents the complex flexural rigidity N⋅m2, in which η denotes the loss factor (–) (loss tangent for a uniform beam) and i=−1 represents the imaginary unit.The general solution of equation (1) can be written as(2)w=expiωtC1coshβ∗x+C2sinhβ∗x+C3cosβ∗x+C4sinβ∗x,(3)β∗4=ρAω2EI1+iη,where ω, Ckk=1,2,3,4, and β∗=βR−iα denote the radial frequency rad⋅s−1, constants (m), and complex wavenumber m−1, respectively. Variables βR and α represent the real wavenumber and attenuation constant (factor), respectively, which were derived from equation (3) and expressed as in [8].


## 3.1.2. Complex Wavenumber Equations

Real wavenumber:(4)βR=ω2ρA8EI1/411+η2+11+η21/21/421+1/1+η21/21/2+11/2.Attenuation constant (factor):(5)α=ω2ρA8EI1/411+η2+11+η21/21/421+1/1+η21/21/2−11/2.


## 3.1.3. Inverse Equations

The following inverse equations were obtained from equations (4) and (5) [9]:(6a)EI=ω2ρA8P,η=Q22−Q22−1,where(6b)Q=1+α/βR21−α/βR2,P=Q422−Q2Q+12βR4,under the conditions 0≤α/βR<2−1,1≤Q<2.Using the inverse equation (equation (6a)), the storage modulus E and loss factor η were obtained for a given complex wavenumber; equation (6b) denotes conditions for the wave solution (Appendix A).


## 3.2. FRF Mobility at Excitation Point

Using the solution of the equation of beam motion, the dimensionless FRF mobilityW can be expressed as [8, 9](7a)W=u+iv=−iρbhωL2Mob∗=ρbhωL2ImMob∗−iReMob∗,where(7b)Mob∗=−iωyd0Fx=L/2=iβ∗ρbhω1+1/coshβ∗L/2cosβ∗L/2tanβ∗L/2+tanhβ∗L/2.Mob∗ denotes the dimensional FRF mobility (velocity amplitude at the excitation point)/(amplitude of the excitation force) ms−1N−1, Fx=L/2=EI∗∂3yd/∂x3x=L/2 represents the shearing force on the rod (N), L denotes the beam length (m), and yd0 denotes the amplitude of rod displacement (m) [9].W defined using equation (7a) can be expressed as a function of the dimensionless complex wavenumber Z [9] as follows.


### 3.2.1. Dimensionless Mobility Equation

(8a)W=Z1+1/coshZcosZtanZ+tanhZ.Equation (8a) is the dimensionless mobility equation, where(8b)Z=β∗L2=X+iY,(8c)X=βRL2,(8d)Y=−αL2,(8e)Ymin≤Y≤0,0≤−YX<2−1.Equations (8b)–(8d) are the dimensionless complex wavenumber, real wavenumber, and attenuation constant ×−1, respectively. Equation (8e) represents the constraints, where the first expression is an additional condition for excluding evanescent waves. Figure 2 portrays the behavior of the complex function W with Y as a parameter, where the following properties of W are portrayed as well [9].Figure 2
Dimensionless mobility in complexZ-plane. (a) Absolute value Wvs.X/π. (b) Phase angle argWvs.X/π. When X>1, the existence condition 0≤−Y/X<2−1 holds for Y<2−1.
(a)(b)


### 3.2.2. Properties of Dimensionless Mobility (W)

(a)
If viscoelasticity is constant, the dimensionless real wavenumberX increases with frequency.(b)
TheX-axis is delimited with resonance points at intervals of approximately π.(c)
The phase angle ofW ranges from 0 to π (180°).(d)
The relationship between the phase angle of dimensional mobilityMob∗ and W is given by argMob∗=argW+π/2.(e)
The phase angle increases from a frequency slightly lower than the antiresonance point and starts to decline from that slightly lower than the resonance point. Although these alterations are stepped at small attenuation, they become gentle as the attenuation increases.When the FRF mobility is measured, the real wavenumberβR and attenuation constant α are obtained for each frequency via the inverse calculation of equation (8a) for Z. The ramp response method developed in a previous study was used for the inverse calculation [9] (Appendix B). In general, properties (b) and (e) are vital for examining the dimensionless mobility data. The real part of Z must be between the dimensionless real wavenumbers for the nR−1th and nRth resonance points. To ensure this, the information regarding the closest resonance point on the high-frequency side of the measurement frequency f designated as the nRth resonance point is used as the table value. Furthermore, viscoelasticity (E, η) was determined for each frequency by applying the inverse formula to the complex wavenumber data.


### 3.2.3. Resonance

The dimensionless real wavenumber of thenRth resonance point is the nRth root of(9)tanXR+tanhXR=0.


### 3.2.4. Antiresonance

The dimensionless real wavenumber of thenth antiresonance point represents the nth root of the following equation. In general, almost no vibration is seen near the exciting rod; however, the free end vibrates vigorously at the antiresonant frequency.(10)1+coshXn′cosXn′=0.


### 3.2.5. Viscoelasticity by Inverse Analysis: For FRF Mobility without Noise

The numerical solution of equations (8a)–(8e) with respect to Z yields the complex wavenumber β∗=βR−iα vs. frequency relationship [9]. Moreover, the viscoelasticity E,η of the beam specimen was determined by applying the complex wavenumber to the inverse formula (equations (6a)–(6b)).We calculated the dimensionless mobilityW of the Oberst beam (Figure 3) based on the viscoelasticity of PVC measured by Wojtowicki et al. [5]. Thereafter, we regarded this as virtual measurement data and performed inverse calculations [9]. The obtained complex wavenumber was applied to equations (6a)–(6b) for obtaining the storage modulus and structural loss factor of the composite beam. As presented in Figure 4, the virtual measurement data (“measured”) for viscoelasticity were reproduced via inverse calculations (“estimated”). In their study [5], a 6 mm-thick and 400 mm-long PVC-based material was glued onto a 1.6 mm-thick aluminum substrate.Figure 3
Dimensionless mobility vs. frequency (Oberst beam). (a) Absolute value. (b) Phase angle.
(a)(b)Figure 4
Viscoelasticity by inverse calculation. (a) Storage modulus. (b) Loss factor.
(a)(b)


## 3.2.1. Dimensionless Mobility Equation

(8a)W=Z1+1/coshZcosZtanZ+tanhZ.Equation (8a) is the dimensionless mobility equation, where(8b)Z=β∗L2=X+iY,(8c)X=βRL2,(8d)Y=−αL2,(8e)Ymin≤Y≤0,0≤−YX<2−1.Equations (8b)–(8d) are the dimensionless complex wavenumber, real wavenumber, and attenuation constant ×−1, respectively. Equation (8e) represents the constraints, where the first expression is an additional condition for excluding evanescent waves. Figure 2 portrays the behavior of the complex function W with Y as a parameter, where the following properties of W are portrayed as well [9].Figure 2
Dimensionless mobility in complexZ-plane. (a) Absolute value Wvs.X/π. (b) Phase angle argWvs.X/π. When X>1, the existence condition 0≤−Y/X<2−1 holds for Y<2−1.
(a)(b)


## 3.2.2. Properties of Dimensionless Mobility (W)

(a)
If viscoelasticity is constant, the dimensionless real wavenumberX increases with frequency.(b)
TheX-axis is delimited with resonance points at intervals of approximately π.(c)
The phase angle ofW ranges from 0 to π (180°).(d)
The relationship between the phase angle of dimensional mobilityMob∗ and W is given by argMob∗=argW+π/2.(e)
The phase angle increases from a frequency slightly lower than the antiresonance point and starts to decline from that slightly lower than the resonance point. Although these alterations are stepped at small attenuation, they become gentle as the attenuation increases.When the FRF mobility is measured, the real wavenumberβR and attenuation constant α are obtained for each frequency via the inverse calculation of equation (8a) for Z. The ramp response method developed in a previous study was used for the inverse calculation [9] (Appendix B). In general, properties (b) and (e) are vital for examining the dimensionless mobility data. The real part of Z must be between the dimensionless real wavenumbers for the nR−1th and nRth resonance points. To ensure this, the information regarding the closest resonance point on the high-frequency side of the measurement frequency f designated as the nRth resonance point is used as the table value. Furthermore, viscoelasticity (E, η) was determined for each frequency by applying the inverse formula to the complex wavenumber data.


## 3.2.3. Resonance

The dimensionless real wavenumber of thenRth resonance point is the nRth root of(9)tanXR+tanhXR=0.


## 3.2.4. Antiresonance

The dimensionless real wavenumber of thenth antiresonance point represents the nth root of the following equation. In general, almost no vibration is seen near the exciting rod; however, the free end vibrates vigorously at the antiresonant frequency.(10)1+coshXn′cosXn′=0.


## 3.2.5. Viscoelasticity by Inverse Analysis: For FRF Mobility without Noise

The numerical solution of equations (8a)–(8e) with respect to Z yields the complex wavenumber β∗=βR−iα vs. frequency relationship [9]. Moreover, the viscoelasticity E,η of the beam specimen was determined by applying the complex wavenumber to the inverse formula (equations (6a)–(6b)).We calculated the dimensionless mobilityW of the Oberst beam (Figure 3) based on the viscoelasticity of PVC measured by Wojtowicki et al. [5]. Thereafter, we regarded this as virtual measurement data and performed inverse calculations [9]. The obtained complex wavenumber was applied to equations (6a)–(6b) for obtaining the storage modulus and structural loss factor of the composite beam. As presented in Figure 4, the virtual measurement data (“measured”) for viscoelasticity were reproduced via inverse calculations (“estimated”). In their study [5], a 6 mm-thick and 400 mm-long PVC-based material was glued onto a 1.6 mm-thick aluminum substrate.Figure 3
Dimensionless mobility vs. frequency (Oberst beam). (a) Absolute value. (b) Phase angle.
(a)(b)Figure 4
Viscoelasticity by inverse calculation. (a) Storage modulus. (b) Loss factor.
(a)(b)


## 4. Results and Discussion: Actual Measurements and Data Processing

The measured FRF mobility (V/F) data contained noise and were not smooth, as illustrated in Figure3. To develop a denoising method as discussed, inverse calculation was applied to these data, which resulted in unrealistic viscoelasticities with large fluctuations of frequency.


### 4.1. Data Processing Flow and Noise Removal

Figure5 depicts the data processing flow that comprises (a) preprocessing, (b) inverse calculation, and (c) postprocessing. In the first step, measured data were normalized and converted to dimensionless mobility values W. Second, the complex wavenumber was obtained using equations (8a)–(8e) with respect to Z. Subsequently, the viscoelasticity was identified from the complex wavenumber using equations (6a) and (6b). If W contained noise, the relationship between the estimated viscoelasticity and frequency becomes unsmooth and unrealistic. Finally, the relationship between the complex wavenumber and frequency was curve fitted with functions of mechanical significance instead of directly eliminating noise from the viscoelastic data. Moreover, equations (8a)–(8e) were applied to the curve-fitted complex wavenumber to obtain the dimensionless mobility Wfit as a smooth frequency function. However, the smooth relationship between viscoelasticity Efit,ηfit and frequency can be found by applying equations (6a) and (6b) to the curve-fitted complex wavenumber.Figure 5
Data processing flow.


### 4.2. Preprocessing (1)

Data points with a coherence of ≤0.99 were first removed from the measured FRF mobility data. In addition, the term “coherence” is defined and summarized in the following section. Furthermore, the phase angle range of the mobilityMob∗ was normalized within 90°–270° [9]. The preprocessed FRF mobility data are presented in Figure 6, wherein the phase angle data include spike-like noise represented by dashed circles. Similar to the properties of W, the phase angle increased from a frequency slightly under the antiresonance point and started declining before the resonance point. A few phase angle data that cannot be normalized within the range of 90°–270° are present, which were subjected to the following equations:(11a)if phase Mob<90°,then phase Mob⟵90°+εL,(11b)if phase Mob>270°,then phase Mob⟵270°−εH.Figure 6
Measured mobilities (materials: ABS and PP). (a) Material: ABS (acrylonitrile-butadiene-styrene copolymer). (b) Material: PP (polypropylene).
(a)(b)In these cases,εL,εH>0 represent small positive numbers and phase Mob denotes the phase angle of dimensional mobility °.In Figure6, εL=εH=4°.


#### 4.2.1. Definition of Coherence

The term “coherence” refers to the output value of the measurement system and is defined ascohf=SVFf/SVVfSFFf1/2, where f denotes the frequency (Hz) and SVVf and SFFf represent the power spectra of the excitation point velocity and excitation force, respectively. Furthermore, SVFf represents the cross spectrum between the excitation point velocity and excitation force. The units of SVVfSFFf1/2 and SVFf are identical to each other such that coherence is a dimensionless number with values in the range of 0-1. Coherence has a meaning similar to the correlation coefficient in the present context. For a small coherence, the velocity of the excitation point is considered to be noise independent of the excitation force.


### 4.3. Preprocessing (2): Gain Correction

The FRF mobility data are influenced by the gain of the measurement system, which is estimated in equation (12). Parameters f0 and kG in equation (12) were determined by employing the least-squares method using the analytical solution and measured mobility data for an aluminum beam.(12)Gestf=ff0−kG,where f denotes the frequency (Hz), f0 (>0) is a constant (Hz), and kG (>0) is a constant (–).Using this gain, the measured dimensional mobilityMobM∗ at f was corrected using(13a)Mob∗=MobM∗Gestf,(13b)20log10Mob∗=20log10MobM∗−20log10Gestf,where equation (13b) represents the absolute mobility in decibels.


#### 4.3.1. Determination of Gain Parameters by Measuring an Aluminum Beam

Assuming that the gain of the measurement system is independent of the beam material, the parameters were estimated using the measured results of the aluminum beam. Consequently,f0=10HzandkG=0.15− were obtained. The impact of gain correction on the absolute value of the dimensionless mobility data of the aluminum beam is presented in Figure 7. For appropriate gain correction, a dimensional mobility of 0 dB corresponds to 1ms−1N−1. In the least-squares method, the residual Sres was minimized:(14)Sres=∑nuCn−uMn/Gestn2+vCn−vMn/Gestn2∑nuCn2+vCn2,where u=ReW,v=ImW, and C, M, and n represent calculated, measured, and n-th frequency, respectively.Figure 7
Dimensionless mobility (absolute value): measured and corrected values compared with analytical solution (aluminum beam).


### 4.4. Preprocessing (3): Rendering FRF Mobility Dimensionless

As mobilities after preprocessing steps (1) and (2) are displayed in decibels, the absolute valueMob∗ms−1N−1, real part ReMob∗, and imaginary part ImMob∗ of the dimensional mobility can be described using(15a)Mob∗=10dBMob/20,(15b)ReMob∗=Mob∗cosphase Mob,(15c)ImMob∗=Mob∗sinphase Mob,where phaseMob denotes the phase angle of the dimensional mobility °.


#### 4.4.1. Dimensionless Mobility

The measured dataMob∗ of the FRF mobility can be converted to dimensionless mobility W using(16)W=−iρbhωL2Mob∗=ρbhωL2ImMob∗−iReMob∗,i=−1,where ρbhL denotes the mass of the beam specimen. The phase angle arg(W) of the dimensionless mobility portrayed in Figure 8 ranges from 0 to π.Figure 8
Dimensionless mobility (materials: ABS and PP, length: 127 mm, and thickness: 1.6 mm). (a) Material: ABS. (b) Material: PP.
(a)(b)


### 4.5. Inverse Calculation

The application of inverse calculations expressed in equations (8a)–(8e) to W in Figure 8 yielded the relationship between the dimensionless complex wavenumber and frequency, as indicated in Figure 9. The thick lines represent the dimensionless real wavenumber X, while the thin lines represent the dimensionless attenuation constant (−Y). A drastically steep rise can be observed in the dimensionless real wavenumber X near 5.5 kHz in Figure 9(a) (ABS) and 6 kHz in Figure 9(b) (PP). However, the dimensionless attenuation constant alternated between small and large values.Figure 9
Dimensionless complex wavenumber (real wavenumber and attenuation constant). (a) Material: ABS. (b) Material: PP.
(a)(b)The storage modulus and loss factor were obtained by applying equations (6a) and (6b) to the real wavenumber βR=2X/L and attenuation constant α=−2Y/L obtained from the results in Figure 9. These results are presented in Figure 10, where the elastic modulus of an elastic body with respect to the resonant frequency is “purely elastic.” In addition, both the storage modulus and loss factor fluctuated significantly with the frequency, which is an anomaly for viscoelastic data. Moreover, the loss factor near the resonant frequency obtained using the half-power bandwidth method was calculated using the measurement system and presented as the “half-power bandwidth.”Figure 10
Viscoelasticity estimated from complex wavenumbers in Figure9. (a) Material: ABS. (b) Material: PP.
(a)(b)


### 4.6. Postprocessing: Noise Removal by Fitting Complex Wavenumber vs. Frequency Relationship to a Theoretical Curve

Instead of directly eliminating the noise from the viscoelasticity vs. frequency data, we developed a method to fit a curve for the relationship between the complex wavenumber and frequency using smooth functions of intermediate data.


#### 4.6.1. Real Wavenumber vs. Frequency Relationship

The generalized Maxwell model is applicable to storage modulusE of materials, which is a monotonically increasing function of frequency. As ω ⟶ ∞ at high frequencies, E approaches the instantaneous modulus Einst, and as ω ⟶ 0 at low frequencies, it approaches the long-term modulus E∞. If η has a small frequency dependence, equation (4) can be considered for the following relationship between the real wavenumber βR and radial frequency ω.(17)βR∼ω1/2asω⟶∞,(18a)∂lnβR∂lnω=∂lnX∂lnf=12−14∂lnE∂lnω=12−14ωE∂E∂ω⟶kf0,0<kf0<12asω⟶+0,0<kf0<∂lnβR∂lnω≤12,(18b)∂E∂ω>0,E∞≤E≤Einst,for0<ω<∞.Therefore, if the relationship oflnβRvs.lnω is plotted, the slope at low frequencies is <1/2 and approaches 1/2 at high frequencies. Thus, we can determine a curve-fitting equation such that the slope satisfies equations (18a) and (18b).


#### 4.6.2. Equation to Fit the Relationship between Real Wavenumber vs. Frequency

Slope from Reference Point. The plot of lnβRvs.lnω has many irregularities. Instead of the gradient ∂lnβR/∂lnω at each frequency, we define the slope kfω of the straight line connecting the points lnω,lnβR and lnω0,lnβ0 from the reference point using the following equations:(19a)lnβR=lnβ0+kfωlnωω0,where(19b)kfω=lnβR−lnβ0lnω−lnω0=lnβR/β0lnω/ω0,forω≠ω0,where ω0 and β0 denote the reference radial frequency s−1 and reference wavenumber (m−1), respectively; furthermore, βR=β0forω=ω0.(19c)β04=ρAω02E0I=12ρω02E0h2,where E0 denotes the storage modulus at ω0.From equations (18a) and (18b), the slope kfω satisfies(20)0<kf0<kfω<12,kfω⟶12,asω⟶∞.Furthermore,ω0 was(A)
Located near the center of the frequency range of interest and(B)
Situated far from the peaks and valleys due to resonance points, antiresonance points, and noise of the dimensionless mobilityWThe slopekfω does not always satisfy equation (20) due to variations in the real wavenumber βR from the inverse calculation on the FRF data. Therefore, the following conditions were used to satisfy equation (20):(21)ifkfω>12,thenkfω⟵12,ifkfω<kmin,thenkfω⟵kminkmindenotes a small positive number.


#### 4.6.3. Curve-Fit Formula for Slopekfω

kfω obtained using the sigmoid function, a smooth function with a value between 0 and 1, was adopted as the curve-fitting equation satisfying equation (20).(22a)kfω=12ysxFd,where(22b)xFd=lnωωI=xF−xFI,(22c)ysxFd=11+exp−xFd,where xFd denotes the difference in the logarithm of frequency (logarithm of frequency ratio) (–), ωI denotes the radial frequency of the inflection point (s−1), and ysxFd represents the sigmoid function of xFd.(22d)xF=lnω,xFI=lnωI.The method for identifyingωI is described in Appendix C. Evaluating the sigmoid function is not necessary once ωI has been calculated. Moreover, equations (22a)–(22d) can be simplified as(23)kfω=12ωω+ωI.Applying the slope of equation (23) to equation (19a), the real wavenumber βR becomes a smooth function of the frequency. Figure 11 presents examples depicting the relationship between logarithms of the dimensionless real wavenumber and frequency. The real wavenumber βR was converted to a dimensionless real wavenumber X=βRL/2. In Figure 11, the points at which the thick line deviated from the regression curve were considered to be affected by noise. Vertical and horizontal axes of Figure 11 were rewritten in Figure 12 in linear forms. Here, for both ABS and PP, the thick and dashed lines were in good agreement for frequencies of 4 kHz or less.Figure 11
Curve fitting of dimensionless real wavenumber using equations (24a) and (24b). (a) Material: ABS. Reference point: f0=2243Hz,E0=2.560GHz,β0=133.4m−1,andX0=8.402−. Curve-fitting parameter: ωI=11.27s−1andfI=ωI/2π=1.794Hz. (b) Material: PP. Reference point: f0=2229Hz,E0=2.741GHz,β0=130.7m−1,andX0=8.234−. Curve-fitting parameter: ωI=13.10s−1andfI=ωI/2π=2.085Hz.
(a)(b)Figure 12
Dimensionless real wavenumber vs. frequency (linear scale). (a) Material: ABS. (b) Material: PP.
(a)(b)


#### 4.6.4. Equation to Fit Relationship between Attenuation Constant and Frequency

In Figure10, the loss factor obtained by the inverse calculation from FRF data fluctuates in the range from 0 to approximately 0.3. However, it must satisfy conditions (A), (B), and (C). In addition, the ratio α/βR of the attenuation constant to real wavenumber was fitted to a theoretical curve based on these conditions. Thereafter, E and η were calculated from the values of βR and α using equations (6a) and (6b).(A) Relationship between Ratioα/βRand Loss Factorη. From the complex wavenumber equations (4) and (5), α/βR can be expressed as a function of only the loss factor η, which is expressed as(24a)αβR=−YX=2/1+1+η2−1/21/2−12/1+1+η2−1/21/2+11/2,(24b)αβR=−YX≈η4,asη⟶0.(B) Frequency at Which Ratioα/βRand Loss FactorηAre Maximized. The following relationship holds from equations (6a) and (6b):(25a)∂η∂ω=dηdQ∂Q∂ω=dηdQ2α/βR1−α/βR22⋅∂∂ωαβR,where(25b)Q=1+α/βR21−α/βR2,dηdQ=Q22−Q22−1−1/2×4Q32−Q3>0,1≤Q<2.Equations (25a) and (25b) maximize α/βR at the frequency corresponding to the maximization of the loss factor η.(C) For Many Viscoelastic Materials, the Loss Factor η Has the Maximum Value at a Single Frequency. From properties (A)–(C), the relationship of the ratio α/βR=−Y/Xvs.ω should be fitted with an appropriate curve to consider the maximum value at a particular ω. In this case, the loss factor in Figure 10 was replaced by a smooth curve with only one maximal peak.


#### 4.6.5. Regression Curve Proportional to Gamma Distribution

Figure13 shows examples of the relationship of α/βR and frequency for both materials, which shows a large fluctuation of α/βR. Therefore, we considered the horizontal axis (frequency (Hz)) as a random variable, and the value on the vertical axis was considered as the number of counts in the histogram. Thus, the curve fit of the data in Figure 13 can be replaced with a case of an appropriate probability distribution applied to the experimental data represented by a histogram.Figure 13
Ratio−Y/X=α/βR vs. frequency. (a) Material: ABS. (b) Material: PP.
(a)(b)Therefore, a function proportional to the density function of the gamma distribution, a continuous distribution similar to the Poisson distribution, was used as the curve-fit function. The proposed equation for the regression curve ofα/βR is presented in equations (26a)–(26e), wherein the normalization by the step size r0 divided into 16 intervals between the minimum and maximum α/βR was adopted in equations (26a)–(26e). The plots obtained using this normalization are presented in Figure 14.(26a)r0nF=−YX=r0n0RfΓpσpgf,where(26b)gf=1Γpσpfp−1exp−fσ,f>0,(26c)Γp=∫0∞xp−1e−xdx,p>0,where gf denotes the probability density of the gamma distribution and Γp denotes the gamma function.(26d)nF=α/βRr0,(26e)r0=rmax−rmin16,where rminandrmax represent the minimum and maximum values of α/βR in the frequency range, respectively; r0 denotes the step width (–) and f=ω/2π denotes the frequency (Hz); σ>0andp>0− are the parameters in the gamma distribution; and n0Rf represents the coefficient for the height of the distribution HZ−p, as defined in(27)r0n0RfσpΓp=r0∫0∞nFdf.Figure 14
Curve fitting ofnF. (a) Material: ABS. (b) Material: PP.
(a)(b)Considering the thin line in Figure14 as a histogram with frequency f as a random variable, the sample mean mf and sample standard deviation sf with respect to frequency f can be defined. Parameters p and σ of the gamma distribution in equation (26b) can be expressed by equations related to mean mf and standard deviation sf. The coefficient n0Rf was determined using equation (27) such that the integral value of the curve-fit equation and plot data is equal.


#### 4.6.6. Determining Gamma Distribution Parameters

For the probability distribution function expressed by equation (26b), the population mean Ef and population variance Vf are given by(28)Ef=pσ,Vf=pσ2.If these are equal to the sample meanmf and sample variance sf2, parameters p and σ are given by(29a)p=mfsf2,σ=sf2mf.The sample meanmf and sample standard deviation sf were calculated using equation (29b). In equation (29b), the region where the convergence points of the inverse calculation (markers in Figure 14) are densely present has large weights.(29b)mf=∑knFkfkN,sf2=∑knFkfk−mf2N,where(29c)N=∑knFk.Here,N is similar to the total count but is generally a noninteger real number. The inverse calculation of equation (8a) was performed for 635 frequency points from 30 to 6400 Hz with 10 Hz increments. The number of points where equation (8a) converged was 490 for ABS and 543 for PP, which are represented by markers in Figures 13 and 14. Plot point numbers k are assigned to these convergence points in the ascending order of frequency.


#### 4.6.7. Curve-Fit for Ratioα/βR vs. Frequency

The thin lines in Figure14 show the normalized plots of ratio α/βR for ABS and PP. The thick lines are regression curves proportional to the gamma distributions, and the curve-fit parameters are shown in Figure 14. The gamma distribution parameters p are 3.4 and 4.1 for ABS and PP, respectively. The anomalous jumps are circled in Figure 14(a), which appear near the spike-like noises of the phase angles in Figures 6 and 8. Although the thick lines did not apparently correspond with the thin lines well, the curve-fit result was evaluated by the degree to which the original dimensionless mobility could be expressed. The reason why the thin lines (value by direct inverse calculation from FRF data) vary is discussed in Appendix D.


#### 4.6.8. Curve-Fitting Parameter

ABS:p=3.430−,σ=918.6Hz,mf=3151Hz,sf=1701Hz,Γp=3.078−,r0=5.534×10−3−,andn0Rf=1.071×10−7Hz−p.PP:p=4.130−,σ=809.3Hz,mf=3342Hz,sf=1645Hz,Γp=7.080−,r0=2.881×10−3−,andn0Rf=4.825×10−9Hz−p.


#### 4.6.9. Viscoelastic Model with Frequency Dependence Whose Loss Factor Is Proportional to the Gamma Distribution

Córtes and Elejabarrieta [10] derived a viscoelastic model represented by equations (30a) and (30b) via a theoretical study of the relaxation function. They experimentally validated this model using a material in which concrete was impregnated with epoxy resin. Equation (30b) is the case where p = 2 in equations (26a)–(26e). The constant η0 is 2×10−2, and therefore, the approximation of equation (24b) holds for such a small loss factor. Equation (30b) is a special case of our curve-fit equations (26a)–(26e).(30a)E′ω=E0=Const,(30b)ηω=η0ωω0exp1−ωω0,where E′ω denotes the storage modulus Pa and ηω denotes the loss factor (–); furthermore, η0=2cω0/E0e, in which c denotes the viscosity Pa⋅s and e represents the base of the natural logarithm.


#### 4.6.10. FRF Mobility Based on Curve-Fitted Complex Wavenumber Data

We applied the curve-fitted data of the dimensionless complex wavenumberZ=X+iY to equations (8a)–(8e) to obtain the dimensionless mobility. The results were then compared with the original FRF data (see Figures 15 and 16). Furthermore, in these figures, the standard deviation StdReg of the difference between the thick and thin gray lines was calculated with respect to the absolute value and phase angle of W. These are shown as (thick line) ± StdReg values via dashed lines. However, (thick line)−StdReg has a small lower limit to prevent the dashed line from being negative.Figure 15
FRF mobility calculated from the curve-fitted dimensionless complex wavenumber (material: ABS). (a) Absolute value vs. frequency. (b) Phase angle vs. frequency.
(a)(b)Figure 16
FRF mobility calculated from the curve-fitted dimensionless complex wavenumber (material: PP). (a) Absolute value vs. frequency. (b) Phase angle vs. frequency.
(a)(b)In these figures, the thin lines are marked with the convergence points of the inverse calculations of equations (8a)–(8e), which are densely present in the frequency range between the antiresonance and resonance points. Here, the dimensionless mobilities based on the measured data (gray lines) and FRF based on the curve-fit of the complex wavenumber (thick lines) agree well. However, in the frequency range that crosses the resonance points, the convergence points of the inverse calculation of equations (8a)–(8e) are sparse, and the deviation between the thick and thin lines is conspicuous. In Figures 15 and 16, the thin gray lines fall within the range of (thick line) ± StdReg in most frequency ranges.


#### 4.6.11. Viscoelasticity with Curve-Fitted Complex Wavenumbers

We converted the curve-fitted dimensionless complex wavenumbers to dimensional complex wavenumbersβ∗=βR−iα,whereβR=2X/Landα=−2Y/L. The results of determining the viscoelasticity by applying equations (6a) and (6b) to these complex wavenumbers are shown in Figures 17 and 18. In these figures, the thin gray lines are the results estimated directly from the FRF data, while the thick lines are based on the curve-fitted complex wavenumbers. Furthermore, the standard deviations of the thick and thin lines for the storage modulus and loss factor were calculated as StdReg. In Figures 17 and 18, the values of (thick line) ± StdReg are drawn as dashed lines. However, (thick line)−StdReg has a small lower limit to prevent negative values. In Figures 17(b) and 18(b), the thick lines have maximum values near 2 kHz and 2.5 kHz, respectively.Figure 17
Storage modulus and loss factor calculated from the curve-fitted dimensionless complex wavenumber (material: ABS). (a) Storage modulus vs. frequency. (b) Loss factor vs. frequency.
(a)(b)Figure 18
Storage modulus and loss factor calculated from the curve-fitted dimensionless complex wavenumber (material: PP). (a) Storage modulus vs. frequency. (b) Loss factor vs. frequency. yReg: calculated from the curve-fitted data (regression). StdReg: standard deviation of the residual.
(a)(b)


#### 4.6.12. Abnormal Data

In Figure17(b), an alternate long and short dashed line is shown near the resonance point at 5.6 kHz. Near this frequency, the absolute value of mobility in Figure 6(a) rose abnormally from the minimum value. The circles near 1 kHz in Figures 17(a) and 17(b) correspond to the spike-like noise of the phase angles in Figure 6(a).


### 4.7. Effect of Resonance/Antiresonance Points on Fluctuations of Storage Modulus and Loss Factor

The resonant and antiresonant frequencies of the beam specimen are indicated by markers on the frequency axes of Figures17 and 18. Here, the estimated storage modulus and loss factor measured via direct inverse calculations from the FRF data (thin lines) are high near the antiresonant frequency and low near the resonant frequency. The method of determining viscoelasticity based on the complex wavenumber is not originally affected by modes such as resonance and antiresonance. However, the inverse calculations for obtaining the complex wavenumbers from the FRF mobility converge well in the frequency range from the antiresonance to resonance point, whereas the number of convergence points becomes sparse in a frequency range exceeding the resonance point. In Figures 6 and 8 for the FRF data, the falling edge of the phase angle near the resonant frequency represents a sudden change close to a step. According to the qualitative study of equations in Appendix D, the inverse calculation converges poorly because dW/dZ∝1/η2 holds true near the resonance point. However, near the antiresonant frequency, noise may occur due to mixing of natural vibrations, which is concerning.


### 4.8. Estimation of the Magnitude of Noise Contained in the Dimensionless Complex Wavenumber and FRFW

#### 4.8.1. Estimating the Magnitude of Noise Contained in the Mobility Measurement Data

The dimensionless mobilityWfit obtained by applying the curve-fitted dimensionless complex wavenumber to equations (8a)–(8e) has a value different from the original measurement data WI as shown in Figures 15 and 16. Figures 19 and 20 show the real and imaginary parts of Wfit=ufit+ivfit and WI=uI+ivI, respectively. Here, we consider the difference ΔW=WI−Wfit as noise in the FRF measurement data.(31)ΔW=WI−Wfit=uI−ufit+ivI−vfit,where uandv denote the real and imaginary parts of the dimensionless mobility W(−), respectively. Subscripts I and fit denote direct inverse calculation value from the measurement and the curve-fitted complex wavenumber, respectively.Figure 19
Comparison of dimensionless mobility (material: ABS). (a) Re(W) vs. frequency. (b) Im(W) vs. frequency.
(a)(b)Figure 20
Comparison of dimensionless mobility (material: PP). (a) Re(W) vs. frequency. (b) Im(W) vs. frequency.
(a)(b)


#### 4.8.2. Noise in Dimensionless Complex Wavenumber

Consider changes in the dimensionless real wavenumberX and dimensionless attenuation constant −Y when the FRF W changes by ∆W. Since W expressed by equation (8a) is differentiable at points other than at the resonance point, the following relationship holds:(32)dWdZ=∂W∂X=1i∂W∂Y.Therefore, the following relationships among the change ∆W of W, that of ∆X of X, and that of ∆Y of Y hold:(33)ΔW=∂W∂XΔX+∂W∂YΔY=dWdZΔX+iΔY.Therefore,(34a)ΔX=ReΔWdW/dZ,ΔY=ImΔWdW/dZ,where(34b)ΔW=Δu+iΔv,Δu=uI−ufit,Δv=vI−vfit.Equations (34a) and (34b) show that the fluctuation of the FRF affects both the real and imaginary parts of the complex wavenumber. On the contrary, the difference between the nondimensional complex wave number calculated directly from the FRF W measurement data and the curve-fitted nondimensional complex wave number can be expressed by(35a)δX=XI−Xfit,δY=YI−Yfit,(35b)ZI=XI+iYI,Zfit=Xfit+iYfit.Let us considerδX or ∆X as the fluctuation of X and δY or ∆Y as the fluctuation of Y.SlopedW/dZ. Figure 21 shows the gradient dW/dZ when Z=Xfit+iYfit is applied to equations (D.7a) and (D.7b). To make the figure clear, the vertical axis of Figure 21(a) is cut off by ±200.Figure 21
SlopedW/dZ vs. frequency. (a) Material: ABS. (b) Material: PP.
(a)(b)Fluctuation of Dimensionless Complex Wavenumber. ΔX/Xfit and δX/Xfit are considered the degree of noise in the data of real wavenumber X. From Figures 22(a) and 23(a), ΔX/Xfit and δX/Xfit are in the range of about 5–7% in most frequency regions. However, δX/Xfit increases sharply above 5.5 kHz in Figures 22(a) and 23(a). This agrees with the fact that the relation of dimensionless real wavenumber X vs. the frequency relationship near 5.5 kHz in Figures 12(a) and 12(b) is not smooth and corresponds to the spike-like noise of the phase angle near 5.5 kHz in Figures 8 and 9. For the fluctuation of the dimensionless attenuation constant −Y, ΔY/Xfit and δY/Xfit were regarded as indicators of the noise magnitude. In Figures 22(b) and 23(b), the ratio of attenuation constant to real wavenumber −Yfit/Xfit is shown by dashed lines, which are the envelopes above ΔY/Xfit and δY/Xfit. For the measurement data of this study, the dimensionless attenuation constant −YI obtained by direct inverse calculation from the FRF measurement data was small −YI≈0, where the frequency was slightly higher than the resonance point. Therefore, δY≈−Yfit≥0 in this frequency range. From Figures 22(b) and 23(b), the noise magnitudes in the attenuation constant ΔY/XfitandδY/Xfit are in the range of approximately 2–4% in most frequency regions below 5.5 kHz.Figure 22
Fluctuation of complex wavenumber (material: ABS). (a) Fluctuation of real wavenumber. (b) Fluctuation of attenuation constant.
(a)(b)Figure 23
Fluctuation of complex wavenumber (material: PP). (a) Fluctuation of real wavenumber. (b) Fluctuation of attenuation constant.
(a)(b)


### 4.9. Effect of Phase Angle Correction on Loss Factor Estimates

The measurement data (thin line) of the phase angle ofW in Figure 16(b) were replaced with data based on the curve-fitted complex wavenumber, represented by a dashed line in Figure 24(a). Nonetheless, the absolute value of W remained on the thick line in Figure 24(a). The loss factor estimated from the inverse calculation at this instant is depicted in Figure 24(b), wherein the fluctuation of the thin line is smaller than that in Figure 18(b). Moreover, the phase angle measurement data of FRF W on the high frequency side of the resonance point significantly influences the fluctuation of the estimated loss factor.Figure 24
Influence of phase angle on loss factor estimates (material: PP). (a) Revision of the phase angle. (b) Effect on the loss factor estimates.
(a)(b)


## 4.1. Data Processing Flow and Noise Removal

Figure5 depicts the data processing flow that comprises (a) preprocessing, (b) inverse calculation, and (c) postprocessing. In the first step, measured data were normalized and converted to dimensionless mobility values W. Second, the complex wavenumber was obtained using equations (8a)–(8e) with respect to Z. Subsequently, the viscoelasticity was identified from the complex wavenumber using equations (6a) and (6b). If W contained noise, the relationship between the estimated viscoelasticity and frequency becomes unsmooth and unrealistic. Finally, the relationship between the complex wavenumber and frequency was curve fitted with functions of mechanical significance instead of directly eliminating noise from the viscoelastic data. Moreover, equations (8a)–(8e) were applied to the curve-fitted complex wavenumber to obtain the dimensionless mobility Wfit as a smooth frequency function. However, the smooth relationship between viscoelasticity Efit,ηfit and frequency can be found by applying equations (6a) and (6b) to the curve-fitted complex wavenumber.Figure 5
Data processing flow.


## 4.2. Preprocessing (1)

Data points with a coherence of ≤0.99 were first removed from the measured FRF mobility data. In addition, the term “coherence” is defined and summarized in the following section. Furthermore, the phase angle range of the mobilityMob∗ was normalized within 90°–270° [9]. The preprocessed FRF mobility data are presented in Figure 6, wherein the phase angle data include spike-like noise represented by dashed circles. Similar to the properties of W, the phase angle increased from a frequency slightly under the antiresonance point and started declining before the resonance point. A few phase angle data that cannot be normalized within the range of 90°–270° are present, which were subjected to the following equations:(11a)if phase Mob<90°,then phase Mob⟵90°+εL,(11b)if phase Mob>270°,then phase Mob⟵270°−εH.Figure 6
Measured mobilities (materials: ABS and PP). (a) Material: ABS (acrylonitrile-butadiene-styrene copolymer). (b) Material: PP (polypropylene).
(a)(b)In these cases,εL,εH>0 represent small positive numbers and phase Mob denotes the phase angle of dimensional mobility °.In Figure6, εL=εH=4°.


### 4.2.1. Definition of Coherence

The term “coherence” refers to the output value of the measurement system and is defined ascohf=SVFf/SVVfSFFf1/2, where f denotes the frequency (Hz) and SVVf and SFFf represent the power spectra of the excitation point velocity and excitation force, respectively. Furthermore, SVFf represents the cross spectrum between the excitation point velocity and excitation force. The units of SVVfSFFf1/2 and SVFf are identical to each other such that coherence is a dimensionless number with values in the range of 0-1. Coherence has a meaning similar to the correlation coefficient in the present context. For a small coherence, the velocity of the excitation point is considered to be noise independent of the excitation force.


## 4.2.1. Definition of Coherence

The term “coherence” refers to the output value of the measurement system and is defined ascohf=SVFf/SVVfSFFf1/2, where f denotes the frequency (Hz) and SVVf and SFFf represent the power spectra of the excitation point velocity and excitation force, respectively. Furthermore, SVFf represents the cross spectrum between the excitation point velocity and excitation force. The units of SVVfSFFf1/2 and SVFf are identical to each other such that coherence is a dimensionless number with values in the range of 0-1. Coherence has a meaning similar to the correlation coefficient in the present context. For a small coherence, the velocity of the excitation point is considered to be noise independent of the excitation force.


## 4.3. Preprocessing (2): Gain Correction

The FRF mobility data are influenced by the gain of the measurement system, which is estimated in equation (12). Parameters f0 and kG in equation (12) were determined by employing the least-squares method using the analytical solution and measured mobility data for an aluminum beam.(12)Gestf=ff0−kG,where f denotes the frequency (Hz), f0 (>0) is a constant (Hz), and kG (>0) is a constant (–).Using this gain, the measured dimensional mobilityMobM∗ at f was corrected using(13a)Mob∗=MobM∗Gestf,(13b)20log10Mob∗=20log10MobM∗−20log10Gestf,where equation (13b) represents the absolute mobility in decibels.


### 4.3.1. Determination of Gain Parameters by Measuring an Aluminum Beam

Assuming that the gain of the measurement system is independent of the beam material, the parameters were estimated using the measured results of the aluminum beam. Consequently,f0=10HzandkG=0.15− were obtained. The impact of gain correction on the absolute value of the dimensionless mobility data of the aluminum beam is presented in Figure 7. For appropriate gain correction, a dimensional mobility of 0 dB corresponds to 1ms−1N−1. In the least-squares method, the residual Sres was minimized:(14)Sres=∑nuCn−uMn/Gestn2+vCn−vMn/Gestn2∑nuCn2+vCn2,where u=ReW,v=ImW, and C, M, and n represent calculated, measured, and n-th frequency, respectively.Figure 7
Dimensionless mobility (absolute value): measured and corrected values compared with analytical solution (aluminum beam).


## 4.3.1. Determination of Gain Parameters by Measuring an Aluminum Beam

Assuming that the gain of the measurement system is independent of the beam material, the parameters were estimated using the measured results of the aluminum beam. Consequently,f0=10HzandkG=0.15− were obtained. The impact of gain correction on the absolute value of the dimensionless mobility data of the aluminum beam is presented in Figure 7. For appropriate gain correction, a dimensional mobility of 0 dB corresponds to 1ms−1N−1. In the least-squares method, the residual Sres was minimized:(14)Sres=∑nuCn−uMn/Gestn2+vCn−vMn/Gestn2∑nuCn2+vCn2,where u=ReW,v=ImW, and C, M, and n represent calculated, measured, and n-th frequency, respectively.Figure 7
Dimensionless mobility (absolute value): measured and corrected values compared with analytical solution (aluminum beam).


## 4.4. Preprocessing (3): Rendering FRF Mobility Dimensionless

As mobilities after preprocessing steps (1) and (2) are displayed in decibels, the absolute valueMob∗ms−1N−1, real part ReMob∗, and imaginary part ImMob∗ of the dimensional mobility can be described using(15a)Mob∗=10dBMob/20,(15b)ReMob∗=Mob∗cosphase Mob,(15c)ImMob∗=Mob∗sinphase Mob,where phaseMob denotes the phase angle of the dimensional mobility °.


### 4.4.1. Dimensionless Mobility

The measured dataMob∗ of the FRF mobility can be converted to dimensionless mobility W using(16)W=−iρbhωL2Mob∗=ρbhωL2ImMob∗−iReMob∗,i=−1,where ρbhL denotes the mass of the beam specimen. The phase angle arg(W) of the dimensionless mobility portrayed in Figure 8 ranges from 0 to π.Figure 8
Dimensionless mobility (materials: ABS and PP, length: 127 mm, and thickness: 1.6 mm). (a) Material: ABS. (b) Material: PP.
(a)(b)


## 4.4.1. Dimensionless Mobility

The measured dataMob∗ of the FRF mobility can be converted to dimensionless mobility W using(16)W=−iρbhωL2Mob∗=ρbhωL2ImMob∗−iReMob∗,i=−1,where ρbhL denotes the mass of the beam specimen. The phase angle arg(W) of the dimensionless mobility portrayed in Figure 8 ranges from 0 to π.Figure 8
Dimensionless mobility (materials: ABS and PP, length: 127 mm, and thickness: 1.6 mm). (a) Material: ABS. (b) Material: PP.
(a)(b)


## 4.5. Inverse Calculation

The application of inverse calculations expressed in equations (8a)–(8e) to W in Figure 8 yielded the relationship between the dimensionless complex wavenumber and frequency, as indicated in Figure 9. The thick lines represent the dimensionless real wavenumber X, while the thin lines represent the dimensionless attenuation constant (−Y). A drastically steep rise can be observed in the dimensionless real wavenumber X near 5.5 kHz in Figure 9(a) (ABS) and 6 kHz in Figure 9(b) (PP). However, the dimensionless attenuation constant alternated between small and large values.Figure 9
Dimensionless complex wavenumber (real wavenumber and attenuation constant). (a) Material: ABS. (b) Material: PP.
(a)(b)The storage modulus and loss factor were obtained by applying equations (6a) and (6b) to the real wavenumber βR=2X/L and attenuation constant α=−2Y/L obtained from the results in Figure 9. These results are presented in Figure 10, where the elastic modulus of an elastic body with respect to the resonant frequency is “purely elastic.” In addition, both the storage modulus and loss factor fluctuated significantly with the frequency, which is an anomaly for viscoelastic data. Moreover, the loss factor near the resonant frequency obtained using the half-power bandwidth method was calculated using the measurement system and presented as the “half-power bandwidth.”Figure 10
Viscoelasticity estimated from complex wavenumbers in Figure9. (a) Material: ABS. (b) Material: PP.
(a)(b)


## 4.6. Postprocessing: Noise Removal by Fitting Complex Wavenumber vs. Frequency Relationship to a Theoretical Curve

Instead of directly eliminating the noise from the viscoelasticity vs. frequency data, we developed a method to fit a curve for the relationship between the complex wavenumber and frequency using smooth functions of intermediate data.


### 4.6.1. Real Wavenumber vs. Frequency Relationship

The generalized Maxwell model is applicable to storage modulusE of materials, which is a monotonically increasing function of frequency. As ω ⟶ ∞ at high frequencies, E approaches the instantaneous modulus Einst, and as ω ⟶ 0 at low frequencies, it approaches the long-term modulus E∞. If η has a small frequency dependence, equation (4) can be considered for the following relationship between the real wavenumber βR and radial frequency ω.(17)βR∼ω1/2asω⟶∞,(18a)∂lnβR∂lnω=∂lnX∂lnf=12−14∂lnE∂lnω=12−14ωE∂E∂ω⟶kf0,0<kf0<12asω⟶+0,0<kf0<∂lnβR∂lnω≤12,(18b)∂E∂ω>0,E∞≤E≤Einst,for0<ω<∞.Therefore, if the relationship oflnβRvs.lnω is plotted, the slope at low frequencies is <1/2 and approaches 1/2 at high frequencies. Thus, we can determine a curve-fitting equation such that the slope satisfies equations (18a) and (18b).


### 4.6.2. Equation to Fit the Relationship between Real Wavenumber vs. Frequency

Slope from Reference Point. The plot of lnβRvs.lnω has many irregularities. Instead of the gradient ∂lnβR/∂lnω at each frequency, we define the slope kfω of the straight line connecting the points lnω,lnβR and lnω0,lnβ0 from the reference point using the following equations:(19a)lnβR=lnβ0+kfωlnωω0,where(19b)kfω=lnβR−lnβ0lnω−lnω0=lnβR/β0lnω/ω0,forω≠ω0,where ω0 and β0 denote the reference radial frequency s−1 and reference wavenumber (m−1), respectively; furthermore, βR=β0forω=ω0.(19c)β04=ρAω02E0I=12ρω02E0h2,where E0 denotes the storage modulus at ω0.From equations (18a) and (18b), the slope kfω satisfies(20)0<kf0<kfω<12,kfω⟶12,asω⟶∞.Furthermore,ω0 was(A)
Located near the center of the frequency range of interest and(B)
Situated far from the peaks and valleys due to resonance points, antiresonance points, and noise of the dimensionless mobilityWThe slopekfω does not always satisfy equation (20) due to variations in the real wavenumber βR from the inverse calculation on the FRF data. Therefore, the following conditions were used to satisfy equation (20):(21)ifkfω>12,thenkfω⟵12,ifkfω<kmin,thenkfω⟵kminkmindenotes a small positive number.


### 4.6.3. Curve-Fit Formula for Slopekfω

kfω obtained using the sigmoid function, a smooth function with a value between 0 and 1, was adopted as the curve-fitting equation satisfying equation (20).(22a)kfω=12ysxFd,where(22b)xFd=lnωωI=xF−xFI,(22c)ysxFd=11+exp−xFd,where xFd denotes the difference in the logarithm of frequency (logarithm of frequency ratio) (–), ωI denotes the radial frequency of the inflection point (s−1), and ysxFd represents the sigmoid function of xFd.(22d)xF=lnω,xFI=lnωI.The method for identifyingωI is described in Appendix C. Evaluating the sigmoid function is not necessary once ωI has been calculated. Moreover, equations (22a)–(22d) can be simplified as(23)kfω=12ωω+ωI.Applying the slope of equation (23) to equation (19a), the real wavenumber βR becomes a smooth function of the frequency. Figure 11 presents examples depicting the relationship between logarithms of the dimensionless real wavenumber and frequency. The real wavenumber βR was converted to a dimensionless real wavenumber X=βRL/2. In Figure 11, the points at which the thick line deviated from the regression curve were considered to be affected by noise. Vertical and horizontal axes of Figure 11 were rewritten in Figure 12 in linear forms. Here, for both ABS and PP, the thick and dashed lines were in good agreement for frequencies of 4 kHz or less.Figure 11
Curve fitting of dimensionless real wavenumber using equations (24a) and (24b). (a) Material: ABS. Reference point: f0=2243Hz,E0=2.560GHz,β0=133.4m−1,andX0=8.402−. Curve-fitting parameter: ωI=11.27s−1andfI=ωI/2π=1.794Hz. (b) Material: PP. Reference point: f0=2229Hz,E0=2.741GHz,β0=130.7m−1,andX0=8.234−. Curve-fitting parameter: ωI=13.10s−1andfI=ωI/2π=2.085Hz.
(a)(b)Figure 12
Dimensionless real wavenumber vs. frequency (linear scale). (a) Material: ABS. (b) Material: PP.
(a)(b)


### 4.6.4. Equation to Fit Relationship between Attenuation Constant and Frequency

In Figure10, the loss factor obtained by the inverse calculation from FRF data fluctuates in the range from 0 to approximately 0.3. However, it must satisfy conditions (A), (B), and (C). In addition, the ratio α/βR of the attenuation constant to real wavenumber was fitted to a theoretical curve based on these conditions. Thereafter, E and η were calculated from the values of βR and α using equations (6a) and (6b).(A) Relationship between Ratioα/βRand Loss Factorη. From the complex wavenumber equations (4) and (5), α/βR can be expressed as a function of only the loss factor η, which is expressed as(24a)αβR=−YX=2/1+1+η2−1/21/2−12/1+1+η2−1/21/2+11/2,(24b)αβR=−YX≈η4,asη⟶0.(B) Frequency at Which Ratioα/βRand Loss FactorηAre Maximized. The following relationship holds from equations (6a) and (6b):(25a)∂η∂ω=dηdQ∂Q∂ω=dηdQ2α/βR1−α/βR22⋅∂∂ωαβR,where(25b)Q=1+α/βR21−α/βR2,dηdQ=Q22−Q22−1−1/2×4Q32−Q3>0,1≤Q<2.Equations (25a) and (25b) maximize α/βR at the frequency corresponding to the maximization of the loss factor η.(C) For Many Viscoelastic Materials, the Loss Factor η Has the Maximum Value at a Single Frequency. From properties (A)–(C), the relationship of the ratio α/βR=−Y/Xvs.ω should be fitted with an appropriate curve to consider the maximum value at a particular ω. In this case, the loss factor in Figure 10 was replaced by a smooth curve with only one maximal peak.


### 4.6.5. Regression Curve Proportional to Gamma Distribution

Figure13 shows examples of the relationship of α/βR and frequency for both materials, which shows a large fluctuation of α/βR. Therefore, we considered the horizontal axis (frequency (Hz)) as a random variable, and the value on the vertical axis was considered as the number of counts in the histogram. Thus, the curve fit of the data in Figure 13 can be replaced with a case of an appropriate probability distribution applied to the experimental data represented by a histogram.Figure 13
Ratio−Y/X=α/βR vs. frequency. (a) Material: ABS. (b) Material: PP.
(a)(b)Therefore, a function proportional to the density function of the gamma distribution, a continuous distribution similar to the Poisson distribution, was used as the curve-fit function. The proposed equation for the regression curve ofα/βR is presented in equations (26a)–(26e), wherein the normalization by the step size r0 divided into 16 intervals between the minimum and maximum α/βR was adopted in equations (26a)–(26e). The plots obtained using this normalization are presented in Figure 14.(26a)r0nF=−YX=r0n0RfΓpσpgf,where(26b)gf=1Γpσpfp−1exp−fσ,f>0,(26c)Γp=∫0∞xp−1e−xdx,p>0,where gf denotes the probability density of the gamma distribution and Γp denotes the gamma function.(26d)nF=α/βRr0,(26e)r0=rmax−rmin16,where rminandrmax represent the minimum and maximum values of α/βR in the frequency range, respectively; r0 denotes the step width (–) and f=ω/2π denotes the frequency (Hz); σ>0andp>0− are the parameters in the gamma distribution; and n0Rf represents the coefficient for the height of the distribution HZ−p, as defined in(27)r0n0RfσpΓp=r0∫0∞nFdf.Figure 14
Curve fitting ofnF. (a) Material: ABS. (b) Material: PP.
(a)(b)Considering the thin line in Figure14 as a histogram with frequency f as a random variable, the sample mean mf and sample standard deviation sf with respect to frequency f can be defined. Parameters p and σ of the gamma distribution in equation (26b) can be expressed by equations related to mean mf and standard deviation sf. The coefficient n0Rf was determined using equation (27) such that the integral value of the curve-fit equation and plot data is equal.


### 4.6.6. Determining Gamma Distribution Parameters

For the probability distribution function expressed by equation (26b), the population mean Ef and population variance Vf are given by(28)Ef=pσ,Vf=pσ2.If these are equal to the sample meanmf and sample variance sf2, parameters p and σ are given by(29a)p=mfsf2,σ=sf2mf.The sample meanmf and sample standard deviation sf were calculated using equation (29b). In equation (29b), the region where the convergence points of the inverse calculation (markers in Figure 14) are densely present has large weights.(29b)mf=∑knFkfkN,sf2=∑knFkfk−mf2N,where(29c)N=∑knFk.Here,N is similar to the total count but is generally a noninteger real number. The inverse calculation of equation (8a) was performed for 635 frequency points from 30 to 6400 Hz with 10 Hz increments. The number of points where equation (8a) converged was 490 for ABS and 543 for PP, which are represented by markers in Figures 13 and 14. Plot point numbers k are assigned to these convergence points in the ascending order of frequency.


### 4.6.7. Curve-Fit for Ratioα/βR vs. Frequency

The thin lines in Figure14 show the normalized plots of ratio α/βR for ABS and PP. The thick lines are regression curves proportional to the gamma distributions, and the curve-fit parameters are shown in Figure 14. The gamma distribution parameters p are 3.4 and 4.1 for ABS and PP, respectively. The anomalous jumps are circled in Figure 14(a), which appear near the spike-like noises of the phase angles in Figures 6 and 8. Although the thick lines did not apparently correspond with the thin lines well, the curve-fit result was evaluated by the degree to which the original dimensionless mobility could be expressed. The reason why the thin lines (value by direct inverse calculation from FRF data) vary is discussed in Appendix D.


### 4.6.8. Curve-Fitting Parameter

ABS:p=3.430−,σ=918.6Hz,mf=3151Hz,sf=1701Hz,Γp=3.078−,r0=5.534×10−3−,andn0Rf=1.071×10−7Hz−p.PP:p=4.130−,σ=809.3Hz,mf=3342Hz,sf=1645Hz,Γp=7.080−,r0=2.881×10−3−,andn0Rf=4.825×10−9Hz−p.


### 4.6.9. Viscoelastic Model with Frequency Dependence Whose Loss Factor Is Proportional to the Gamma Distribution

Córtes and Elejabarrieta [10] derived a viscoelastic model represented by equations (30a) and (30b) via a theoretical study of the relaxation function. They experimentally validated this model using a material in which concrete was impregnated with epoxy resin. Equation (30b) is the case where p = 2 in equations (26a)–(26e). The constant η0 is 2×10−2, and therefore, the approximation of equation (24b) holds for such a small loss factor. Equation (30b) is a special case of our curve-fit equations (26a)–(26e).(30a)E′ω=E0=Const,(30b)ηω=η0ωω0exp1−ωω0,where E′ω denotes the storage modulus Pa and ηω denotes the loss factor (–); furthermore, η0=2cω0/E0e, in which c denotes the viscosity Pa⋅s and e represents the base of the natural logarithm.


### 4.6.10. FRF Mobility Based on Curve-Fitted Complex Wavenumber Data

We applied the curve-fitted data of the dimensionless complex wavenumberZ=X+iY to equations (8a)–(8e) to obtain the dimensionless mobility. The results were then compared with the original FRF data (see Figures 15 and 16). Furthermore, in these figures, the standard deviation StdReg of the difference between the thick and thin gray lines was calculated with respect to the absolute value and phase angle of W. These are shown as (thick line) ± StdReg values via dashed lines. However, (thick line)−StdReg has a small lower limit to prevent the dashed line from being negative.Figure 15
FRF mobility calculated from the curve-fitted dimensionless complex wavenumber (material: ABS). (a) Absolute value vs. frequency. (b) Phase angle vs. frequency.
(a)(b)Figure 16
FRF mobility calculated from the curve-fitted dimensionless complex wavenumber (material: PP). (a) Absolute value vs. frequency. (b) Phase angle vs. frequency.
(a)(b)In these figures, the thin lines are marked with the convergence points of the inverse calculations of equations (8a)–(8e), which are densely present in the frequency range between the antiresonance and resonance points. Here, the dimensionless mobilities based on the measured data (gray lines) and FRF based on the curve-fit of the complex wavenumber (thick lines) agree well. However, in the frequency range that crosses the resonance points, the convergence points of the inverse calculation of equations (8a)–(8e) are sparse, and the deviation between the thick and thin lines is conspicuous. In Figures 15 and 16, the thin gray lines fall within the range of (thick line) ± StdReg in most frequency ranges.


### 4.6.11. Viscoelasticity with Curve-Fitted Complex Wavenumbers

We converted the curve-fitted dimensionless complex wavenumbers to dimensional complex wavenumbersβ∗=βR−iα,whereβR=2X/Landα=−2Y/L. The results of determining the viscoelasticity by applying equations (6a) and (6b) to these complex wavenumbers are shown in Figures 17 and 18. In these figures, the thin gray lines are the results estimated directly from the FRF data, while the thick lines are based on the curve-fitted complex wavenumbers. Furthermore, the standard deviations of the thick and thin lines for the storage modulus and loss factor were calculated as StdReg. In Figures 17 and 18, the values of (thick line) ± StdReg are drawn as dashed lines. However, (thick line)−StdReg has a small lower limit to prevent negative values. In Figures 17(b) and 18(b), the thick lines have maximum values near 2 kHz and 2.5 kHz, respectively.Figure 17
Storage modulus and loss factor calculated from the curve-fitted dimensionless complex wavenumber (material: ABS). (a) Storage modulus vs. frequency. (b) Loss factor vs. frequency.
(a)(b)Figure 18
Storage modulus and loss factor calculated from the curve-fitted dimensionless complex wavenumber (material: PP). (a) Storage modulus vs. frequency. (b) Loss factor vs. frequency. yReg: calculated from the curve-fitted data (regression). StdReg: standard deviation of the residual.
(a)(b)


### 4.6.12. Abnormal Data

In Figure17(b), an alternate long and short dashed line is shown near the resonance point at 5.6 kHz. Near this frequency, the absolute value of mobility in Figure 6(a) rose abnormally from the minimum value. The circles near 1 kHz in Figures 17(a) and 17(b) correspond to the spike-like noise of the phase angles in Figure 6(a).


## 4.6.1. Real Wavenumber vs. Frequency Relationship

The generalized Maxwell model is applicable to storage modulusE of materials, which is a monotonically increasing function of frequency. As ω ⟶ ∞ at high frequencies, E approaches the instantaneous modulus Einst, and as ω ⟶ 0 at low frequencies, it approaches the long-term modulus E∞. If η has a small frequency dependence, equation (4) can be considered for the following relationship between the real wavenumber βR and radial frequency ω.(17)βR∼ω1/2asω⟶∞,(18a)∂lnβR∂lnω=∂lnX∂lnf=12−14∂lnE∂lnω=12−14ωE∂E∂ω⟶kf0,0<kf0<12asω⟶+0,0<kf0<∂lnβR∂lnω≤12,(18b)∂E∂ω>0,E∞≤E≤Einst,for0<ω<∞.Therefore, if the relationship oflnβRvs.lnω is plotted, the slope at low frequencies is <1/2 and approaches 1/2 at high frequencies. Thus, we can determine a curve-fitting equation such that the slope satisfies equations (18a) and (18b).


## 4.6.2. Equation to Fit the Relationship between Real Wavenumber vs. Frequency

Slope from Reference Point. The plot of lnβRvs.lnω has many irregularities. Instead of the gradient ∂lnβR/∂lnω at each frequency, we define the slope kfω of the straight line connecting the points lnω,lnβR and lnω0,lnβ0 from the reference point using the following equations:(19a)lnβR=lnβ0+kfωlnωω0,where(19b)kfω=lnβR−lnβ0lnω−lnω0=lnβR/β0lnω/ω0,forω≠ω0,where ω0 and β0 denote the reference radial frequency s−1 and reference wavenumber (m−1), respectively; furthermore, βR=β0forω=ω0.(19c)β04=ρAω02E0I=12ρω02E0h2,where E0 denotes the storage modulus at ω0.From equations (18a) and (18b), the slope kfω satisfies(20)0<kf0<kfω<12,kfω⟶12,asω⟶∞.Furthermore,ω0 was(A)
Located near the center of the frequency range of interest and(B)
Situated far from the peaks and valleys due to resonance points, antiresonance points, and noise of the dimensionless mobilityWThe slopekfω does not always satisfy equation (20) due to variations in the real wavenumber βR from the inverse calculation on the FRF data. Therefore, the following conditions were used to satisfy equation (20):(21)ifkfω>12,thenkfω⟵12,ifkfω<kmin,thenkfω⟵kminkmindenotes a small positive number.


## 4.6.3. Curve-Fit Formula for Slopekfω

kfω obtained using the sigmoid function, a smooth function with a value between 0 and 1, was adopted as the curve-fitting equation satisfying equation (20).(22a)kfω=12ysxFd,where(22b)xFd=lnωωI=xF−xFI,(22c)ysxFd=11+exp−xFd,where xFd denotes the difference in the logarithm of frequency (logarithm of frequency ratio) (–), ωI denotes the radial frequency of the inflection point (s−1), and ysxFd represents the sigmoid function of xFd.(22d)xF=lnω,xFI=lnωI.The method for identifyingωI is described in Appendix C. Evaluating the sigmoid function is not necessary once ωI has been calculated. Moreover, equations (22a)–(22d) can be simplified as(23)kfω=12ωω+ωI.Applying the slope of equation (23) to equation (19a), the real wavenumber βR becomes a smooth function of the frequency. Figure 11 presents examples depicting the relationship between logarithms of the dimensionless real wavenumber and frequency. The real wavenumber βR was converted to a dimensionless real wavenumber X=βRL/2. In Figure 11, the points at which the thick line deviated from the regression curve were considered to be affected by noise. Vertical and horizontal axes of Figure 11 were rewritten in Figure 12 in linear forms. Here, for both ABS and PP, the thick and dashed lines were in good agreement for frequencies of 4 kHz or less.Figure 11
Curve fitting of dimensionless real wavenumber using equations (24a) and (24b). (a) Material: ABS. Reference point: f0=2243Hz,E0=2.560GHz,β0=133.4m−1,andX0=8.402−. Curve-fitting parameter: ωI=11.27s−1andfI=ωI/2π=1.794Hz. (b) Material: PP. Reference point: f0=2229Hz,E0=2.741GHz,β0=130.7m−1,andX0=8.234−. Curve-fitting parameter: ωI=13.10s−1andfI=ωI/2π=2.085Hz.
(a)(b)Figure 12
Dimensionless real wavenumber vs. frequency (linear scale). (a) Material: ABS. (b) Material: PP.
(a)(b)


## 4.6.4. Equation to Fit Relationship between Attenuation Constant and Frequency

In Figure10, the loss factor obtained by the inverse calculation from FRF data fluctuates in the range from 0 to approximately 0.3. However, it must satisfy conditions (A), (B), and (C). In addition, the ratio α/βR of the attenuation constant to real wavenumber was fitted to a theoretical curve based on these conditions. Thereafter, E and η were calculated from the values of βR and α using equations (6a) and (6b).(A) Relationship between Ratioα/βRand Loss Factorη. From the complex wavenumber equations (4) and (5), α/βR can be expressed as a function of only the loss factor η, which is expressed as(24a)αβR=−YX=2/1+1+η2−1/21/2−12/1+1+η2−1/21/2+11/2,(24b)αβR=−YX≈η4,asη⟶0.(B) Frequency at Which Ratioα/βRand Loss FactorηAre Maximized. The following relationship holds from equations (6a) and (6b):(25a)∂η∂ω=dηdQ∂Q∂ω=dηdQ2α/βR1−α/βR22⋅∂∂ωαβR,where(25b)Q=1+α/βR21−α/βR2,dηdQ=Q22−Q22−1−1/2×4Q32−Q3>0,1≤Q<2.Equations (25a) and (25b) maximize α/βR at the frequency corresponding to the maximization of the loss factor η.(C) For Many Viscoelastic Materials, the Loss Factor η Has the Maximum Value at a Single Frequency. From properties (A)–(C), the relationship of the ratio α/βR=−Y/Xvs.ω should be fitted with an appropriate curve to consider the maximum value at a particular ω. In this case, the loss factor in Figure 10 was replaced by a smooth curve with only one maximal peak.


## 4.6.5. Regression Curve Proportional to Gamma Distribution

Figure13 shows examples of the relationship of α/βR and frequency for both materials, which shows a large fluctuation of α/βR. Therefore, we considered the horizontal axis (frequency (Hz)) as a random variable, and the value on the vertical axis was considered as the number of counts in the histogram. Thus, the curve fit of the data in Figure 13 can be replaced with a case of an appropriate probability distribution applied to the experimental data represented by a histogram.Figure 13
Ratio−Y/X=α/βR vs. frequency. (a) Material: ABS. (b) Material: PP.
(a)(b)Therefore, a function proportional to the density function of the gamma distribution, a continuous distribution similar to the Poisson distribution, was used as the curve-fit function. The proposed equation for the regression curve ofα/βR is presented in equations (26a)–(26e), wherein the normalization by the step size r0 divided into 16 intervals between the minimum and maximum α/βR was adopted in equations (26a)–(26e). The plots obtained using this normalization are presented in Figure 14.(26a)r0nF=−YX=r0n0RfΓpσpgf,where(26b)gf=1Γpσpfp−1exp−fσ,f>0,(26c)Γp=∫0∞xp−1e−xdx,p>0,where gf denotes the probability density of the gamma distribution and Γp denotes the gamma function.(26d)nF=α/βRr0,(26e)r0=rmax−rmin16,where rminandrmax represent the minimum and maximum values of α/βR in the frequency range, respectively; r0 denotes the step width (–) and f=ω/2π denotes the frequency (Hz); σ>0andp>0− are the parameters in the gamma distribution; and n0Rf represents the coefficient for the height of the distribution HZ−p, as defined in(27)r0n0RfσpΓp=r0∫0∞nFdf.Figure 14
Curve fitting ofnF. (a) Material: ABS. (b) Material: PP.
(a)(b)Considering the thin line in Figure14 as a histogram with frequency f as a random variable, the sample mean mf and sample standard deviation sf with respect to frequency f can be defined. Parameters p and σ of the gamma distribution in equation (26b) can be expressed by equations related to mean mf and standard deviation sf. The coefficient n0Rf was determined using equation (27) such that the integral value of the curve-fit equation and plot data is equal.


## 4.6.6. Determining Gamma Distribution Parameters

For the probability distribution function expressed by equation (26b), the population mean Ef and population variance Vf are given by(28)Ef=pσ,Vf=pσ2.If these are equal to the sample meanmf and sample variance sf2, parameters p and σ are given by(29a)p=mfsf2,σ=sf2mf.The sample meanmf and sample standard deviation sf were calculated using equation (29b). In equation (29b), the region where the convergence points of the inverse calculation (markers in Figure 14) are densely present has large weights.(29b)mf=∑knFkfkN,sf2=∑knFkfk−mf2N,where(29c)N=∑knFk.Here,N is similar to the total count but is generally a noninteger real number. The inverse calculation of equation (8a) was performed for 635 frequency points from 30 to 6400 Hz with 10 Hz increments. The number of points where equation (8a) converged was 490 for ABS and 543 for PP, which are represented by markers in Figures 13 and 14. Plot point numbers k are assigned to these convergence points in the ascending order of frequency.


## 4.6.7. Curve-Fit for Ratioα/βR vs. Frequency

The thin lines in Figure14 show the normalized plots of ratio α/βR for ABS and PP. The thick lines are regression curves proportional to the gamma distributions, and the curve-fit parameters are shown in Figure 14. The gamma distribution parameters p are 3.4 and 4.1 for ABS and PP, respectively. The anomalous jumps are circled in Figure 14(a), which appear near the spike-like noises of the phase angles in Figures 6 and 8. Although the thick lines did not apparently correspond with the thin lines well, the curve-fit result was evaluated by the degree to which the original dimensionless mobility could be expressed. The reason why the thin lines (value by direct inverse calculation from FRF data) vary is discussed in Appendix D.


## 4.6.8. Curve-Fitting Parameter

ABS:p=3.430−,σ=918.6Hz,mf=3151Hz,sf=1701Hz,Γp=3.078−,r0=5.534×10−3−,andn0Rf=1.071×10−7Hz−p.PP:p=4.130−,σ=809.3Hz,mf=3342Hz,sf=1645Hz,Γp=7.080−,r0=2.881×10−3−,andn0Rf=4.825×10−9Hz−p.


## 4.6.9. Viscoelastic Model with Frequency Dependence Whose Loss Factor Is Proportional to the Gamma Distribution

Córtes and Elejabarrieta [10] derived a viscoelastic model represented by equations (30a) and (30b) via a theoretical study of the relaxation function. They experimentally validated this model using a material in which concrete was impregnated with epoxy resin. Equation (30b) is the case where p = 2 in equations (26a)–(26e). The constant η0 is 2×10−2, and therefore, the approximation of equation (24b) holds for such a small loss factor. Equation (30b) is a special case of our curve-fit equations (26a)–(26e).(30a)E′ω=E0=Const,(30b)ηω=η0ωω0exp1−ωω0,where E′ω denotes the storage modulus Pa and ηω denotes the loss factor (–); furthermore, η0=2cω0/E0e, in which c denotes the viscosity Pa⋅s and e represents the base of the natural logarithm.


## 4.6.10. FRF Mobility Based on Curve-Fitted Complex Wavenumber Data

We applied the curve-fitted data of the dimensionless complex wavenumberZ=X+iY to equations (8a)–(8e) to obtain the dimensionless mobility. The results were then compared with the original FRF data (see Figures 15 and 16). Furthermore, in these figures, the standard deviation StdReg of the difference between the thick and thin gray lines was calculated with respect to the absolute value and phase angle of W. These are shown as (thick line) ± StdReg values via dashed lines. However, (thick line)−StdReg has a small lower limit to prevent the dashed line from being negative.Figure 15
FRF mobility calculated from the curve-fitted dimensionless complex wavenumber (material: ABS). (a) Absolute value vs. frequency. (b) Phase angle vs. frequency.
(a)(b)Figure 16
FRF mobility calculated from the curve-fitted dimensionless complex wavenumber (material: PP). (a) Absolute value vs. frequency. (b) Phase angle vs. frequency.
(a)(b)In these figures, the thin lines are marked with the convergence points of the inverse calculations of equations (8a)–(8e), which are densely present in the frequency range between the antiresonance and resonance points. Here, the dimensionless mobilities based on the measured data (gray lines) and FRF based on the curve-fit of the complex wavenumber (thick lines) agree well. However, in the frequency range that crosses the resonance points, the convergence points of the inverse calculation of equations (8a)–(8e) are sparse, and the deviation between the thick and thin lines is conspicuous. In Figures 15 and 16, the thin gray lines fall within the range of (thick line) ± StdReg in most frequency ranges.


## 4.6.11. Viscoelasticity with Curve-Fitted Complex Wavenumbers

We converted the curve-fitted dimensionless complex wavenumbers to dimensional complex wavenumbersβ∗=βR−iα,whereβR=2X/Landα=−2Y/L. The results of determining the viscoelasticity by applying equations (6a) and (6b) to these complex wavenumbers are shown in Figures 17 and 18. In these figures, the thin gray lines are the results estimated directly from the FRF data, while the thick lines are based on the curve-fitted complex wavenumbers. Furthermore, the standard deviations of the thick and thin lines for the storage modulus and loss factor were calculated as StdReg. In Figures 17 and 18, the values of (thick line) ± StdReg are drawn as dashed lines. However, (thick line)−StdReg has a small lower limit to prevent negative values. In Figures 17(b) and 18(b), the thick lines have maximum values near 2 kHz and 2.5 kHz, respectively.Figure 17
Storage modulus and loss factor calculated from the curve-fitted dimensionless complex wavenumber (material: ABS). (a) Storage modulus vs. frequency. (b) Loss factor vs. frequency.
(a)(b)Figure 18
Storage modulus and loss factor calculated from the curve-fitted dimensionless complex wavenumber (material: PP). (a) Storage modulus vs. frequency. (b) Loss factor vs. frequency. yReg: calculated from the curve-fitted data (regression). StdReg: standard deviation of the residual.
(a)(b)


## 4.6.12. Abnormal Data

In Figure17(b), an alternate long and short dashed line is shown near the resonance point at 5.6 kHz. Near this frequency, the absolute value of mobility in Figure 6(a) rose abnormally from the minimum value. The circles near 1 kHz in Figures 17(a) and 17(b) correspond to the spike-like noise of the phase angles in Figure 6(a).


## 4.7. Effect of Resonance/Antiresonance Points on Fluctuations of Storage Modulus and Loss Factor

The resonant and antiresonant frequencies of the beam specimen are indicated by markers on the frequency axes of Figures17 and 18. Here, the estimated storage modulus and loss factor measured via direct inverse calculations from the FRF data (thin lines) are high near the antiresonant frequency and low near the resonant frequency. The method of determining viscoelasticity based on the complex wavenumber is not originally affected by modes such as resonance and antiresonance. However, the inverse calculations for obtaining the complex wavenumbers from the FRF mobility converge well in the frequency range from the antiresonance to resonance point, whereas the number of convergence points becomes sparse in a frequency range exceeding the resonance point. In Figures 6 and 8 for the FRF data, the falling edge of the phase angle near the resonant frequency represents a sudden change close to a step. According to the qualitative study of equations in Appendix D, the inverse calculation converges poorly because dW/dZ∝1/η2 holds true near the resonance point. However, near the antiresonant frequency, noise may occur due to mixing of natural vibrations, which is concerning.


## 4.8. Estimation of the Magnitude of Noise Contained in the Dimensionless Complex Wavenumber and FRFW

### 4.8.1. Estimating the Magnitude of Noise Contained in the Mobility Measurement Data

The dimensionless mobilityWfit obtained by applying the curve-fitted dimensionless complex wavenumber to equations (8a)–(8e) has a value different from the original measurement data WI as shown in Figures 15 and 16. Figures 19 and 20 show the real and imaginary parts of Wfit=ufit+ivfit and WI=uI+ivI, respectively. Here, we consider the difference ΔW=WI−Wfit as noise in the FRF measurement data.(31)ΔW=WI−Wfit=uI−ufit+ivI−vfit,where uandv denote the real and imaginary parts of the dimensionless mobility W(−), respectively. Subscripts I and fit denote direct inverse calculation value from the measurement and the curve-fitted complex wavenumber, respectively.Figure 19
Comparison of dimensionless mobility (material: ABS). (a) Re(W) vs. frequency. (b) Im(W) vs. frequency.
(a)(b)Figure 20
Comparison of dimensionless mobility (material: PP). (a) Re(W) vs. frequency. (b) Im(W) vs. frequency.
(a)(b)


### 4.8.2. Noise in Dimensionless Complex Wavenumber

Consider changes in the dimensionless real wavenumberX and dimensionless attenuation constant −Y when the FRF W changes by ∆W. Since W expressed by equation (8a) is differentiable at points other than at the resonance point, the following relationship holds:(32)dWdZ=∂W∂X=1i∂W∂Y.Therefore, the following relationships among the change ∆W of W, that of ∆X of X, and that of ∆Y of Y hold:(33)ΔW=∂W∂XΔX+∂W∂YΔY=dWdZΔX+iΔY.Therefore,(34a)ΔX=ReΔWdW/dZ,ΔY=ImΔWdW/dZ,where(34b)ΔW=Δu+iΔv,Δu=uI−ufit,Δv=vI−vfit.Equations (34a) and (34b) show that the fluctuation of the FRF affects both the real and imaginary parts of the complex wavenumber. On the contrary, the difference between the nondimensional complex wave number calculated directly from the FRF W measurement data and the curve-fitted nondimensional complex wave number can be expressed by(35a)δX=XI−Xfit,δY=YI−Yfit,(35b)ZI=XI+iYI,Zfit=Xfit+iYfit.Let us considerδX or ∆X as the fluctuation of X and δY or ∆Y as the fluctuation of Y.SlopedW/dZ. Figure 21 shows the gradient dW/dZ when Z=Xfit+iYfit is applied to equations (D.7a) and (D.7b). To make the figure clear, the vertical axis of Figure 21(a) is cut off by ±200.Figure 21
SlopedW/dZ vs. frequency. (a) Material: ABS. (b) Material: PP.
(a)(b)Fluctuation of Dimensionless Complex Wavenumber. ΔX/Xfit and δX/Xfit are considered the degree of noise in the data of real wavenumber X. From Figures 22(a) and 23(a), ΔX/Xfit and δX/Xfit are in the range of about 5–7% in most frequency regions. However, δX/Xfit increases sharply above 5.5 kHz in Figures 22(a) and 23(a). This agrees with the fact that the relation of dimensionless real wavenumber X vs. the frequency relationship near 5.5 kHz in Figures 12(a) and 12(b) is not smooth and corresponds to the spike-like noise of the phase angle near 5.5 kHz in Figures 8 and 9. For the fluctuation of the dimensionless attenuation constant −Y, ΔY/Xfit and δY/Xfit were regarded as indicators of the noise magnitude. In Figures 22(b) and 23(b), the ratio of attenuation constant to real wavenumber −Yfit/Xfit is shown by dashed lines, which are the envelopes above ΔY/Xfit and δY/Xfit. For the measurement data of this study, the dimensionless attenuation constant −YI obtained by direct inverse calculation from the FRF measurement data was small −YI≈0, where the frequency was slightly higher than the resonance point. Therefore, δY≈−Yfit≥0 in this frequency range. From Figures 22(b) and 23(b), the noise magnitudes in the attenuation constant ΔY/XfitandδY/Xfit are in the range of approximately 2–4% in most frequency regions below 5.5 kHz.Figure 22
Fluctuation of complex wavenumber (material: ABS). (a) Fluctuation of real wavenumber. (b) Fluctuation of attenuation constant.
(a)(b)Figure 23
Fluctuation of complex wavenumber (material: PP). (a) Fluctuation of real wavenumber. (b) Fluctuation of attenuation constant.
(a)(b)


## 4.8.1. Estimating the Magnitude of Noise Contained in the Mobility Measurement Data

The dimensionless mobilityWfit obtained by applying the curve-fitted dimensionless complex wavenumber to equations (8a)–(8e) has a value different from the original measurement data WI as shown in Figures 15 and 16. Figures 19 and 20 show the real and imaginary parts of Wfit=ufit+ivfit and WI=uI+ivI, respectively. Here, we consider the difference ΔW=WI−Wfit as noise in the FRF measurement data.(31)ΔW=WI−Wfit=uI−ufit+ivI−vfit,where uandv denote the real and imaginary parts of the dimensionless mobility W(−), respectively. Subscripts I and fit denote direct inverse calculation value from the measurement and the curve-fitted complex wavenumber, respectively.Figure 19
Comparison of dimensionless mobility (material: ABS). (a) Re(W) vs. frequency. (b) Im(W) vs. frequency.
(a)(b)Figure 20
Comparison of dimensionless mobility (material: PP). (a) Re(W) vs. frequency. (b) Im(W) vs. frequency.
(a)(b)


## 4.8.2. Noise in Dimensionless Complex Wavenumber

Consider changes in the dimensionless real wavenumberX and dimensionless attenuation constant −Y when the FRF W changes by ∆W. Since W expressed by equation (8a) is differentiable at points other than at the resonance point, the following relationship holds:(32)dWdZ=∂W∂X=1i∂W∂Y.Therefore, the following relationships among the change ∆W of W, that of ∆X of X, and that of ∆Y of Y hold:(33)ΔW=∂W∂XΔX+∂W∂YΔY=dWdZΔX+iΔY.Therefore,(34a)ΔX=ReΔWdW/dZ,ΔY=ImΔWdW/dZ,where(34b)ΔW=Δu+iΔv,Δu=uI−ufit,Δv=vI−vfit.Equations (34a) and (34b) show that the fluctuation of the FRF affects both the real and imaginary parts of the complex wavenumber. On the contrary, the difference between the nondimensional complex wave number calculated directly from the FRF W measurement data and the curve-fitted nondimensional complex wave number can be expressed by(35a)δX=XI−Xfit,δY=YI−Yfit,(35b)ZI=XI+iYI,Zfit=Xfit+iYfit.Let us considerδX or ∆X as the fluctuation of X and δY or ∆Y as the fluctuation of Y.SlopedW/dZ. Figure 21 shows the gradient dW/dZ when Z=Xfit+iYfit is applied to equations (D.7a) and (D.7b). To make the figure clear, the vertical axis of Figure 21(a) is cut off by ±200.Figure 21
SlopedW/dZ vs. frequency. (a) Material: ABS. (b) Material: PP.
(a)(b)Fluctuation of Dimensionless Complex Wavenumber. ΔX/Xfit and δX/Xfit are considered the degree of noise in the data of real wavenumber X. From Figures 22(a) and 23(a), ΔX/Xfit and δX/Xfit are in the range of about 5–7% in most frequency regions. However, δX/Xfit increases sharply above 5.5 kHz in Figures 22(a) and 23(a). This agrees with the fact that the relation of dimensionless real wavenumber X vs. the frequency relationship near 5.5 kHz in Figures 12(a) and 12(b) is not smooth and corresponds to the spike-like noise of the phase angle near 5.5 kHz in Figures 8 and 9. For the fluctuation of the dimensionless attenuation constant −Y, ΔY/Xfit and δY/Xfit were regarded as indicators of the noise magnitude. In Figures 22(b) and 23(b), the ratio of attenuation constant to real wavenumber −Yfit/Xfit is shown by dashed lines, which are the envelopes above ΔY/Xfit and δY/Xfit. For the measurement data of this study, the dimensionless attenuation constant −YI obtained by direct inverse calculation from the FRF measurement data was small −YI≈0, where the frequency was slightly higher than the resonance point. Therefore, δY≈−Yfit≥0 in this frequency range. From Figures 22(b) and 23(b), the noise magnitudes in the attenuation constant ΔY/XfitandδY/Xfit are in the range of approximately 2–4% in most frequency regions below 5.5 kHz.Figure 22
Fluctuation of complex wavenumber (material: ABS). (a) Fluctuation of real wavenumber. (b) Fluctuation of attenuation constant.
(a)(b)Figure 23
Fluctuation of complex wavenumber (material: PP). (a) Fluctuation of real wavenumber. (b) Fluctuation of attenuation constant.
(a)(b)


## 4.9. Effect of Phase Angle Correction on Loss Factor Estimates

The measurement data (thin line) of the phase angle ofW in Figure 16(b) were replaced with data based on the curve-fitted complex wavenumber, represented by a dashed line in Figure 24(a). Nonetheless, the absolute value of W remained on the thick line in Figure 24(a). The loss factor estimated from the inverse calculation at this instant is depicted in Figure 24(b), wherein the fluctuation of the thin line is smaller than that in Figure 18(b). Moreover, the phase angle measurement data of FRF W on the high frequency side of the resonance point significantly influences the fluctuation of the estimated loss factor.Figure 24
Influence of phase angle on loss factor estimates (material: PP). (a) Revision of the phase angle. (b) Effect on the loss factor estimates.
(a)(b)


## 5. Conclusion

In this study, we presented a method for estimating the frequency dispersion of viscoelasticity from the measured FRF mobility data using the center impedance method for a beam specimen. To remove the noise in the measured data, the relationship between the complex wavenumber and frequency was curve fitted with smooth functions having mechanical meanings. Furthermore, the viscoelasticity was determined from the curve-fitted complex wavenumber. Using the proposed method, viscoelasticity could be determined over a wide frequency range using only the data obtained from FRF measured for several seconds. This proposed method can be utilized to materials comprising polymer alloys, fibers, etc. as it is difficult to expand the estimated frequency range by simply using TTS for these materials. This method can be applied to materials with large loss factors of 1 or more.Only FRF mobility from existing measurement systems is utilized here. In the measured FRF data, the phase angle drops in steps on the high frequency side of the resonance point. On the contrary, for the FRF analytical solution, this drop has a smooth frequency dependence. However, the analytical solution only deals with bending vibrations, and vibration modes such as torsion can also cause noise if the beam test specimen is not installed properly. We confirmed that our data contains noise near the resonance points and hypothesized that one of the causes is random excitation by white noise, which is obtained by imposing impulses at random time intervals following a Poisson distribution. By adjusting the average and variance of the number of impulses per unit time, we anticipate that the noise in FRF data may be reduced. We would like to discuss and cooperate with measuring instrument manufacturers to test this hypothesis and optimize the excitation force containing white noise.



---

*Source: 1644823-2021-09-27.xml*