A Amor

A Fernandes
email: amancio.fernandes@sorbonne-universite.fr

J Pouget

C Maurini

Nonlinear dynamics and snap-through regimes of a bistable buckled beam excited by an electromagnetic Laplace force

Keywords: Bistable beam, buckling, snap-through, electromagnetic actuation, bifurcation, nonlinear dynamics, Laplace force

Nonlinear dynamics and snap-through regimes of a bistable buckled beam excited by an electromagnetic

Laplace force

Introduction

Bistable mechanisms are at the basis of many advances engineering applications [START_REF] Cao | Bistable structures for advanced functional systems[END_REF][START_REF] Zhang | Bistable morphing composite structures: A review[END_REF], including Micro-Electro-Mechanical-Systems [Charlot et al., 2008, Pane andAsano, 2008], micro-robotics, non-volatile memories, medical endoscopy (guiding mechanism in small diameter blood vessels) [START_REF] Schomburg | Design optimization of bistable microdiaphragm valves[END_REF]], tactile surfaces for visually impaired people [START_REF] Chouvardas | Tactile displays: Overview and recent advances[END_REF][START_REF] Vitushinsky | Bistable thin-film shape memory actuators for applications in tactile displays[END_REF]. Their nonlinear dynamical properties are widely exploited for energy harvesting [Cottone et al., 2009a, Emam and Inman, 2015[START_REF] Van Blarigan | Dynamic analysis of a buckled asymmetric piezoelectric beam for energy harvesting[END_REF][START_REF] Zhu | Enhanced buckled-beam piezoelectric energy harvesting using midpoint magnetic force[END_REF] and shape control applications [START_REF] Fernandes | Multiparameter actuation for shape control of bistable composite plates[END_REF], Hagood et al., Schoeftner et al., 2015].

Buckled beams constitute a simple and effective way to obtain bistable devices at different length-scales. The bistability is obtained by applying an axial pre-stress to an initial straight beam.

If the pre-stress is sufficiently large, the beam shows two stable equilibria, where the beam bends upwards or downwards. When the beam is actuated by a transverse force, it can snap between the two states. In the case of a dynamic actuation, different regimes are possible, including intrawell oscillations around each of the two stable equilibrium configurations, or chaotic cross-well 1 dynamics with persistent snapping. Characterizing the conditions leading to the different regimes is fundamental for shape-control or energy harvesting applications. In this work, we investigate the nonlinear dynamics of a bistable buckled beam through numerical and experimental approaches.

Considering the case of harmonic load, we determine the conditions on the forcing amplitude and frequency leading to the different dynamic regimes.

Our original experimental setup exploits electromagnetic Laplace forces to produce a finely tunable actuation of the mechanical structure. Extending to the dynamic context the ideas presented in [START_REF] Amor | Snap-through of elastic bistable beam under contactless magnetic actuation[END_REF][START_REF] Amor | Numerical and experimental investigations of bistable beam snapping using distributed laplace force[END_REF], we control the forcing frequency and amplitude by generating an electric current inside a bistable metallic beam placed between two permanent magnets. The resulting Laplace force is a follower force constantly normal to the beam's axis in the deformed configuration. It is well-known that the Laplace force can induce subtle instability for a wire in a magnetic field [START_REF] Healey | Large rotating states of a conducting elastic wire in a magnetic field: subtle symmetry and multiparameter bifurcation[END_REF][START_REF] Valverde | Magnetically-induced buckling of a whirling conducting rod with applications to electrodynamic space tethers[END_REF][START_REF] Wolfe | Equilibrium states of an elastic conductor in a magnetic field: A paradigm of bifurcation theory[END_REF]]. Yet, we were not able to find other recent studies using this effect for shape control. Other contactless actuation devices proposed in the literature exploit different kinds of actuating forces such as electrostatic forces which are effective only at microscopic scale [START_REF] Chen | Snap-through buckling of initially curved microbeam subject to an electrostatic force[END_REF][START_REF] Krylov | Bistability of curved microbeams actuated by fringing electrostatic fields[END_REF][START_REF] Park | Pre-shaped buckled-beam actuators: Theory and experiments[END_REF][START_REF] Younis | Nonlinear dynamics of mems arches under harmonic electrostatic actuation[END_REF] with extension to dynamics of arch-shaped beams [START_REF] Krylov | Dynamic stability of electrostatically actuated initially curved shallow micro beams[END_REF][START_REF] Ouakad | Static response and natural frequencies of microbeams actuated by out-of-plane electrostatic fringing-fields[END_REF][START_REF] Ramini | Theoretical and experimental investigation of the nonlinear behavior of an electrostatically actuated in-plane mems arch[END_REF][START_REF] Wu | Excitation-Induced Stability in a Bistable Duffing Oscillator: Analysis and Experiments[END_REF]. Micro-beams actuated by electro-thermo-mechanical effect or by shape memory alloys heated by laser are reported in [START_REF] Barth | A bistable shape memory microswitch with high energy density[END_REF][START_REF] Zaidi | Contactless and selective energy transfer to a bistable micro-actuator using laser heated shape memory alloy[END_REF]. Our method based on Laplace forces provides an excellent control of the beam dynamics for a beam of 20cm length.

Our modeling approach is based on the one-dimensional elastic beam theory [START_REF] Chen | Static snapping load of a hinged extensible elastica[END_REF][START_REF] Ouakad | Static response and natural frequencies of microbeams actuated by out-of-plane electrostatic fringing-fields[END_REF][START_REF] Goss | The History of the Planar Elastica: Insights into Mechanics and Scientific Method[END_REF][START_REF] Magnusson | Behaviour of the extensible elastica solution[END_REF][START_REF] Patrício | An elastica problem: instabilities of an elastic arch[END_REF]] within the hypothesis of moderate rotation [START_REF] Neukirch | Comparison of the Von Kármán and Kirchhoff models for the post-buckling and vibrations of elastic beams[END_REF]. To correctly describe the dynamics of bistable beam snap-through process, we take into account the effect of the beam extensibility [START_REF] Camescasse | Bistable buckled beam: Elastica modeling and analysis of static actuation[END_REF][START_REF] Neukirch | Vibrations of post-buckled rods: The singular inextensible limit[END_REF]. Hence, we adopt a Galërkin approximation on the beam buckling modes for the space-discretisation and we integrate the resulting system of the ordinary differential equation with standard time-stepping schemes, as in [START_REF] Cazottes | Bistable Buckled Beam: Modeling of Actuating Force and Experimental Validations[END_REF][START_REF] Qiu | A curved-beam bistable mechanism[END_REF][START_REF] Vangbo | An analytical analysis of a compressed bistable buckled beam[END_REF]]. While most of the available works retain a single buckling mode in the Galërkin approximation [START_REF] Chandra | A numerical investigation of snap-through in a shallow arch-like model[END_REF], Cottone et al., 2009b[START_REF] Erturk | A piezomagnetoelastic structure for broadband vibration energy harvesting[END_REF][START_REF] Medina | Dynamic release condition in latched curved micro beams[END_REF], we show that at least three modes should be included in the reduced model to correctly approximate the nonlinear snapping dynamics.

The paper proceeds as follows. Section 2 formulates the weakly nonlinear model for the beam dynamics. We progressively introduce several simplifying assumptions starting from a variational formulation of a nonlinear extensible beam model. Section 3 reports the numerical results for the non-linear forced dynamics based on the reduced modal model. The results obtained from the numerical simulations include the identification of the possible snapping regimes and the analysis of their dependence on the forcing parameters. Our results clearly show the classical transition from the regular to a chaotic behavior. Section 4 presents the experimental results and reports a fine comparison between the experimental finding and the numerical result obtained with the reduced models. Section 5 draws the conclusions.

2 Beam model

Kinematics

The bistable mechanism considered in the present study consists of an isotropic and homogeneous elastic simply-supported beam subject to an end-shortening ∆L, see Figure 1. In its natural configuration the beam is straight, with length L, and has a rectangular cross-section of width b and thickness h. When the end-shorting overtakes a critical value ∆L c , the straight configuration becomes unstable. The beam buckles and two possible stable equilibrium configurations can be observed, with upwards or downwards bending. We formulate here the beam model used in the rest of this paper. A fixed Cartesian reference frame R 0 : {A; e 1 , e 2 , e 3 } is attached to the straight reference configuration, with the unit vector e 1 in the direction of the beam axis. The beam is regarded as an elastic curve that coincides with the neutral axis of the beam. A material point G 0 of the beam in the reference configuration is given by the vector

-→ AG 0 = q 0 (s) = s e 1 with s ∈ [0, L],
where s is the abscissa along the beam axis. The current placement G of the point G 0 in the reference configuration is

-→ AG(s) = q(s) = x(s) e 1 + y(s) e 2 = (s + u(s)) e 1 + v(s) e 2 , s ∈ [0, L] , (1)
where u(s) = u(s) e 1 + v(s) e 2 is the displacement vector.

We denote by τ (s) = d q(s)/ds d q(s)ds = cos θ(s) e 1 + sin θ(s) e 2 , n(s) =sin θ(s) e 1 + cos θ(s) e 2 ,

the tangent and normal vectors to the current configuration at the point of abscissa s in the reference configuration, where θ is the angle formed by the tangent with e 1 . Hence, adopting an unshearable beam model, we characterize the beam deformation by the axial stretch ε and curvature κ defined by

ε(s) = ds -ds ds = d q(s) ds -1, κ(s) = dθ(s) ds , (3)
where we denote by s the arc-length in the current configuration. Expanding the expression d q(s)/ds = (1 + ε(s)) τ (s) for small angles θ, we get

u (s) = x (s) -1 = (1 + ε(s)) cos θ(s) -1 = ε(s) + (1 + ε(s)) - θ(s) 2 2 + θ(s) 4 24 + o(θ(s) 5 ) (4) v (s) = y (s) = (1 + ε(s)) sin θ(s) = (1 + ε(s)) θ(s) - θ(s) 3 6 + o(θ(s) 4 ). (5)
Here and henceforth, we denote by a f the derivative of a function f with respect to the abscissa s. Under the approximation of moderate rotations and infinitesimal but not vanishing extension (ε ∼ θ 2 ), we keep the following approximation for the displacement vector as a function of the rotation and the extension

u (s) = ε(s) -(1 + ε(s)) θ(s) 2 2 + θ(s) 4 24 , v (s) = θ(s). (6)
In particular, we neglect the third order terms for the transverse displacement v (s). This greatly simplifies the formulation in the dynamical settings without significantly compromising the geometrical accuracy for slender beams [START_REF] Neukirch | Comparison of the Von Kármán and Kirchhoff models for the post-buckling and vibrations of elastic beams[END_REF].

Variational formulation of the equilibrium

For a linear elastic and isotropic material behavior, the elastic energy density per unit line reads as

ϕ(ε, κ) = EA 2 ε 2 + EI 2 κ 2 , ( 7
)
where E is the Young modulus of the material, A and I are the cross-sectional area and moment of inertia. We will consider a beam with a rectangular cross-section of thickness h and width b, for which A = hb and I = bh 3 /12. We define the normal force and the bending moment by

N := ∂ϕ(ε, κ) ∂ε = EA ε, M := ∂ϕ(ε, κ) ∂κ = EI κ. (8)
The total elastic energy can be written as a function of the transversal displacement v and the extension ε as follows, where we use that κ(s) = θ (s) = v (s):

E(ε, v) = L 0 EA 2 ε(s) 2 + EI 2 v (s) 2 ds. (9)
Adopting the principle of the minimum of the potential energy, the equilibrium configurations are found as the solution of the following constrained minimization problem:

min{E(ε, v), (ε, v) : v(0) = 0, v(L) = 0, u(0) = 0, u(L) = -∆L}. (10)
Integrating (6), one gets

u(L) = L 0 ε(s) -(1 + ε(s)) v (s) 2 2 + v (s) 4 24 ds. (11)
Using the Lagrange multiplier method to impose the nonlinear constraint on the end-shortening given by L -∆L = L 0

x (s) ds yielding u(L) = -∆L, and denoting by P the corresponding multiplier, (10) can be reformulated as:

min v∈C,ε∈L 2 (R 3 ) max P ∈R Ẽ(ε, v, P ) = E(ε, v) + P ∆L + L 0 ε(s) -(1 + ε(s)) v (s) 2 2 + v (s) 4 24 ds . ( 12
)
where

C ≡ {v ∈ H 1 (0, L) : v(0) = 0, v(L) = 0}, (13)
and L 2 and H 1 denote the spaces of square integrable functions and of functions with square integrable first derivatives, respectively, whilst R is the set of real numbers. The stationarity condition with respect to ε gives:

ε(s) = - P EA 1 - v (s) 2 2 . ( 14
)
Substituting this solution in (12), we can finally recast the minimum principle (10) as

min v∈C max P ∈R Ẽ(v, P ), (15)
with

Ẽ(v, P ) = L 0 EI 2 v (s) 2 - P 2 2EA 1 - v (s) 2 2 2 -P v (s) 2 2 - v (s) 4 24 ds + P ∆L. (16)
In the following, we will rewrite the energy functional (16) in the form

Ẽ * (v * , P * ) = 1 2 1 0 v * (s * ) 2 ds * - P * 2 (1 -kP * ) 1 0 v * (s * ) 2 ds * + P * (1 -3kP * ) 24 1 0 v * (s * ) 4 ds * + P * ∆L * -k P * 2 2 , (17)
where we introduced the dimensionless extensibility parameter

k = I AL 2 = h 2 12L 2 , ( 18
)
and we rescaled the variables as follows

s * = s L , v * (s * ) = v(s) L , P * = P L 2 EI , Ẽ * = L EI Ẽ, ∆L * = ∆L L . (19)
The weak formulation of the equilibrium conditions for the beam is given by the following stationarity conditions with respect to v * and P * :

Ẽ * (v * , P * )(v * ) = 0, ∀v * ∈ C, (20a)
∂ Ẽ * ∂P * (v * , P * ) = 0, (20b)
where Ẽ * (v * , P * )(v * ) is the directional derivative of the functional E * with respect to v * in the direction v * :

Ẽ * (v * , P * )(v * ) = 1 0 v * v * ds * -P * (1 -kP * ) 1 0 v * v * ds + P * (1 -3kP * ) 6 1 0 v * 3 v * ds * . (21)
Equation (20b) gives a further scalar equation to calculate the force P * to get the imposed end-shortening

P * = 1 k ∆L * - 1 2 1 0 v * (s * ) 2 ds * + 1 24 1 0 v * (s * ) 4 ds * 1 - 1 0 v * (s * ) 2 ds * + 1 4 1 0 v * (s * ) 4 ds * . ( 22
)

Reduced dynamical model

We will study the nonlinear dynamics of the beam with a reduced model based on the following hypotheses:

• We assume that the dynamics is dominated by the transverse inertia. Neglecting the terms due to the longitudinal motion, we write the kinetic energy as

T dv dt = 1 2 L 0 ρS dv dt 2 ds. (23)
• We introduce a Galërkin approximation based on a modal expansion of the transverse displacement as

v * (s * , t) = N i=1 a i (t) v i (s * ), v i (s * ) = √ 2 sin(iπs * ). (24)
Within these approximations, we get the discretized elastic energy functional

V (a, P * ) = Ẽ * (a, P * ) = N i=1 (i 2 π 2 -P * (1 -kP * )) i 2 π 2 a 2 i 2 + P * (1 -3kP * ) N i,j,m,n=1 G ijmn a i a j a m a n 4 + P * ∆L * -k P * 2 2 , (25)
where a = [a 1 , . . . , a N ] T and

G ijmn = 2 i j m n π 4 3 1 0 cos(iπs) cos(jπs) cos(mπs) cos(nπs) ds. (26)
The discrete dimensionless expression of the kinetic energy reads as

T ( ȧ) := L EI T dv dt = 1 2 N i=1 ȧi 2 , (27)
where we denote by ȧi := dai dt * = Ω 0 dai dt the derivative with respect to the dimensionless time by t * = Ω 0 t scaled with

Ω 0 = 1 L 2 EI ρS . ( 28
)
Finally, using the standard Lagrangian formalism, the equations of motion for the modal coefficients read as

d dt * ∂T ∂ ȧi + ∂V ∂a i = - ∂R ∂ ȧi + Q i , (29)
∂V ∂P * = 0, (30)
where R( ȧ) and Q i model the effect of viscous damping and external forces, respectively. Here and henceforth, we consider only the effect of a transverse distribute loading per unit length

f * (s * , t * ) = f * (s * , t * ) e 2 = L EI f (s, t) e 2 , (31)
and a modal damping, setting

R( ȧ) = 1 2 N i=1 C i ȧ2 i , Q i (t * ) = 1 0 √ 2 sin (iπs * )f * (s * , t * ) ds * . (32)
Finally, the dynamical equations ( 29) in terms of the modal coefficients read as

äi (t * ) + (i 2 π 2 -P * (1 -kP * )) i 2 π 2 a i (t * ) + P * (1 -3kP * ) N j,l,m=1 G ijlm a j (t * )a l (t * )a m (t * ) + C i ȧi (t * ) = Q i (t * ), (33)
where the end-force P * is determined solving the constraint equation ( 30):

P * = 1 k ∆L * - N i=1 i 2 π 2 a 2 i 2 + N i,j,m,n=1 G ijmn a i a j a m a n 4 1 - N i=1 i 2 π 2 a 2 i + 3 2 N i,j,m,n=1 G ijmn a i a j a m a n . ( 34
)
3 Numerical simulations

The dynamic analysis of the bistable beam is here and henceforth limited to the modal model including only the three buckling modes (N = 3) in ( 24). Accordingly, the dynamics of the bistable beam is governed by the set of nonlinear differential equations ( 33) for the modal coordinates (a 1 , a 2 , a 3 ), where the Lagrange multiplier P * is given explicitly by ( 34). The free parameters of the model are the extensibility k, the dimensionless modal damping factor C i , and the end-shortening ∆L * . The numerical simulations of this section are for k = 10 -4 and C i = 2. We consider the external forcing in the form of a dimensionless transverse point force

f * = F p δ(s * -s * F ) sin(ω F t * ) e 2
, where δ is the Dirac delta function. In the following, we investigate the dynamical behavior of the beam as a function of the end-shorting ∆L * , the forcing frequency ω F , amplitude F p and application point s * F . To simplify the notation, we omit the * in the dimensionless frequency ω F and amplitude F p . The corresponding scaling factors are given in ( 28) and (31).

Static equilibria and natural frequency of the buckled beam

As a preliminary step, we calculate with our modal model the static equilibria of the beam as a function of the end shortening. Figure 2-left reports the beam midpoint deflection. We observe the classical pitchfork bifurcation diagram, where the beam becomes bistable (v m > 0) for

∆L > ∆L c = L 2 1 -1 -4π 2 k .
After the buckling threshold, we show only the stable solution with v m > 0, not representing the symmetric stable equilibrium with v m < 0. Figure 2-right plots as a function of the endshorting the natural frequencies ω 1 , ω 2 , ω 3 for small oscillations around the stable equilibrium.

The sharp jump of the derivatives corresponds to the buckling threshold. In the buckled regime, the frequencies ω 2 and ω 3 associated to the shape of the type sin(2π/L) and sin(3π/L) are almost independent of the end-shortening, whilst the frequency ω 1 associated to the modal shape of the type sin(π/L) increases almost quadratically.

Numerical method for solving the equations of the switching bistable beam

The dynamical response of the bistable beam subject to a harmonic actuation is obtained by solving in time the differential-algebraic system of equations including the three equations for the dynamics of the modal coefficients (33) and the expression (34) for the Lagrange multiplier P * that assure the prescribed end-shortening. Substituting directly (34) in ( 33) results in a very stiff dynamics for k → 0 and possible numerical issues. Hence, we apply a staggered procedure, where at each time step t i , given the values of the modal coefficients and velocities at the previous time-steps, we (i) update the value of the Lagrange multiplier P * that impose the end-shorting using ( 34) and

(ii) update the modal coefficients and velocities by a classical Runge-Kutta scheme.

We report below numerical results illustrating the observed dynamical behavior. 3. Inverted snap-through: the bistable undergoes some snap-through events as in the previous case, but the bistable oscillates around the other stable state, as in Figure 3(c).

Global long-time bistable response

4. Persistent snap-through: the bistable switches permanently between the two stable states, as in Figure 3(d).

We perform a parametric analysis to determine the qualitative properties of the response as a function of the forcing frequency and amplitude. Figure 4 reports the results for a point load applied at s * F = 0.4 and an end-shortening ∆L * = 1%. We color each point of the (ω F , F p )-plane according to the type of response observed. To assess the influence of the number of modes N retained in the modal model in the description of the beam dynamics, we report the results obtained This force was determined in a previous study dealing with the problem of static actuation of the bistable beam by a punctual force [START_REF] Camescasse | Bistable buckled beam and force actuation: Experimental validations[END_REF]. We note that the regions of the possible snapping of the bistable beam become wider as the number of modes accounted for increases. Indeed, for N = 1, the region of the snapping (either inverted or reversed or persistent snap) is narrower than those for N = 2 or N = 3. These results suggest that the accurate modelling of the nonlinear dynamics requires at least 3 modes in the considered post-buckling regime. For N = 1, we observe a minimal snapping force F p = 2.7 for ω F = 36.2, while for N = 2 or N = 3 the minimal force is F p = 1.6 for ω F = 32, close to second natural frequency ω 2 . These observations show that, in the vicinity of the second eigenmode a dynamical snapping is obtained with a lower actuating force. Such an observation has been conjectured by [START_REF] Casals-Terre | Dynamic analysis of a snap-action micromechanism[END_REF] and is therefore, here, confirmed by our model. We can compare our results with those of similar studies, in particular the work of [START_REF] Chandra | A numerical investigation of snap-through in a shallow arch-like model[END_REF]. For N = 1 the amplitude of the force necessary to trigger the snapping is very close to that of the static case. In addition, beyond the circular frequency ω F = 75, the model with only one mode does not show any snapping, in contrast to the model with 2 or 3 modes.

Figure 5 reports similar results obtained for a force applied at the midpoint of the beam (s * F = 0.5). In this case, the second skew-symmetric mode 2 is not excited. The static snapping force is F stat = 8.78. Moreover, the introduction of the third mode has a minor impact on the minimal snapping force in the vicinity of ω 2 . On the other hand, in the vicinity of ω 3 the inverted and reversed snaps are more probable. Finally, we examine the influence of the end-shortening on the snapping behavior. Figure 6 shows the results of the numerical simulations for N = 3. We report results for the end-shortenings ∆L * = 0.25 %, ∆L * = 0.5 %, and ∆L * = 1 %, for which the non-dimensional circular frequencies of the first mode are respectively ω 1 = 17.28, ω 1 = 28.14, and ω 1 = 42.2. As the end-shortening and consequently the value of ω 1 increases, the snapping zone is getting wider. Thus, we have an enrichment of the snapping possibilities of the bistable beam. Figure 6(d) shows the boundary that separates the possible snapping zones and the no-snapping zones, normalizing the force with respect to the force necessary for the snapping in the corresponding static case. We also note that for a small end-shortening ∆L * = 0.25% switching is not possible for non-dimensional frequencies ω * greater than 25.

Period doubling and transition to chaos

We investigate in further detail the transition from a regular to a chaotic dynamics. As before, we consider the case of a beam excited by a transversal harmonic force of amplitude F p and frequency ω F applied at the curvilinear abscissa s F = 0.4L, and we observe the vertical displacement v m = v(L/2) at the same point. Since we are interested in the long-time behavior of the bistable beam, we perform the time integration for a sufficiently long time interval and report the results only after that the effect of the initial transient can be neglected. Figure 8 shows the effect of increasing the excitation amplitude at a fixed frequency, reporting the Poincaré section (left) and frequency spectrum (right) of the responses. For F p = 5.6 (Figure 8a) the response is quasi-periodic. The Poincaré section has two compact orbits and the frequency spectrum is discrete. Increasing the amplitude, the Poincaré section stretches and forms strange attractor as shown in Figure 8b-8c, which is the signature of a chaotic response of the bistable beam.

Figure 9 reports the one-dimensional projections of the Poincaré sections on v * m -axis as a function of the forcing amplitude F p for two different frequencies. We observe a cascade of perioddoubling bifurcations as F p increases. In Figure 9 for ω F = 6π, the first bifurcation appears at F p = 4.3. For larger forces, we observe a dense network of bifurcations, interrupted by an interesting zone of regular responses for F p ∈ (8.5 -9). Figure 9 for ω F = 20π shows a cascade of bifurcations leading to chaos with a persistent snap-through. Our results are coherent with the observations of [START_REF] Guo | Numerical and experimental and experimental studies on nonlinear dynamics and performance of a bistable piezoelectric cantilever generator[END_REF] in the context of a bistable beam excited by an electro-magnetic system. In this Section, we compare the numerical results against the experimental findings. We consider the case of a simply-supported beam. Obtaining accurate and reproducible experimental results on the nonlinear dynamics and snapping properties implies to overcome two main difficulties:

• to control precisely the initial shape of the beam by limiting the imperfection to a tight tolerance;

• to design non-invasive forcing and measurement devices that do not perturb the intrinsic dynamical properties of the beam.

Here we report experimental results based on specific designed experimental set-up. In particular, we leverage an electromagnetic forcing device that exploits the action of Laplace forces [see also [START_REF] Amor | Snap-through of elastic bistable beam under contactless magnetic actuation[END_REF]Amor et al., , 2022, for static applications], for static applications]. A metallic beam hosting a controlled current flow is placed in a magnetic field. The modulation of the electric current allows us to easily control the excitation.

Description of the experimental set up

We perform experiments on a thin elastic beam of thickness h = 0.2mm, length L = 200mm and width b = 10mm obtained from a sheet of conductive amagnetic stainless steel A316L, with Young modulus E = 207GPa and volume mass density ρ = 7932kg/m 3 . We carefully selected the material and manufacturing process of the original sheet to ensure the planarity of the beam and to limit the residual stresses. Two identical magnets placed on either side of the beam generate a background magnetic induction in the covered region, as shown in Figure 10. Each of the two magnet is a planar plate made of Neodyne magnetic material (NdFeB alloy, magnetic grade N45, adhesive force 32Kg) with length 80mm, width 17mm, thickness 5mm. They generate an almost uniform magnetic field B = B e 3 in the region x C = 0 ≤ x ≤ x D = 80mm and vanishing elsewhere, with B = 0.2 Tesla.

We generate a current flow I = i τ (s) through the beam by connecting the two ends of the beam to the series connection of a resistance R and a power amplifier driven by a suitable function generator (BK ® Precision Model 4055). The power amplifier (Brüel & Kjaer ® , Type 2718, gain 40 dB, maximum power 75W, maximum current of approximately 7A) controls the amplitude of the current flowing through the beam. An oscilloscope and a precision (±1mA) ammeter provides the measure of the electric current. A laser high-speed laser profilometer (Keyence ® LJ-V7060) is used to measure the beam displacement at selected points.

Modelling of the Laplace force

When an electric charge q moves in a magnetic field B with a velocity v, it experiences the electromagnic Lorentz force f q = q v × B [START_REF] Jackson | Classical electrodynamics[END_REF]. An electric current flowing through a metallic beam inside a magnetic field induces a force per unit line

f (s, t) = I(t) × B(s) = I(t) B(s) n(s, t), (35)
called Laplace force, where n(s) is the normal to the beam axis in its current configuration, see (2).

The Laplace force is the macroscopic resultant of the Lorentz forces acting on each moving electron when an electric current flows into the beam. Viceversa, the mechanical velocity of the beam as a continuum in the magnetic field B does not induce any Lorentz force, because the net charge of the metallic beam is null. The force is a follower loading because its orientation depends on beam normal director in the current configuration n [START_REF] Healey | Large rotating states of a conducting elastic wire in a magnetic field: subtle symmetry and multiparameter bifurcation[END_REF][START_REF] Valverde | Magnetically-induced buckling of a whirling conducting rod with applications to electrodynamic space tethers[END_REF][START_REF] Wolfe | Equilibrium states of an elastic conductor in a magnetic field: A paradigm of bifurcation theory[END_REF]. Follower forces can induce peculiar instabilities, including Höpf-type bifurcation.

This is especially important for beam buckling instabilities when the force is purely tangential and several classical and recent studies are devoted to this aspect [START_REF] Bigoni | Nonlinear Solid Mechanics: Bifurcation Theory and Material Instability[END_REF][START_REF] Bigoni | Experimental evidence of flutter and divergence instabilities induced by dry friction[END_REF][START_REF] Ziegler | Linear elastic stability[END_REF]. In our specific application, the force is purely normal and the beam slope is relatively small (less than 10 degrees). In this context, the effect of the dependence of the force on the orientation of the beam normal can safely be neglected in a first approximation. The validity of this hypothesis will be confirmed in the following sections by the good agreement between the numerical results obtained in this way and the experimental findings. Hence, we model the effect of the Laplace force only as a dead load f (s, t) = I(t) B(s) e 2 and use the following expression to determine the loading coefficients (32) representing the effect of the Laplace force in our modal model ( 33):

Q i (t * ) = 1 0 √ 2 sin (iπs * )f * (s * , t * ) ds * = √ 2 iπ (cos (iπs * C ) -cos (iπs * D )) B I(t * ) L EI (36)
where we assumed that the magnetic induction B(s) is a constant B for s ∈ (x C , x D ) and vanishing elsewhere.

Modal properties of the buckled beam

We consider the bistable system obtained when imposing end shorting ∆L = 0.1mm (∆L * = 0.05%) or ∆L = 0.2mm (∆L * = 0.1%) to the initial straight beam of length L = 200mm and thickness h = 0.2mm. The aspect ratio of the beam used in the experimental setup is h/L = 10 -4 , Before investigating the nonlinear dynamics, we performed a classical experimental linear modal analysis on the buckled beams. The Table 1 reports the natural frequencies and damping ratio of the buckled beam associated to the modes v i (1 ≤ i ≤ 3) considered in our modal model ( 24), where the modes are ordered in function of their wavelength, and not their frequency. We can observe a very good agreement between the predictions of the simplified beam model and the experiments, with relative errors smaller than 4%. The modal damping ratio in the last column of Table 1 are used to estimate the damping coefficients C i (32) of the numerical model. The extensibility parameter ( 18) is very small for the experimental case k 8.3 × 10 -8 . In this context, the extensibility of our model should be regarded merely as a regularization for the purely inextensible model. The precise value of k does not have an important influence on the result.

Numerical and experimental results on the nonlinear dynamics

We compare here the qualitative and quantitive properties of the nonlinear forced dynamics predicted by our numerical models with the experimental results. For the applications, we identify the values of the amplitude and frequencies of the harmonic forcing for which the snapping be- 33) are calculated from the frequencies and damping ratios as

C i = 2ω i ζ i .
tween the two bistable state is possible. Figure 12 current required to induce a snapping behavior is similar for the two considered end-shorting. It happens for a forcing frequency close to 30Hz, slightly lower than the natural frequency of the second mode of the buckled beam, see Table 1. The permanent snap-through regions are wider than those obtained with the point force in Figure 4-6.

Conclusions

We studied the nonlinear dynamics of a bistable buckled beam. We presented an original setup for the dynamic non-invasive actuation of the buckled beam. The setup exploits the Laplace force generated by controlling the flow of an electric current through the beam placed in a magnetic induction generated by permanent magnets. We proved that this system can provide a sufficient high actuation to induce intra-well dynamic snapping between the two stable equilibria of the buckled beam. Comparisons between experimental and numerical results show that a weakly nonlinear simplified modal beam model is able quantitatively predict the complex dynamics and the regular or chaotic snapping regimes. Our results show that the modal model must include at least three modes to retrieve acceptable results [START_REF] Cazottes | Bistable Buckled Beam: Modeling of Actuating Force and Experimental Validations[END_REF][START_REF] Qiu | A curved-beam bistable mechanism[END_REF]. Even if the Laplace force at the basis of the actuating mechanisms is a follower force, we approximated it as a dead load. The good agreement between experimental and numerical results confirms that this approximation is acceptable in the considered regime. We proved that the proposed actuating system can be effective to control the dynamics of the bistable beam. Further work should be aimed at optimizing the device for engineering applications. From the fundamental point of view, it would be interesting to investigate possible static or dynamical instabilities due to the follower character of the Laplace force. Laplace forces can constitute an effective mean also for the experimental investigation of flutter instabilities, as an alternative to friction-induced devices proposed by Bigoni and coworkers Bigoni and Misseroni [2020].

Figure 1 :

 1 Figure 1: Elastic simply supported beam with an end-shorting ∆L and a transverse actuating force f : (a) the straight reference configuration of the beam, (b) the buckled configuration with the actuating force. The force P = -P e 1 denotes the constraint reaction required to impose the end-shorting ∆L.

Figure 2 :

 2 Figure 2: Left: the midpoint deflection of the beam v m (L/2) and Right: first three natural frequencies of the buckled beam, as a function of the end-shorting for k = 10 -4 .

Figure 3

 3 Figure3reports typical patterns for time-dependence of the vertical displacement of the beam midspan for a buckled beam with an end-shortening ∆L * = 1%. They can be classified in four scenarios with respect to the snap-through behavior:1. No snap-through: the bistable oscillates around one of both stable states (first buckling mode), Figure3(a).2. Reversed snap-through: the bistable experiences a number of snap-through movements and it oscillates around the starting stable state, as in Figure3(b).

Figure 3 :

 3 Figure 3: Time evolution of the vertical movement of the beam midspan according to the four kinds of dynamical behavior of the bistable beam: (a) no snap-through, (b) reversed snap-through, (c) inverted snap-through and (d) permanent snap-through

Figure 4 :

 4 Figure 4: Maps of the classification of the bistable beam switching response in the (ω F , F p ) space according to the number of buckling modes involved in the reduced model for actuation located at s * F = 0.4 (∆L * = ∆L/L = 1 %, k = 10 -4 and C 1 = C 2 = C 3 = 1).

Figure 5 :

 5 Figure 5: Maps of the classification of the bistable beam switching response in the (ω F , F p ) space according to the number of buckling modes involved in the reduced model for actuation located at s * F = 0.5 (∆L * = ∆L/L = 1 %, k = 10 -4 and C 1 = C 2 = C 3 = 1).

Figure 6 :

 6 Figure 6: Maps of the classification of the bistable beam switching response in the (ω F , F p ) space according to the end-shortening for N = 3 (k = 10 -4 , C 1 = C 2 = C 3 = 1 and s * F = 0.4 ).

Figure 7

 7 Figure 7 reports the dynamic response for different values of the forcing amplitude and frequency. We plot the time response (left), the projection of the phase-portrait in the v * mv * m plane (center), and corresponding frequency spectrum (right). The latter is obtained by performing a Fast Fourier Transform (FFT) of the time signals. The red-dots in the graphs in the center columns are the stroboscopic view of the state of the system at time intervals corresponding to the forcing period t F = 2π/ω F . They represent the projection of the Poincaré's section of the dynamical system in the v * mv * m plane [see e.g. Ott, 2002]. For these simulations the end-shortening and the damping factors are set to ∆L * = 1% and C 1 = C 2 = C 3 = 1. The results of Figures 7a-7b for (F p , ω F ) = (7, 85) and (F p , ω F ) = (3, 80) show two typical responses without snap-through. Even when the beam oscillates around the original equilibrium, nonlinear effect can be important, including possible super-harmonics and period-doubling phenomena, as show clearly by the red dots in the Poincaré maps. Figures 7c-7f report responses including the oscillation of the beam between the two stable equilibrium conditions, with the case of inverted (Figure 7c), reverted (Figure 7d), and persistent snap-through (Figure 7e-7f), with possible period-doubling (Figure 7d and 7f).

Figure 8

 8 Figure8illustrates the route from periodic or quasi-periodic responses to a chaotic response.

  (a) No snap-through for Fp = 7 and ω F = 85 (b) No snap-through for Fp = 3 and ω F = 80 (c) inverted snap-through for Fp = 3 and ω F = 34 (d) snap-through for Fp = 4.75 and ω F = 14 (e) persistent snap-through for Fp = 7 and ω F = 21 (f) persistent snap-through for Fp = 7 and ω F = 50

Figure 7 :Figure 8 :

 78 Figure 7: Time evolution of the vertical displacement of the beam midpoint v * m (left), phase portrait (center) including the Poincaré section (red dots) and frequency spectrum (right) of the Fast Fourier Transform |FRF(v * m )| of v * m . Numerical parameters: k = 10 -4 , ∆L * = ∆L/L = 1%, C 1 = C 2 = C 3 = 1 and s * F = 0.4 .

  (a) ω F = 6π (b) ω F = 20π

Figure 9 :

 9 Figure 9: Projection of the Poincaré section on the v * m -axis as a function of the forcing amplitude F p for ω F = 6π (a) and ω F = 20π (b), showing the transition to chaos. To better highlight the initial period doublings, the diagrams on the right column report a zoom inside the red box of the graphics of the left column. The simulations are for k = 10 -4 , ∆L * = ∆L/L = 1%, C 1 = C 2 = C 3 = 1 and s * F = 0.4

Figures 10 -

 10 Figures 10-11 show the schematics and different views of the experimental set-up. The beam is simply-supported at both ends, leaving free the rotation of the cross-section. The supports are made of a low friction material (PolyEtherEtherKetone) to limit the dissipation. The displacement of one end (point A) is blocked, whilst the axial displacement of the other end (point B) is imposed by micro-metric translator stage. A digital micrometer with an accuracy of 1 µm is used to control the end-shortening.

Figure 10 :

 10 Figure 10: Experimental setup. Left: Schematics of the measurement and actuating system, showing the laser profilometers and the electric circuit imposing a current flow through the beam. Right: Top and side view of the geometry of the beam and the permanent magnets generating the magnetic induction B.

Figure 11 :

 11 Figure 11: Experimental set-up for the simply-supported bistable beam with an electromagnetic actuation. The metallic conductive amagnetic beam (200 × 10 × 0.2mm) is connected to two wire imposing a current flow through the beam. The two transparent plates on the left includes permanent magnets inducing a constant magnetic induction in the corresponding region. The two ends of the beam are clamped to axes that can rotates with negligible friction to approximate simple-support conditions.

  resumes our main results. As in Figures4-7, we classify each point in the amplitude (current intensity I in Ampere) and frequency (in Hertz) plane according to the regimes observed in the numerical simulation: no-snap (white), reverted (red), inverted (green) or permanent (blue) snap. The pink solid lines report the boundary between the region with snap and no-snap observed experimentally. The agreement between experimental and numerical results is qualitatively and quantitatively good, both for ∆L * = 0.05% (left) and ∆L * = 0.1% (right). Obviously in the experiments the distinction between the different regimes is more difficult and generating data is time-consuming. The discrepancies become more important at higher frequencies, where higher modes, neglected in the model, can have an influence on the results.

Figure 12 :

 12 Figure 12: Maps of the classifications of the bistable beam switching response actuated by the Laplace force in the (f r, I) space according to the end-shortening : (left) ∆L * = 0.05% and (right) ∆L * = 0.1%.

Figure 13

 13 Figure 13 compares the numerical and results for the minimal force amplitude required to obtain the snap as a function of the frequency. Figure 14 reports explicitly the numerical and experimental time responses corresponding to the points (a)-(f) in Figure 13 and Table 2.

Figure 15

 15 Figure15shows the corresponding frequency spectrums. The agreement is generally good. The noisy frequency spectrums of Figures15(d)-(f) are a clear signature of the chaotic behavior.

Figure 13 :

 13 Figure 13: Boundaries between the regions of no-snap-through and snap-through of the bistable response in (f r, I) space, model and experimental results (∆L * = 0.1%). cases f r [Hz] I [A] snap-through (a) 2 5 persistent (b) 10 2 no (c) 20 4 persistent (d) 40 3 persistent (e) 60 3 reverted (f) 70 4 inverted Table2: Dynamical response of the bistable according to the frequency and the amplitude of the harmonic electric current for an actuation by means of the Laplace force.

Figure 14 :

 14 Figure14: The time evolution of the midpoint of the bistable beam actuated by the Laplace force for ∆L * = 0.1% (see Table2 and Figure13for information about the corresponding kind of the bistable response). (i) numerical simulations : curves in blue and (ii) experimental results : curves in red.

Figure 15 :

 15 Figure 15: Numerical (left column, in blue) and experimental (right column, in red) frequency spectra associated to the temporal responses displayed in Figure 14 for ∆L * = 0.1%).

Table 1 :

 1 Natural frequencies and damping ratios of the buckled beam used in the experimental setup for two different end-shortening. The damping coefficients in equation (

	∆L	∆L *	Deflection (v m ) Mode	Model	Frequencies Exp	Error	Modal damping (ζ)
				1	403.7Hz 402Hz 0.42%	0.5%
	0.1mm 0.05%	2.84mm	2	40.1Hz	40Hz	0.25%	4%
				3	98.2Hz	102Hz 3.73%	3%
				1	571.4Hz 586Hz 1.82%	0.4%
	0.2mm	1%	4.03mm	2	40.1Hz 39.3Hz 2.04%	4%
				3	98.2Hz 97.9Hz 0.31%	4%

  us remarquons que l'écart entre les résultats expérimentaux et ceux du modèle augm nt vers les hautes fréquences. Cette observation peut être expliquée par les erreur ation des amortissements modaux. En e et, le rôle de ces derniers est plus impo s cas de grandes vitesses de la poutre. De plus, le faible jeu au niveau des pivots eng uvements de corps rigide de translation. Lors de l'analyse modale expérimentale, repéré les pics de fréquences dûs à ces déplacements aux alentours de 87Hz. Lorsq est excitée avec des fréquences proches de cette valeur, les mouvements de corps is en résonance et les mesures sont impactées ayant une moindre qualité. Nous n' as pu obtenir des basculements avec des fréquences d'excitation supérieures à 80Hz alle du courant applicable (le courant maximal est 6.5 A).nt donné que le système étudié est un système chaotique, le moindre écart entre le es expérimentaux et ceux du modèle peut engendrer des réponses totalement di ére e nous attendons donc pas à obtenir une concordance parfaite entre les réponses ux procédures. L'objectif est de mettre en évidence un comportement global équiv i érentes valeurs d'amplitudes et de fréquences. Les déplacements du centre de la p tion du temps pour le modèle et l'expérience sont illustrés en Figure (4.44).

te étude a mis en évidence di érents comportements de la poutre. Les Figures (4 et(4.44d) correspondent à un comportement en basculement persistant. Cependant s courbes sont distinctes. Dans le cas (a), aucun comportement transitoire est ob tre suit l'allure du signal du courant et bascule régulièrement en suivant la faibl de l'excitation (2Hz). Un temps de stabilisation a été nécessaire avant d'atteind

Acknowledgements

The research work reported in the paper has been supported by the research project BISCELTECH funded by Fonds Unique Interministériel (FUI-APP21).


Florian Rigaud-Minet
email: florian.rigaud-minet@cea.fr

Julien Buckley

William Vandendaele

Stéphane Bécu

Matthew Charles

Jérôme Biscarrat

Romain Gwoziecki

Charlotte Gillot

Veronique Sousa

Hervé Morel
email: herve.morel@insa-lyon.fr

Dominique Planson
email: dominique.planson@insa-lyon.fr

Vertical Current Temperature Analysis of GaN-on-Si Epitaxy through Analytical Modelling

The current temperature dependence of a 3.9 µm GaN-on-Silicon epitaxy at temperatures ranging from 25 °C to 200 °C and for voltages beneath the plateau limited by the carrier generation in the p-doped substrate is analyzed. Based upon a fit with analytical transport models, the vertical current is reported for the first time to be the superposition of three different transport mechanisms: Recombination limited at low and middle voltages for low temperatures, Ohmic at low and middle voltages for high temperatures and Nearest Neighbor Hopping (NNH) at high voltages. Thanks to their temperature dependences and therefore to their energy, the transports can be related to either point defects clustered around threading dislocations or to interstitial nitrogen.

Introduction

For about ten years, tremendous efforts have been made to tackle with thermal and lattice mismatch of GaN-on-Si epitaxy to produce high quality GaN active layers [1] for power diodes and transistors in the 100 V to 1000 V range.

Many transport mechanisms were identified as governing the vertical leakage that may contribute to the device leakage. Near the breakdown, the study of Borga [2] reveals that the current depends on the nucleation layer/Substrate interface with a plateau in the I-V characteristics limited by the thermal carrier generation preceding an avalanche at the AlN/Si interface. At lower voltage, several mechanisms such as Space Charge Limited Current or Poole Frenkel were identified by means of TCAD [3], or analytical fit [4]. However, no analytical fitting was done for a wide range of temperature and for voltage below the plateau limited by thermal carrier generation to identify the possible origin of the leakage.

Experiments and models

The 3.9 µm AlGaN/GaN-on-Si epitaxy analyzed in this study was grown on a 1 mm thick p-doped 200 mm (111)-oriented silicon substrate with a resistivity between 3 and 20 Ω.cm. It consists of four main parts as illustrated in Fig. 1(a) deposited by Metal Organic Chemical Vapor Deposition (MOCVD). There is an AlN nucleation layer followed by an AlxGa1-xN Strain Relief Layer (SRL), and a 1.3 µm thick GaN:C with carbon concentration of roughly 5x10 18 cm -3 . The last hundred nanometers are composed of GaN layers. The top contact is a Ti/Al bilayer circular Ohmic contact with a diameter of 357 µm.

The measurements were done between 0 V and 800 V with a sweep rate of 20 V/s to prevent trapping from occurring and an integration over two Power Line Cycles (PLC). White light illumination was used to have a wider voltage range without being limited by the generation of electrons in the substrate. This phenomenon can still be identified in the temperature curves of the vertical current in Fig. 1(b) above 750 V. Verification was made to enhance that this illumination does not affect the leakage level studied here.

𝐼𝐼 𝛼𝛼 𝐴𝐴 1 × 𝑽𝑽 Recombination-limited model 𝐼𝐼 = 𝑞𝑞 × 𝐴𝐴 2 × 𝐺𝐺 × �µ 𝑛𝑛 × 𝜏𝜏 𝑛𝑛 × √𝑽𝑽 (ℎ𝑦𝑦𝑦𝑦. : µ 𝑛𝑛 ≫ µ 𝑝𝑝 )
Nearest Neighbor Hopping (NNH)

𝐼𝐼 = 𝐴𝐴 3 × 𝐶𝐶 1 × 𝑞𝑞 × 𝐾𝐾 𝐵𝐵 × 𝑇𝑇 × 𝑎𝑎 × 𝑒𝑒 𝑞𝑞×𝑎𝑎×𝑬𝑬 𝐾𝐾 𝐵𝐵 ×𝑇𝑇 - 𝜙𝜙 𝑡𝑡 𝐾𝐾 𝐵𝐵 ×𝑇𝑇 -2× 𝑎𝑎 𝛼𝛼 (ℎ𝑦𝑦𝑦𝑦. : 𝑞𝑞 × 𝑎𝑎 × 𝑬𝑬 ≈ 10 × 𝐾𝐾 𝐵𝐵 × 𝑇𝑇) I = Current, A i =

Results and discussion

A clear temperature and bias dependence is visible in Fig. 1(b) and is analyzed with 1D analytical models reported in TABLE I assuming a temperature independent effective surface. Transport models such as Variable Range Hopping have not been used because either their temperature or bias dependence could not fit the data. Noting that the uncertainty levels that are reported were determined by visual graphic method when fitting data.

The first analytical model that was used for this study is the ohmic model as described in TABLE I with an Arrhenius temperature dependence as in [5]. Indeed, a linear dependence in the log-log plot is visible as already exposed in the literature [6], [7]. An activation energy of 0.76 ± 0.03 eV was obtained for this transport mechanism.

Nevertheless, at lower temperature, the log-log slope is smaller. This is due to the recombination-limited transport contribution that vanishes with temperature. In fact, its square root dependence with voltage, as shown in TABLE I, leads to a slope of one-half in the log-log plot and allows fitting the experimental at lower temperatures but still below 400 V. The Arrhenius temperature dependence gives an activation energy of 0.27 ± 0.02 eV.

The third analytical model that has been used to fit the experimental curves especially between 575 V and 750 V is the Nearest Neighbor Hopping (NNH) with its equation given in TABLE I. Thanks to the ln(I) vs. E graph in Fig. 2, plotted by assuming that the electric field was uniformed in the III-V buffer, the trap energy level is 0.75 ± 0.01 eV.

To be complete, the authors proved that the Poole Frenkel (PF) model could not explain the experimental data because the relative permittivity extracted with the formula in [4] from the fit of 1.59 is inconsistent with the expected value of 8.9. The addition of all of these three different models in parallel as illustrated in Fig. 3(a) with the parameters exposed before gives the fitting plotted in Fig. 3(b).

Having identified three different energies related to all the contributions of the vertical leakages with two that are similar, their possible origins will be discussed now.

The energy range: 0.27 ± 0.02 eV related to the recombination-limited model embeds energy levels that were reported in the bulk GaN or AlGaN Deep Level Trap Spectroscopy (DLTS) literature but there is no consensus concerning its precise origin. The common signature of these energy levels is the logarithmic capture kinetics [8]- [10], revealing an extended defect such as a dislocation. Furthermore, this is more likely to be related to a cluster of punctual defects around the dislocation indiscriminately of its kind (screw, edge or mixed) rather than a dislocation alone according to Soh et al. [8]. According to the authors, the most relevant hypothesis for the cluster is the gallium vacancies VGa or a V N -V Ga pairs as presented by Fang [11]. In fact, abinitio simulations that in an N-rich environment (the case of the MOCVD deposition used in this study), the formation of V Ga near the dislocation has quite a low formation energy and even lower than the V N [12]. The Ohmic and hopping models revealed energies of 0.76 ± 0.03 eV and 0.75 ± 0.01 eV respectively. This can be related to interstitial nitrogen as reported in the paper of Asghar [13] on MOCVD grown vertical p-n diodes. This hypothesis is consistent since the GaN epitaxy studied here was grown in an N-rich environment.

Conclusions

The vertical current temperature dependence of a GaNon-Si epitaxy has been studied from 25 °C to 200 °C for voltages below 800 V where the current is limited by carrier thermal generation. Under 400 V, at low temperature, the current is expected to be limited by the recombination of carriers through edge, screw or mixed dislocations clustered with Gallium vacancies or Nitrogen vacancies or both in pair because of its activation energy of 0.27 ± 0.02 eV. At higher temperature and still below 400 V, the main current contribution is ohmic. The activation energy of the mechanism is equal to 0.76 ± 0.03 eV that could be due to the present of nitrogen in interstitial sites. Above 400 V and whatever the temperature, the transport is enhanced with a NNH mechanism possibly related to interstitial nitrogen because its energy is equal to 0.75 ± 0.01 eV. This is consistent with the N-rich environment during growth.

  Effective surface, V = Voltage, q = Elementary charge, G = Generation rate, µ n = Electron mobility, 𝜏𝜏 𝑛𝑛 = Electron lifetime, C 1 = Arbitrary constant, Φ T = trap level, E = Electric field, K B = Boltzmann constant, T = Temperature, a = Mean hopping distance, α = Decay length

Fig. 1 .

 1 Fig. 1. (a) Schematic of the vertical epitaxy (b) Experimental vertical current characteristics at temperature ranging from 25 °C to 200 °C

Fig. 2 .

 2 Fig. 2. Fitting between 575 V and 750 V with the nearest neighbor hopping model (fitting parameters: Φ T = 0.75 eV & a = 2.42 nm & α = 0.24 nm) with respect to the experimental values

Fig. 3 .

 3 Fig. 3. (a) The different contributions modeling for the vertical current (b) Final fitting with the experimental data

TABLE I :

 I TRANSPORT MECHANISMS

	Transport model	Analytical model
	Ohmic conduction

Acknowledgements

This study was funded by the French national program "Programme d'Investissements d'Avenir PSPC n°9 -G-Mobility".