Dynamics and Control of a Magnetic Transducer Array Using Multi-Physics Models and Artificial Neural Networks

A linear mechanical oscillator is non-linearly coupled with an electromagnet and its driving circuit through a magnetic field. The resulting non-linear dynamics are investigated using magnetic circuit approximations without major loss of accuracy and in the interest of brevity. Different computational approaches to simulate the setup in terms of dynamical system response and design parameters optimization are pursued. A current source operating in baseband without modulation directly feeds the electromagnet, which consists commonly of a solenoid and a horseshoe-shaped core. The electromagnet is then magnetically coupled to a mass made of soft magnetic material and attached to a spring with damping. The non-linear system is described by a linearized steady-space representation while is examined for controllability and observability. A controller using a pole placement approach is built to stabilize the element. Drawing upon the fact that coupling works both ways, enabling estimation of the mass position and velocity (state variables) by processing the induced voltage across the electromagnet, a state observer is constructed. Accurate and fast tracking of the state variables, along with the possibility of driving more than one module from the same source using modulation, proves the applicability of the electro-magneto-mechanical transducer for sensor applications. Next, a three-layer feed-forward artificial neural network (ANN) system equivalent was trained using the non-linear plant-linear controller-linear observer configuration. Simulations to investigate the robustness of the system with respect to different equilibrium points and input currents were carried out. The ANN proved robust with respect to position accuracy.


Introduction
Advancements in the field of magnetic materials in terms of better efficiency and energy densities [1], together with the integration of mechanics, electromagnetics, power and control electronics into the system, have enabled intelligent electromagnetic actuators and sensors. Electromechanical devices are fundamentally based on quasi-static magnetic and electric fields, where force is generated [2]. Widely spread nowadays due to their simple design and system integration, they can measure changes in the magnetic reluctance of the circuit, which can be correspondingly related to a change in speed or position [3]. Precise operation can be realized using control techniques such as feedback control, which enables them both as sensors and actuators [4]. Typical applications can be found in aerospace [5], automotive [6], ocean [7] and biomedical devices [8] as well as flywheel and transmission speed [9], crank and cam shaft position for engine timing, throttle valve position for air intake, steering wheel position, pedal position, fluid level, chassis height, and in electronic door locks [6]. Given the joint sensing-actuation functionality, electro-magneto-mechanical devices are a promising solution for multi-agent architectures. Such architectures have a vast application field, such as mechatronics and robotics [10], smart grids [11] as well  Figure 1a,b show the electromechanical oscillator investigated in this work and in [25] respectively.

Lagrangian Formulation of the System
Analysis of dynamic systems requires the equations of motion for such systems to be derived. The ability of computation and analytical approaches depends on the mathematical description used to formulate system dynamics. Often Newtonian mechanics is replaced by Lagrangian dynamics, which solely depends on the energy balance within a non-dissipative system: Lagrange's equation can be augmented, using the so-called "power" function, to include dissipative forces as well as other non-dissipative, non-conservative external forces: Let us now consider the electromechanical system shown in Figure 1a. In specific, the system's circuit consists of an electric current source, a conductance, a capacitor and the Sensors 2021, 21, 6788 3 of 28 electromagnet's coil all connected in parallel. The Lagrangian of the system can then be written as: By substituting (4)- (8) into (3): (9) For the mechanical subsystem, the Lagrange equation can be written as: From the right-hand side of the above, the electromagnetic force applied to the mass by the electromagnet is obtained as follows.
∂x (11) The system has two degrees of freedom, thus for the electromagnetic part: (12) Note that for C = 0: i c = C dV dt = 0 equivalent to an open circuit.

Electromagnetic Subsystem
In this work the magnetic circuit solution is adopted, based on Maxwell's equations and the correspondence of the differential equations for linear magnetization problems with those for steady current problems. This is only true in the special case where boundary conditions for the two cases are identical [26]. In order to simplify the analysis of magnetic circuits and make them suitable and accurate enough for engineering applications, assumptions must be made. Initially, the terms of the Maxwell equations containing timevarying electric fields are neglected and the system is regarded as magneto quasi-static. The second assumption refers to the magnetic circuit, which is thought of as a structure of high-permeability materials. This results in the magnetic flux to flow entirely through the path defined by the structure.
An electromagnet featuring a solenoid around a horseshoe-shaped core, an air gap and a moving mass is shown in Figure 2.  An electromagnet featuring a solenoid around a horseshoe-shaped core, an air and a moving mass is shown in Figure 2. The magnetic reluctances ℜ of the circuit are defined as: The magnetomotive force (mmf) drives a magnetic flux through the magnetic luctances of the core, air gap and moving mass. From (13), (14) and (15) it is clear that magnetic reluctance depends on the material permeability. High magnetic permeabi can result in small reluctances, in fact much smaller than that of the air gap. The result magnetic field stores energy in the air gap equal to: The force in the gap acting on the mass is given by the rate of change of energy w gap length, thus: It is known that: thus we can write (17) as: The magnetic reluctances of the circuit are defined as: The magnetomotive force (mmf) drives a magnetic flux Φ through the magnetic reluctances of the core, air gap and moving mass. From (13), (14) and (15) it is clear that the magnetic reluctance depends on the material permeability. High magnetic permeability can result in small reluctances, in fact much smaller than that of the air gap. The resulting magnetic field stores energy in the air gap equal to: The force in the gap acting on the mass is given by the rate of change of energy with gap length, thus: It is known that: thus we can write (17) as: Sensors 2021, 21, 6788 5 of 28 or alternatively: The force, acceleration, velocity and displacement are vectors. The electromagnetic force is acting on the mass minimizing the gap length.
The dynamic equations for the electrical subsystem must be derived as well. Knowledge of the induced voltage and variable inductance are key for both sensory and actuating applications.
It is known that: and: We assume: and (27) can now be written as:

Electromechanical System Simulation
The dynamical electromechanical system described in this work is simulated using Matlab [27], in particular: 2.3.1. Simscape Implementation Matlab Simscape uses physical blocks to model complex systems (Figure 3). Each block models an item which is interconnected with others, forming a diagram equivalent to the mathematical model of the system. The interaction of the blocks is accomplished through signal ports permitting non-directional energy exchange. Additional signal ports are available to the blocks in order to specify initial conditions and enhance the model's computational speed. Blocks are available for both active and passive elements depending on their ability to deliver energy to the system, store it, or dissipate it. Results of the simulated system are shown in Figure 4.
• Matlab Simulink; 2.3.1. Simscape Implementation Matlab Simscape uses physical blocks to model complex systems (Figure 3). Each block models an item which is interconnected with others, forming a diagram equivalent to the mathematical model of the system. The interaction of the blocks is accomplished through signal ports permitting non-directional energy exchange. Additional signal ports are available to the blocks in order to specify initial conditions and enhance the model's computational speed. Blocks are available for both active and passive elements depending on their ability to deliver energy to the system, store it, or dissipate it. Results of the simulated system are shown in Figure 4.

Simulink Implementation
Matlab Simulink allows dynamical modelling using a graphic environment throughout signal blocks. These blocks are complete only when input and output signals are fully defined. The magnetic force function is built and integrated into the Simulink model using the approximation described in (29) for the inductance.

Simscape Implementation
Matlab Simscape uses physical blocks to model complex systems (Figure 3). Each block models an item which is interconnected with others, forming a diagram equivalent to the mathematical model of the system. The interaction of the blocks is accomplished through signal ports permitting non-directional energy exchange. Additional signal ports are available to the blocks in order to specify initial conditions and enhance the model's computational speed. Blocks are available for both active and passive elements depending on their ability to deliver energy to the system, store it, or dissipate it. Results of the simulated system are shown in Figure 4.

Simulink Implementation
Matlab Simulink allows dynamical modelling using a graphic environment throughout signal blocks. These blocks are complete only when input and output signals are fully defined. The magnetic force function is built and integrated into the Simulink model using the approximation described in (29) for the inductance.

Simulink Implementation
Matlab Simulink allows dynamical modelling using a graphic environment throughout signal blocks. These blocks are complete only when input and output signals are fully defined. The magnetic force function is built and integrated into the Simulink model using the approximation described in (29) for the inductance. Figures 5 and 6 shows the air gap length and the Simulink coil voltage over a simulation time of 10 s, for the implementation in Figure 5. Initially the output voltage is equal to zero, since the system is in equilibrium (no spring deformation). When a current i is applied to the circuit an attractive force is pulling the mass toward the electromagnet (Figure 2), changing the air gap length thus inducing voltage on the inductance. When the spring force overcomes the magnetic force the mass moves away from the electromagnet.   Figure 5. Initially the output voltage is equal to zero, since the system is in equilibrium (no spring deformation). When a current i is applied to the circuit an attractive force is pulling the mass toward the electromagnet ( Figure  2), changing the air gap length thus inducing voltage on the inductance. When the spring force overcomes the magnetic force the mass moves away from the electromagnet. Motion on the positive x direction induces a positive voltage (negative x direction motion results in negative induced voltage).

Methods Comparison
Next the two methods described above to perform the dynamic simulation of the system are compared. Figure 7 confirms that the results are comparable hence both methods are suitable. lation time of 10 s, for the implementation in Figure 5. Initially the output voltage is equal to zero, since the system is in equilibrium (no spring deformation). When a current i is applied to the circuit an attractive force is pulling the mass toward the electromagnet ( Figure  2), changing the air gap length thus inducing voltage on the inductance. When the spring force overcomes the magnetic force the mass moves away from the electromagnet. Motion on the positive x direction induces a positive voltage (negative x direction motion results in negative induced voltage).

Methods Comparison
Next the two methods described above to perform the dynamic simulation of the system are compared. Figure 7 confirms that the results are comparable hence both methods are suitable.

Methods Comparison
Next the two methods described above to perform the dynamic simulation of the system are compared. Figure 7 confirms that the results are comparable hence both methods are suitable.

State-Space Problem Formulation
The state-space model is defined in terms of the derivatives of the states and its representation is given below: In order to formulate the state equations of the electro-magneto-mechanical oscillator shown in Figure 1, constraints limiting the system according to its physical and performance characteristics must be defined. For the inductance, we assume the expression described by (29), derived using the magnetic circuit approximation. In favor of simplicity: Figure 1 shows the origin of the mechanical coordinate (x), assumed equivalent to the neutral position of the spring while that of the electromagnetic system (z) at the position where the mass is marginally in contact with the magnet armature, such that: thus: Based on the expression derived for the electromagnetic force earlier where the magnetic flux and consequently the current appears squared, one can easily derive that the force will be unidirectional, and always attractive, regardless of the electric current's sign: Additionally, given that we want to avoid any physical contact in the system, the velocity of a moving part must be controlled and constrained as follows [28]:

State-Space Problem Formulation
The state-space model is defined in terms of the derivatives of the states and its representation is given below: In order to formulate the state equations of the electro-magneto-mechanical oscillator shown in Figure 1, constraints limiting the system according to its physical and performance characteristics must be defined. For the inductance, we assume the expression described by (29), derived using the magnetic circuit approximation. In favor of simplicity: where Figure 1 shows the origin of the mechanical coordinate (x), assumed equivalent to the neutral position of the spring while that of the electromagnetic system (z) at the position where the mass is marginally in contact with the magnet armature, such that: thus: Based on the expression derived for the electromagnetic force earlier where the magnetic flux and consequently the current appears squared, one can easily derive that the force will be unidirectional, and always attractive, regardless of the electric current's sign: Additionally, given that we want to avoid any physical contact in the system, the velocity of a moving part must be controlled and constrained as follows [28]: ..
with states: the mass potion x, the mass velocity .
x and the magnetic flux ψ. Assuming the voltage across the inductance to be the systems' output: (39)-(41) represent the dynamics of a nonlinear system thus linearization is applied around some equilibrium point x 0 , .
and for the equilibrium point .
We define as state variables the displacement x, the velocity .
x, the magnetic flux ψ and adopt the variable (30) is now written as: and: The partial derivatives described in (43)-(45) to linearize the system are calculated using both the analytical and numerical methods around an equilibrium point point x 0 = 0.015 m (Table 1). Table 1. Analytical and numerical values for the linearization partial derivatives.

Partial Derivative
Analytical Numerical  Table 2 summarizes the state variables, inputs and outputs of the linearized system: Table 2. State variables, inputs and outputs summary for the linearized system.

System Controllability and Observability
A linear time invariant system . x = Ax + Bu, x(0) = x 0 is controllable if over the interval 0, t f a control input u(t) ∀ t ∈ 0, t f steers the state from x init to x t f according to: or in other words an LTI system is controllable at time t f > 0 if for any initial state and for any target state x t f , a control input u(t) exists that can steer the system from x(0) to x t f over the defined interval. It can be shown that The observability matrix O has full rank since rank (O) is equal to number of state variables.
Let's now consider the linearized system as described in Section 3 with state variables the spring displacement x and the velocity .
x, and output the inductance voltage V L . As stated earlier, the system contains an unstable pole with a corresponding damping ratio of −1. A negative damping ratio, called driving force, increases the system oscillation rather than driving the system to stability (damping force).

Linear Controller Design Using Pole Placement
To stabilize the system around a specific position, a feedback controller is designed using poles placement (Figure 8). In order to implement full-state feedback all state variables must be known to the controller at all times.

variables.
Let's now consider the linearized system as described in Section 3 with state variables the spring displacement and the velocity , and output the inductance voltage . As stated earlier, the system contains an unstable pole with a corresponding damping ratio of -1. A negative damping ratio, called driving force, increases the system oscillation rather than driving the system to stability (damping force).

Linear Controller Design Using Pole Placement
To stabilize the system around a specific position, a feedback controller is designed using poles placement (Figure 8). In order to implement full-state feedback all state variables must be known to the controller at all times The system specific closed loop matrix is now equal to ( − ) and the gain matrix is a 1 × 2 matrix hence, The controller equation is:

Linear Controller
where u is the system input, K gain the negative state-feedback gain matrix and x the state variables, and: where: The system specific closed loop A cl matrix is now equal to A − BK gain and the gain K gain matrix is a 1 × 2 matrix hence, The eigenvalues of A cl are given by det(A cl − λI): Let us now assume that we would like to place two poles p 1 and p 2 such that the characteristic equation is: The gain matrix K gain can now be calculated by setting the coefficients equal to each other. The location of the poles p 1 and p 2 characterize the stability and time-domain performance of the system.
The controller design and the selection of the gain matrix K gain are based on two criteria: settling time < 0.5 sec For an overshoot less than 5% we have: and for a settling time (3% error) of less than 0.5 s: For a second order system with ζ < 1 the poles are given by: The closed loop system showed in Figure 8 is simulated using the following three alternative pole sets:  respectively. The controller is based on the linearized system while the plant is the non-linear system described by (39)-(41). The current applied to the plant features a small delta with respect to the equilibrium current, δi s0 = 0.01 A, to remain in the region where the linearization is valid. Figures 9-11 show the controlled mass displacement for different equilibrium points using poles set as p 31 , p 32 and p 33 .

Observer Design
The pole placement method described earlier requires knowledge of the system states at all times. In practice, this is very demanding due to unavailable or expensive sensors often driven by high noise as well. Moreover, some types of application require the system to operate as a sensor rather than as an actuator. Given the proven observability of the studied system, a state observer is designed to overcome the issue. A state observer is capable of estimating and observing the system states regardless of whether some state variables are available for direct measurement or not. The only measurable output for the hereby-described electromechanical oscillator is the inductance voltage. The observer mimics the model of the plant, except for an additional term, introduced to compensate for inaccuracies in the matrices A and B. This new term, called estimation error, is equal to the difference between the estimated output value and the measured value. The observer equations are: wherex the estimated states,ŷ the estimated output and L the observer gain matrix. The system's output is compared to the estimated output. The resulting error is corrected to help the observer converge to the system state values. The observer's error dynamic is given by: The input to the open-loop plant is equal to: Thus: and: . x (88) shows that the state dynamics is decoupled from the observer's error e (separation principle) thus the gain K can be computed independently of the observe gain L. The gain matrix K gain is computed, for p 31 = −30 + i30, p 32 = −30 − i30, p 33 = −150 and is equal to K gain = [262.01 8.23 − 351.22]. For the observer, we select eigenvalues located at least five times further left on the complex plane than the plant eigenvalues, since we want the dynamics of the observer to be much faster than the system itself. In other words, the observer must converge to the plants' states before they converge to zero. Hence we choose: p obs1 = −150 + 6i, p obs2 = −150 − 6i and p obs3 = −750. Poles placed further apart than that could result in complex systems with large bandwidth, making the system vulnerable to noise effects. The analytical procedure to derive the observer matrix is identical to the one described for the gain matrix K gain . Alternatively, one can use the place command in Matlab and: L = −9.405 × 10 3 1.914 × 10 3 − 1.29 3.23 × 10 3 . The selection of matrix L is a tradeoff between a fast dynamic response and noise reduction. Figure 12 shows the non-linear plant-linear controller-linear observer block diagram. Figure 13 shows theMatlab Simulink system implementation. The states are now unknown to the controller, hence the observer is used to generate estimates by comparing the plant output to the observer output. These estimates of the states are fed into the linear controller to stabilize the system around an equilibrium point. The observer is a virtual system, designed to approximate the state variables not available for sensing. A properly designed observer has all initial conditions specified. Despite there isn't a standard rule, for this study are chosen such that the initial error is equal to e = [0.004, −0.04, −0.07]. Responses of the state variables, for different equilibrium points are depicted in Figures 17-19: it is clear that the estimate states converge to the actual state variables while tracking reasonably well the equilibrium value. rs 2021, 21, x FOR PEER REVIEW 16 the plant output to the observer output. These estimates of the states are fed into the l controller to stabilize the system around an equilibrium point.          The observer is a virtual system, designed to approximate the state variables not available for sensing. A properly designed observer has all initial conditions specified. Despite there isn't a standard rule, for this study are chosen such that the initial error is equal to = [0.004, −0.04, −0.07]. Responses of the state variables, for different equilibrium points are depicted in Figures 17-19: it is clear that the estimate states converge to The observer is a virtual system, designed to approximate the state variables available for sensing. A properly designed observer has all initial conditions speci Despite there isn't a standard rule, for this study are chosen such that the initial err equal to = [0.004, −0.04, −0.07]. Responses of the state variables, for different equ rium points are depicted in Figures 17-19: it is clear that the estimate states converg the actual state variables while tracking reasonably well the equilibrium value.

Artificial Neural Network (ANN) Representation of the System
An ANN is an adaptive, often non-linear system, that learns to perform a function (an input/output map) when an input and the corresponding desired or target response set is presented to the untrained ANN [29].
An ANN is structured in three layers ( Figure 20): • the input layer, where there is no real processing done, is essentially a "fan-out" layer where the input vector is distributed to the hidden layer; • the hidden layer, being the computational core of the ANN; • the output layer, which combines all the "votes" of the hidden layer.

Artificial Neural Network (ANN) Representation of the System
An ANN is an adaptive, often non-linear system, that learns to perform a function (an input/output map) when an input and the corresponding desired or target response set is presented to the untrained ANN [29].
An ANN is structured in three layers ( Figure 20): • the input layer, where there is no real processing done, is essentially a "fan-out" layer where the input vector is distributed to the hidden layer; • the hidden layer, being the computational core of the ANN; • the output layer, which combines all the "votes" of the hidden layer. The general mathematical expression that describes an ANN is [29]: is multiplied by the weight w_n. Second, the weighted input _ is added to the scalar bias b to form the net input n. The bias is similar to a weight, except that it has a constant input of 1. Finally, the net input is passed through the The general mathematical expression that describes an ANN is [29]: First, the input x ANN is multiplied by the weight w _n . Second, the weighted input xw _ n is added to the scalar bias b ANN to form the net input n. The bias is similar to a weight, except that it has a constant input of 1. Finally, the net input is passed through the transfer function f ANN , which produces the output y ANN .
The equation can be more detailed and further expanded to [29]: with p the number of input nodes, M the number of hidden nodes and K the number of output nodes.
The activation function f ANN is commonly chosen to be one of the following: • the logistic sigmoid function (Figure 6), commonly abbreviated as logsig, • the hyperbolic tangent function (Figure 6), commonly abbreviated as tansig, The change in error due to output layers weight can be found by Least Mean Square (LSM) algorithms. The output error is given by [29]: The back-propagating error corresponding to the hidden layer output can be expressed as the partial derivative [29]: Successful training of the network requires calibration of the weights to counteract the error gradient to minimize the output error according to [29]: w_ new ji = w_ old ji + ∆w _ji (100)

ANN System Configuration
In this section, a feed-forward network organized in layers is proposed to approximate the electro-magneto-mechanical system. The input layer has three neurons, the hidden layer 10 neurons while the output layer has three neurons. The workflow design has the following steps: •

. Generation of ANN Data
The training data is generated using the system described in Section 4. In order to provide some variance to the system, the input current is modelled as: The training data set for the ANN correspond to the source current (input), the inductance voltage (input), the simulation time (input) and the state variables (output).

ANN Implementation and Training
The Matlab Neural Network Toolbox Mathworks (The MathWorks, Inc.) is used to construct the network. The selected training algorithm is the Levenberg-Marquardt with the following samples set aside for: • training, 70% of the data; • validation, 15% of the data; • and the remaining 15% for testing.

ANN Simulation of the Non-Linear System-Linear Controller-Observer
The system described in Section 4.2 is now been simulated by the trained ANN analogue. Figure            The dashed line (Figures 22, 24 and 26) represents the ideal case where the ANN simulated results are an exact replica of the targets. The solid line represents the best-fit linear regression line between trained net outputs and targets, which in our case is adjacent to the ideal case. All cases indicate a good fit, with the scatter manifesting the data points exhibiting a poor fit. R correlates the targets to the simulated net outputs. If R = 1, there is an exact linear relationship between outputs and targets. If R = 0, no linear rela- The dashed line (Figures 22, 24 and 26) represents the ideal case where the ANN simulated results are an exact replica of the targets. The solid line represents the best-fit linear regression line between trained net outputs and targets, which in our case is adjacent to the ideal case. All cases indicate a good fit, with the scatter manifesting the data points exhibiting a poor fit. R correlates the targets to the simulated net outputs. If R = 1, there is an exact linear relationship between outputs and targets. If R = 0, no linear relationship between outputs and targets can be inferred. All three cases, for different equilibrium points, show a good approximation of the non-linear plant-linear controller-linear observer system output by the ANN (Figures 23, 25 and 27).

Conclusions
The dynamics of an electro-magneto-mechanical system was investigated in this work. The coupled system simulation was developed and comparatively assessed using different modelling techniques in Matalb/Simulink and Matlab/Simscape. Next, the nonlinear state-space problem formulation and the linearization method using partial derivatives were presented. The resulting linearized system proved both controllable and observable, yet marginally stable. To stabilize and control the system, pole-placement with full-state feedback is employed. To implement the full-state feedback a state observer is introduced and designed. Using the models developed and the simulation results, an ANN equivalent has been configured and tuned for the non-linear plant-linear controller-linear observer system. The model is flexible and accurate; specifically, it is able to adapt to different equilibrium conditions (equilibrium position was varied in the range 0.005 m and 0.015 m) and compensates satisfactorily for changes in current source input (sinusoidal reference input). The system may operate as a contactless displacement sensor or actuator, or in combined mode. In ongoing work, the connection of more than one transducer to the same source will be demonstrated using amplitude modulation. The tuned frequency of each transducer will be defined by assigning different values to the capacitor as shown in Figure 1b. Preliminary analysis shows that such transducer arrays can be used in maritime and aerospace applications, where avoiding physical contact is often a sine qua non requirement, to achieve coordinated sensing and actuation of e.g., space structures like in an orbital station to detect and mitigate the vibrational effects of hits by space debris; or in the maritime, for energy harvesting without violating hull integrity with penetrations to implement physical contact. Other applications are possible too, especially if promising theoretical tools like feedback linearization (instead of the local linearization demonstrated here) are employed.