Add missing Matlab function

figs/inkscape/convert_svg.sh

@@ -0,0 +1,18 @@
+#!/bin/bash
+
+# Directory containing SVG files
+INPUT_DIR="."
+
+# Loop through all SVG files in the directory
+for svg_file in "$INPUT_DIR"/*.svg; do
+    # Check if there are SVG files in the directory
+    if [ -f "$svg_file" ]; then
+        # Output PDF file name
+        pdf_file="../${svg_file%.svg}.pdf"
+        png_file="../${svg_file%.svg}"
+        
+        # Convert SVG to PDF using Inkscape
+        inkscape "$svg_file" --export-filename="$pdf_file" && \
+            pdftocairo -png -singlefile -cropbox "$pdf_file" "$png_file"
+    fi
+done

figs/inkscape/detail_fem_joints_wire.svg

@@ -23,13 +23,13 @@
      inkscape:pagecheckerboard="0"
      inkscape:deskcolor="#d1d1d1"
      inkscape:document-units="mm"
-     inkscape:zoom="3.0060173"
-     inkscape:cx="-10.977981"
-     inkscape:cy="86.659515"
+     inkscape:zoom="6.0120346"
+     inkscape:cx="13.971975"
+     inkscape:cy="76.596366"
      inkscape:window-width="2534"
-     inkscape:window-height="1367"
+     inkscape:window-height="1387"
      inkscape:window-x="11"
-     inkscape:window-y="60"
+     inkscape:window-y="38"
      inkscape:window-maximized="1"
      inkscape:current-layer="layer1" /><defs
      id="defs1" /><g
@@ -59,5 +59,5 @@
        sodipodi:nodetypes="sccsccs" /><path
        id="path1-6"
        style="font-variation-settings:normal;opacity:1;vector-effect:none;fill:#ccebff;fill-opacity:1;stroke:#000000;stroke-width:0.264583;stroke-linecap:round;stroke-linejoin:round;stroke-miterlimit:4;stroke-dasharray:none;stroke-dashoffset:0;stroke-opacity:1;-inkscape-stroke:none;stop-color:#000000;stop-opacity:1"
-       d="m 54.316011,85.45368 c -1.379215,4.832397 4.285358,6.124171 5.992665,6.094811 -6e-6,1.148557 -2.957169,2.079657 -6.605013,2.079657 -3.647844,0 -6.605007,-0.9311 -6.605014,-2.079663 1.707306,0.02936 7.371879,-1.262421 5.992664,-6.094818 z"
-       sodipodi:nodetypes="ccscc" /></g></svg>
+       d="m 54.316011,85.45368 c -0.887898,3.110957 2.373522,3.991648 4.443957,4.756737 1.040155,0.38437 1.282585,0.436539 1.548708,1.338074 -6e-6,1.148557 -2.957169,2.079657 -6.605013,2.079657 -3.647844,0 -6.605007,-0.9311 -6.605014,-2.079663 0.293655,-0.938556 0.931993,-1.036998 1.983564,-1.485812 2.025194,-0.864359 4.789205,-1.875725 4.0091,-4.609006 z"
+       sodipodi:nodetypes="cscscscc" /></g></svg>

matlab/detail_fem_1_flexible_body.m

@@ -0,0 +1,307 @@
+%% Clear Workspace and Close figures
+clear; close all; clc;
+
+%% Intialize Laplace variable
+s = zpk('s');
+
+%% Path for functions, data and scripts
+addpath('./src/'); % Path for scripts
+addpath('./mat/'); % Path for data
+addpath('./STEPS/'); % Path for Simscape Model
+addpath('./subsystems/'); % Path for Subsystems Simulink files
+
+%% Linearization options
+opts = linearizeOptions;
+opts.SampleTime = 0;
+
+%% Open Simscape Model
+mdl = 'detail_fem_super_element'; % Name of the Simulink File
+open(mdl); % Open Simscape Model
+
+%% Colors for the figures
+colors = colororder;
+freqs = logspace(1,3,500); % Frequency vector [Hz]
+
+%% Estimate "Sensor Constant" - (THP5H)
+d33 = 680e-12;        % Strain constant [m/V]
+n   = 160;            % Number of layers per stack
+eT  = 4500*8.854e-12; % Permittivity under constant stress [F/m]
+sD  = 21e-12;         % Compliance under constant electric displacement [m2/N]
+
+gs = d33/(eT*sD*n);   % Sensor Constant [V/m]
+
+%% Estimate "Actuator Constant" - (THP5H)
+d33 = 680e-12;   % Strain constant [m/V]
+n   = 320;       % Number of layers
+
+cE  = 1/sD;      % Youngs modulus [N/m^2]
+A   = (10e-3)^2; % Area of the stacks [m^2]
+L   = 40e-3;     % Length of the two stacks [m]
+ka  = cE*A/L;    % Stiffness of the two stacks [N/m]
+
+ga = d33*n*ka;   % Actuator Constant [N/V]
+
+%% Load reduced order model
+K = readmatrix('APA95ML_K.CSV'); % order: 48
+M = readmatrix('APA95ML_M.CSV');
+[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('APA95ML_out_nodes_3D.txt');
+
+%% Stiffness estimation
+m = 0.0001; % block-free condition, no payload
+k_support = 1e9;
+c_support = 1e3;
+
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Fd'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '/y'],  1, 'openoutput'); io_i = io_i + 1;
+
+G = linearize(mdl, io);
+
+% The inverse of the DC gain of the transfer function
+% from vertical force to vertical displacement is the axial stiffness of the APA
+k_est = 1/dcgain(G); % [N/m]
+
+%% Estimated compliance of the APA95ML
+freqs = logspace(2, log10(5000), 1000);
+
+% Get first resonance indice
+i_max = find(abs(squeeze(freqresp(G, freqs(2:end), 'Hz'))) - abs(squeeze(freqresp(G, freqs(1:end-1), 'Hz'))) < 0, 1);
+
+figure;
+hold on;
+plot(freqs, abs(squeeze(freqresp(G, freqs, 'Hz'))), 'DisplayName', 'Compliance');
+plot([freqs(1), freqs(end)], [1/k_est, 1/k_est], 'k--', 'DisplayName', sprintf('$1/k$ ($k = %.0f N/\\mu m$)', 1e-6*k_est))
+xline(freqs(i_max), '--', 'linewidth', 1, 'color', [0,0,0], 'DisplayName', sprintf('$f_0 = %.0f$ Hz', freqs(i_max)))
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+xlabel('Frequency [Hz]'); ylabel('Amplitude [m/N]');
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+xlim([100, 5000]);
+
+%% Estimation of the amplification factor and Stroke
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Fa'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '/y'],  1, 'openoutput'); io_i = io_i + 1;
+io(io_i) = linio([mdl, '/d'],  1, 'openoutput'); io_i = io_i + 1;
+
+G = linearize(mdl, io);
+
+% Estimated amplification factor
+ampl_factor = abs(dcgain(G(1,1))./dcgain(G(2,1)));
+
+% Estimated stroke
+apa_stroke = ampl_factor * 3 * 20e-6; % [m]
+
+%% Experimental plant identification
+% with PD200 amplifier (gain of 20) - 2 stacks as an actuator, 1 as a sensor
+load('apa95ml_5kg_2a_1s.mat')
+
+Va = 20*u; % Voltage amplifier gain: 20
+
+% Spectral Analysis parameters
+Ts = t(end)/(length(t)-1);
+Nfft = floor(1/Ts);
+win = hanning(Nfft);
+Noverlap = floor(Nfft/2);
+
+% Identification of the transfer function from Va to di
+[G_y,  f]  = tfestimate(detrend(Va, 0), detrend(y, 0), win, Noverlap, Nfft, 1/Ts);
+[G_Vs, ~]  = tfestimate(detrend(Va, 0), detrend(v, 0), win, Noverlap, Nfft, 1/Ts);
+
+%% Plant Identification from Multi-Body model
+% Load Reduced Order Matrices
+K = readmatrix('APA95ML_K.CSV'); % order: 48
+M = readmatrix('APA95ML_M.CSV');
+[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('APA95ML_out_nodes_3D.txt');
+
+m = 5.5; % Mass of the suspended granite [kg]
+k_support = 4e7;
+c_support = 3e2;
+
+% Compute transfer functions
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Va'], 1, 'openinput');  io_i = io_i + 1; % Voltage accros piezo stacks [V]
+io(io_i) = linio([mdl, '/y'],  1, 'openoutput'); io_i = io_i + 1; % Vertical Displacement [m]
+io(io_i) = linio([mdl, '/Vs'], 1, 'openoutput'); io_i = io_i + 1; % Sensor stack voltage [V]
+
+Gm = linearize(mdl, io);
+Gm.InputName  = {'Va'};
+Gm.OutputName = {'y', 'Vs'};
+
+%% Comparison of the identified transfer function from Va to di to the multi-body model
+freqs = logspace(1, 3, 500);
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+plot(f, abs(G_y), '-', 'color', [colors(2,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'Measured FRF');
+plot(freqs, abs(squeeze(freqresp(Gm('y', 'Va'), freqs, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', 'Model')
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude $y/V_a$ [m/V]'); set(gca, 'XTickLabel',[]);
+hold off;
+ylim([1e-8, 1e-5]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+ax2 = nexttile;
+hold on;
+plot(f, 180/pi*angle(G_y), '-', 'color' , [colors(2,:), 0.5], 'linewidth', 2.5);
+plot(freqs, 180/pi*angle(squeeze(freqresp(Gm('y', 'Va'), freqs, 'Hz'))), '--', 'color', colors(2,:))
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
+hold off;
+yticks(-360:45:360);
+ylim([-45, 180]);
+
+linkaxes([ax1,ax2],'x');
+xlim([10, 1e3]);
+
+%% Comparison of the identified transfer function from Va to Vs to the multi-body model
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+plot(f, abs(G_Vs), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'Measured FRF');
+plot(freqs, abs(squeeze(freqresp(Gm('Vs', 'Va'), freqs, 'Hz'))), '--', 'color', colors(1,:), 'DisplayName', 'Model')
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude $V_s/V_a$ [V/V]'); set(gca, 'XTickLabel',[]);
+hold off;
+ylim([1e-3, 1e1]);
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+ax2 = nexttile;
+hold on;
+plot(f, 180/pi*angle(G_Vs), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5);
+plot(freqs, 180/pi*angle(squeeze(freqresp(Gm('Vs', 'Va'), freqs, 'Hz'))), '--', 'color', colors(1,:))
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
+hold off;
+yticks(-360:90:360); ylim([-180, 180]);
+
+linkaxes([ax1,ax2],'x');
+xlim([10, 1e3]);
+
+%% Integral Force Feedback Controller
+K_iff = (1/(s + 2*2*pi))*(s/(s + 0.5*2*pi));
+K_iff.inputname = {'Vs'};
+K_iff.outputname = {'u_iff'};
+
+% New damped plant input
+S1 = sumblk("u = u_iff + u_damp");
+
+% Voltage amplifier with gain of 20
+voltage_amplifier = tf(20);
+voltage_amplifier.inputname = {'u'};
+voltage_amplifier.outputname = {'Va'};
+
+%% Load experimental data with IFF implemented with different gains
+load('apa95ml_iff_test.mat', 'results');
+
+% Tested gains
+g_iff = [0, 10, 50, 100, 500, 1000];
+
+% Spectral Analysis parameters
+Ts = t(end)/(length(t)-1);
+Nfft = floor(1/Ts);
+win = hanning(Nfft);
+Noverlap = floor(Nfft/2);
+
+%% Computed the identified FRF of the damped plants
+tf_iff = {zeros(1, length(g_iff))};
+for i=1:length(g_iff)
+    [tf_est, f] = tfestimate(results{i}.u, results{i}.y, win, Noverlap, Nfft, 1/Ts);
+    tf_iff(i) = {tf_est};
+end
+
+%% Estimate the damped plants from the multi-body model
+Gm_iff = {zeros(1, length(g_iff))};
+
+for i=1:length(g_iff)
+    K_iff_g = -K_iff*g_iff(i); K_iff_g.inputname = {'Vs'}; K_iff_g.outputname = {'u_iff'};
+    Gm_iff(i) = {connect(Gm, K_iff_g, S1, voltage_amplifier, {'u_damp'}, {'y'})};
+end
+
+%% Identify second order plants from the experimental data
+% This is mandatory to estimate the experimental "poles"
+% an place them in the root-locus plot
+G_id = {zeros(1,length(results))};
+
+f_start = 70; % [Hz]
+f_end = 500; % [Hz]
+
+for i = 1:length(results)
+    tf_id = tf_iff{i}(sum(f<f_start):length(f)-sum(f>f_end));
+    f_id = f(sum(f<f_start):length(f)-sum(f>f_end));
+
+    gfr = idfrd(tf_id, 2*pi*f_id, Ts);
+    G_id(i) = {procest(gfr,'P2UDZ')};
+end
+
+%% Comparison of the Root-Locus computed from the multi-body model and the identified closed-loop poles
+gains = logspace(0, 5, 1000);
+figure;
+hold on;
+plot(real( pole(Gm('Vs', 'Va'))), imag( pole(Gm('Vs', 'Va'))), 'kx', 'HandleVisibility', 'off');
+plot(real(tzero(Gm('Vs', 'Va'))), imag(tzero(Gm('Vs', 'Va'))), 'ko', 'HandleVisibility', 'off');
+for i = 1:length(gains)
+    cl_poles = pole(feedback(Gm('Vs', 'Va'), gains(i)*K_iff));
+    plot(real(cl_poles(imag(cl_poles)>100)), imag(cl_poles(imag(cl_poles)>100)), 'k.', 'HandleVisibility', 'off');
+end
+for i = 1:length(g_iff)
+    cl_poles = pole(Gm_iff{i});
+    plot(real(cl_poles(imag(cl_poles)>100)), imag(cl_poles(imag(cl_poles)>100)), '.', 'MarkerSize', 20, 'color', colors(i,:), 'HandleVisibility', 'off');
+    plot(real(pole(G_id{i})), imag(pole(G_id{i})), 'x', 'color', colors(i,:), 'DisplayName', sprintf('g = %0.f', g_iff(i)), 'DisplayName', sprintf('$g = %.0f$', g_iff(i)));
+end
+xlabel('Real Part');
+ylabel('Imaginary Part');
+axis equal;
+ylim([-100, 2100]);
+xlim([-2100,100]);
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+%% Experimental damped plant for several IFF gains and estimated damped plants from the model
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2, 1]);
+hold on;
+plot(f, abs(tf_iff{1}), '-', 'DisplayName', '$g = 0$', 'color', [0,0,0, 0.5], 'linewidth', 2.5)
+plot(f, abs(squeeze(freqresp(Gm_iff{1}, f, 'Hz'))), 'k--', 'HandleVisibility', 'off')
+for i = 2:length(results)
+    plot(f, abs(tf_iff{i}), '-', 'DisplayName', sprintf('g = %0.f', g_iff(i)), 'color', [colors(i-1,:), 0.5], 'linewidth', 2.5)
+end
+for i = 2:length(results)
+    plot(f, abs(squeeze(freqresp(Gm_iff{i}, f, 'Hz'))), '--', 'color', colors(i-1,:), 'HandleVisibility', 'off')
+end
+set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
+ylabel('Amplitude $y/V_a$ [m/N]'); set(gca, 'XTickLabel',[]);
+hold off;
+ylim([1e-6, 2e-4]);
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+ax2 = nexttile;
+hold on;
+plot(f, 180/pi*angle(tf_iff{1}./squeeze(freqresp(exp(-s*2e-4), f, 'Hz'))), '-', 'color', [0,0,0, 0.5], 'linewidth', 2.5)
+plot(f, 180/pi*angle(squeeze(freqresp(Gm_iff{1}, f, 'Hz'))), 'k--')
+for i = 2:length(results)
+    plot(f, 180/pi*angle(tf_iff{i}./squeeze(freqresp(exp(-s*2e-4), f, 'Hz'))), '-', 'color', [colors(i-1,:), 0.5], 'linewidth', 2.5)
+    plot(f, 180/pi*angle(squeeze(freqresp(Gm_iff{i}, f, 'Hz'))), '--', 'color', colors(i-1,:))
+end
+set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+hold off;
+yticks(-360:45:360);
+ylim([-10, 190]);
+
+linkaxes([ax1,ax2], 'x');
+xlim([150, 500]);

matlab/detail_fem_2_actuators.m

@@ -0,0 +1,318 @@
+%% Clear Workspace and Close figures
+clear; close all; clc;
+
+%% Intialize Laplace variable
+s = zpk('s');
+
+%% Path for functions, data and scripts
+addpath('./src/'); % Path for scripts
+addpath('./mat/'); % Path for data
+addpath('./STEPS/'); % Path for Simscape Model
+addpath('./subsystems/'); % Path for Subsystems Simulink files
+
+%% Linearization options
+opts = linearizeOptions;
+opts.SampleTime = 0;
+
+%% Open Simscape Model
+mdl = 'detail_fem_APA300ML'; % Name of the Simulink File
+open(mdl); % Open Simscape Model
+
+% Piezoelectric parameters
+ga = -25.9; % [N/V]
+gs = -5.08e6; % [V/m]
+
+%% Colors for the figures
+colors = colororder;
+freqs = logspace(1,3,500); % Frequency vector [Hz]
+
+%% Identify dynamics with "Reduced Order Flexible Body"
+K = readmatrix('APA300ML_mat_K.CSV');
+M = readmatrix('APA300ML_mat_M.CSV');
+[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('APA300ML_out_nodes_3D.txt');
+
+m = 5; % [kg]
+ga = 25.9; % [N/V]
+gs = 5.08e6; % [V/m]
+
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Va'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '/Fd'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '/z'], 1, 'openoutput'); io_i = io_i + 1;
+io(io_i) = linio([mdl, '/Vs'], 1, 'openoutput'); io_i = io_i + 1;
+
+G_fem = linearize(mdl, io);
+G_fem_z = G_fem('z','Va');
+G_fem_Vs = G_fem('Vs', 'Va');
+G_fem_comp = G_fem('z', 'Fd');
+
+%% Determine c1 and k1 from the zero
+G_zeros = zero(minreal(G_fem_Vs));
+G_zeros = G_zeros(imag(G_zeros)>0);
+[~, i_sort] = sort(imag(G_zeros));
+G_zeros = G_zeros(i_sort);
+G_zero = G_zeros(1);
+
+% Solving 2nd order equations
+c1 = -2*m*real(G_zero);
+k1 = m*(imag(G_zero)^2 + real(G_zero)^2);
+
+%% Determine ka, ke, ca, ce from the first pole
+G_poles = pole(minreal(G_fem_z));
+G_poles = G_poles(imag(G_poles)>0);
+[~, i_sort] = sort(imag(G_poles));
+G_poles = G_poles(i_sort);
+G_pole = G_poles(1);
+
+% Solving 2nd order equations
+ce = 3*(-2*m*real(G_pole(1)) - c1);
+ca = 1/2*ce;
+
+ke = 3*(m*(imag(G_pole)^2 + real(G_pole)^2) - k1);
+ka = 1/2*ke;
+
+%% Matching sensor/actuator constants
+% ga = dcgain(G_fem_z) / (1/(ka + k1*ke/(k1 + ke)));
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '_2dof', '/Fa'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '_2dof', '/z'], 1, 'openoutput'); io_i = io_i + 1;
+ga = dcgain(G_fem_z)/dcgain(linearize([mdl, '_2dof'], io));
+
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '_2dof', '/Va'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '_2dof', '/dL'], 1, 'openoutput'); io_i = io_i + 1;
+gs = dcgain(G_fem_Vs)/dcgain(linearize([mdl, '_2dof'], io));
+
+%% Identify dynamics with tuned 2DoF model
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '_2dof', '/Va'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '_2dof', '/Fd'], 1, 'openinput');  io_i = io_i + 1;
+io(io_i) = linio([mdl, '_2dof', '/z'], 1, 'openoutput'); io_i = io_i + 1;
+io(io_i) = linio([mdl, '_2dof', '/Vs'], 1, 'openoutput'); io_i = io_i + 1;
+
+G_2dof = linearize([mdl, '_2dof'], io);
+G_2dof_z = G_2dof('z','Va');
+G_2dof_Vs = G_2dof('Vs', 'Va');
+G_2dof_comp = G_2dof('z', 'Fd');
+
+%% Comparison of the transfer functions from Va to vertical motion - FEM vs 2DoF
+freqs = logspace(1, 3, 500);
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+plot(freqs, abs(squeeze(freqresp(G_fem_z, freqs, 'Hz'))), '-', 'color', [colors(2,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'FEM')
+plot(freqs, abs(squeeze(freqresp(G_2dof_z, freqs, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', '2DoF Model')
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude $y/V_a$ [m/V]'); set(gca, 'XTickLabel',[]);
+hold off;
+ylim([1e-8, 2e-4]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+ax2 = nexttile;
+hold on;
+plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_fem_z, freqs, 'Hz')))), '-', 'color', [colors(2,:), 0.5], 'linewidth', 2.5);
+plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_2dof_z, freqs, 'Hz')))), '--', 'color', colors(2,:))
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
+hold off;
+yticks(-360:45:360); ylim([-20, 200]);
+
+linkaxes([ax1,ax2],'x');
+xlim([10, 1e3]);
+
+%% Comparison of the transfer functions from Va to Vs - FEM vs 2DoF
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+plot(freqs, abs(squeeze(freqresp(G_fem_Vs, freqs, 'Hz'))), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'FEM');
+plot(freqs, abs(squeeze(freqresp(G_2dof_Vs, freqs, 'Hz'))), '--', 'color', colors(1,:), 'DisplayName', '2DoF Model')
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude $V_s/V_a$ [V/V]'); set(gca, 'XTickLabel',[]);
+hold off;
+ylim([6e-4, 3e1]);
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+ax2 = nexttile;
+hold on;
+plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_fem_Vs, freqs, 'Hz')))), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5);
+plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_2dof_Vs, freqs, 'Hz')))), '--', 'color', colors(1,:))
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
+hold off;
+yticks(-360:45:360); ylim([-20, 200]);
+
+linkaxes([ax1,ax2],'x');
+xlim([10, 1e3]);
+
+%% Effect of electrical boundaries on the
+oc = load('detail_fem_apa95ml_open_circuit.mat',  't', 'encoder', 'u');
+sc = load('detail_fem_apa95ml_short_circuit.mat', 't', 'encoder', 'u');
+
+% Spectral Analysis parameters
+Ts = sc.t(end)/(length(sc.t)-1);
+Nfft = floor(2/Ts);
+win = hanning(Nfft);
+Noverlap = floor(Nfft/2);
+
+% Identification of the transfer function from Va to di
+[G_oc,  f]  = tfestimate(detrend(oc.u, 0), detrend(oc.encoder, 0), win, Noverlap, Nfft, 1/Ts);
+[G_sc,  f]  = tfestimate(detrend(sc.u, 0), detrend(sc.encoder, 0), win, Noverlap, Nfft, 1/Ts);
+
+% Find resonance frequencies
+[~, i_oc] = max(abs(G_oc(f<300)));
+[~, i_sc] = max(abs(G_sc(f<300)));
+
+%% Effect of the electrical bondaries of the force sensor stack on the APA95ML resonance frequency
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+plot(f, abs(G_oc), '-', 'DisplayName', sprintf('Open-Circuit - $f_0 = %.1f Hz$', f(i_oc)))
+plot(f, abs(G_sc), '-', 'DisplayName', sprintf('Short-Circuit - $f_0 = %.1f Hz$', f(i_sc)))
+set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
+ylabel('Amplitude [m/N]'); set(gca, 'XTickLabel',[]);
+hold off;
+ylim([1e-6, 1e-4]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+ax2 = nexttile;
+hold on;
+plot(f, 180/pi*angle(G_oc), '-')
+plot(f, 180/pi*angle(G_sc), '-')
+set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
+ylabel('Phase'); xlabel('Frequency [Hz]');
+hold off;
+yticks(-360:45:360);
+ylim([-20, 200]);
+axis padded 'auto x'
+
+linkaxes([ax1,ax2], 'x');
+xlim([100, 300]);
+
+%% Compare Dynamics between "Reduced Order" flexible joints and "2-dof and 3-dof" joints
+% Let's initialize all the stages with default parameters.
+initializeGround('type', 'rigid');
+initializeGranite('type', 'rigid');
+initializeTy('type', 'rigid');
+initializeRy('type', 'rigid');
+initializeRz('type', 'rigid');
+initializeMicroHexapod('type', 'rigid');
+initializeSample('m', 50);
+
+initializeSimscapeConfiguration();
+initializeDisturbances('enable', false);
+initializeLoggingConfiguration('log', 'none');
+initializeController('type', 'open-loop');
+initializeReferences();
+
+mdl = 'detail_fem_nass';
+
+% Input/Output definition
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Controller'],     1, 'openinput');              io_i = io_i + 1; % Actuator Inputs
+io(io_i) = linio([mdl, '/Tracking Error'], 1, 'openoutput', [], 'EdL');  io_i = io_i + 1; % Errors in the frame of the struts
+io(io_i) = linio([mdl, '/NASS'],       3, 'openoutput', [], 'fn');  io_i = io_i + 1; % Force Sensors
+
+% Flexible actuators
+initializeSimplifiedNanoHexapod('actuator_type', 'flexible', ...
+            'flex_type_F', '2dof', ...
+            'flex_type_M', '3dof');
+
+G_flex = linearize(mdl, io);
+G_flex.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+G_flex.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+
+% Actuators modeled as 2DoF system
+initializeSimplifiedNanoHexapod('actuator_type', 'apa300ml', ...
+            'flex_type_F', '2dof', ...
+            'flex_type_M', '3dof');
+
+G_ideal = linearize(mdl, io);
+G_ideal.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+G_ideal.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+
+%% Comparison of the dynamics for actuators modeled using "reduced order flexible body" and using 2DoF system - HAC plant
+freqs = logspace(1, 4, 1000);
+
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for j = 1:5
+    for k = j+1:6
+        plot(freqs, abs(squeeze(freqresp(G_flex("l"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(1,:), 0.1], ...
+            'HandleVisibility', 'off');
+        plot(freqs, abs(squeeze(freqresp(G_ideal("l"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(2,:), 0.1], ...
+            'HandleVisibility', 'off');
+    end
+end
+plot(freqs, abs(squeeze(freqresp(G_flex("l1","f1"), freqs, 'Hz'))), 'color', colors(1,:), 'DisplayName', 'FEM');
+plot(freqs, abs(squeeze(freqresp(G_ideal("l1","f1"), freqs, 'Hz'))), 'color', colors(2,:), 'DisplayName', '2-DoF');
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [m/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-10, 1e-4]);
+
+ax2 = nexttile;
+hold on;
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex("l1","f1"), freqs, 'Hz'))));
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_ideal("l1","f1"), freqs, 'Hz'))));
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-20, 200]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');
+
+%% Comparison of the dynamics for actuators modeled using "reduced order flexible body" and using 2DoF system - IFF plant
+freqs = logspace(0, 3, 1000);
+
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for j = 1:5
+    for k = j+1:6
+        plot(freqs, abs(squeeze(freqresp(G_flex("fm"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(1,:), 0.1], ...
+            'HandleVisibility', 'off');
+        plot(freqs, abs(squeeze(freqresp(G_ideal("fm"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(2,:), 0.1], ...
+            'HandleVisibility', 'off');
+    end
+end
+plot(freqs, abs(squeeze(freqresp(G_flex("fm1","f1"), freqs, 'Hz'))), 'color', colors(1,:), 'DisplayName', 'FEM');
+plot(freqs, abs(squeeze(freqresp(G_ideal("fm1","f1"), freqs, 'Hz'))), 'color', colors(2,:), 'DisplayName', '2-DoF');
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [N/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-5, 1e1]);
+
+ax2 = nexttile;
+hold on;
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex("fm1","f1"), freqs, 'Hz'))));
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_ideal("fm1","f1"), freqs, 'Hz'))));
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-20, 200]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');

matlab/detail_fem_3_flexible_joints.m

@@ -0,0 +1,623 @@
+%% Clear Workspace and Close figures
+clear; close all; clc;
+
+%% Intialize Laplace variable
+s = zpk('s');
+
+%% Path for functions, data and scripts
+addpath('./src/'); % Path for scripts
+addpath('./mat/'); % Path for data
+addpath('./STEPS/'); % Path for Simscape Model
+addpath('./subsystems/'); % Path for Subsystems Simulink files
+
+%% Linearization options
+opts = linearizeOptions;
+opts.SampleTime = 0;
+
+%% Open Simscape Model
+mdl = 'detail_fem_nass'; % Name of the Simulink File
+open(mdl); % Open Simscape Model
+
+%% Colors for the figures
+colors = colororder;
+freqs = logspace(1,3,500); % Frequency vector [Hz]
+
+%% Identify the dynamics for several considered bending stiffnesses
+% Let's initialize all the stages with default parameters.
+initializeGround('type', 'rigid');
+initializeGranite('type', 'rigid');
+initializeTy('type', 'rigid');
+initializeRy('type', 'rigid');
+initializeRz('type', 'rigid');
+initializeMicroHexapod('type', 'rigid');
+initializeSample('m', 50);
+
+initializeSimscapeConfiguration();
+initializeDisturbances('enable', false);
+initializeLoggingConfiguration('log', 'none');
+initializeController('type', 'open-loop');
+initializeReferences();
+
+% Input/Output definition
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Controller'],     1, 'openinput');              io_i = io_i + 1; % Actuator Inputs
+io(io_i) = linio([mdl, '/Tracking Error'], 1, 'openoutput', [], 'EdL');  io_i = io_i + 1; % Errors in the frame of the struts
+io(io_i) = linio([mdl, '/NASS'],       3, 'openoutput', [], 'fn');  io_i = io_i + 1; % Force Sensors
+
+% Effect of bending stiffness
+Kf = [0, 50, 100, 500]; % [Nm/rad]
+G_Kf = {zeros(length(Kf), 1)};
+
+for i = 1:length(Kf)
+    % Limited joint axial compliance
+    initializeSimplifiedNanoHexapod('actuator_type', '1dof', ...
+                          'flex_type_F', '2dof', ...
+                          'flex_type_M', '3dof', ...
+                          'actuator_k', 1e6, ...
+                          'actuator_c', 1e1, ...
+                          'actuator_kp', 0, ...
+                          'actuator_cp', 0, ...
+                          'Fsm', 56e-3, ... % APA300ML weight 112g
+                          'Msm', 56e-3, ...
+                          'Kf_F', Kf(i), ...
+                          'Kf_M', Kf(i));
+
+    G_Kf(i) = {linearize(mdl, io)};
+    G_Kf{i}.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+    G_Kf{i}.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+end
+
+freqs = logspace(0, 3, 1000);
+
+%% Effect of the flexible joint bending stiffness on the HAC-plant
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for i = 1:length(Kf)
+    for j = 1:5
+        for k = j+1:6
+            plot(freqs, abs(squeeze(freqresp(G_Kf{i}("l"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(i,:), 0.1], ...
+                 'HandleVisibility', 'off');
+        end
+    end
+    plot(freqs, abs(squeeze(freqresp(G_Kf{i}("l1","f1"), freqs, 'Hz'))), 'color', colors(i,:), 'DisplayName', sprintf('$k_f = %.0f $ [Nm/rad]', Kf(i)));
+    for j = 2:6
+        plot(freqs, abs(squeeze(freqresp(G_Kf{i}("l"+j,"f"+j), freqs, 'Hz'))), 'color', colors(i,:), 'HandleVisibility', 'off');
+    end
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [m/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-10, 1e-4]);
+
+ax2 = nexttile;
+hold on;
+for i = 1:length(Kf)
+    plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_Kf{i}(1, 1), freqs, 'Hz')))));
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-200, 20]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');
+
+%% Effect of the flexible joint bending stiffness on the IFF plant
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for i = 1:length(Kf)
+    for j = 1:5
+        for k = j+1:6
+            plot(freqs, abs(squeeze(freqresp(G_Kf{i}("fm"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(i,:), 0.1], ...
+                 'HandleVisibility', 'off');
+        end
+    end
+    plot(freqs, abs(squeeze(freqresp(G_Kf{i}("fm1","f1"), freqs, 'Hz'))), 'color', colors(i,:), 'DisplayName', sprintf('$k_f = %.0f $ [Nm/rad]', Kf(i)));
+    for j = 2:6
+        plot(freqs, abs(squeeze(freqresp(G_Kf{i}("fm"+j,"f"+j), freqs, 'Hz'))), 'color', colors(i,:), 'HandleVisibility', 'off');
+    end
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [N/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-4, 1e2]);
+
+ax2 = nexttile();
+hold on;
+for i = 1:length(Kf)
+    plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_Kf{i}("fm1", "f1"), freqs, 'Hz')))));
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-20, 200]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');
+
+%% Decentalized IFF
+Kiff = -200 * ...              % Gain
+       1/s * ... % LPF: provides integral action
+       eye(6);                 % Diagonal 6x6 controller (i.e. decentralized)
+
+Kiff.InputName = {'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+Kiff.OutputName = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+
+%% Root Locus for decentralized IFF - 1dof actuator - Effect of joint bending stiffness
+gains = logspace(-1, 2, 400);
+
+figure;
+tiledlayout(1, 1, 'TileSpacing', 'compact', 'Padding', 'None');
+nexttile();
+hold on;
+
+for i = 1:length(Kf)
+    plot(real(pole(G_Kf{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))),  imag(pole(G_Kf{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))),  'x', 'color', colors(i,:), ...
+        'DisplayName', sprintf('$k_f = %.0f$ Nm/rad', Kf(i)));
+    plot(real(tzero(G_Kf{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))), imag(tzero(G_Kf{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))), 'o', 'color', colors(i,:), ...
+        'HandleVisibility', 'off');
+
+    for g = gains
+        clpoles = pole(feedback(G_Kf{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}), g*Kiff, +1));
+        plot(real(clpoles), imag(clpoles), '.', 'color', colors(i,:), ...
+            'HandleVisibility', 'off');
+    end
+
+end
+
+xline(0, 'HandleVisibility', 'off'); yline(0, 'HandleVisibility', 'off');
+hold off;
+axis equal;
+xlim(1.1*[-900, 100]); ylim(1.1*[-100, 900]);
+xticks(1.1*[-900:100:0]);
+yticks(1.1*[0:100:900]);
+set(gca, 'XTickLabel',[]); set(gca, 'YTickLabel',[]);
+xlabel('Real part'); ylabel('Imaginary part');
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+%% Identify the dynamics for several considered bending stiffnesses - APA300ML
+G_Kf_apa300ml = {zeros(length(Kf), 1)};
+
+for i = 1:length(Kf)
+    % Limited joint axial compliance
+    initializeSimplifiedNanoHexapod('actuator_type', 'apa300ml', ...
+                          'flex_type_F', '2dof', ...
+                          'flex_type_M', '3dof', ...
+                          'Fsm', 56e-3, ... % APA300ML weight 112g
+                          'Msm', 56e-3, ...
+                          'Kf_F', Kf(i), ...
+                          'Kf_M', Kf(i));
+
+    G_Kf_apa300ml(i) = {linearize(mdl, io)};
+    G_Kf_apa300ml{i}.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+    G_Kf_apa300ml{i}.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+end
+
+Kiff = -1000 * ...              % Gain
+       1/(s) * ... % LPF: provides integral action
+       eye(6);                 % Diagonal 6x6 controller (i.e. decentralized)
+
+Kiff.InputName = {'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+Kiff.OutputName = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+
+%% Root Locus for decentralized IFF - APA300ML actuator - Effect of joint bending stiffness
+gains = logspace(-1, 2, 300);
+
+figure;
+tiledlayout(1, 1, 'TileSpacing', 'compact', 'Padding', 'None');
+nexttile();
+hold on;
+
+for i = 1:length(Kf)
+    plot(real(pole(G_Kf_apa300ml{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))),  imag(pole(G_Kf_apa300ml{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))),  'x', 'color', colors(i,:), ...
+        'DisplayName', sprintf('$k_f = %.0f$ [Nm/rad]', Kf(i)));
+    plot(real(tzero(G_Kf_apa300ml{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))), imag(tzero(G_Kf_apa300ml{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))), 'o', 'color', colors(i,:), ...
+        'HandleVisibility', 'off');
+
+    for g = gains
+        clpoles = pole(feedback(G_Kf_apa300ml{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}), g*Kiff, +1));
+        plot(real(clpoles), imag(clpoles), '.', 'color', colors(i,:), ...
+            'HandleVisibility', 'off');
+    end
+
+end
+
+xline(0, 'HandleVisibility', 'off'); yline(0, 'HandleVisibility', 'off');
+hold off;
+axis equal;
+xlim(1.4*[-900, 100]); ylim(1.4*[-100, 900]);
+xticks(1.4*[-900:100:0]);
+yticks(1.4*[0:100:900]);
+set(gca, 'XTickLabel',[]); set(gca, 'YTickLabel',[]);
+xlabel('Real part'); ylabel('Imaginary part');
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+%% Identify the dynamics for several considered axial stiffnesses
+% Let's initialize all the stages with default parameters.
+initializeGround('type', 'rigid');
+initializeGranite('type', 'rigid');
+initializeTy('type', 'rigid');
+initializeRy('type', 'rigid');
+initializeRz('type', 'rigid');
+initializeMicroHexapod('type', 'rigid');
+initializeSample('m', 50);
+
+initializeSimscapeConfiguration();
+initializeDisturbances('enable', false);
+initializeLoggingConfiguration('log', 'none');
+initializeController('type', 'open-loop');
+initializeReferences();
+
+% Input/Output definition
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Controller'],     1, 'openinput');              io_i = io_i + 1; % Actuator Inputs
+io(io_i) = linio([mdl, '/Tracking Error'], 1, 'openoutput', [], 'EdL');  io_i = io_i + 1; % Errors in the frame of the struts
+io(io_i) = linio([mdl, '/NASS'],       3, 'openoutput', [], 'fn');  io_i = io_i + 1; % Force Sensors
+
+% Effect of bending stiffness
+Ka = 1e6*[1000, 100, 10, 1]; % [Nm/rad]
+G_Ka = {zeros(length(Ka), 1)};
+
+for i = 1:length(Ka)
+    % Limited joint axial compliance
+    initializeSimplifiedNanoHexapod('actuator_type', '1dof', ...
+                          'flex_type_F', '2dof_axial', ...
+                          'flex_type_M', '4dof', ...
+                          'actuator_k', 1e6, ...
+                          'actuator_c', 1e1, ...
+                          'actuator_kp', 0, ...
+                          'actuator_cp', 0, ...
+                          'Fsm', 56e-3, ... % APA300ML weight 112g
+                          'Msm', 56e-3, ...
+                          'Ca_F', 1, ...
+                          'Ca_M', 1, ...
+                          'Ka_F', Ka(i), ...
+                          'Ka_M', Ka(i));
+
+    G_Ka(i) = {linearize(mdl, io)};
+    G_Ka{i}.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+    G_Ka{i}.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+end
+
+freqs = logspace(1, 4, 1000);
+
+%% Effect of the flexible joint axial stiffness on the HAC-plant
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for i = 1:length(Ka)
+    for j = 1:5
+        for k = j+1:6
+            plot(freqs, abs(squeeze(freqresp(G_Ka{i}("l"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(i,:), 0.1], ...
+                 'HandleVisibility', 'off');
+        end
+    end
+end
+for i = 1:length(Ka)
+    plot(freqs, abs(squeeze(freqresp(G_Ka{i}("l1","f1"), freqs, 'Hz'))), 'color', colors(i,:), 'DisplayName', sprintf('$k_a = %.0f$ [N/$\\mu$m]', 1e-6*Ka(i)));
+    % for j = 2:6
+    %     plot(freqs, abs(squeeze(freqresp(G_Ka{i}("l"+j,"f"+j), freqs, 'Hz'))), 'color', colors(i,:), 'HandleVisibility', 'off');
+    % end
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [m/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-10, 1e-4]);
+
+ax2 = nexttile;
+hold on;
+for i = 1:length(Ka)
+    plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_Ka{i}(1, 1), freqs, 'Hz')))));
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-200, 20]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');
+
+%% Effect of the flexible joint axial stiffness on the IFF plant
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for i = 1:length(Ka)
+    for j = 1:5
+        for k = j+1:6
+            plot(freqs, abs(squeeze(freqresp(G_Ka{i}("fm"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(i,:), 0.1], ...
+                 'HandleVisibility', 'off');
+        end
+    end
+end
+for i = 1:length(Ka)
+    plot(freqs, abs(squeeze(freqresp(G_Ka{i}("fm1","f1"), freqs, 'Hz'))), 'color', colors(i,:), 'DisplayName', sprintf('$k_a = %.0f$ [N/$\\mu$m]', 1e-6*Ka(i)));
+    % for j = 2:6
+    %     plot(freqs, abs(squeeze(freqresp(G_Ka{i}("fm"+j,"f"+j), freqs, 'Hz'))), 'color', colors(i,:), 'HandleVisibility', 'off');
+    % end
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [N/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-4, 1e2]);
+
+ax2 = nexttile();
+hold on;
+for i = 1:length(Ka)
+    plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(G_Ka{i}("fm1", "f1"), freqs, 'Hz')))));
+end
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-20, 200]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');
+
+%% Decentalized IFF
+Kiff = -200 * ...              % Gain
+       1/s * ... % LPF: provides integral action
+       eye(6);                 % Diagonal 6x6 controller (i.e. decentralized)
+
+Kiff.InputName = {'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+Kiff.OutputName = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+
+%% Root Locus for decentralized IFF - 1dof actuator - Effect of joint bending stiffness
+gains = logspace(-1, 2, 400);
+
+figure;
+tiledlayout(1, 1, 'TileSpacing', 'compact', 'Padding', 'None');
+nexttile();
+hold on;
+
+for i = 1:length(Ka)
+    plot(real(pole(G_Ka{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))),  imag(pole(G_Ka{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))),  'x', 'color', colors(i,:), ...
+        'DisplayName', sprintf('$k_a = %.0f$ N/$\\mu$m', 1e-6*Ka(i)));
+    plot(real(tzero(G_Ka{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))), imag(tzero(G_Ka{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}))), 'o', 'color', colors(i,:), ...
+        'HandleVisibility', 'off');
+
+    for g = gains
+        clpoles = pole(feedback(G_Ka{i}({"fm1", "fm2", "fm3", "fm4", "fm5", "fm6"}, {"f1", "f2", "f3", "f4", "f5", "f6"}), g*Kiff, +1));
+        plot(real(clpoles), imag(clpoles), '.', 'color', colors(i,:), ...
+            'HandleVisibility', 'off');
+    end
+
+end
+
+xline(0, 'HandleVisibility', 'off'); yline(0, 'HandleVisibility', 'off');
+hold off;
+axis equal;
+xlim(1.1*[-900, 100]); ylim(1.1*[-100, 900]);
+xticks(1.1*[-900:100:0]);
+yticks(1.1*[0:100:900]);
+set(gca, 'XTickLabel',[]); set(gca, 'YTickLabel',[]);
+xlabel('Real part'); ylabel('Imaginary part');
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+%% Compute the damped plants
+Kiff = -500 * ...              % Gain
+       1/(s + 2*pi*0.1) * ... % LPF: provides integral action
+       eye(6);                 % Diagonal 6x6 controller (i.e. decentralized)
+
+Kiff.InputName = {'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+Kiff.OutputName = {'u1iff', 'u2iff', 'u3iff', 'u4iff', 'u5iff', 'u6iff'};
+
+% New damped plant input
+S1 = sumblk("f1 = u1iff + u1");
+S2 = sumblk("f2 = u2iff + u2");
+S3 = sumblk("f3 = u3iff + u3");
+S4 = sumblk("f4 = u4iff + u4");
+S5 = sumblk("f5 = u5iff + u5");
+S6 = sumblk("f6 = u6iff + u6");
+
+G_Ka_iff = {zeros(1,length(Ka))};
+for i=1:length(Ka)
+    G_Ka_iff(i) = {connect(G_Ka{i}, Kiff, S1, S2, S3, S4, S5, S6, {'u1', 'u2', 'u3', 'u4', 'u5', 'u6'}, {'l1', 'l2', 'l3', 'l4', 'l5', 'l6'})};
+end
+
+%% Interaction Analysis - RGA Number
+rga = zeros(length(Ka), length(freqs));
+for i = 1:length(Ka)
+    for j = 1:length(freqs)
+        rga(i,j) = sum(sum(abs(inv(evalfr(G_Ka_iff{i}({"l1", "l2", "l3", "l4", "l5", "l6"}, {"u1", "u2", "u3", "u4", "u5", "u6"}), 1j*2*pi*freqs(j)).').*evalfr(G_Ka_iff{i}({"l1", "l2", "l3", "l4", "l5", "l6"}, {"u1", "u2", "u3", "u4", "u5", "u6"}), 1j*2*pi*freqs(j)) - eye(6))));
+    end
+end
+
+%% RGA number for the damped plants - Effect of the flexible joint axial stiffness
+figure;
+hold on;
+for i = 1:length(Ka)
+    plot(freqs, rga(i,:), 'DisplayName', sprintf('$k_a = %.0f$ N/$\\mu$m', 1e-6*Ka(i)))
+end
+hold off;
+xlabel('Frequency [Hz]'); ylabel('RGA number');
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylim([0, 10]); xlim([10, 5e3]);
+leg = legend('location', 'northwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+
+%% Extract stiffness of the joint from the reduced order model
+% We first extract the stiffness and mass matrices.
+K = readmatrix('flex025_mat_K.CSV');
+M = readmatrix('flex025_mat_M.CSV');
+% Then, we extract the coordinates of the interface nodes.
+[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('flex025_out_nodes_3D.txt');
+
+m = 1;
+
+%% Name of the Simulink File
+mdl = 'detail_fem_joint';
+
+%% Input/Output definition
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/T'], 1, 'openinput');  io_i = io_i + 1; % Forces and Torques
+io(io_i) = linio([mdl, '/D'], 1, 'openoutput'); io_i = io_i + 1; % Translations and Rotations
+
+G = linearize(mdl, io);
+
+% Stiffness extracted from the Simscape model
+k_a = 1/dcgain(G(3,3)); % Axial stiffness [N/m]
+k_f = 1/dcgain(G(4,4)); % Bending stiffness [N/m]
+k_t = 1/dcgain(G(6,6)); % Torsion stiffness [N/m]
+
+
+% Stiffness extracted from the Stiffness matrix
+k_s = K(1,1); % shear [N/m]
+% k_s = K(2,2); % shear [N/m]
+k_a = K(3,3); % axial [N/m]
+k_f = K(4,4); % bending [Nm/rad]
+% k_f = K(5,5); % bending [Nm/rad]
+k_t =  K(6,6); % torsion [Nm/rad]
+
+%% Compare Dynamics between "Reduced Order" flexible joints and "2-dof and 3-dof" joints
+% Let's initialize all the stages with default parameters.
+initializeGround('type', 'rigid');
+initializeGranite('type', 'rigid');
+initializeTy('type', 'rigid');
+initializeRy('type', 'rigid');
+initializeRz('type', 'rigid');
+initializeMicroHexapod('type', 'rigid');
+initializeSample('m', 50);
+
+initializeSimscapeConfiguration();
+initializeDisturbances('enable', false);
+initializeLoggingConfiguration('log', 'none');
+initializeController('type', 'open-loop');
+initializeReferences();
+
+mdl = 'detail_fem_nass';
+
+% Input/Output definition
+clear io; io_i = 1;
+io(io_i) = linio([mdl, '/Controller'],     1, 'openinput');              io_i = io_i + 1; % Actuator Inputs
+io(io_i) = linio([mdl, '/Tracking Error'], 1, 'openoutput', [], 'EdL');  io_i = io_i + 1; % Errors in the frame of the struts
+io(io_i) = linio([mdl, '/NASS'],       3, 'openoutput', [], 'fn');  io_i = io_i + 1; % Force Sensors
+
+% Fully flexible joints
+initializeSimplifiedNanoHexapod('actuator_type', 'apa300ml', ...
+            'flex_type_F', 'flexible', ...
+            'flex_type_M', 'flexible', ...
+            'Fsm', 56e-3, ... % APA300ML weight 112g
+            'Msm', 56e-3);
+
+G_flex = linearize(mdl, io);
+G_flex.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+G_flex.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+
+% Flexible joints modelled by 2DoF and 3DoF joints
+initializeSimplifiedNanoHexapod('actuator_type', 'apa300ml', ...
+            'flex_type_F', '2dof_axial', ...
+            'flex_type_M', '4dof', ...
+            'Kf_F', k_f, ...
+            'Kt_F', k_t, ...
+            'Ka_F', k_a, ...
+            'Kf_M', k_f, ...
+            'Kt_M', k_t, ...
+            'Ka_M', k_a, ...
+            'Cf_F', 1e-2, ...
+            'Ct_F', 1e-2, ...
+            'Ca_F', 1e-2, ...
+            'Cf_M', 1e-2, ...
+            'Ct_M', 1e-2, ...
+            'Ca_M', 1e-2, ...
+            'Fsm', 56e-3, ... % APA300ML weight 112g
+            'Msm', 56e-3);
+
+G_ideal = linearize(mdl, io);
+G_ideal.InputName  = {'f1', 'f2', 'f3', 'f4', 'f5', 'f6'};
+G_ideal.OutputName = {'l1', 'l2', 'l3', 'l4', 'l5', 'l6', 'fm1', 'fm2', 'fm3', 'fm4', 'fm5', 'fm6'};
+
+%% Comparison of the dynamics with joints modelled with FEM and modelled with "ideal joints" - HAC plant
+freqs = logspace(1, 4, 1000);
+
+%% Effect of the flexible joint axial stiffness on the HAC-plant
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for j = 1:5
+    for k = j+1:6
+        plot(freqs, abs(squeeze(freqresp(G_flex("l"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(1,:), 0.1], ...
+            'HandleVisibility', 'off');
+        plot(freqs, abs(squeeze(freqresp(G_ideal("l"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(2,:), 0.1], ...
+            'HandleVisibility', 'off');
+    end
+end
+plot(freqs, abs(squeeze(freqresp(G_flex("l1","f1"), freqs, 'Hz'))), 'color', colors(1,:), 'DisplayName', 'Reduced Order Flexible Joints');
+plot(freqs, abs(squeeze(freqresp(G_ideal("l1","f1"), freqs, 'Hz'))), 'color', colors(2,:), 'DisplayName', 'Bot: $k_f$, $k_a$, Top: $k_f$, $k_t$, $k_a$');
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [m/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-10, 1e-4]);
+
+ax2 = nexttile;
+hold on;
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex("l1","f1"), freqs, 'Hz'))));
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_ideal("l1","f1"), freqs, 'Hz'))));
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-20, 200]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');
+
+freqs = logspace(0, 3, 1000);
+
+%% Effect of the flexible joint axial stiffness on the HAC-plant
+figure;
+tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+
+ax1 = nexttile([2,1]);
+hold on;
+for j = 1:5
+    for k = j+1:6
+        plot(freqs, abs(squeeze(freqresp(G_flex("fm"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(1,:), 0.1], ...
+            'HandleVisibility', 'off');
+        plot(freqs, abs(squeeze(freqresp(G_ideal("fm"+k,"f"+j), freqs, 'Hz'))), 'color', [colors(2,:), 0.1], ...
+            'HandleVisibility', 'off');
+    end
+end
+plot(freqs, abs(squeeze(freqresp(G_flex("fm1","f1"), freqs, 'Hz'))), 'color', colors(1,:), 'DisplayName', 'Reduced Order Flexible Joints');
+plot(freqs, abs(squeeze(freqresp(G_ideal("fm1","f1"), freqs, 'Hz'))), 'color', colors(2,:), 'DisplayName', 'Bot: $k_f$, $k_a$, Top: $k_f$, $k_t$, $k_a$');
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
+ylabel('Amplitude [N/N]'); set(gca, 'XTickLabel',[]);
+leg = legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
+ylim([1e-5, 1e1]);
+
+ax2 = nexttile;
+hold on;
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex("fm1","f1"), freqs, 'Hz'))));
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_ideal("fm1","f1"), freqs, 'Hz'))));
+hold off;
+set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
+ylabel('Phase [deg]'); xlabel('Frequency [Hz]');
+ylim([-20, 200]);
+yticks([-360:45:360]);
+
+linkaxes([ax1,ax2],'x');

matlab/src/extractNodes.m

@@ -0,0 +1,91 @@
+function [int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes(filename)
+% extractNodes -
+%
+% Syntax: [n_xyz, nodes] = extractNodes(filename)
+%
+% Inputs:
+%    - filename - relative or absolute path of the file that contains the Matrix
+%
+% Outputs:
+%    - n_xyz -
+%    - nodes - table containing the node numbers and corresponding dof of the interfaced DoFs
+
+arguments
+    filename
+end
+
+fid = fopen(filename,'rt');
+
+if fid == -1
+    error('Error opening the file');
+end
+
+n_xyz = []; % Contains nodes coordinates
+n_i = []; % Contains nodes indices
+
+n_num = []; % Contains node numbers
+n_dof = {}; % Contains node directions
+
+while 1
+    % Read a line
+    nextline = fgetl(fid);
+
+    % End of the file
+    if ~isstr(nextline), break, end
+
+    % Line just before the list of nodes coordinates
+    if contains(nextline, 'NODE') && ...
+            contains(nextline, 'X') && ...
+            contains(nextline, 'Y') && ...
+            contains(nextline, 'Z')
+
+        while 1
+            nextline = fgetl(fid);
+
+            if nextline < 0, break, end
+
+            c = sscanf(nextline, ' %f');
+
+            if isempty(c), break, end
+
+            n_xyz = [n_xyz; c(2:4)'];
+            n_i = [n_i; c(1)];
+        end
+    end
+
+    if nextline < 0, break, end
+
+    % Line just before the list of node DOF
+  if contains(nextline, 'NODE') && ...
+         contains(nextline, 'LABEL')
+
+        while 1
+            nextline = fgetl(fid);
+
+            if nextline < 0, break, end
+
+            c = sscanf(nextline, ' %d %s');
+
+            if isempty(c), break, end
+
+            n_num = [n_num; c(1)];
+
+            n_dof{length(n_dof)+1} = char(c(2:end)');
+        end
+
+        nodes = table(n_num, string(n_dof'), 'VariableNames', {'node_i', 'node_dof'});
+    end
+
+    if nextline < 0, break, end
+end
+
+fclose(fid);
+
+int_i = unique(nodes.('node_i')); % indices of interface nodes
+
+% Extract XYZ coordinates of only the interface nodes
+if length(n_xyz) > 0 && length(n_i) > 0
+    int_xyz = n_xyz(logical(sum(n_i.*ones(1, length(int_i)) == int_i', 2)), :);
+else
+    int_xyz = n_xyz;
+end

nass-fem.bib

@@ -1,15 +1,15 @@
-@article{souleille18_concep_activ_mount_space_applic,
-  author          = {Souleille, Adrien and Lampert, Thibault and Lafarga, V and
-                  Hellegouarch, Sylvain and Rondineau, Alan and Rodrigues,
-                  Gon{\c{c}}alo and Collette, Christophe},
-  title           = {A Concept of Active Mount for Space Applications},
-  journal         = {CEAS Space Journal},
-  volume          = 10,
-  number          = 2,
-  pages           = {157--165},
-  year            = 2018,
-  publisher       = {Springer},
-  keywords        = {parallel robot, iff},
+@article{mcinroy02_model_desig_flexur_joint_stewar,
+  author          = {J.E. McInroy},
+  title           = {Modeling and Design of Flexure Jointed Stewart Platforms
+                  for Control Purposes},
+  journal         = {IEEE/ASME Transactions on Mechatronics},
+  volume          = 7,
+  number          = 1,
+  pages           = {95-99},
+  year            = 2002,
+  doi             = {10.1109/3516.990892},
+  url             = {https://doi.org/10.1109/3516.990892},
+  keywords        = {parallel robot, flexure},
 }
 
 
@@ -66,6 +66,8 @@
 }
 
 
+
+
 @article{claeyssen07_amplif_piezoel_actuat,
   author          = {Frank Claeyssen and R. Le Letty and F. Barillot and O.
                   Sosnicki},
@@ -109,6 +111,53 @@
 
 
 
+@book{pintelon12_system_ident,
+  author          = {Rik Pintelon and Johan Schoukens},
+  title           = {System Identification : a Frequency Domain Approach},
+  year            = 2012,
+  publisher       = {Wiley IEEE Press},
+  url             = {https://doi.org/10.1002/9781118287422},
+  address         = {Hoboken, N.J. Piscataway, NJ},
+  doi             = {10.1002/9781118287422},
+  isbn            = 9780470640371,
+}
+
+
+
+@article{souleille18_concep_activ_mount_space_applic,
+  author          = {Souleille, Adrien and Lampert, Thibault and Lafarga, V and
+                  Hellegouarch, Sylvain and Rondineau, Alan and Rodrigues,
+                  Gon{\c{c}}alo and Collette, Christophe},
+  title           = {A Concept of Active Mount for Space Applications},
+  journal         = {CEAS Space Journal},
+  volume          = 10,
+  number          = 2,
+  pages           = {157--165},
+  year            = 2018,
+  publisher       = {Springer},
+  keywords        = {parallel robot, iff},
+}
+
+
+
+@article{verma20_dynam_stabil_thin_apert_light,
+  author          = {Mohit Verma and Adrien Pece and Sylvain Hellegouarch and
+                  Jennifer Watchi and Gilles Durand and Simon Chesn{\'e} and
+                  Christophe Collette},
+  title           = {Dynamic Stabilization of Thin Aperture Light Collector
+                  Space Telescope Using Active Rods},
+  journal         = {Journal of Astronomical Telescopes, Instruments, and
+                  Systems},
+  volume          = 6,
+  number          = 01,
+  pages           = 1,
+  year            = 2020,
+  doi             = {10.1117/1.jatis.6.1.014002},
+  url             = {http://dx.doi.org/10.1117/1.JATIS.6.1.014002},
+  DATE_ADDED      = {Thu Apr 3 21:25:20 2025},
+}
+
+
 @phdthesis{hanieh03_activ_stewar,
   author          = {Hanieh, Ahmed Abu},
   keywords        = {parallel robot},
@@ -120,37 +169,13 @@
 
 
 
-@article{mcinroy02_model_desig_flexur_joint_stewar,
-  author          = {J.E. McInroy},
-  title           = {Modeling and Design of Flexure Jointed Stewart Platforms
-                  for Control Purposes},
-  journal         = {IEEE/ASME Transactions on Mechatronics},
-  volume          = 7,
-  number          = 1,
-  pages           = {95-99},
-  year            = 2002,
-  doi             = {10.1109/3516.990892},
-  url             = {https://doi.org/10.1109/3516.990892},
-  keywords        = {parallel robot, flexure},
-}
-
-
-
-@article{yang19_dynam_model_decoup_contr_flexib,
-  author          = {Yang, XiaoLong and Wu, HongTao and Chen, Bai and Kang,
-                  ShengZheng and Cheng, ShiLi},
-  title           = {Dynamic Modeling and Decoupled Control of a Flexible
-                  Stewart Platform for Vibration Isolation},
-  journal         = {Journal of Sound and Vibration},
-  volume          = 439,
-  pages           = {398-412},
-  year            = 2019,
-  doi             = {10.1016/j.jsv.2018.10.007},
-  url             = {https://doi.org/10.1016/j.jsv.2018.10.007},
-  issn            = {0022-460X},
-  keywords        = {parallel robot, flexure, decoupled control},
-  month           = {Jan},
-  publisher       = {Elsevier BV},
+@book{schmidt20_desig_high_perfor_mechat_third_revis_edition,
+  author          = {Schmidt, R Munnig and Schitter, Georg and Rankers, Adrian},
+  title           = {The Design of High Performance Mechatronics - Third Revised
+                  Edition},
+  year            = 2020,
+  publisher       = {Ios Press},
+  keywords        = {favorite},
 }
 
 
@@ -173,6 +198,25 @@
 
 
 
+@article{yang19_dynam_model_decoup_contr_flexib,
+  author          = {Yang, XiaoLong and Wu, HongTao and Chen, Bai and Kang,
+                  ShengZheng and Cheng, ShiLi},
+  title           = {Dynamic Modeling and Decoupled Control of a Flexible
+                  Stewart Platform for Vibration Isolation},
+  journal         = {Journal of Sound and Vibration},
+  volume          = 439,
+  pages           = {398-412},
+  year            = 2019,
+  doi             = {10.1016/j.jsv.2018.10.007},
+  url             = {https://doi.org/10.1016/j.jsv.2018.10.007},
+  issn            = {0022-460X},
+  keywords        = {parallel robot, flexure, decoupled control},
+  month           = 1,
+  publisher       = {Elsevier BV},
+}
+
+
+
 @article{du14_piezo_actuat_high_precis_flexib,
   author          = {Zhijiang Du and Ruochong Shi and Wei Dong},
   title           = {A Piezo-Actuated High-Precision Flexible Parallel Pointing
@@ -187,3 +231,16 @@
   keywords        = {parallel robot},
 }
 
+
+
+@book{preumont18_vibrat_contr_activ_struc_fourt_edition,
+  author          = {Andre Preumont},
+  title           = {Vibration Control of Active Structures - Fourth Edition},
+  year            = 2018,
+  publisher       = {Springer International Publishing},
+  url             = {https://doi.org/10.1007/978-3-319-72296-2},
+  doi             = {10.1007/978-3-319-72296-2},
+  keywords        = {favorite, parallel robot},
+  series          = {Solid Mechanics and Its Applications},
+}
+

nass-fem.tex

@@ -1,4 +1,4 @@
-% Created 2025-02-27 Thu 01:24
+% Created 2025-04-03 Thu 22:01
 % Intended LaTeX compiler: pdflatex
 \documentclass[a4paper, 10pt, DIV=12, parskip=full, bibliography=totoc]{scrreprt}
 
@@ -8,13 +8,6 @@
 \author{Dehaeze Thomas}
 \date{\today}
 \title{Optimization using Finite Element Models}
-\hypersetup{
- pdfauthor={Dehaeze Thomas},
- pdftitle={Optimization using Finite Element Models},
- pdfkeywords={},
- pdfsubject={},
- pdfcreator={Emacs 29.4 (Org mode 9.6)}, 
- pdflang={English}}
 \usepackage{biblatex}
 
 \begin{document}
@@ -23,34 +16,26 @@
 \tableofcontents
 
 \clearpage
+During the nano-hexapod's detailed design phase, a hybrid modeling approach combining finite element analysis with multi-body dynamics was developed.
+This methodology, utilizing reduced-order flexible bodies, was created to enable both detailed component optimization and efficient system-level simulation, addressing the impracticality of a full FEM for real-time control scenarios.
 
-During the detailed design phase of the nano-hexapod, optimizing individual components while ensuring their dynamic compatibility with the complete system presents significant challenges.
-While Finite Element Analysis (FEA) serves as a powerful tool for component-level optimization, understanding how the dynamics of each element interacts within the complete nano-active stabilization system (NASS) becomes crucial.
-A full finite element model of the assembled system, while theoretically possible, would prove impractical for simulating real-time control scenarios due to its computational complexity.
-
-This chapter presents a hybrid modeling approach that combines finite element analysis with multi-body dynamics, enabling both detailed component optimization and efficient system-level simulation.
-The methodology employs reduced-order flexible bodies, whereby components whose dynamic properties are determined through FEA can be effectively integrated into the multi-body framework.
-The theoretical foundations and practical implementation of this approach are presented in Section \ref{sec:detail_fem_super_element}, where experimental validation using an amplified piezoelectric actuator demonstrates the method's accuracy in predicting both open and closed-loop dynamic behavior.
-
-This validated modeling framework is then applied to optimize two critical elements of the nano-hexapod: the actuators and the flexible joints.
-Section \ref{sec:detail_fem_actuator} examines the selection and characterization of amplified piezoelectric actuators, developing both high-fidelity and computationally efficient models that capture essential dynamic characteristics.
-Section \ref{sec:detail_fem_joint} addresses the design of flexible joints, where precise control of directional stiffness proves crucial for system performance.
-In both cases, the hybrid modeling approach enables detailed component optimization while maintaining the ability to predict system-level dynamic behavior, particularly under closed-loop control conditions.
-
+The theoretical foundations and implementation are presented in Section \ref{sec:detail_fem_super_element}, where experimental validation was performed using an Amplified Piezoelectric Actuator.
+The framework was then applied to optimize two critical nano-hexapod elements: the actuators (Section \ref{sec:detail_fem_actuator}) and the flexible joints (Section \ref{sec:detail_fem_joint}).
+Through this approach, system-level dynamic behavior under closed-loop control conditions could be successfully predicted while detailed component-level optimization was facilitated.
 \chapter{Reduced order flexible bodies}
-\label{sec:orgfd22661}
+\label{sec:org9c118d2}
 \label{sec:detail_fem_super_element}
 Components exhibiting complex dynamical behavior are frequently found to be unsuitable for direct implementation within multi-body models.
 These components are traditionally analyzed using Finite Element Analysis (FEA) software.
 However, a methodological bridge between these two analytical approaches has been established, whereby components whose dynamical properties have been determined through FEA can be successfully integrated into multi-body models \cite{hatch00_vibrat_matlab_ansys}.
-This combined multibody-FEA modeling approach presents significant advantages, as it enables the selective application of FEA modeling to specific elements while maintaining the computational efficiency of multi-body analysis for the broader system \cite{rankers98_machin}.
+This combined multibody-FEA modeling approach presents significant advantages, as it enables the accurate FE modeling to specific elements while maintaining the computational efficiency of multi-body analysis for the broader system \cite{rankers98_machin}.
 
 The investigation of this hybrid modeling approach is structured in three sections.
 First, the fundamental principles and methodological approaches of this modeling framework are introduced (Section \ref{ssec:detail_fem_super_element_theory}).
 It is then illustrated through its practical application to the modelling of an Amplified Piezoelectric Actuator (APA) (Section \ref{ssec:detail_fem_super_element_example}).
 Finally, the validity of this modeling approach is demonstrated through experimental validation, wherein the obtained dynamics from the hybrid modelling approach is compared with measurements (Section \ref{ssec:detail_fem_super_element_validation}).
 \section{Procedure}
-\label{sec:org93ab665}
+\label{sec:orgd6255e7}
 \label{ssec:detail_fem_super_element_theory}
 
 In this modeling approach, some components within the multi-body framework are represented as \emph{reduced-order flexible bodies}, wherein their modal behavior is characterized through reduced mass and stiffness matrices derived from finite element analysis (FEA) models.
@@ -64,23 +49,22 @@ Initially, the component is modeled in a finite element software with appropriat
 Subsequently, interface frames are defined at locations where the multi-body model will establish connections with the component.
 These frames serve multiple functions, including connecting to other parts, applying forces and torques, and measuring relative motion between defined frames.
 
-Following the establishment of these interface parameters, modal reduction is performed using the Craig-Bampton method \cite{craig68_coupl_subst_dynam_analy} (also known as the ``fixed-interface method''), a technique that transforms the extensive FEA degrees of freedom into a significantly reduced set of retained degrees of freedom.
+Following the establishment of these interface parameters, modal reduction is performed using the Craig-Bampton method \cite{craig68_coupl_subst_dynam_analy} (also known as the ``fixed-interface method''), a technique that significantly reduces the number of DoF while while still presenting the main dynamical characteristics.
 This transformation typically reduces the model complexity from hundreds of thousands to fewer than 100 DoF.
 The number of degrees of freedom in the reduced model is determined by \eqref{eq:detail_fem_model_order} where \(n\) represents the number of defined frames and \(p\) denotes the number of additional modes to be modeled.
-The outcome of this procedure is an \(m \times m\) set of reduced mass and stiffness matrices, which can subsequently be incorporated into the multi-body model to represent the component's dynamic behavior.
+The outcome of this procedure is an \(m \times m\) set of reduced mass and stiffness matrices, \(m\) being the total retained number of degrees of freedom, which can subsequently be incorporated into the multi-body model to represent the component's dynamic behavior.
 
 \begin{equation}\label{eq:detail_fem_model_order}
 m = 6 \times n + p
 \end{equation}
-
 \section{Example with an Amplified Piezoelectric Actuator}
-\label{sec:org1e66a5f}
+\label{sec:org99a5f8a}
 \label{ssec:detail_fem_super_element_example}
 The presented modeling framework was first applied to an Amplified Piezoelectric Actuator (APA) for several reasons.
 Primarily, this actuator represents an excellent candidate for implementation within the nano-hexapod, as will be elaborated in Section \ref{sec:detail_fem_actuator}.
 Additionally, an Amplified Piezoelectric Actuator (the APA95ML shown in Figure \ref{fig:detail_fem_apa95ml_picture}) was available in the laboratory for experimental testing.
 
-The APA consists of multiple piezoelectric stacks arranged horizontally (depicted in blue in Figure \ref{fig:detail_fem_apa95ml_picture}) and of an amplifying shell structure (shown in red) that serves two purposes: the application of pre-stress to the piezoelectric elements and the amplification of their displacement into the vertical direction \cite{claeyssen07_amplif_piezoel_actuat}.
+The APA consists of multiple piezoelectric stacks arranged horizontally (depicted in blue in Figure \ref{fig:detail_fem_apa95ml_picture}) and of an amplifying shell structure (shown in red) that serves two purposes: the application of pre-stress to the piezoelectric elements and the amplification of their displacement in the vertical direction \cite{claeyssen07_amplif_piezoel_actuat}.
 The selection of the APA for validation purposes was further justified by its capacity to simultaneously demonstrate multiple aspects of the modeling framework.
 The specific design of the APA allows for the simultaneous modeling of a mechanical structure analogous to a flexible joint, piezoelectric actuation, and piezoelectric sensing, thereby encompassing the principal elements requiring validation.
 
@@ -104,10 +88,10 @@ Stiffness & \(21\,N/\mu m\)\\
 \end{tabularx}
 \captionof{table}{\label{tab:detail_fem_apa95ml_specs}APA95ML specifications}
 \end{minipage}
-\paragraph{Finite Element Model}
-\label{sec:orgce5afdb}
+\subsubsection{Finite Element Model}
+\label{sec:org9ffb40b}
 
-The development of the finite element model for the APA95ML necessitated the specification of appropriate material properties, as summarized in Table \ref{tab:detail_fem_material_properties}.
+The development of the finite element model for the APA95ML required the knowledge of the material properties, as summarized in Table \ref{tab:detail_fem_material_properties}.
 The finite element mesh, shown in Figure \ref{fig:detail_fem_apa95ml_mesh}, was then generated.
 
 \begin{table}[htbp]
@@ -123,7 +107,7 @@ Piezoelectric Ceramics (PZT) & \(49.5\,GPa\) & \(0.31\) & \(7800\,\text{kg}/m^3\
 \end{tabularx}
 \end{table}
 
-The definition of interface frames, or ``remote points'', constitute a critical aspect of the model preparation.
+The definition of interface frames constitutes a critical aspect of the model preparation.
 Seven frames were established: one frame at the two ends of each piezoelectric stack to facilitate strain measurement and force application, and additional frames at the top and bottom of the structure to enable connection with external elements in the multi-body simulation.
 
 Six additional modes were considered, resulting in total model order of \(48\).
@@ -134,19 +118,18 @@ The modal reduction procedure was then executed, yielding the reduced mass and s
 \begin{center}
 \includegraphics[scale=1,scale=1]{figs/detail_fem_apa95ml_mesh.png}
 \end{center}
-\subcaption{\label{fig:detail_fem_apa95ml_mesh}Obtained mesh and "remote points"}
+\subcaption{\label{fig:detail_fem_apa95ml_mesh} }
 \end{subfigure}
 \begin{subfigure}{0.48\textwidth}
 \begin{center}
 \includegraphics[scale=1,scale=1]{figs/detail_fem_apa_modal_schematic.png}
 \end{center}
-\subcaption{\label{fig:detail_fem_apa_model_schematic}Inclusion in multi-body model}
+\subcaption{\label{fig:detail_fem_apa_model_schematic} }
 \end{subfigure}
-\caption{\label{fig:detail_fem_apa95ml_model}Obtained mesh and defined interface frames (or ``remote points'') in the finite element model of the APA95ML (\subref{fig:detail_fem_apa95ml_mesh}). Interface with the multi-body model is shown in (\subref{fig:detail_fem_apa_modal_schematic}).}
+\caption{\label{fig:detail_fem_apa95ml_model}Obtained mesh and defined interface frames (or ``remote points'') in the finite element model of the APA95ML (\subref{fig:detail_fem_apa95ml_mesh}). Interface with the multi-body model is shown in (\subref{fig:detail_fem_apa_model_schematic}).}
 \end{figure}
-
-\paragraph{Super Element in the Multi-Body Model}
-\label{sec:org809b5e7}
+\subsubsection{Super Element in the Multi-Body Model}
+\label{sec:orgf1583cd}
 
 Previously computed reduced order mass and stiffness matrices were imported in a multi-body model block called ``Reduced Order Flexible Solid''.
 This block has several interface frames corresponding to the ones defined in the FEA software.
@@ -156,9 +139,8 @@ Therefore, a force source \(F_a\) operating between frames \(\{3\}\) and \(\{2\}
 This is illustrated in Figure \ref{fig:detail_fem_apa_model_schematic}.
 
 However, to have access to the physical voltage input of the actuators stacks \(V_a\) and to the generated voltage by the force sensor \(V_s\), conversion between the electrical and mechanical domains need to be determined.
-
-\paragraph{Sensor and Actuator ``constants''}
-\label{sec:orgb3a075f}
+\subsubsection{Sensor and Actuator ``constants''}
+\label{sec:orga72c195}
 
 To link the electrical domain to the mechanical domain, an ``actuator constant'' \(g_a\) and a ``sensor constant'' \(g_s\) were introduced as shown in Figure \ref{fig:detail_fem_apa_model_schematic}.
 
@@ -174,8 +156,8 @@ From \cite{fleming10_integ_strain_force_feedb_high} the relation between the for
   F_a = g_a \cdot V_a, \quad g_a = d_{33} n k_a, \quad k_a = \frac{c^{E} A}{L}
 \end{equation}
 
-Unfortunately, it is difficult to know exactly which material is used in the amplified piezoelectric actuator\footnote{The manufacturer of the APA95ML was not willing to share the piezoelectric material properties of the stack.}.
-However, based on the available properties of the stacks in the data-sheet (summarized in Table \ref{tab:detail_fem_stack_parameters}), the soft Lead Zirconate Titanate ``THP5H'' from Thorlabs seemed to match quite well the observed properties.
+Unfortunately, it is difficult to know exactly which material is used for the piezoelectric stacks\footnote{The manufacturer of the APA95ML was not willing to share the piezoelectric material properties of the stack.}.
+Yet, based on the available properties of the stacks in the data-sheet (summarized in Table \ref{tab:detail_fem_stack_parameters}), the soft Lead Zirconate Titanate ``THP5H'' from Thorlabs seemed to match quite well the observed properties.
 
 \begin{table}[htbp]
 \caption{\label{tab:detail_fem_stack_parameters}Stack Parameters}
@@ -195,11 +177,11 @@ Stack Area & \(mm^2\) & 10x10\\
 \end{tabularx}
 \end{table}
 
-The properties of this ``THP5H'' material used to compute \(g_a\) and \(g_s\) are listed in Table \ref{tab:test_apa_piezo_properties}.
+The properties of this ``THP5H'' material used to compute \(g_a\) and \(g_s\) are listed in Table \ref{tab:detail_fem_piezo_properties}.
 From these parameters, \(g_s = 5.1\,V/\mu m\) and \(g_a = 26\,N/V\) were obtained.
 
 \begin{table}[htbp]
-\caption{\label{tab:test_apa_piezo_properties}Piezoelectric properties used for the estimation of the sensor and actuators sensitivities}
+\caption{\label{tab:detail_fem_piezo_properties}Piezoelectric properties used for the estimation of the sensor and actuators sensitivities}
 \centering
 \begin{tabularx}{1\linewidth}{ccX}
 \toprule
@@ -215,9 +197,8 @@ From these parameters, \(g_s = 5.1\,V/\mu m\) and \(g_a = 26\,N/V\) were obtaine
 \bottomrule
 \end{tabularx}
 \end{table}
-
-\paragraph{Identification of the APA Characteristics}
-\label{sec:orge041867}
+\subsubsection{Identification of the APA Characteristics}
+\label{sec:org30104dc}
 
 Initial validation of the finite element model and its integration as a reduced-order flexible model within the multi-body model was accomplished through comparative analysis of key actuator characteristics against manufacturer specifications.
 
@@ -225,7 +206,7 @@ The stiffness of the APA95ML was estimated from the multi-body model by computin
 The inverse of the DC gain this transfer function corresponds to the axial stiffness of the APA95ML.
 A value of \(23\,N/\mu m\) was found which is close to the specified stiffness in the datasheet of \(k = 21\,N/\mu m\).
 
-The multi-body model predicted a resonant frequency under block-free conditions of \(2024\,\text{Hz}\) (Figure \ref{fig:detail_fem_apa95ml_compliance}), which is in agreement with the nominal specification of \(2000\,\text{Hz}\).
+The multi-body model predicted a resonant frequency under block-free conditions of \(\approx 2\,\text{kHz}\) (Figure \ref{fig:detail_fem_apa95ml_compliance}), which is in agreement with the nominal specification.
 
 \begin{figure}[htbp]
 \centering
@@ -233,7 +214,7 @@ The multi-body model predicted a resonant frequency under block-free conditions
 \caption{\label{fig:detail_fem_apa95ml_compliance}Estimated compliance of the APA95ML}
 \end{figure}
 
-In order to estimate the stroke of the APA95ML, first the mechanical amplification factor, defined as the ratio between vertical displacement and horizontal stack displacement, needs to be determined.
+In order to estimate the stroke of the APA95ML, the mechanical amplification factor, defined as the ratio between vertical displacement and horizontal stack displacement, was first determined.
 This characteristic was quantified through analysis of the transfer function relating horizontal stack motion to vertical actuator displacement, from which an amplification factor of \(1.5\) was derived.
 
 The piezoelectric stacks, exhibiting a typical strain response of \(0.1\,\%\) relative to their length (here equal to \(20\,mm\)), produce an individual nominal stroke of \(20\,\mu m\) (see data-sheet of the piezoelectric stacks on Table \ref{tab:detail_fem_stack_parameters}, page \pageref{tab:detail_fem_stack_parameters}).
@@ -241,44 +222,32 @@ As three stacks are used, the horizontal displacement is \(60\,\mu m\).
 Through the established amplification factor of 1.5, this translates to a predicted vertical stroke of \(90\,\mu m\) which falls within the manufacturer-specified range of \(80\,\mu m\) and \(120\,\mu m\).
 
 The high degree of concordance observed across multiple performance metrics provides a first validation of the ability to include FEM into multi-body model.
-
 \section{Experimental Validation}
-\label{sec:org354cea4}
+\label{sec:orga0c1c4b}
 \label{ssec:detail_fem_super_element_validation}
 Further validation of the reduced-order flexible body methodology was undertaken through experimental investigation.
-The goal is to measure the dynamics of the APA95ML and compared it with predictions derived from the multi-body model incorporating the actuator as a flexible element.
+The goal was to measure the dynamics of the APA95ML and to compare it with predictions derived from the multi-body model incorporating the actuator as a flexible element.
 
-The test bench illustrated in Figure \ref{fig:detail_fem_apa95ml_bench} was used, which consists of a \(5.7\,kg\) granite suspended on top of the APA95ML.
+The test bench illustrated in Figure \ref{fig:detail_fem_apa95ml_bench_schematic} was used, which consists of a \(5.7\,kg\) granite suspended on top of the APA95ML.
 The granite's motion was vertically guided with an air bearing system, and a fibered interferometer was used to measured its vertical displacement \(y\).
-A digital-to-analog converter (DAC) was used to generate the control signal \(u\), which was subsequently conditioned through a voltage amplification stage providing a gain factor of \(20\), ultimately yielding the effective voltage \(V_a\) across the two piezoelectric stacks.
+A digital-to-analog converter (DAC) was used to generate the control signal \(u\), which was subsequently conditioned through a voltage amplifier with a gain of \(20\), ultimately yielding the effective voltage \(V_a\) across the two piezoelectric stacks.
 Measurement of the sensor stack voltage \(V_s\) was performed using an analog-to-digital converter (ADC).
 
 \begin{figure}[htbp]
-\begin{subfigure}{0.34\textwidth}
-\begin{center}
-\includegraphics[scale=1,width=0.95\linewidth]{figs/detail_fem_apa95ml_bench_picture.png}
-\end{center}
-\subcaption{\label{fig:detail_fem_apa95ml_bench_picture}Picture of the test bench}
-\end{subfigure}
-\begin{subfigure}{0.72\textwidth}
-\begin{center}
-\includegraphics[scale=1,width=0.95\linewidth]{figs/detail_fem_apa95ml_bench_schematic.png}
-\end{center}
-\subcaption{\label{fig:detail_fem_apa95ml_bench_schematic}Schematic with signals}
-\end{subfigure}
-\caption{\label{fig:detail_fem_apa95ml_bench}Test bench used to validate ``reduced order solid bodies'' using an APA95ML. Picture of the bench is shown in (\subref{fig:detail_fem_apa95ml_bench_picture}). Schematic is shown in (\subref{fig:detail_fem_apa95ml_bench_schematic}).}
+\centering
+\includegraphics[scale=1,width=\linewidth]{figs/detail_fem_apa95ml_bench_schematic.png}
+\caption{\label{fig:detail_fem_apa95ml_bench_schematic}Test bench used to validate ``reduced order solid bodies'' using an APA95ML.}
 \end{figure}
-\paragraph{Comparison of the dynamics}
-\label{sec:orgc0ac08b}
+\subsubsection{Comparison of the dynamics}
+\label{sec:orgd50cd8c}
 
 Frequency domain system identification techniques were used to characterize the dynamic behavior of the APA95ML.
-The identification procedure necessitated careful choice of the excitation signal \cite[, chap. 5]{pintelon12_system_ident}.
-The most used ones are impulses (particularly suited to modal analysis), steps, random noise signals, and multi-sine excitations.
+The identification procedure required careful choice of the excitation signal \cite[, chap. 5]{pintelon12_system_ident}.
 During all this experimental work, random noise excitation was predominantly employed.
 
 The designed excitation signal is then generated and both input and output signals are synchronously acquired.
 From the obtained input and output data, the frequency response functions were derived.
-To improve the quality of the obtained frequency domain data, averaging and windowing were used \cite[, chap. 13]{pintelon12_system_ident}..
+To improve the quality of the obtained frequency domain data, averaging and windowing were used \cite[, chap. 13]{pintelon12_system_ident}.
 
 The obtained frequency response functions from \(V_a\) to \(V_s\) and to \(y\) are compared with the theoretical predictions derived from the multi-body model in Figure \ref{fig:detail_fem_apa95ml_comp_plant}.
 
@@ -300,11 +269,10 @@ Regarding the amplitude characteristics, the constants \(g_a\) and \(g_s\) could
 \end{center}
 \subcaption{\label{fig:detail_fem_apa95ml_comp_plant_sensor}from $V_a$ to $V_s$}
 \end{subfigure}
-\caption{\label{fig:detail_fem_apa95ml_comp_plant}Comparison of the measured frequency response functions and the identified dynamics from the finite element model of the APA95ML. Both for the dynamics from \(V_a\) to \(y\) (\subref{fig:detail_fem_apa95ml_comp_plant_actuator}) and from \(V_a\) to \(V_s\) (\subref{fig:detail_fem_apa95ml_comp_plant_sensor})}
+\caption{\label{fig:detail_fem_apa95ml_comp_plant}Comparison of the measured frequency response functions and the finite element model of the APA95ML. Both for the dynamics from \(V_a\) to \(y\) (\subref{fig:detail_fem_apa95ml_comp_plant_actuator}) and from \(V_a\) to \(V_s\) (\subref{fig:detail_fem_apa95ml_comp_plant_sensor})}
 \end{figure}
-
-\paragraph{Integral Force Feedback with APA}
-\label{sec:org4b65d74}
+\subsubsection{Integral Force Feedback with APA}
+\label{sec:org3d5f71f}
 
 To further validate this modeling methodology, its ability to predict closed-loop behavior was verified experimentally.
 Integral Force Feedback (IFF) was implemented using the force sensor stack, and the measured dynamics of the damped system were compared with model predictions across multiple feedback gains.
@@ -315,11 +283,11 @@ The IFF controller implementation, defined in equation \ref{eq:detail_fem_iff_co
   K_{\text{IFF}}(s) = \frac{g}{s + 2\cdot 2\pi} \cdot \frac{s}{s + 0.5 \cdot 2\pi}
 \end{equation}
 
-The theoretical damped dynamics of the closed-loop system was analyzed through using the model by computed the root locus plot shown in Figure \ref{fig:detail_fem_apa95ml_iff_root_locus}.
+The theoretical damped dynamics of the closed-loop system was estimated using the model by computed the root locus plot shown in Figure \ref{fig:detail_fem_apa95ml_iff_root_locus}.
 For experimental validation, six gain values were tested: \(g = [0,\,10,\,50,\,100,\,500,\,1000]\).
 The measured frequency responses for each gain configuration were compared with model predictions, as presented in Figure \ref{fig:detail_fem_apa95ml_damped_plants}.
 
-The close agreement between experimental measurements and theoretical predictions across all gain configurations demonstrates the model's capability to accurately predict both open-loop and closed-loop system dynamics, thereby validating its utility for control system design and analysis.
+The close agreement between experimental measurements and theoretical predictions across all gain configurations demonstrates the model's capability to accurately predict both open-loop and closed-loop system dynamics.
 
 \begin{figure}[htbp]
 \begin{subfigure}{0.48\textwidth}
@@ -334,44 +302,35 @@ The close agreement between experimental measurements and theoretical prediction
 \end{center}
 \subcaption{\label{fig:detail_fem_apa95ml_damped_plants}Damped plants}
 \end{subfigure}
-\caption{\label{fig:detail_fem_apa95ml_iff_results}Obtained results using Integral Force Feedback with the APA95ML. Obtained closed-loop poles as a function of the controller gain \(g\) are prediction by root Locus plot (\subref{fig:detail_fem_apa95ml_iff_root_locus}). Circles are predictions from the model while crosses are poles estimated from the experimental data. Damped plants estimated from the model (dashed curves) and measured ones (solid curves) are compared in (\subref{fig:detail_fem_apa95ml_damped_plants}) for all tested controller gains.}
+\caption{\label{fig:detail_fem_apa95ml_iff_results}Results using Integral Force Feedback with the APA95ML. Closed-loop poles as a function of the controller gain \(g\) are predicted by root Locus plot (\subref{fig:detail_fem_apa95ml_iff_root_locus}). Circles are predictions from the model while crosses are poles estimated from the experimental data. Damped plants estimated from the model (dashed curves) and measured ones (solid curves) are compared in (\subref{fig:detail_fem_apa95ml_damped_plants}) for all tested controller gains.}
 \end{figure}
-
 \section*{Conclusion}
-\label{sec:org105aef7}
-The modeling procedure presented in this section will demonstrate significant utility for the optimization of complex mechanical components within multi-body systems, particularly in the design of actuators (Section \ref{sec:detail_fem_actuator}) and flexible joints (Section \ref{sec:detail_fem_joint}).
-
-Through experimental validation using an Amplified Piezoelectric Actuator, the methodology has been shown to accurately predict both open-loop and closed-loop dynamic behavior, thereby establishing its reliability for component design and system analysis.
-
-While this modeling approach provides accurate predictions of component behavior, the resulting model order can become prohibitively high for practical time-domain simulations.
-This is exemplified by the nano-hexapod configuration, where the implementation of six Amplified Piezoelectric Actuators, each modeled with 48 degrees of freedom, yields 288 degrees of freedom only for the actuators.
-However, the methodology remains valuable for system analysis, as the extraction of frequency domain characteristics can be efficiently performed even with such high-order models.
+\label{sec:orgc9b9f6b}
+The experimental validation with an Amplified Piezoelectric Actuator confirms that this methodology accurately predicts both open-loop and closed-loop dynamic behaviors.
+This verification establishes its effectiveness for component design and system analysis applications.
 
+The approach will be especially beneficial for optimizing actuators (Section \ref{sec:detail_fem_actuator}) and flexible joints (Section \ref{sec:detail_fem_joint}) for the nano-hexapod.
 \chapter{Actuator Selection}
-\label{sec:org3ec4809}
+\label{sec:orgdbca6b6}
 \label{sec:detail_fem_actuator}
-The selection and modeling of actuators constitutes a critical step in the development of the nano-hexapod.
-This chapter presents the approach to actuator selection and modeling.
-First, specifications for the nano-hexapod actuators are derived from previous analyses, leading to the selection of the actuator type and ultimately to a specific model (Section \ref{ssec:detail_fem_actuator_specifications}).
-Then, the chosen actuator is modeled using the reduced-order flexible body approach developed in the previous section, enabling validation of this selection through detailed dynamical analysis (Section \ref{ssec:detail_fem_actuator_apa300ml}).
-Finally, a simplified two-degree-of-freedom model is developed to facilitate time-domain simulations while maintaining accurate representation of the actuator's essential characteristics (Section \ref{ssec:detail_fem_actuator_apa300ml_2dof}).
 \section{Choice of the Actuator based on Specifications}
-\label{sec:org929a34b}
+\label{sec:orgc816a2f}
 \label{ssec:detail_fem_actuator_specifications}
 
 The actuator selection process was driven by several critical requirements derived from previous dynamic analyses.
 A primary consideration is the actuator stiffness, which significantly impacts system dynamics through multiple mechanisms.
 The spindle rotation induces gyroscopic effects that modify plant dynamics and increase coupling, necessitating sufficient stiffness.
-Conversely, the actuator stiffness must be carefully limited to ensure the nano-hexapod's suspension modes remain below the problematic modes of the micro-stations to limit the coupling between the two structures.
+Conversely, the actuator stiffness must be carefully limited to ensure the nano-hexapod's suspension modes remain below the problematic modes of the micro-station to limit the coupling between the two structures.
 These competing requirements suggest an optimal stiffness of approximately \(1\,N/\mu m\).
 
 Additional specifications arise from the control strategy and physical constraints.
-The implementation of a HAC-LAC (High Authority Control-Low Authority Control) architecture necessitates integrated force sensing capability.
+The implementation of the decentralized Integral Force Feedback (IFF) architecture necessitates force sensors to be collocated with each actuator.
 The system's geometric constraints limit the actuator height to 50mm, given the nano-hexapod's maximum height of 95mm and the presence of flexible joints at each strut extremity.
-Furthermore, the actuator stroke must exceed the micro-station positioning errors while providing additional margin for mounting adjustments and operational flexibility, which is estimated at \(\approx 100\,\mu m\).
+Furthermore, the actuator stroke must exceed the micro-station positioning errors while providing additional margin for mounting adjustments and operational flexibility.
+An actuator stroke of \(\approx 100\,\mu m\) is therefore required.
 
-Three actuator technologies were evaluated (examples are shown in Figure \ref{fig:detail_fem_actuator_pictures}): voice coil actuators, piezoelectric stack actuators, and amplified piezoelectric actuators.
-Variable reluctance actuators were not considered despite their superior efficiency compared to voice coil actuators, as their inherent nonlinearity would introduce unnecessary control complexity.
+Three actuator technologies were evaluated (examples of such actuators are shown in Figure \ref{fig:detail_fem_actuator_pictures}): voice coil actuators, piezoelectric stack actuators, and amplified piezoelectric actuators.
+Variable reluctance actuators were not considered despite their superior efficiency compared to voice coil actuators, as their inherent nonlinearity would introduce control complexity.
 
 \begin{figure}[htbp]
 \begin{subfigure}{0.25\textwidth}
@@ -395,17 +354,17 @@ Variable reluctance actuators were not considered despite their superior efficie
 \caption{\label{fig:detail_fem_actuator_pictures}Example of actuators considered for the nano-hexapod. Voice coil from Sensata Technologies (\subref{fig:detail_fem_voice_coil_picture}). Piezoelectric stack actuator from Physik Instrumente (\subref{fig:detail_fem_piezo_picture}). Amplified Piezoelectric Actuator from DSM (\subref{fig:detail_fem_fpa_picture}).}
 \end{figure}
 
-Voice coil actuators (shown in Figure \ref{fig:detail_fem_voice_coil_picture}), when combined with flexure guides of wanted stiffness \(\approx 1\,N/\mu m\), would require forces above \(100\,N\) to achieve the specified \(100\,\mu m\) displacement.
-While these actuators offer excellent linearity and long strokes, the constant force requirement would result in significant steady-state current, leading to thermal loads that could compromise system stability.
-Their advantages were not considered adapted for this application, diminishing their benefits relative to piezoelectric solutions.
+Voice coil actuators (shown in Figure \ref{fig:detail_fem_voice_coil_picture}), when combined with flexure guides of wanted stiffness (\(\approx 1\,N/\mu m\)), would require forces in the order of \(100\,N\) to achieve the specified \(100\,\mu m\) displacement.
+While these actuators offer excellent linearity and long strokes capabilities, the constant force requirement would result in significant steady-state current, leading to thermal loads that could compromise system stability.
+Their advantages (linearity and long stroke) were not considered adapted for this application, diminishing their benefits relative to piezoelectric solutions.
 
 Conventional piezoelectric stack actuators (shown in Figure \ref{fig:detail_fem_piezo_picture}) present two significant limitations for the current application.
 Their stroke is inherently limited to approximately \(0.1\,\%\) of their length, meaning that even with the maximum allowable height of \(50\,mm\), the achievable stroke would only be \(50\,\mu m\), insufficient for the application.
 Additionally, their extremely high stiffness, typically around \(100\,N/\mu m\), exceeds the desired specifications by two orders of magnitude.
 
-Amplified Piezoelectric Actuators (APAs) emerged as the optimal solution by addressing these limitations through an specific mechanical design.
+Amplified Piezoelectric Actuators (APAs) emerged as the optimal solution by addressing these limitations through a specific mechanical design.
 The incorporation of a shell structure serves multiple purposes: it provides mechanical amplification of the piezoelectric displacement, reduces the effective axial stiffness to more suitable levels for the application, and creates a compact vertical profile.
-Furthermore, the multi-stack configuration enables one stack to be dedicated to force sensing, ensuring excellent collocation with the actuator stacks, a critical feature for implementing robust decentralized control strategies.
+Furthermore, the multi-stack configuration enables one stack to be dedicated to force sensing, ensuring excellent collocation with the actuator stacks, a critical feature for implementing robust decentralized IFF \cite{souleille18_concep_activ_mount_space_applic,verma20_dynam_stabil_thin_apert_light}.
 Moreover, using APA for active damping has been successfully demonstrated in similar applications \cite{hanieh03_activ_stewar}.
 
 Several specific APA models were evaluated against the established specifications (Table \ref{tab:detail_fem_piezo_act_models}).
@@ -429,9 +388,8 @@ Height  \(< 50\, [mm]\) & 22 & 30 & 24 & 27 & 16\\
 \bottomrule
 \end{tabularx}
 \end{table}
-
 \section{APA300ML - Reduced Order Flexible Body}
-\label{sec:orgeff9b1b}
+\label{sec:org3980880}
 \label{ssec:detail_fem_actuator_apa300ml}
 
 The validation of the APA300ML started by incorporating a ``reduced order flexible body'' into the multi-body model as explained in Section \ref{sec:detail_fem_super_element}.
@@ -457,12 +415,11 @@ While this high order provides excellent accuracy for validation purposes, it pr
 \end{figure}
 
 The sensor and actuator ``constants'' (\(g_s\) and \(g_a\)) derived in Section \ref{ssec:detail_fem_super_element_example} for the APA95ML were used for the APA300ML model, as both actuators employ identical piezoelectric stacks.
-
 \section{Simpler 2DoF Model of the APA300ML}
-\label{sec:org120c274}
+\label{sec:org8816367}
 \label{ssec:detail_fem_actuator_apa300ml_2dof}
 
-To facilitate efficient time-domain simulations while maintaining essential dynamic characteristics, a simplified two-degree-of-freedom model was developed, adapted from \cite{souleille18_concep_activ_mount_space_applic}.
+To facilitate efficient time-domain simulations while maintaining essential dynamic characteristics, a simplified two-degree-of-freedom model, adapted from \cite{souleille18_concep_activ_mount_space_applic}, was developed.
 
 This model, illustrated in Figure \ref{fig:detail_fem_apa_2dof_model}, comprises three components.
 The mechanical shell is characterized by its axial stiffness \(k_1\) and damping \(c_1\).
@@ -485,15 +442,15 @@ The model requires tuning of 8 parameters (\(k_1\), \(c_1\), \(k_e\), \(c_e\), \
 
 The shell parameters \(k_1\) and \(c_1\) were determined first through analysis of the zero in the \(V_a\) to \(V_s\) transfer function.
 The physical interpretation of this zero can be understood through Root Locus analysis: as controller gain increases, the poles of a closed-loop system converge to the open-loop zeros.
-In this context, the zero corresponds to the poles of the system with a theoretical infinite-gain controller that ensures zero force in the sensor stack.
-This condition effectively represents the dynamics of an APA without the force sensor stack.
+The open-loop zero therefore corresponds to the poles of the system with a theoretical infinite-gain controller that ensures zero force in the sensor stack.
+This condition effectively represents the dynamics of an APA without the force sensor stack (i.e. an APA with only the shell).
 This physical interpretation enables straightforward parameter tuning: \(k_1\) determines the frequency of the zero, while \(c_1\) defines its damping characteristic.
 
 The stack parameters (\(k_a\), \(c_a\), \(k_e\), \(c_e\)) were then derived from the first pole of the \(V_a\) to \(y\) response.
 Given that identical piezoelectric stacks are used for both sensing and actuation, the relationships \(k_e = 2k_a\) and \(c_e = 2c_a\) were enforced, reflecting the series configuration of the dual actuator stacks.
 Finally, the sensitivities \(g_s\) and \(g_a\) were adjusted to match the DC gains of the respective transfer functions.
 
-The resulting parameters, documented in Table \ref{tab:detail_fem_apa300ml_2dof_parameters}, yield dynamic behavior that closely matches the high-order finite element model, as demonstrated in Figure \ref{fig:detail_fem_apa300ml_comp_fem_2dof_fem_2dof}.
+The resulting parameters, listed in Table \ref{tab:detail_fem_apa300ml_2dof_parameters}, yield dynamic behavior that closely matches the high-order finite element model, as demonstrated in Figure \ref{fig:detail_fem_apa300ml_comp_fem_2dof_fem_2dof}.
 While higher-order modes and non-axial flexibility are not captured, the model accurately represents the fundamental dynamics within the operational frequency range.
 
 \begin{table}[htbp]
@@ -530,9 +487,8 @@ While higher-order modes and non-axial flexibility are not captured, the model a
 \end{subfigure}
 \caption{\label{fig:detail_fem_apa300ml_comp_fem_2dof_fem_2dof}Comparison of the transfer functions extracted from the finite element model of the APA300ML and of the 2DoF model. Both for the dynamics from \(V_a\) to \(d_i\) (\subref{fig:detail_fem_apa300ml_comp_fem_2dof_actuator}) and from \(V_a\) to \(V_s\) (\subref{fig:detail_fem_apa300ml_comp_fem_2dof_force_sensor})}
 \end{figure}
-
 \section{Electrical characteristics of the APA}
-\label{sec:orga7af5a1}
+\label{sec:org6794c26}
 \label{ssec:detail_fem_actuator_apa300ml_electrical}
 
 The behavior of piezoelectric actuators is characterized by coupled constitutive equations that establish relationships between electrical properties (charges, voltages) and mechanical properties (stress, strain) \cite[, chapter 5.5]{schmidt20_desig_high_perfor_mechat_third_revis_edition}.
@@ -540,7 +496,7 @@ The behavior of piezoelectric actuators is characterized by coupled constitutive
 To evaluate the impact of electrical boundary conditions on the system dynamics, experimental measurements were conducted using the APA95ML, comparing the transfer function from \(V_a\) to \(y\) under two distinct configurations.
 With the force sensor stack in open-circuit condition (analogous to voltage measurement with high input impedance) and in short-circuit condition (similar to charge measurement with low output impedance).
 As demonstrated in Figure \ref{fig:detail_fem_apa95ml_effect_electrical_boundaries}, short-circuiting the force sensor stack results in a minor decrease in resonance frequency.
-This relatively modest effect validates the simplifying assumption made in the model of the APA.
+The developed models of the APA do not represent such behavior, but as this effect is quite small, this validates the simplifying assumption made in the models.
 
 \begin{figure}[htbp]
 \centering
@@ -550,21 +506,20 @@ This relatively modest effect validates the simplifying assumption made in the m
 
 However, the electrical characteristics of the APA remain crucial for instrumentation design.
 Proper consideration must be given to voltage amplifier specifications and force sensor signal conditioning requirements.
-These aspects, being fundamental to system implementation, will be addressed in the instrumentation chapter.
-
+These aspects will be addressed in the instrumentation chapter.
 \section{Validation with the Nano-Hexapod}
-\label{sec:orgba951c9}
+\label{sec:org0a27bad}
 \label{ssec:detail_fem_actuator_apa300ml_validation}
 
 The integration of the APA300ML model within the nano-hexapod simulation framework served two validation objectives: to validate the APA300ML choice through analysis of system dynamics with APA modelled as flexible bodies, and to validate the simplified 2DoF model through comparative analysis with the full FEM implementation.
 
-The dynamic characteristics predicted using the flexible body model align well with the design requirements established during the conceptual phase.
+The dynamics predicted using the flexible body model align well with the design requirements established during the conceptual phase.
 The dynamics from \(\bm{u}\) to \(\bm{V}_s\) exhibits the desired alternating pole-zero pattern (Figure \ref{fig:detail_fem_actuator_fem_vs_perfect_hac_plant}), a critical characteristic for implementing robust decentralized Integral Force Feedback.
 Additionally, the model predicts no problematic high-frequency modes in the dynamics from \(\bm{u}\) to \(\bm{\epsilon}_{\mathcal{L}}\) (Figure \ref{fig:detail_fem_actuator_fem_vs_perfect_iff_plant}), maintaining consistency with earlier conceptual simulations.
 These findings suggest that the control performance targets established during the conceptual phase remain achievable with the selected actuator.
 
 Comparative analysis between the high-order FEM implementation and the simplified 2DoF model (Figure \ref{fig:detail_fem_actuator_fem_vs_perfect_plants}) demonstrates remarkable agreement in the frequency range of interest.
-This validates the use of the simplified model for time-domain simulations, where computational efficiency is paramount.
+This validates the use of the simplified model for time-domain simulations.
 The reduction in model order is substantial: while the FEM implementation results in approximately 300 states (36 states per actuator plus 12 additional states), the 2DoF model requires only 24 states for the complete nano-hexapod.
 
 These results validate both the selection of the APA300ML and the effectiveness of the simplified modeling approach for the nano-hexapod.
@@ -584,9 +539,8 @@ These results validate both the selection of the APA300ML and the effectiveness
 \end{subfigure}
 \caption{\label{fig:detail_fem_actuator_fem_vs_perfect_plants}Comparison of the dynamics obtained between a nano-hexpod having the actuators modeled with FEM and a nano-hexapod having actuators modelled a 2DoF system. Both from actuator force \(\bm{f}\) to strut motion measured by external metrology \(\bm{\epsilon}_{\mathcal{L}}\) (\subref{fig:detail_fem_actuator_fem_vs_perfect_iff_plant}) and to the force sensors \(\bm{f}_m\) (\subref{fig:detail_fem_actuator_fem_vs_perfect_hac_plant}).}
 \end{figure}
-
 \chapter{Flexible Joint Design}
-\label{sec:org93c9b2c}
+\label{sec:orgee992e3}
 \label{sec:detail_fem_joint}
 High-precision position control at the nanometer scale requires systems to be free from friction and backlash, as these nonlinear phenomena severely limit achievable positioning accuracy.
 This fundamental requirement prevents the use of conventional joints, necessitating instead the implementation of flexible joints that achieve motion through elastic deformation.
@@ -594,6 +548,12 @@ For Stewart platforms requiring nanometric precision, numerous flexible joint de
 For design simplicity and component standardization, identical joints are employed at both ends of the nano-hexapod struts.
 
 \begin{figure}[htbp]
+\begin{subfigure}{0.3\textwidth}
+\begin{center}
+\includegraphics[scale=1,height=5cm]{figs/detail_fem_joints_preumont.png}
+\end{center}
+\subcaption{\label{fig:detail_fem_joints_preumont}}
+\end{subfigure}
 \begin{subfigure}{0.35\textwidth}
 \begin{center}
 \includegraphics[scale=1,height=5cm]{figs/detail_fem_joints_yang.png}
@@ -602,17 +562,11 @@ For design simplicity and component standardization, identical joints are employ
 \end{subfigure}
 \begin{subfigure}{0.3\textwidth}
 \begin{center}
-\includegraphics[scale=1,height=5cm]{figs/detail_fem_joints_preumont.png}
-\end{center}
-\subcaption{\label{fig:detail_fem_joints_preumont}}
-\end{subfigure}
-\begin{subfigure}{0.3\textwidth}
-\begin{center}
 \includegraphics[scale=1,height=5cm]{figs/detail_fem_joints_wire.png}
 \end{center}
 \subcaption{\label{fig:detail_fem_joints_wire}}
 \end{subfigure}
-\caption{\label{fig:detail_fem_joints_examples}Example of different flexible joints geometry used for Stewart platforms. (\subref{fig:detail_fem_joints_yang}) \cite{yang19_dynam_model_decoup_contr_flexib}. (\subref{fig:detail_fem_joints_preumont}) \cite{preumont07_six_axis_singl_stage_activ}. (\subref{fig:detail_fem_joints_wire}) \cite{du14_piezo_actuat_high_precis_flexib}.}
+\caption{\label{fig:detail_fem_joints_examples}Example of different flexible joints geometry used for Stewart platforms. (\subref{fig:detail_fem_joints_preumont}) Typical ``universal'' flexible joint used in \cite{preumont07_six_axis_singl_stage_activ}. (\subref{fig:detail_fem_joints_yang}) Torsional stiffness can be explicitely specified as done in \cite{yang19_dynam_model_decoup_contr_flexib}. (\subref{fig:detail_fem_joints_wire}) ``Thin'' flexible joints having ``notch curves'' are also used \cite{du14_piezo_actuat_high_precis_flexib}.}
 \end{figure}
 
 While ideally these joints would permit free rotation about defined axes while maintaining infinite rigidity in other degrees of freedom, practical implementations exhibit parasitic stiffness that can impact control performance \cite{mcinroy02_model_desig_flexur_joint_stewar}.
@@ -620,16 +574,15 @@ This section examines how these non-ideal characteristics affect system behavior
 
 The analysis of bending and axial stiffness effects enables the establishment of comprehensive specifications for the flexible joints.
 These specifications guide the development and optimization of a flexible joint design through finite element analysis (Section \ref{ssec:detail_fem_joint_specs}).
-The validation process, detailed in Section \ref{ssec:detail_fem_joint_validation}, begins with the integration of the joints as ``reduced order flexible bodies'' in the nano-hexapod model, followed by the development of computationally efficient lower-order models that preserve the essential dynamic characteristics.
+The validation process, detailed in Section \ref{ssec:detail_fem_joint_validation}, begins with the integration of the joints as ``reduced order flexible bodies'' in the nano-hexapod model, followed by the development of computationally efficient lower-order models that preserve the essential dynamic characteristics of the flexible joints.
 \section{Bending and Torsional Stiffness}
-\label{sec:org582c93a}
+\label{sec:org76084a2}
 \label{ssec:detail_fem_joint_bending}
 
-The presence of bending stiffness in flexible joints causes the forces applied by the struts to deviate from the strut direction.
-This additional spring forces can affect system dynamics.
+The presence of bending stiffness in flexible joints causes the forces applied by the struts to deviate from the strut direction \cite{mcinroy02_model_desig_flexur_joint_stewar} and can affect system dynamics.
 
-To isolate and quantify these effects, simulations were conducted with the micro-station considered rigid and using simplified single-degree-of-freedom actuators (stiffness \(1\,N/\mu m\)) without parallel force sensor stiffness.
-Joint bending stiffness was varied from 0 (ideal case) to 500 Nm/rad.
+To quantify these effects, simulations were conducted with the micro-station considered rigid and using simplified 1DoF actuators (stiffness of \(1\,N/\mu m\)) without parallel stiffness to the force sensors.
+Flexible joint bending stiffness was varied from 0 (ideal case) to \(500\,Nm/\text{rad}\).
 
 Analysis of the plant dynamics reveals two significant effects.
 For the transfer function from \(\bm{f}\) to \(\bm{\epsilon}_{\mathcal{L}}\), bending stiffness increases low-frequency coupling, though this remains small for realistic stiffness values (Figure \ref{fig:detail_fem_joints_bending_stiffness_hac_plant}).
@@ -644,7 +597,7 @@ This characteristic impacts the achievable damping using decentralized Integral
 This is confirmed by the Root Locus plot in Figure \ref{fig:detail_fem_joints_bending_stiffness_iff_locus_1dof}.
 This effect becomes less significant when using the selected APA300ML actuators (Figure \ref{fig:detail_fem_joints_bending_stiffness_iff_locus_apa300ml}), which already incorporate parallel stiffness by design which is higher than the one induced by flexible joint stiffness.
 
-A parallel analysis of torsional stiffness revealed similar dynamic effects, though these proved less critical for system performance.
+A parallel analysis of torsional stiffness revealed similar effects, though these proved less critical for system performance.
 
 \begin{figure}[h!tbp]
 \begin{subfigure}{0.48\textwidth}
@@ -677,9 +630,8 @@ A parallel analysis of torsional stiffness revealed similar dynamic effects, tho
 \end{subfigure}
 \caption{\label{fig:detail_fem_joints_bending_stiffness_iff_locus}Effect of bending stiffness of the flexible joints on the attainable damping with decentralized IFF. When having an actuator modelled as 1DoF without parallel stiffness to the force sensor (\subref{fig:detail_fem_joints_bending_stiffness_iff_locus_1dof}), and with the 2DoF model of the APA300ML (\subref{fig:detail_fem_joints_bending_stiffness_iff_locus_apa300ml})}
 \end{figure}
-
 \section{Axial Stiffness}
-\label{sec:org05206a1}
+\label{sec:org4f7778b}
 \label{ssec:detail_fem_joint_axial}
 
 The limited axial stiffness (\(k_a\)) of flexible joints introduces an additional compliance between the actuation point and the measurement point.
@@ -692,7 +644,7 @@ The resulting frequency responses (Figure \ref{fig:detail_fem_joints_axial_stiff
 
 The force-sensor (IFF) plant exhibits minimal sensitivity to axial compliance, as evidenced by both frequency response data (Figure \ref{fig:detail_fem_joints_axial_stiffness_iff_plant}) and root locus analysis (Figure \ref{fig:detail_fem_joints_axial_stiffness_iff_locus}).
 
-However, the externally measured (HAC) plant demonstrates significant effects: internal strut modes appear at high frequencies, introducing substantial cross-coupling between axes.
+However, the transfer function from \(\bm{f}\) to \(\bm{\epsilon}_{\mathcal{L}}\) demonstrates significant effects: internal strut modes appear at high frequencies, introducing substantial cross-coupling between axes.
 This coupling is quantified through RGA analysis of the damped system (Figure \ref{fig:detail_fem_joints_axial_stiffness_rga_hac_plant}), which confirms increasing interaction between control channels at frequencies above the joint-induced resonance.
 
 Above this resonance frequency, two critical limitations emerge.
@@ -733,9 +685,8 @@ Based on this analysis, an axial stiffness specification of \(100\,N/\mu m\) was
 \end{subfigure}
 \caption{\label{fig:detail_fem_joints_axial_stiffness_iff_results}Effect of axial stiffness of the flexible joints on the attainable damping with decentralized IFF (\subref{fig:detail_fem_joints_axial_stiffness_iff_locus}). Estimation of the coupling of the damped plants using the RGA-number (\subref{fig:detail_fem_joints_axial_stiffness_rga_hac_plant})}
 \end{figure}
-
 \section{Specifications and Design flexible joints}
-\label{sec:org2310ef0}
+\label{sec:org1531b2d}
 \label{ssec:detail_fem_joint_specs}
 
 The design of flexible joints for precision applications requires careful consideration of multiple mechanical characteristics.
@@ -759,13 +710,13 @@ Bending Stroke & \(> 1\,\text{mrad}\) & 24.5\\
 \end{table}
 
 Among various possible flexible joint architectures, the design shown in Figure \ref{fig:detail_fem_joints_design} was selected for three key advantages.
-First, the geometry creates coincident x and y rotation axes, ensuring well-defined kinematic behavior through precise definition of the system's Jacobian matrix.
+First, the geometry creates coincident \(x\) and \(y\) rotation axes, ensuring well-defined kinematic behavior, important for the precise definition of the nano-hexapod Jacobian matrix.
 Second, the design allows easy tuning of different directional stiffnesses through a limited number of geometric parameters.
 Third, the architecture inherently provides high axial stiffness while maintaining the required compliance in rotational degrees of freedom.
 
 The joint geometry was optimized through parametric finite element analysis.
 The optimization process revealed an inherent trade-off between maximizing axial stiffness and achieving sufficiently low bending/torsional stiffness, while maintaining material stresses within acceptable limits.
-The final design, featuring a neck dimension of 0.25mm, achieves mechanical properties closely matching the target specifications, as verified through finite element analysis and documented in Table \ref{tab:detail_fem_joints_specs}.
+The final design, featuring a neck dimension of 0.25mm, achieves mechanical properties closely matching the target specifications, as verified through finite element analysis and summarized in Table \ref{tab:detail_fem_joints_specs}.
 
 \begin{figure}[htbp]
 \begin{subfigure}{0.48\textwidth}
@@ -782,9 +733,8 @@ The final design, featuring a neck dimension of 0.25mm, achieves mechanical prop
 \end{subfigure}
 \caption{\label{fig:detail_fem_joints_design}Designed flexible joints.}
 \end{figure}
-
 \section{Validation with the Nano-Hexapod}
-\label{sec:org4e972fc}
+\label{sec:orgdb5365d}
 \label{ssec:detail_fem_joint_validation}
 
 The designed flexible joint was first validated through integration into the nano-hexapod model using reduced-order flexible bodies derived from finite element analysis.
@@ -803,7 +753,7 @@ To improve computational efficiency, a low order representation was developed us
 
 After evaluating various configurations, a compromise was achieved by modeling bottom joints with bending and axial stiffness (\(k_f\) and \(k_a\)), and top joints with bending, torsional, and axial stiffness (\(k_f\), \(k_t\) and \(k_a\)).
 This simplification reduces the total model order to 48 states: 12 for the payload, 12 for the struts, and 24 for the joints (12 each for bottom and top joints).
-While additional degrees of freedom could potentially capture more dynamic features, the selected configuration preserves essential behavioral characteristics while minimizing computational complexity.
+While additional degrees of freedom could potentially capture more dynamic features, the selected configuration preserves essential system characteristics while minimizing computational complexity.
 
 \begin{figure}[htbp]
 \begin{subfigure}{0.48\textwidth}
@@ -820,14 +770,13 @@ While additional degrees of freedom could potentially capture more dynamic featu
 \end{subfigure}
 \caption{\label{fig:detail_fem_joints_fem_vs_perfect_plants}Comparison of the dynamics obtained between a nano-hexpod including joints modelled with FEM and a nano-hexapod having bottom joint modelled by bending stiffness \(k_f\) and axial stiffness \(k_a\) and top joints modelled by bending stiffness \(k_f\), torsion stiffness \(k_t\) and axial stiffness \(k_a\). Both from actuator force \(\bm{f}\) to strut motion measured by external metrology \(\bm{\epsilon}_{\mathcal{L}}\) (\subref{fig:detail_fem_joints_fem_vs_perfect_iff_plant}) and to the force sensors \(\bm{f}_m\) (\subref{fig:detail_fem_joints_fem_vs_perfect_hac_plant}).}
 \end{figure}
-
 \chapter*{Conclusion}
-\label{sec:org5382fe7}
+\label{sec:orgbb92db4}
 \label{sec:detail_fem_conclusion}
 
-In this chapter, the methodology of combining finite element analysis with multi-body modeling has been demonstrated and validated, proving particularly valuable for the detailed design phase of the nano-hexapod.
+In this chapter, the methodology of combining finite element analysis with multi-body modeling has been demonstrated and validated, proving particularly valuable for the detailed design of nano-hexapod components.
 The approach was first validated using an amplified piezoelectric actuator, where predicted dynamics showed excellent agreement with experimental measurements for both open and closed-loop behavior.
-This validation established confidence in the method's ability to accurately predict component behavior within the broader system context.
+This validation established confidence in the method's ability to accurately predict component behavior within a larger system.
 
 The methodology was then successfully applied to optimize two critical components.
 For the actuators, it enabled validation of the APA300ML selection while providing both high-fidelity and computationally efficient models for system simulation.
@@ -836,6 +785,5 @@ In both cases, the ability to seamlessly integrate finite element models into th
 
 A key outcome of this work is the development of reduced-order models that maintain prediction accuracy while enabling efficient time-domain simulation.
 Such model reduction, guided by detailed understanding of component behavior, provides the foundation for subsequent control system design and optimization.
-
 \printbibliography[heading=bibintoc,title={Bibliography}]
 \end{document}