Dual-axis MEMS micromirror

Choose the clamp

Every boundary surface of the mesh is a facet group surface_<id>. List the surfaces to clamp in clamped; their numbers are in FIG · 01. The default, surface 8, is bonded to the handle wafer.

clamped = [8] # surface numbers, read in ParaView

grid = SVK.FerriteGmsh.togrid(joinpath(@__DIR__, "micromirror_device_layer.msh"))
Ferrite.addfacetset!(grid, "Dirichlet",
    union((Ferrite.getfacetset(grid, "surface_$s") for s in clamped)...))

pv_dir = joinpath(@__DIR__, "results", "paraview")
mkpath(pv_dir)
write_surfaces(pv_dir, grid, clamped) # print the surface numbers

State the case

mechanical_model fixes every displacement on the clamped facets. The silicon is a cubic crystal, CubicCrystal, with the die edges along ⟨110⟩. There is no damping.

case = SVK.mechanical_model(grid;
    material = SVK.CubicCrystal(c11 = 165.7e3, c12 = 63.9e3, c44 = 79.6e3, ρ = 2.329e-9,
        rotation = π / 4),
    damping = SVK.RayleighDamping(α = 0.0, β = 0.0),
    dirichlet = "Dirichlet")
case.info.n_dofs # print the number of free degrees of freedom

Look at the modes

SVK.spectrum computes the modes once, and build_model reuses them.

nev = 10
sp = SVK.spectrum(case; nev);

The table lists each mode's frequency and how much it tilts the mirror (step 06).

coords = Ferrite.get_node_coordinate.(Ferrite.getnodes(grid))
mirror = findall(x -> abs(x[3]) < 1e-6 && hypot(x[1], x[2]) < 480, coords)
fit = pinv([ones(length(mirror)) getindex.(coords[mirror], 1) getindex.(coords[mirror], 2)])
rows = [SVK.probe_dof(case, n, 3) for n in mirror] # free-DOF rows of u_z on the mirror face
functional(c) = (l = zeros(case.info.n_dofs); l[rows] .= c; l)
θx, θy = functional(fit[3, :]), functional(-fit[2, :])

modes = map(1:nev) do p
    ϕ = mode_shape(sp, p)
    (; p, frequency = abs(imag(sp.eigenvalues[2p - 1])) / 2π,
        tilt_x = 500abs(θx ⋅ ϕ), tilt_y = 500abs(θy ⋅ ϕ))
end
println(" p   frequency   mirror tilt x   mirror tilt y")
foreach(m -> @printf("%2d  %10.4f  %14.3f  %14.3f\n", m...), modes)
write_modes(pv_dir, case, sp)

Build a physics-independent model for MORFE

master = [p] selects mode 3, the lowest mode that tilts the mirror about $y$. The rest is as in From a mesh to a ROM.

p = 3 # master mode, read off the table and ParaView
order = 5
(; model, spectral, meta) = build_model(case; master = [p], spectrum = sp, expansion_order = order)
meta.spectrum.eigenvalues[meta.master_indices] # print master eigenvalues

Compute the invariant manifold

As in From a mesh to a ROM, with tol_relative = 0.05: a resonance threshold relative to each eigenvalue.

W, R = parametrise(model, spectral, order;
    resonance = ResonanceConfig(style = :complex_normal_form, tol_relative = 0.05))
R # print reduced dynamics

Obtain the backbone

normal_form_branch reads the backbone off R, as in From a mesh to a ROM.

The observable is the mirror tilt $\theta$. A plane fitted to $u_z$ on the mirror face has slope $\sin\theta$ about the rotation axis (FIG · 03). observable_polynomial projects W onto that slope, cycle_amplitude takes its amplitude, and asind turns it into degrees.

θ = observable_polynomial(W, modes[p].tilt_x ≥ modes[p].tilt_y ? θx : θy) # sine of the mirror tilt
θ₁ = cycle_amplitude(MORFE.Polynomials.restrict_polynomial_to_degree(θ, 1), 1.0) # its linear part per unit ρ
ρ = range(0, sind(10) / θ₁; length = 301) # modal amplitude, up to a linear tilt of 10°
ω₀ = abs(imag(meta.spectrum.eigenvalues[meta.master_indices[1]])) # linear eigenfrequency

curves = map(3:2:order) do N
    b = normal_form_branch(restrict_ReducedDynamics_to_degree(R, N); amplitudes = ρ)
    P = MORFE.Polynomials.restrict_polynomial_to_degree(θ, N)
    (; N, b.amplitude, b.frequency, ratio = b.frequency ./ ω₀,
        tilt = asind.(cycle_amplitude.(Ref(P), b.amplitude))) # mirror tilt in degrees
end

Plot it

The mirror stiffens: at a tilt of 10°, Δω ≈ 0.085 %. Orders 3 and 5 nearly coincide up to 5°.

using CairoMakie
fig = Figure(size = (520, 380))
ax = Axis(fig[1, 1]; xlabel = "Δω [%]", ylabel = "mirror tilt [°]")
foreach(c -> lines!(ax, 100 .* (c.ratio .- 1), c.tilt; label = "order $(c.N)"), curves)
axislegend(ax; position = :lt)
fig

Save the ROM and the backbone

save_rom writes W, R and a summary; two CSVs hold the backbone and the tilt row of W.

MORFE.save_rom(joinpath(@__DIR__, "results"), W, R; drop_below = 0.0)

open(joinpath(@__DIR__, "results", "data", "backbone.csv"), "w") do io
    println(io, "order,r,omega,omega_ratio,tilt_deg")
    for c in curves, i in eachindex(c.amplitude)
        println(io, join((c.N, c.amplitude[i], c.frequency[i], c.ratio[i], c.tilt[i]), ","))
    end
end

open(joinpath(@__DIR__, "results", "data", "W_tilt_coefficients.csv"), "w") do io
    println(io, "exp_1,exp_2,W_tilt_re,W_tilt_im")
    for (e, c) in zip(MORFE.Polynomials.multiindex_set(θ).exponents,
        MORFE.Polynomials.coefficients(θ))
        println(io, join((e[1], e[2], real(c), imag(c)), ","))
    end
end