A Comparison of Methods to Measure the Coupling Coefficient of Electromagnetic Vibration Energy Harvesters

Vibration energy harvesters transform environmental vibration energy into usable electrical energy. The transformation is only possible because of a coupling between the mechanical part of the energy harvester and the electric circuit. This paper compares several measurement methods to determine the electromagnetic coupling coefficient. These methods are applied to various implementations of an energy harvester and the results are compared with one another and with simulation data by analyzing the magnetic flux. The average deviation between the measurement methods and the simulation data in our study was 5%. This good agreement validates the methods. Based on this, we recommend determination of the coupling coefficient and the optimum load resistance for maximum power harvesting on the basis of simulations and the open circuit method, because this procedure leads to the shortest measurement times.


Introduction
The deployment of many wireless sensor networks only makes sense commercially and environmentally if the sensors can be powered by energy harvesting (i.e., by converting energy freely available in the environment to usable electrical energy) [1]. The energy sources can be vibration, light, temperature gradients, airflow, or electromagnetic waves [2,3]. Vibration energy can be converted by electromagnetic, piezoelectric, or electrostatic principles [4]. Electromagnetic conversion is achieved by the relative movement of a permanent magnet and a coil, in which a voltage is induced.
Vibration harvesters are designed for a specific environment. The energy output is maximized when the resonance frequency matches the vibration frequency. In time-variant environments, the average energy output can be increased by adapting the resonance frequency of a narrowband harvester to the instantaneous vibration frequency or by using a broadband harvester [5][6][7]. An example of the latter case are nonlinear transducers [8,9].
An essential aspect of vibration energy harvesters is the coupling between the mechanical and electrical part of harvester. In piezoelectric harvesters, this coupling is described by the piezoelectric coefficient [10,11], while in electromagnetic harvesters it is described by the electromagnetic coupling coefficient [12,13].
Several authors have dealt with the simulation or measurement of the latter. Glynne-Jones et al. [14] experimentally investigated the optimum load resistance leading to the maximum generated power. They also derive the magnetic flux density from the induction voltage, but did not deal with the coupling. Stephen [12] derived the coupling coefficient and its influence on the optimum load 2 of 14 resistance in detail, but did not provide experimental validation. O'Donnell et al. [13] investigated the coupling by finite element analysis (FEA), but again no measurements were performed to corroborate the simulation results.
Cheng et al. [15] derived a lumped element equivalent circuit of electromagnetic energy harvesters, built a prototype, and measured its parameters. The coupling was determined by applying a current to the pickup coil and measuring the resulting magnetic force; the coupling coefficient then is the slope of the linear fit function between the force and the current. The damping coefficient was measured by an open circuit resonance curve, as well as by the impulse response. The optimum load resistance was estimated by applying different load resistances and measuring the load voltage. This paper is very detailed, but lacks a comparison between simulated and measured coupling coefficients.
Spreemann et al. [16] demonstrated an interesting numerical optimization approach for the parameters of electromagnetic harvesters, such as the resting position of the mass and geometrical magnet and coil parameters. They calculated the magnetic field analytically and reported the normalized gradient of the magnetic flux density. The optimum load resistance of a prototype was measured by applying different load resistances and searching for maximum power. A concrete coupling coefficient value was neither given by calculation nor by measurement.
Mann and Owens [17] investigated a nonlinear electromagnetic energy harvester with a bistable potential well. They reported a measured value of the coupling, but did not clearly state how the value was measured. Overall, the physical dimension of the parameter appears to be wrong. Cepnik et al. [18] described the time-optimized simulation of the electromagnetic coupling by analytical calculation, but did not report simulated or measured values. Szarka et al. [19] examined power conditioning systems for electromagnetic transducers and stated the numerical value of the electromagnetic coupling coefficient of a prototype, but did not reveal how the value was obtained.
It is obvious that the coupling coefficient in electromagnetic energy harvesters has not received sufficient attention up to now. Many values cited in the literature cannot be considered reliable and can hardly be compared to one another for lack of detail. This paper aims at shedding some light on the issue. It deals with the simulation and measurement of the electromagnetic coupling coefficient. Section 2 explains the electromagnetic coupling theory and how the coupling coefficient can be extracted from simulation and measurement data. Section 3 is devoted to the simulation (by FEA), realization, and experimental characterization of a prototype harvester. The results of the measurement methods are discussed in Section 4. Finally, Section 5 gives a summary.

Theory
Electromagnetic vibration energy harvesters mainly consist of permanent magnets and induction coils. When exposed to vibration, the permanent magnets move relative to the coils, in which a voltage is induced by Faraday's law. The voltage can be used to power a load. The mechanical part of the harvester is a distributed system and must be described by partial differential equations. In many instances, however, it suffices to track the movement of a representative element or reference point on the harvester, which can be described by ordinary differential equations. In this sense, a light cantilever beam with heavy permanent magnets and iron keepers near its free end can be modeled by an equivalent spring-mass-damper system, with the mass of the magnets and the keepers as effective mass m, an effective spring constant k, and a mechanical damping coefficient c m (Figure 1a). The mass and the spring constant define the natural resonance frequency ω n = √ k/m. The electrical part consists of the induction coil with the induction voltage u ind as the voltage source, L C and R C as coil inductance and resistance, respectively, and R L as load resistance (Figure 1b).
The equivalent lumped element model is excited by a base acceleration a(t) = − ..
y(t) with the base position y. With the relative displacement z(t) = x(t) − y(t) (the absolute displacement x(t) may be taken to coincide with the center of mass of the permanent magnet in the case of the tip-loaded Micromachines 2019, 10, 826 3 of 14 cantilever beam) and the coupling coefficient K, the system is described by the following coupled differential equations [15,18] Equation (1), with K = 0, is the differential equation of an unloaded mechanical oscillator [20]. The extra term is due to the back action of the induction current in the linear approximation [16].
Micromachines 2019, 10, 826 3 of 14 loaded cantilever beam) and the coupling coefficient K , the system is described by the following coupled differential equations [15,18] Equation (1), with = 0 K , is the differential equation of an unloaded mechanical oscillator [20]. The extra term is due to the back action of the induction current in the linear approximation [16].
The left part of Equation (2) is derived from Faraday's law by with N the number of coil turns and ϕ the magnetic flux through a single coil turn [18]. K depends on the coil geometry and the magnetic circuit details. For example, for a given coil volume, K can be increased by using thinner wires and more coil turns. Neglecting the coil inductance C L in Equation (2) at the commonly low frequencies [13], solving the current ( ) i t and inserting the result into Equation (1) shows that the back action of the induction current on the mechanical part of the harvester is equivalent to an additional electrical damping coefficient e c : Here, e c depends on both the harvester and the load resistance L R . With small L R , the electrical damping coefficient e c , and consequently the total damping coefficient = + high. This reduces the displacement ( ) z t , the velocity ( ) z t , the induction voltage ind u , and the load power L P . With high L R , e c is small, and therefore ( ) z t and ( ) z t are high, but the load power is again small because the load current is small.
Stephen [12] showed that care must be exercised when looking for the optimum load resistance L,opt R providing the highest load power. The maximum power transfer theorem would lead one to The left part of Equation (2) is derived from Faraday's law by with N the number of coil turns and ϕ the magnetic flux through a single coil turn [18]. K depends on the coil geometry and the magnetic circuit details. For example, for a given coil volume, K can be increased by using thinner wires and more coil turns. Neglecting the coil inductance L C in Equation (2) at the commonly low frequencies [13], solving the current i(t) and inserting the result into Equation (1) shows that the back action of the induction current on the mechanical part of the harvester is equivalent to an additional electrical damping coefficient c e : Here, c e depends on both the harvester and the load resistance R L . With small R L , the electrical damping coefficient c e , and consequently the total damping coefficient c = c e + c m , are high. This reduces the displacement z(t), the velocity . z(t), the induction voltage u ind , and the load power P L . With high R L , c e is small, and therefore z(t) and . z(t) are high, but the load power is again small because the load current is small.
Stephen [12] showed that care must be exercised when looking for the optimum load resistance R L,opt providing the highest load power. The maximum power transfer theorem would lead one to set R L,opt = R C , but this is not true here because of the back action of the load current. It rather follows from dP L /dR L = 0 that [12]

Four Methods of Measuring the Electromagnetic Coupling Coefficient
The above theoretical considerations lead to several methods of extracting the coupling coefficient from simulation or measurement data. First, K can be found by an electromagnetic FEA. Modeling the magnetic circuit and the induction coil shown in Figure 2a in a FE program (e.g. Ansys Electromagnetics) and making the parts move relative to each other yields the total magnetic flux Φ(z) = Nϕ(z). Figure 2b shows an exemplary curve with a third-order polynomial fitting function. The linear gradient of the curve at z = 0 is the simulated coupling coefficient K sim used for comparison in this work. K sim is approximately valid for small oscillation amplitudes. The geometry and the simulation details behind

Four Methods of Measuring the Electromagnetic Coupling Coefficient
The above theoretical considerations lead to several methods of extracting the coupling coefficient from simulation or measurement data. First, K can be found by an electromagnetic FEA. Modeling the magnetic circuit and the induction coil shown in Figure 2a in a FE program (e.g. Ansys Electromagnetics) and making the parts move relative to each other yields the total magnetic flux Φ = ϕ ( ) ( ) z N z . Figure 2b shows an exemplary curve with a third-order polynomial fitting function.
with the mechanical damping factor ζ = m m /(2 ) c km . Fitting measured data with functions of this type and choosing the parameter ζ m for optimum agreement between measured data and fit curve allows one to identify the coupling coefficient K ( Figure 3). The corresponding estimate of K will be called oc K in the following.
The resonance curve of the absolute position ( ) x t in the case of a short circuit (  To extract coupling coefficient values from measurements, the base of an energy harvester is excited harmonically with constant acceleration amplitudeâ and varying frequencies. With R L → ∞ , the measured open circuit voltage is the induction voltage. The transfer function, defined as the ratio of the induction voltage phasor U ind to the acceleration phasor A, is derived from Equation (1) (underlined symbols are understood to represent complex-valued quantities.) Its magnitude is with the mechanical damping factor ζ m = c m /(2 √ km). Fitting measured data with functions of this type and choosing the parameter ζ m for optimum agreement between measured data and fit curve allows one to identify the coupling coefficient K (Figure 3). The corresponding estimate of K will be called K oc in the following.
Here, X and Y are the phasors of ( ) x t and ( ) y t . Fitting such a curve to measured data respectively yields the short-circuit damping factor ζ e,sc and the electrical damping coefficient e,sc c . From this, an estimate sc K of the coupling coefficient is obtained from Equation (4) by This method requires knowledge of the mechanical damping coefficient ζ m , and therefore serves as additional validation of the open circuit method. A third way of estimating the coupling coefficient is by evaluating the measured optimum load resistance L,opt R . The load voltage spectrum near resonance is measured with different load resistances, and the resistance resulting in the highest power at resonance is L,opt R [14][15][16]. An estimate of the coupling coefficient then follows from Equation (5): The procedure is time consuming because one has to determine the peaks of resonance curves for various load resistances. It does not suffice to just measure the output power at a fixed frequency (e.g., the open circuit resonance frequency) because the resonance frequency varies with the damping [20].
The fourth measurement method for the coupling coefficient is based on a linear variation of the current i through the coil [15]. The resulting magnetic force mag F on the permanent magnets is measured by the cantilever beam deflection x via the effective spring constant k . One obtains the following estimate of the coupling coefficient i K : This can be exploited fully automatically with a slowly rising and falling current. The slow change reduces measurement errors because the mass is slowly swinging. The spring constant can The resonance curve of the absolute position x(t) in the case of a short circuit (R L = 0) yields the transfer function magnitude.
Here, X and Y are the phasors of x(t) and y(t). Fitting such a curve to measured data respectively yields the short-circuit damping factor ζ e,sc and the electrical damping coefficient c e,sc . From this, an estimate K sc of the coupling coefficient is obtained from Equation (4) by This method requires knowledge of the mechanical damping coefficient ζ m , and therefore serves as additional validation of the open circuit method.
A third way of estimating the coupling coefficient is by evaluating the measured optimum load resistance R L,opt . The load voltage spectrum near resonance is measured with different load resistances, and the resistance resulting in the highest power at resonance is R L,opt [14][15][16]. An estimate of the coupling coefficient then follows from Equation (5): The procedure is time consuming because one has to determine the peaks of resonance curves for various load resistances. It does not suffice to just measure the output power at a fixed frequency (e.g., the open circuit resonance frequency) because the resonance frequency varies with the damping [20].
The fourth measurement method for the coupling coefficient is based on a linear variation of the current i through the coil [15]. The resulting magnetic force F mag on the permanent magnets is measured by the cantilever beam deflection x via the effective spring constant k. One obtains the following estimate of the coupling coefficient K i : This can be exploited fully automatically with a slowly rising and falling current. The slow change reduces measurement errors because the mass is slowly swinging. The spring constant can be taken from k = mω n 2 , which can be measured with a force testing system or can be extracted from simulations.
When the mechanical damping coefficient and the coupling coefficient are known, the optimum load resistance can be predicted by Equation (5). A comparison between the measured load resistance R L,opt and the calculated one, R L,opt,c , provides a measure of goodness of the measurement and the simulation procedures.

Energy Harvester Implementation
To compare the respective characteristics of the various measurement methods for the coupling coefficient, we applied the methods to a typical vibration energy harvester type. This type consists of a cantilever made from copper, of four NdFeB magnets, and of two iron legs to guide the magnetic flux ( Figure 4). The seismic mass at the end of the cantilever weighed m = 4.07g. The magnet grade was N50, with a remnant flux density of B r = 1.43T and a coercive field strength of H c = 955kA/m. The respective dimensions of the magnets and the iron legs were 10 mm × 5 mm × 1 mm and 10.5 mm × 10 mm × 1 mm. The copper beam was 0.5 mm thick and comprised a rectangular slot 14 mm long and 3 mm wide, in which an induction coil could freely move up and down. The shortest distances of the cantilever end and the clamping point to the slot were 3 mm and , respectively. Hence, the full cantilever length amounted to + 14mm + 3mm.

Energy Harvester Implementation
To compare the respective characteristics of the various measurement methods for the coupling coefficient, we applied the methods to a typical vibration energy harvester type. This type consists of a cantilever made from copper, of four NdFeB magnets, and of two iron legs to guide the magnetic flux ( Figure 4). The seismic mass at the end of the cantilever weighed  Table 1.
The absolute uncertainties of the magnet and cantilever geometries were about 0.1 mm, the outer coil dimension uncertainty amounted to 0.2 mm according to the producer, and the remaining coil dimension uncertainties were negligible. The residual flux density uncertainty of the magnets was 30 mT according to the data sheet. Finally, the relative uncertainties of the mass and of the coil resistance were 1%, caused by the scales and multimeters used.  Various implementations of the harvester type described were used in the experiments. The implementations differed in the free clamping lengths and the coils. Three different coils were tested. They all had identical outer dimensions (inner and outer diameters D i = 2.9mm and D o = 9.7mm, respectively, and height h = 1.49mm), but different wire diameters. Enameled copper wires with diameters D w between 30 and 50 µm were used, and the value of D w determined the number of coil turns N and the coil resistance R C . Further values are listed in Table 1. The absolute uncertainties of the magnet and cantilever geometries were about 0.1 mm, the outer coil dimension uncertainty amounted to 0.2 mm according to the producer, and the remaining coil dimension uncertainties were negligible. The residual flux density uncertainty of the magnets was 30 mT according to the data sheet. Finally, the relative uncertainties of the mass and of the coil resistance were 1%, caused by the scales and multimeters used.

Finite Element Simulation
The magnetic circuit with the four magnets and the coil was modeled by the FE program Ansys Electromagnetics (Figure 2; dimensions and magnet data as in Section 3.1). The FE mesh comprised over 10,000 tetrahedron elements. The magnetic flux Φ through the coils vanishes in the quiescent (equilibrium) state z = 0, where z denotes the relative displacement of the magnets and the coils, because the magnets generate a quadrupole field. In the simulation, the dependence of the magnetic flux on the displacement, Φ(z) was easily obtained by varying z. The effects of dimensional and magnet parameter deviations were studied in the same manner.

Measurement Details
On the one hand, we characterized harvesters with identical clamping lengths , different coils, and an identical harmonic acceleration amplitude ofâ = 0.75m/s 2 . On the other hand, we tested harvesters with identical coils, different clamping lengths, and an identical acceleration amplitude of a = 1m/s 2 .
The energy harvesting test bench consisted of a waveform generator connected to a permanent magnet shaker (B&K LDS V406, Naerum, Denmark), an accelerometer (B&K 4534-B-001) for measuring the vibration strength, and a voltage measurement card (NI PIXe-6341, Austin, TX, USA). The measurement card captured the load or open circuit voltages and the voltage from the acceleration sensor. The absolute movement was tracked with a triangulation laser sensor (Micro-Epsilon optoNCDT 2300, Ortenburg, Germany). The load resistance was set by a resistance decade box. The whole bench was controlled by a PC, so the user could vary the frequency and the acceleration amplitude. The center frequency and span of the measurement window in the frequency domain were respectively adjusted to the resonance frequency and a value larger than the resonance-curve bandwidth of the harvester to be characterized. The step width was set to 50 mHz. The voltage and movement resonance curves of all setups are contained in Supplementary Materials.
In the experiments, the actual acceleration amplitudes deviated from the preset values. If ignored, this would introduce errors, as both the load voltage and the cantilever tip displacement are proportional to the acceleration. Therefore, we converted the measured voltages and displacements to the values that would have resulted at the preset reference acceleration. The relative uncertainties of the measured voltages, displacements, and accelerations did not exceed 1%, and were neglected with respect to model and dimensional errors.
The load power was calculated by P L = U 2 /R L , with U being the load voltage. The optimum load resistance R L,opt is the resistance maximizing the power. R L was changed in steps by way of the resistance decade box.

Propagation of Uncertainty
The electromagnetic coupling coefficient and the optimum load resistance are deduced from other parameters with uncertainties as discussed above. We assume that these parameter uncertainties result from independent and normally distributed random errors. We may then write as where u K sim and u K sim are the absolute and relative uncertainties of the simulated coupling coefficient, respectively, x i is the uncertain input parameters, and S K sim x i = ∂K sim /∂x i is the sensitivity of K sim with respect to x i [23]. The remnant flux density B r , the gap width w g , and the outer coil dimension D o are likely to be the most uncertain influence quantities. We, therefore, use Equation (11) with these three parameters only (x 1 = B r etc.) and neglect the remaining parameters: The numbers for setups 1 to 3 are u B r = 0.03T, The uncertainty of the predicted optimum load resistance R L,opt,c is calculated with the simulated coupling coefficient K sim . The uncertainties of K oc and the open-and closed-circuit damping factor values ζ m and ζ 0 follow from the fit of Equations (6) and (7). The uncertainties for R L,opt,c , K sc , K R , and K i are derived from Equations. (5), (8), (9), and (10), respectively, using Equation (11) and c = 4π f n ζm: Figure 5 shows the influence of the load resistance in a harvester setup 1L at a vibration amplitude ofâ = 0.75m/s 2 . The voltage resonance curves for all measured load resistances are contained in Supplementary Materials. The load voltage amplitude U, the electrical damping factor ζ e , and the power P L delivered to the load at the resonance frequency all depend on the load resistance R L (Figure 5a). The expectations from Section 2.1 are met: for a low load resistance (e.g., R L = R C ), the damping factor ζ e is large and the displacement is small, thus limiting the load power. A high load resistance leads to a reduced damping factor, but this is not enough to offset the current reduction, which again limits the load power. The optimum load resistance R L,opt lies between the limiting cases of small and high values. Surprisingly, R L,opt is an order of magnitude higher than the coil resistance R C (Table 2).

Damping Influence and Optimum Load
which again limits the load power. The optimum load resistance L,opt R lies between the limiting cases of small and high values. Surprisingly, L,opt R is an order of magnitude higher than the coil resistance C R ( Table 2).
(a) (b)  The frequency dependence of the load power is in agreement with these observations (Figure 5b). Any load different from L,opt R leads to a lower load power at resonance, max P . The much higher damping associated with small L R can be used for frequency adaptivity: at vibration frequencies near (but not equal to) the harvester resonance frequency, the flat resonance curve for small L R effectively leads to more load power than the narrowband resonance curve for L,opt R .
The maximum load powers of our harvesters with the long free clamping lengths (setups 1L, 2, and 3) amounted to between 223 and 257 µW at = 2 0.75m/s a . In Section 2.1, it was mentioned that a higher number of coil turns N increases the coupling coefficient K . From Equation (3), we may write = K K C N with a constant K C independent of N . However, the load power does not depend on N when the coil dimensions remain unchanged, the load resistance is optimum, and the system is operated at resonance. To understand this, first note that the coil resistance is approximately  The frequency dependence of the load power is in agreement with these observations (Figure 5b). Any load different from R L,opt leads to a lower load power at resonance, P max . The much higher damping associated with small R L can be used for frequency adaptivity: at vibration frequencies near (but not equal to) the harvester resonance frequency, the flat resonance curve for small R L effectively leads to more load power than the narrowband resonance curve for R L,opt .
The maximum load powers of our harvesters with the long free clamping lengths (setups 1L, 2, and 3) amounted to between 223 and 257 µW atâ = 0.75m/s 2 . In Section 2.1, it was mentioned that a higher number of coil turns N increases the coupling coefficient K. From Equation (3), we may write K = C K N with a constant C K independent of N. However, the load power does not depend on N when the coil dimensions remain unchanged, the load resistance is optimum, and the system is operated at resonance. To understand this, first note that the coil resistance is approximately with a constant C R independent of N (cf. derivation in Appendix A). The time-averaged load power then becomes [12] Up until non-idealities and measurement errors, our measurement results demonstrate this independence of the power on the number of coil turns (Table 2). In contrast, the free clamping length clearly has an influence.

4.2.
Measuring K i Figure 6 shows the measured magnetic force F mag as a function of the coil current i for setup 2. The slope of the linear fit is the coupling coefficient K i from Equation (10). The measurement results are in excellent agreement with Equation (10). They also indicate a linear elastic behavior, or else the slope of the measured curve would not be constant.
with a constant R C independent of N (cf. derivation in Appendix A). The time-averaged load power then becomes [12] − Up until non-idealities and measurement errors, our measurement results demonstrate this independence of the power on the number of coil turns (Table 2). In contrast, the free clamping length  clearly has an influence. Figure 6 shows the measured magnetic force mag F as a function of the coil current i for setup 2. The slope of the linear fit is the coupling coefficient i K from Equation (10). The measurement results are in excellent agreement with Equation (10). They also indicate a linear elastic behavior, or else the slope of the measured curve would not be constant.

Measuring
(a) (b) 95% [23]. The scatter in the repeated measurements is much stronger than i K u would suggest, so there must be additional influences.
In fact, the large sample standard deviation is a consequence of the sensitivity of the deflection measurement to laser triangulation. The deflection of a cantilever element depends on the position of the element along the cantilever. Inadvertently focusing the laser beam on the wrong spot results in a deflection value error, which will translate to a coupling coefficient error. The same may be said about an inadvertently oblique (instead of normal) incidence of the laser beam on the cantilever surface. In our series of measurements, we removed and then re-installed the triangulation sensor before each measurement to get an idea of the magnitude of these effects. It is obvious from However, the measurement is more difficult for shorter clamping lengths, as Figure 6b shows for five single measurements for each coil setup. The black dots with error bars u K i , computed by Equation (10), mark single measurements, while the red error bars mark the overall uncertainty of the result of five subsequent measurements according to u K i,Σ = s K i ·t α,n−1 / √ n. Here, s K i is the sample standard deviation, n = 5 is the number of measurements, and t α,n−1 denotes the quantile of Student's probability distribution function with n − 1 degrees of freedom at a confidence level of 1 − α = 95% [23]. The scatter in the repeated measurements is much stronger than u K i would suggest, so there must be additional influences.
In fact, the large sample standard deviation is a consequence of the sensitivity of the deflection measurement to laser triangulation. The deflection of a cantilever element depends on the position of the element along the cantilever. Inadvertently focusing the laser beam on the wrong spot results in a deflection value error, which will translate to a coupling coefficient error. The same may be said about an inadvertently oblique (instead of normal) incidence of the laser beam on the cantilever surface. In our series of measurements, we removed and then re-installed the triangulation sensor before each measurement to get an idea of the magnitude of these effects. It is obvious from Figure 6b that the error increases with decreasing clamping lengths. It follows that Equation (10) substantially underestimates the error of the identified coupling coefficient at small clamping lengths.

Comparison and Discussion
The results of determining the electromagnetic coupling coefficient K by different methods are shown in Figure 7 and Table 3. For K i , the sample mean of five measurements was used at a confidence level of 95%, as explained in Section 4.2.
substantially underestimates the error of the identified coupling coefficient at small clamping lengths.

Comparison and Discussion
The results of determining the electromagnetic coupling coefficient K by different methods are shown in Figure 7 and Table 3. For i K , the sample mean of five measurements was used at a confidence level of 95%, as explained in Section 4.2.  . This clearly shows the appropriateness of all approaches. The maximum observed deviation was between 5% and 8%.
As to the uncertainties of the experimentally identified coupling coefficients, oc K and sc K showed the best performance with a relative uncertainty of 3% (average for all setups), followed by R K with about 4%. The large errors in i K for short clamping lengths were already discussed in Section 4.2.
As the different coupling-coefficient measurement methods agree well, one can select one's preferential method using other criteria. The open circuit method requires only one resonance curve to provide the coupling coefficient K , the mechanical damping factor ζ m , and the resonance frequency r f , which can be equated to the natural resonance frequency n f for small ζ m . Hence, the analysis of a harvester only takes a short time. This can be used to predict the optimum load resistance depending on the vibration strength. The mechanical damping increases with stronger vibration amplitudes [24,25], and therefore the optimum load decreases with the amplitude by Equation (5). Determining the optimum load for the full range of vibration amplitudes would take Figure 7. Comparison of the simulated and measured coupling coefficients of five harvester implementations and the measured and predicted optimum load resistance. Table 3. Simulated and measured coupling coefficients. The data are visualized in Figure 7.

Setup No.
Coupling Coefficient in Wb/m On average, the experimentally determined coupling coefficients deviate by about 3% from the simulated value K sim . This clearly shows the appropriateness of all approaches. The maximum observed deviation was between 5% and 8%.
As to the uncertainties of the experimentally identified coupling coefficients, K oc and K sc showed the best performance with a relative uncertainty of 3% (average for all setups), followed by K R with about 4%. The large errors in K i for short clamping lengths were already discussed in Section 4.2.
As the different coupling-coefficient measurement methods agree well, one can select one's preferential method using other criteria. The open circuit method requires only one resonance curve to provide the coupling coefficient K, the mechanical damping factor ζ m , and the resonance frequency f r , which can be equated to the natural resonance frequency f n for small ζ m . Hence, the analysis of a harvester only takes a short time. This can be used to predict the optimum load resistance depending on the vibration strength. The mechanical damping increases with stronger vibration amplitudes [24,25], and therefore the optimum load decreases with the amplitude by Equation (5). Determining the optimum load for the full range of vibration amplitudes would take much longer than the open circuit method. Knowing the mechanical damping coefficient for different vibrations amplitudes, one can predict R L,opt .
The short circuit method is perfect in theory but less suitable in practice because it requires knowledge of the mechanical damping coefficient and because the high damping under short circuit conditions is associated with a damped oscillation, which is not easy to measure.
The constant current method produces precise results and is easy to implement, but the calculation of the force from the cantilever beam deflection is error-prone, as discussed in Section 4.2. A direct force measurement would be preferable from this point of view, but, of course, is also more complex to implement.
As a result of the good agreement of the coupling coefficients estimated by various measurement methods or simulations, the prediction of the optimum load resistance by Equation (5) also agrees well with the measured optimum load resistance R L,opt , as can be seen in Table 2. The average and maximum deviations were 5% and 14%, respectively, in our experiments.

Summary
An important design parameter of electromagnetic vibration energy harvesters is their electromagnetic coupling coefficient, which influences the output power. There are various ways of experimentally identifying the coefficient, but the agreement of the respective methods and their advantages or disadvantages have not been considered until now in the open literature. In this work, we have experimentally tested four methods for the direct or indirect measurement of the coupling coefficient and have compared the results with data from FE simulations. The extracted coupling coefficients agree to within 3% among all methods tested. This validates both the measurement methods and the theory of electromagnetic coupling and optimum load resistance. This means that the coupling behavior can be predicted on the basis of simulation data alone and the optimum load resistance can be predicted on the basis of estimated (measured) mechanical damping coefficients. In this situation, the method of choice is the one that involves the least measurement effort. This is the open circuit method, which provides all necessary parameters within a short measurement time.

Conflicts of Interest:
The authors declared no potential conflicts of interest with respect to the research, authorship, or publication of this article. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A
An air coil (see Figure A1) consisting of a wire with a total length of tot , a resistivity of ρ el , and a constant cross-sectional area of A wire has the resistance The total length is calculated by the number of coil turns N and the mean length of one coil turn . With a wire diameter much smaller than the outer coil dimensions, and therefore with many coil turns, the total length is As measured by the outer coil dimensions, a single winding encloses the area The ratio between the winding area and the wire cross-sectional area is expressed by the filling factor F := A wire /A winding . One always has F < 1 because of imperfect wiring and the wire coating. Applying these equations to Equation (19), we arrive at the following expression for the coil resistance as a function of the coil turn number N: factor = wire winding : / F A A . One always has < 1 F because of imperfect wiring and the wire coating. Applying these equations to Equation (19), we arrive at the following expression for the coil resistance as a function of the coil turn number N : π + = ρ −