Tangle matlab files and test them

matlab/test_joints_1_dim_meas.m

@@ -11,16 +11,17 @@ addpath('./mat/'); % Path for data
 colors = colororder;
 
 % Measurement Results
-% # - Strange shape: 5
+% The specified flexible beam thickness (gap) is $250\,\mu m$.
+% Four gaps are measured for each flexible joints (2 in the $x$ direction and 2 in the $y$ direction).
+% The "beam thickness" is then estimated to be the mean between the gaps measured on opposite sides.
 
-% The expected flexible beam thickness is $250\,\mu m$.
+% An histogram of the measured beam thicknesses is shown in Figure ref:fig:test_joints_size_hist.
-% However, it is more important that the thickness of all beams are close to each other.
+% The measured thickness is less thant the specified value of $250\,\mu m$, but this optical method may not be very accurate as the estimated gap can depend on the lighting of the part and of its proper alignment.
 
-% The dimension of the beams are been measured at each end to be able to estimate the mean of the beam thickness.
+% However, what is more important than the true value of the thickness is the consistency between all the flexible joints.
-
-% All the measured dimensions are summarized in Table ref:tab:test_joints_flex_dim.
 
 
+%% Measured gap for the 16 flexible joints
 meas_flex = [[223, 226, 224, 214];
              [229, 231, 237, 224];
              [234, 230, 239, 231];
@@ -38,38 +39,13 @@ meas_flex = [[223, 226, 224, 214];
              [209, 214, 220, 221];
              [213, 210, 230, 229]];
 
-
+%% Histogram of the measured gap
-
-% #+name: tab:test_joints_flex_dim
-% #+caption: Measured Dimensions of the flexible beams in $\mu m$
-% #+attr_latex: :environment tabularx :width 0.4\linewidth :align Xcccc
-% #+attr_latex: :center t :booktabs t :float t
-% #+RESULTS:
-% |    |  Y1 |  Y2 |  X1 |  X2 |
-% |----+-----+-----+-----+-----|
-% |  1 | 223 | 226 | 224 | 214 |
-% |  2 | 229 | 231 | 237 | 224 |
-% |  3 | 234 | 230 | 239 | 231 |
-% |  4 | 233 | 227 | 229 | 232 |
-% |  5 | 225 | 212 | 228 | 228 |
-% |  6 | 220 | 221 | 224 | 220 |
-% |  7 | 206 | 207 | 228 | 226 |
-% |  8 | 230 | 224 | 224 | 223 |
-% |  9 | 223 | 231 | 228 | 233 |
-% | 10 | 228 | 230 | 235 | 231 |
-% | 11 | 197 | 207 | 211 | 204 |
-% | 12 | 227 | 226 | 225 | 226 |
-% | 13 | 215 | 228 | 231 | 220 |
-% | 14 | 216 | 224 | 224 | 221 |
-% | 15 | 209 | 214 | 220 | 221 |
-% | 16 | 213 | 210 | 230 | 229 |
-
-% An histogram of these measured dimensions is shown in Figure ref:fig:test_joints_size_hist.
-
-
 figure;
+hold on;
 histogram([(meas_flex(:,1)+meas_flex(:,2))/2,(meas_flex(:,3)+meas_flex(:,4))/2], 7)
-xlabel("Beam's Thickness [$\mu m$]");
+hold off;
+xlabel("Measured beam thickness [$\mu m$]");
+xticks([200, 205, 210, 215, 220, 225, 230, 235])
 
 %% Save beam sizes
 save('./mat/flex_meas_dim.mat', 'meas_flex');

matlab/test_joints_2_bench_dimensioning.m

@@ -10,50 +10,123 @@ addpath('./mat/'); % Path for data
 %% Colors for the figures
 colors = colororder;
 
-% Flexible joint Geometry
-% The flexible joint used for the Nano-Hexapod is shown in Figure ref:fig:test_joints_bend_geometry.
-% Its bending stiffness is foreseen to be $k_{R_y}\approx 5\,\frac{Nm}{rad}$ and its stroke $\theta_{y,\text{max}}\approx 25\,mrad$.
-
-% #+name: fig:test_joints_bend_geometry
-% #+caption: Geometry of the flexible joint
-% [[file:figs/test_joints_bend_geometry.png]]
-
-% The height between the flexible point (center of the joint) and the point where external forces are applied is $h = 20\,mm$.
-
-% Let's define the parameters on Matlab.
-
-kRx = 5; % Bending Stiffness [Nm/rad]
-Rxmax = 25e-3; % Bending Stroke [rad]
-h = 20e-3; % Height [m]
-
 % Required external applied force
 
-% The bending $\theta_y$ of the flexible joint due to the force $F_x$ is:
+% The bending stiffness is foreseen to be $k_{R_y} \approx k_{R_x} \approx 5\,\frac{Nm}{rad}$ and its stroke $\theta_{y,\text{max}}\approx \theta_{x,\text{max}}\approx 25\,mrad$.
-% \begin{equation}
+% The height between the flexible point (center of the joint) and the point where external forces are applied is $h = 22.5\,mm$ (see Figure ref:fig:test_joints_bench_working_principle).
-%   \theta_y = \frac{M_y}{k_{R_y}} = \frac{F_x h}{k_{R_y}}
-% \end{equation}
-
-% Therefore, the applied force to test the full range of the flexible joint is:
-% \begin{equation}
-%   F_{x,\text{max}} = \frac{k_{R_y} \theta_{y,\text{max}}}{h}
-% \end{equation}
 
 
+%% Parameters for study
+kRx = 5;       % Bending Stiffness [Nm/rad]
+Rxmax = 25e-3; % Bending Stroke [rad]
+h = 22.5e-3;     % Height [m]
+
+%% Estimation of the force to test the full stroke
 Fxmax = kRx*Rxmax/h; % Force to induce maximum stroke [N]
 
-
-
-% And we obtain:
-
-sprintf('\\begin{equation} F_{x,max} = %.1f\\, [N] \\end{equation}', Fxmax)
-
 % Required actuator stroke and sensors range
 
 % The flexible joint is designed to allow a bending motion of $\pm 25\,mrad$.
-% The corresponding stroke at the location of the force sensor is:
+% The corresponding stroke at the location of the force sensor is given by eqref:eq:test_joints_max_stroke.
-% \[ d_{x,\text{max}} = h \tan(R_{x,\text{max}}) \]
+% In order to test the full range of the flexible joint, the means of applying a force (explained in the next section) should allow a motion of at least $0.5\,mm$.
+% Similarly, the measurement range of the displacement sensor should also be higher than $0.5\,mm$.
+
+% \begin{equation}\label{eq:test_joints_max_stroke}
+% d_{x,\text{max}} = h \tan(R_{x,\text{max}}) \approx 0.5\,mm
+% \end{equation}
 
 
+%% Estimated maximum stroke [m]
 dxmax = h*tan(Rxmax);
 
-sprintf('\\begin{equation} d_{max} = %.1f\\, [mm] \\end{equation}', 1e3*dxmax)
+% Introduction                                                  :ignore:
+% In order to estimate the accuracy of the measured bending stiffness that can be obtained using this measurement bench, an error budget is performed.
+
+% Based on equation eqref:eq:test_joints_stiff_displ_force, several errors can impact the accuracy of the measured bending stiffness:
+% - Errors in the measured torque $M_x, M_y$: this is mainly due to inaccuracies of the load cell and of the height estimation $h$
+% - Errors in the measured bending motion of the flexible joints $\theta_x, \theta_y$: errors from limited shear stiffness, from the deflection of the load cell itself, and from inaccuracy of the height estimation $h$
+
+% Let's first estimate the displacement induced only by the bending stiffness.
+
+% \begin{equation}\label{eq:test_joints_dbx}
+%   d_{x,b} = h \tan(\theta_y) = h \tan\left( \frac{F_x \cdot h}{k_{R_y}} \right)
+% \end{equation}
+
+% # From equation eqref:eq:test_joints_dbx, we can compute the bending stiffness:
+% # \begin{equation}
+% #   k_{R_y} = \frac{F_x \cdot h}{\tan^{-1}\left( \frac{D_b}{h} \right)}
+% # \end{equation}
+
+% # For small displacement, we have
+% # \begin{equation}
+% #   \boxed{k_{R_y} \approx h^2 \frac{F_x}{d_x}}
+% # \end{equation}
+
+% # And therefore, to precisely measure $k_{R_y}$, we need to:
+% # - precisely measure the motion $d_x$
+% # - precisely measure the applied force $F_x$
+% # - precisely now the height of the force application point $h$
+
+
+%% Stiffness
+ka = 94e6; % Axial Stiffness [N/m]
+ks = 13e6; % Shear Stiffness [N/m]
+kb = 5;    % Bending Stiffness [Nm/rad]
+kt = 260;  % Torsional Stiffness [Nm/rad]
+
+%% Maximum force
+Fa = 469;    % Axial Force before yield [N]
+Fs = 242;    % Shear Force before yield [N]
+Fb = 0.118;  % Bending Force before yield [Nm]
+Ft = 1.508; % Torsional Force before yield [Nm]
+
+%% Compute the corresponding stroke
+Xa = Fa/ka; % Axial Stroke before yield [m]
+Xs = Fs/ks; % Shear Stroke before yield [m]
+Xb = Fb/kb; % Bending Stroke before yield [rad]
+Xt = Ft/kt; % Torsional Stroke before yield [rad]
+
+%% Height between the joint's center and the force application point
+h = 22.5e-3; % [m]
+
+% Effect of Shear
+% The applied force $F_x$ will induce some shear $d_{x,s}$ which is described by eqref:eq:test_joints_shear_displ with $k_s$ the shear stiffness of the flexible joint.
+
+% \begin{equation}\label{eq:test_joints_shear_displ}
+%   d_{x,s} = \frac{F_x}{k_s}
+% \end{equation}
+
+% The measured displacement $d_x$ is affected by the shear as shown in equation eqref:eq:test_joints_displ_shear.
+
+% \begin{equation}\label{eq:test_joints_displ_shear}
+%   d_x = d_{x,b} + d_{x,s} = h \tan\left( \frac{F_x \cdot h}{k_{R_y}} \right) + \frac{F_x}{k_s} \approx F_x \left( \frac{h^2}{k_{R_y}} + \frac{1}{k_s} \right)
+% \end{equation}
+
+% The estimated bending stiffness $k_{\text{est}}$ then depends on the shear stiffness eqref:eq:test_joints_error_shear.
+
+% \begin{equation}\label{eq:test_joints_error_shear}
+%   k_{R_y,\text{est}} = h^2 \frac{F_x}{d_x} \approx k_{R_y} \frac{1}{1 + \frac{k_{R_y}}{k_s h^2}} \approx k_{R_y} \Bigl( 1 - \underbrace{\frac{k_{R_y}}{k_s h^2}}_{\epsilon_{s}} \Bigl)
+% \end{equation}
+
+% With an estimated shear stiffness $k_s = 13\,N/\mu m$ from the finite element model and an height $h=25\,mm$, the estimation errors of the bending stiffness due to shear is $\epsilon_s < 0.1\,\%$
+
+
+%% Estimated error due to shear
+epsilon_s = 100*abs(1-1/(1 + kb/(ks*h^2))); % Error in %
+
+% Effect of load cell limited stiffness
+% As explained in the previous section, because the measurement of the flexible joint deflection is indirectly made with the encoder, errors will be made if the load cell experiences some compression.
+
+% Suppose the load cell has an internal stiffness $k_f$, the same reasoning that was made for the effect of shear can be applied here.
+% The estimation error of the bending stiffness due to the limited stiffness of the load cell is then described by eqref:eq:test_joints_error_load_cell_stiffness.
+
+% \begin{equation}\label{eq:test_joints_error_load_cell_stiffness}
+%   k_{R_y,\text{est}} = h^2 \frac{F_x}{d_x} \approx k_{R_y} \frac{1}{1 + \frac{k_{R_y}}{k_F h^2}} \approx k_{R_y} \Bigl( 1 - \underbrace{\frac{k_{R_y}}{k_F h^2}}_{\epsilon_f} \Bigl)
+% \end{equation}
+
+% With an estimated load cell stiffness of $k_f \approx 1\,N/\mu m$ (from the documentation), the errors due to the load cell limited stiffness is around $\epsilon_f = 1\,\%$.
+
+
+%% Estimated error due to limited load cell stiffness
+kF = 50/0.05e-3; % Estimated load cell stiffness [N/m]
+epsilon_f = 100*abs(1-1/(1 + kb/(kF*h^2))); % Error in %

matlab/test_joints_3_bending_stiff_meas.m

@@ -0,0 +1,336 @@
+%% Clear Workspace and Close figures
+clear; close all; clc;
+
+%% Intialize Laplace variable
+s = zpk('s');
+
+%% Path for functions, data and scripts
+addpath('./mat/'); % Path for data
+
+%% Colors for the figures
+colors = colororder;
+
+% Load Cell Calibration
+% In order to estimate the measured errors of the load cell "FC2231", it is compared against another load cell[fn:5].
+% The two load cells are measured simultaneously while they are pushed against each other (see Figure ref:fig:test_joints_force_sensor_calib_picture).
+% The contact between the two load cells is well defined as one has a spherical interface while the other has a flat surface.
+
+% The measured forces are compared in Figure ref:fig:test_joints_force_sensor_calib_fit.
+% The gain mismatch between the two load cells is approximately $4\,\%$ which is higher than what was specified in the data-sheets.
+% However, the estimated non-linearity is bellow $1\,\%$ for forces between $0.5\,N$ and $5\,N$.
+
+
+%% Force Sensor Calibration
+% Load measurement data
+load('calibration_force_sensor.mat', 't', 'F', 'Fc')
+
+% Remove any offset such that they are both measuring no force when not in contact.
+F  = F  - mean(F( t > 0.5 & t < 1.0)); % FC2231 [N]
+Fc = Fc - mean(Fc(t > 0.5 & t < 1.0)); % XFL212R [N]
+
+% Only get useful stroke
+F  = F( t > 1.55 & t < 2.4);
+Fc = Fc(t > 1.55 & t < 2.4);
+
+% Make a line fit
+fit_F = polyfit(Fc, F, 1);
+
+% Estimate the gain mismatch
+F_gain_mismatch = 100*(1 - fit_F(1)); % in %
+
+% Estimate non-linearity of the sensors
+F_non_linearity = 100*(F - (Fc*fit_F(1) + fit_F(2)))./F; % in %
+
+%% Measured two forces and linear fit
+figure;
+hold on;
+plot(Fc, F, 'k-', 'DisplayName', 'Raw Data');
+plot(Fc([1,end]), Fc([1,end])*fit_F(1) + fit_F(2), '--', 'DisplayName', 'Line Fit');
+hold off;
+xlabel('XFL212R [N]'); ylabel('FC2231 [N]');
+xlim([0.5,5.5]); ylim([0.5,5.5]);
+legend('location', 'southeast', 'FontSize', 8);
+
+% Load Cell Stiffness
+% The objective of this measurement is to estimate the stiffness $k_F$ of the force sensor.
+% To do so, a stiff element (much stiffer than the estimated $k_F \approx 1\,N/\mu m$) is fixed in front of the force sensor as shown in Figure ref:fig:test_joints_meas_force_sensor_stiffness_picture.
+% Then, the force sensor is pushed again this stiff element while the force and the encoder displacement are measured.
+% Measured displacement as a function of the force is shown in Figure ref:fig:test_joints_force_sensor_stiffness_fit.
+% The load cell stiffness can then be estimated by computing a linear fit, and is found to be $k_F \approx 0.75\,N/\mu m$.
+
+
+%% Estimaetd load cell stiffness
+% Load measurement data
+load('force_sensor_stiffness_meas.mat', 't', 'F', 'd')
+
+% Remove offset
+F  = F - mean(F(t > 0.5 & t < 1.0));
+
+% Select important part of data
+F  = F( t > 4.55 & t < 7.24);
+d  = d( t > 4.55 & t < 7.24); d = d - d(1);
+t  = t( t > 4.55 & t < 7.24);
+
+% Linear fit
+fit_k = polyfit(F, d, 1);
+
+%% Displacement as a function of the measured force
+figure;
+hold on;
+plot(F, 1e6*d, 'k-', 'DisplayName', 'Raw Data');
+plot(F([1,end]), 1e6*(F([1,end])*fit_k(1) + fit_k(2)), '--', 'DisplayName', sprintf('Fit, $k_F \\approx %.2f N/\\mu m$', 1e-6/fit_k(1)));
+hold off;
+xlabel('Force [$N$]'); ylabel('Displacement [$\mu m$]');
+xlim([0,45]); ylim([0,60]);
+legend('location', 'southeast', 'FontSize', 8);
+
+% Bending Stiffness estimation
+% The actual stiffness measurement in now performed by manually moving the translation stage from a start position where the force sensor is not yet in contact with the flexible joint to a position where flexible joint is on its mechanical stop.
+
+% The measured force and displacement as a function of time are shown in Figure ref:fig:test_joints_meas_bend_time.
+% Three regions can be observed: first the force sensor tip is not in contact with the flexible joint and the measured force is zero, then the flexible joint deforms linearly, finally the flexible joint comes in contact with the mechanical stops.
+
+% The angular motion $\theta_{y}$ computed from the displacement $d_x$ is displacement as function of the measured torque $T_{y}$ in displayed in Figure ref:fig:test_joints_meas_F_d_lin_fit.
+% The bending stiffness of the flexible joint can be estimated by computing the slope of the curve in the linear regime (red dashed line) and is found to be $k_{R_y} = 4.4\,Nm/\text{rad}$.
+% The bending stroke can also be estimated as shown in Figure ref:fig:test_joints_meas_F_d_lin_fit and is found to be $\theta_{y,\text{max}} = 20.9\,\text{mrad}$.
+
+
+%% Estimate the bending stiffness and stroke from the measurement - First Flexible joint
+% Load Measured Data for first flexible joint
+load('meas_stiff_flex_1_x.mat', 't', 'F', 'd');
+
+% Start measurement at t = 0.2 s
+d = d(t > 0.2);
+F = F(t > 0.2);
+t = t(t > 0.2); t = t - t(1);
+
+% Zero the force
+F = F - mean(F(t < 0.2));
+
+% Find when the force sensor touches the flexible joint
+i_l_start = find(F > 0.3, 1, 'first');
+
+% Compute torque and angular displacement
+h = 22.5e-3; % Height [m]
+Tx = h * F; % Applied torque in [Nm]
+thetax = atan2(d - d(i_l_start), h); % Measured angle in [rad]
+
+% Find then the maximum torque is applied
+[~, i_s_stop] = max(Tx);
+% Linear region stops ~ when 90% of the stroke is reached
+i_l_stop = find(thetax > 0.9*thetax(i_s_stop), 1, 'first');
+% Mechanical "Stop" region start ~20Nmm before maximum torque is applied
+i_s_start = find(Tx > max(Tx)-20e-3, 1, 'first');
+
+% Linear fit in the "linear" region
+fit_l = polyfit(Tx(i_l_start:i_l_stop), thetax(i_l_start:i_l_stop), 1);
+
+% Linear fit in the "mechanical stop" region
+fit_s = polyfit(Tx(i_s_start:i_s_stop), thetax(i_s_start:i_s_stop), 1);
+
+% Reset displacement more precisely based on fit
+thetax = thetax - fit_l(2);
+fit_s(2) = fit_s(2) - fit_l(2);
+fit_l(2) = 0;
+
+%% Estimation of the bending stiffness
+kRx_l = 1/fit_l(1); % Bending Stiffness [Nm/rad]
+kRx_s = 1/fit_s(1); % Mechanical "Stop" Stiffness [Nm/rad]
+
+%% Estimation of the bending stroke
+% This is done by finding the intersection of the two linear fits
+theta_max = fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1)); % Maximum angular stroke [rad]
+Tx_at_theta_max = (fit_s(2) - fit_l(2))/(fit_l(1) - fit_s(1));
+
+%% Measured force and displacement as a function of time
+figure;
+yyaxis left
+hold on;
+plot(t, F);
+plot(0.8*cos(0:0.01:2*pi)+0.8, ...
+     2.8*sin(0:0.01:2*pi)-1, 'k--', 'HandleVisibility', 'off');
+text(1.8, -1.2, sprintf('Not in\ncontact'), 'horizontalalignment', 'left');
+plot(0.4*cos(0:0.01:2*pi)+3, ...
+     1.1*sin(0:0.01:2*pi)+4.8, 'k--');
+text(3.5, 4.8, sprintf('Mechanical\nStop'), 'horizontalalignment', 'left');
+hold off;
+ylabel('Force $F_y$ [N]');
+ylim([-4, 6])
+ylimr = get(gca,'Ylim');
+
+yyaxis right
+plot(t, 1e3*d);
+xlabel('Time [s]');
+ylabel('Displacement $d_y$ [mm]');
+xlim([0,5]);
+% Make the force and displacement superimpose
+ylim([0.364 - 4*(0.8315-0.364)/4.095, 0.364 + 6*(0.8315-0.364)/4.095])
+
+% Regions where to plot the fitted data
+Tx_l_fit = [0, Tx_at_theta_max];
+Tx_s_fit = [Tx_at_theta_max, Tx(i_s_stop)];
+
+figure;
+hold on;
+plot(Tx(1:i_s_stop), 1e3*thetax(1:i_s_stop), '-k', 'DisplayName', 'Raw data')
+plot(Tx_l_fit, 1e3*(Tx_l_fit*fit_l(1) + fit_l(2)), '--', 'DisplayName', sprintf('$k_{R_x} = %.1f$ Nm/rad', kRx_l))
+plot(Tx_s_fit, 1e3*(Tx_s_fit*fit_s(1) + fit_s(2)), '--', 'DisplayName', sprintf('$k_{R_x,stop} = %.0f$ Nm/rad', kRx_s))
+plot([0, Tx_at_theta_max], [1e3*theta_max, 1e3*theta_max], 'k--', 'HandleVisibility', 'off')
+plot([0, Tx_at_theta_max], [0, 0], 'k--', 'HandleVisibility', 'off')
+anArrow = annotation('doublearrow', 'LineWidth', 0.5);
+anArrow.Parent = gca;
+anArrow.Position = [0.05, 0, 0, 1e3*fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1))];
+text(0.052, 0.4*1e3*fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1)), sprintf('$\\theta_{x,\\max} = %.1f$ mrad', 1e3*theta_max), 'horizontalalignment', 'left');
+hold off;
+xlabel('Torque $T_x$ [Nm]');
+ylabel('Angle $\theta_x$ [mrad]');
+xlim([0, 0.15]);
+ylim([-5,25]);
+leg = legend('location', 'southeast', 'FontSize', 8);
+leg.ItemTokenSize(1) = 15;
+
+% Measured flexible joint stiffnesses
+% The same measurement is performed for all the 16 flexible joints, both in the $x$ and $y$ directions.
+% The measured angular motion as a function of the applied torque are shown in Figure ref:fig:test_joints_meas_bending_all_raw_data for all the 16 flexible joints.
+% This gives a first idea of the dispersion of the measured bending stiffnesses (i.e. slope of the linear region) and of the angular stroke.
+
+% An histogram of the measured bending stiffnesses is show in Figure ref:fig:test_joints_bend_stiff_hist.
+% Most of the bending stiffnesses are between $4.6\,Nm/rad$ and $5.0\,Nm/rad$.
+
+
+%% Measure the bending stiffness and Stroke for all the flexible joints
+figure;
+hold on;
+
+%% Start with the X-bending
+% Initialize variables
+Rx = zeros(1,16); % Bending stiffnesses [Nm/rad]
+kSx = zeros(1,16); % Bending stiffnesses at "stop" [Nm/rad]
+Rmx = zeros(1,16); % Bending stroke [rad]
+
+for i = 1:16
+    % Load the data
+    load(sprintf('meas_stiff_flex_%i_x.mat', i), 't', 'F', 'd');
+
+    % Start measurement at t = 0.2 s
+    d = d(t > 0.2);
+    F = F(t > 0.2);
+    t = t(t > 0.2); t = t - t(1);
+
+    % Zero the force
+    F = F - mean(F(t < 0.2));
+
+    % Find when the force sensor touches the flexible joint
+    i_l_start = find(F > 0.3, 1, 'first');
+
+    % Zero the displacement when it comes in contact
+    d = d - d(i_l_start);
+
+    % Compute torque and angular displacement
+    h = 22.5e-3; % Height [m]
+    Tx = h * F; % Applied torque in [Nm]
+    thetax = atan2(d, h); % Measured angle in [rad]
+
+    % Find then the maximum torque is applied
+    [~, i_s_stop] = max(Tx);
+    % Linear region stops ~ when 90% of the stroke is reached
+    i_l_stop = find(thetax > 0.9*thetax(i_s_stop), 1, 'first');
+    % Mechanical "Stop" region start ~20Nmm before maximum torque is applied
+    i_s_start = find(Tx > max(Tx)-20e-3, 1, 'first');
+
+    % Linear fit in the "linear" region
+    fit_l = polyfit(Tx(i_l_start:i_l_stop), thetax(i_l_start:i_l_stop), 1);
+
+    % Linear fit in the "mechanical stop" region
+    fit_s = polyfit(Tx(i_s_start:i_s_stop), thetax(i_s_start:i_s_stop), 1);
+
+    % Reset displacement more precisely based on fit
+    thetax = thetax - fit_l(2);
+    fit_s(2) = fit_s(2) - fit_l(2);
+    fit_l(2) = 0;
+
+    % Estimation of the bending stiffness and bending stroke
+    kRx(i) = 1/fit_l(1); % Bending Stiffness [Nm/rad]
+    kSx(i) = 1/fit_s(1); % Mechanical "Stop" Stiffness [Nm/rad]
+    Rmx(i) = fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1)); % Maximum angular stroke [rad]
+
+    if i == 1
+        plot(Tx(1:10:i_s_stop), 1e3*thetax(1:10:i_s_stop), '-', 'color', [colors(1,:), 0.4], ...
+            'DisplayName', '$k_{R_x}$')
+    else
+        plot(Tx(1:10:i_s_stop), 1e3*thetax(1:10:i_s_stop), '-', 'color', [colors(1,:), 0.4], ...
+            'HandleVisibility', 'off')
+    end
+end
+
+%% Continue with the Y-bending
+% Initialize variables
+kRy = zeros(1,16); % Bending stiffnesses [Nm/rad]
+kSy = zeros(1,16); % Bending stiffnesses at "stop" [Nm/rad]
+Rmy = zeros(1,16); % Bending stroke [rad]
+
+for i = 1:16
+    % Load the data
+    load(sprintf('meas_stiff_flex_%i_y.mat', i), 't', 'F', 'd');
+
+    % Start measurement at t = 0.2 s
+    d = d(t > 0.2);
+    F = F(t > 0.2);
+    t = t(t > 0.2); t = t - t(1);
+
+    % Zero the force
+    F = F - mean(F(t < 0.2));
+
+    % Find when the force sensor touches the flexible joint
+    i_l_start = find(F > 0.3, 1, 'first');
+
+    % Zero the displacement when it comes in contact
+    d = d - d(i_l_start);
+
+    % Compute torque and angular displacement
+    h = 22.5e-3; % Height [m]
+    Ty = h * F; % Applied torque in [Nm]
+    thetay = atan2(d, h); % Measured angle in [rad]
+
+    % Find then the maximum torque is applied
+    [~, i_s_stop] = max(Ty);
+    % Linear region stops ~ when 90% of the stroke is reached
+    i_l_stop = find(thetay > 0.9*thetay(i_s_stop), 1, 'first');
+    % Mechanical "Stop" region start ~20Nmm before maximum torque is applied
+    i_s_start = find(Ty > max(Ty)-20e-3, 1, 'first');
+
+    % Linear fit in the "linear" region
+    fit_l = polyfit(Ty(i_l_start:i_l_stop), thetay(i_l_start:i_l_stop), 1);
+
+    % Linear fit in the "mechanical stop" region
+    fit_s = polyfit(Ty(i_s_start:i_s_stop), thetay(i_s_start:i_s_stop), 1);
+
+    % Reset displacement more precisely based on fit
+    thetay = thetay - fit_l(2);
+    fit_s(2) = fit_s(2) - fit_l(2);
+    fit_l(2) = 0;
+
+    % Estimation of the bending stiffness and bending stroke
+    kRy(i) = 1/fit_l(1); % Bending Stiffness [Nm/rad]
+    kSy(i) = 1/fit_s(1); % Mechanical "Stop" Stiffness [Nm/rad]
+    Rmy(i) = fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1)); % Maximum angular stroke [rad]
+
+    if i == 1
+        plot(Ty(1:10:i_s_stop), 1e3*thetay(1:10:i_s_stop), '-', 'color', [colors(2,:), 0.4], ...
+            'DisplayName', '$k_{R_y}$')
+    else
+        plot(Ty(1:10:i_s_stop), 1e3*thetay(1:10:i_s_stop), '-', 'color', [colors(2,:), 0.4], ...
+            'HandleVisibility', 'off')
+    end
+end
+
+xlabel('Torque $T$ [Nm]');
+ylabel('Angle $\theta$ [mrad]');
+xlim([0, 0.15]);
+ylim([-5,25]);
+legend('location', 'southeast', 'FontSize', 8);
+
+figure;
+histogram([kRx, kRy], [4:0.2:5])
+xlabel('Bending stiffness [Nm/rad]')
+xticks([4:0.2:5])

matlab/test_joints_3_error_budget.m

@@ -1,73 +0,0 @@
-%% Clear Workspace and Close figures
-clear; close all; clc;
-
-%% Intialize Laplace variable
-s = zpk('s');
-
-%% Path for functions, data and scripts
-addpath('./mat/'); % Path for data
-
-%% Colors for the figures
-colors = colororder;
-
-% Finite Element Model
-% From the Finite Element Model, the stiffness and stroke of the flexible joint have been computed and summarized in Tables ref:tab:test_joints_axial_shear_prop and ref:tab:test_joints_bending_torsion_prop.
-
-
-%% Stiffness
-ka = 94e6; % Axial Stiffness [N/m]
-ks = 13e6; % Shear Stiffness [N/m]
-kb = 5;    % Bending Stiffness [Nm/rad]
-kt = 260;  % Torsional Stiffness [Nm/rad]
-
-%% Maximum force
-Fa = 469;    % Axial Force before yield [N]
-Fs = 242;    % Shear Force before yield [N]
-Fb = 0.118;  % Bending Force before yield [Nm]
-Ft = 1.508; % Torsional Force before yield [Nm]
-
-%% Compute the corresponding stroke
-Xa = Fa/ka; % Axial Stroke before yield [m]
-Xs = Fs/ks; % Shear Stroke before yield [m]
-Xb = Fb/kb; % Bending Stroke before yield [rad]
-Xt = Ft/kt; % Torsional Stroke before yield [rad]
-
-% Setup
-
-% The setup is schematically represented in Figure ref:fig:test_joints_bench_side_bis.
-
-% The force is applied on top of the flexible joint with a distance $h$ with the joint's center.
-% The displacement of the flexible joint is also measured at the same height.
-
-% The height between the joint's center and the force application point is:
-
-h = 25e-3; % Height [m]
-
-% Estimation error due to force sensor compression
-% The measured displacement is not done directly at the joint's location.
-% The force sensor compression will then induce an error on the joint's stiffness.
-
-% The force sensor stiffness $k_F$ is estimated to be around:
-
-kF = 50/0.05e-3; % [N/m]
-
-sprintf('k_F = %.1e [N/m]', kF)
-
-% Estimation error due to height estimation error
-% Let's consider an error in the estimation of the height from the application of the force to the joint's center:
-% \begin{equation}
-% h_{\text{est}} = h (1 + \epsilon)
-% \end{equation}
-
-% The computed bending stiffness will be:
-% \begin{equation}
-%   k_\text{est} \approx h_{\text{est}}^2 \frac{F_x}{d_x}
-% \end{equation}
-
-% And the stiffness estimation error is:
-% \begin{equation}
-% \frac{k_{\text{est}}}{k_{R_y}} = (1 + \epsilon)^2
-% \end{equation}
-
-
-h_err = 0.2e-3; % Height estimation error [m]

matlab/test_joints_4_first_meas.m

@@ -1,175 +0,0 @@
-%% Clear Workspace and Close figures
-clear; close all; clc;
-
-%% Intialize Laplace variable
-s = zpk('s');
-
-%% Path for functions, data and scripts
-addpath('./mat/'); % Path for data
-
-%% Colors for the figures
-colors = colororder;
-
-% Force Sensor Calibration
-
-% #+begin_note
-% *Load Cells*:
-% - [[file:doc/A700000007147087.pdf][FC2231-0000-0010-L]]
-% - [[file:doc/FRE_DS_XFL212R_FR_A3.pdf][XFL212R]]
-% #+end_note
-
-% There are both specified to have $\pm 1 \%$ of non-linearity over the full range.
-
-% The XFL212R has a spherical interface while the FC2231 has a flat surface.
-% Therefore, we should have a nice point contact when using the two force sensors as shown in Figure ref:fig:test_joints_force_sensor_calib.
-
-% #+name: fig:test_joints_force_sensor_calib
-% #+caption: Zoom on the two force sensors in contact
-% #+attr_latex: :width 0.8\linewidth
-% [[file:figs/test_joints_force_sensor_calib.jpg]]
-
-% The two force sensors are therefore measuring the exact same force, and we can compare the two measurements.
-
-% Let's load the measured force of both sensors.
-
-%% Load measurement data
-load('calibration_force_sensor.mat', 't', 'F', 'Fc')
-
-
-
-% We remove any offset such that they are both measuring no force when not in contact.
-
-%% Remove offset
-F  = F  - mean(F( t > 0.5 & t < 1.0));
-Fc = Fc - mean(Fc(t > 0.5 & t < 1.0));
-
-figure;
-hold on;
-plot(t, F, 'DisplayName', 'FC2231');
-plot(t, Fc, 'DisplayName', 'XFL212R');
-hold off;
-xlabel('Time [s]'); ylabel('Measured Force [N]');
-xlim([0,15]); ylim([0,55]);
-legend('location', 'southeast');
-
-
-
-% #+name: fig:test_joints_force_sensor_calib_time
-% #+caption: Measured force using both sensors as a function of time
-% #+RESULTS:
-% [[file:figs/test_joints_force_sensor_calib_time.png]]
-
-% Let's select only the first part from the moment they are in contact until the maximum force is reached.
-
-
-%% Only get the first part until maximum force
-F  = F( t > 1.55 & t < 4.65);
-Fc = Fc(t > 1.55 & t < 4.65);
-
-
-
-% Then, let's make a linear fit between the two measured forces.
-
-
-%% Make a line fit
-fit_F = polyfit(Fc, F, 1);
-
-
-
-% The two forces are plotted against each other as well as the linear fit in Figure ref:fig:test_joints_force_sensor_calib_fit.
-
-
-figure;
-hold on;
-plot(Fc, F, '-', 'DisplayName', 'Raw Data');
-plot(Fc([1,end]), Fc([1,end])*fit_F(1) + fit_F(2), '--', 'DisplayName', 'Line Fit');
-hold off;
-xlabel('XFL212R [N]'); ylabel('FC2231 [N]');
-xlim([0,50]); ylim([0,50]);
-legend('location', 'southeast');
-
-
-
-% #+name: fig:test_joints_force_sensor_calib_fit
-% #+caption: Measured two forces and linear fit
-% #+RESULTS:
-% [[file:figs/test_joints_force_sensor_calib_fit.png]]
-
-% The measurement error between the two sensors is shown in Figure ref:fig:test_joints_force_sensor_calib_error.
-% It is below 0.1N for the full measurement range.
-
-
-figure;
-hold on;
-plot(Fc, F - (Fc*fit_F(1) + fit_F(2)), 'k-');
-hold off;
-xlim([0,50]); ylim([-0.12, 0.12]);
-xlabel('Measured Force [N]');
-ylabel('Error [N]')
-
-% Force Sensor Stiffness
-
-% The objective of this measurement is to estimate the stiffness of the force sensor [[file:doc/A700000007147087.pdf][FC2231-0000-0010-L]].
-
-% To do so, a very stiff element is fixed in front of the force sensor as shown in Figure ref:fig:test_joints_meas_force_sensor_stiffness.
-
-% Then, we apply a force on the stiff element through the force sensor.
-% We measure the deflection of the force sensor using an encoder.
-
-% Then, having the force and the deflection, we should be able to estimate the stiffness of the force sensor supposing the stiffness of the other elements are much larger.
-
-% #+name: fig:test_joints_meas_force_sensor_stiffness
-% #+caption: Bench used to measured the stiffness of the force sensor
-% #+attr_latex: :width 0.6\linewidth
-% [[file:figs/test_joints_meas_force_sensor_stiffness.jpg]]
-
-% From the documentation, the deflection of the sensor at the maximum load (50N) is 0.05mm, the stiffness is therefore foreseen to be around $1\,N/\mu m$.
-
-% Let's load the measured force as well as the measured displacement.
-
-%% Load measurement data
-load('force_sensor_stiffness_meas.mat', 't', 'F', 'd')
-
-
-
-% Some pre-processing is applied on the data.
-
-%% Remove offset
-F  = F - mean(F(t > 0.5 & t < 1.0));
-
-%% Select important part of data
-F  = F( t > 4.55 & t < 7.24);
-d  = d( t > 4.55 & t < 7.24); d = d - d(1);
-t  = t( t > 4.55 & t < 7.24);
-
-
-
-% The linear fit is performed.
-
-%% Linear fit
-fit_k = polyfit(F, d, 1);
-
-
-
-% The displacement as a function of the force as well as the linear fit are shown in Figure ref:fig:test_joints_force_sensor_stiffness_fit.
-
-figure;
-hold on;
-plot(F, 1e6*d, '-', 'DisplayName', 'Raw Data');
-plot(F([1,end]), 1e6*(F([1,end])*fit_k(1) + fit_k(2)), '--', 'DisplayName', 'Line Fit');
-hold off;
-xlabel('Force [$N$]'); ylabel('Displacement [$\mu m$]');
-xlim([0,45]); ylim([0,60]);
-legend('location', 'southeast');
-
-
-
-% #+name: fig:test_joints_force_sensor_stiffness_fit
-% #+caption: Displacement as a function of the measured force
-% #+RESULTS:
-% [[file:figs/test_joints_force_sensor_stiffness_fit.png]]
-
-% And we obtain the following stiffness:
-
-%% Force Sensor Stiffness
-sprintf('k = %.2f [N/um]', 1e-6*1/fit_k(1));

matlab/test_joints_5_bending_stiff_meas.m

@@ -1,285 +0,0 @@
-%% Clear Workspace and Close figures
-clear; close all; clc;
-
-%% Intialize Laplace variable
-s = zpk('s');
-
-%% Path for functions, data and scripts
-addpath('./mat/'); % Path for data
-
-%% Colors for the figures
-colors = colororder;
-
-% Analysis of one measurement
-
-% In this section is shown how the data are analysis in order to measured:
-% - the bending stiffness
-% - the bending stroke
-% - the stiffness once the mechanical stops are in contact
-
-
-% The height from the flexible joint's center and the point of application force $h$ is defined below:
-
-h = 25e-3; % [m]
-
-%% Load Data
-load('meas_stiff_flex_1_x.mat', 't', 'F', 'd');
-
-%% Zero the force
-F = F - mean(F(t > 0.1 & t < 0.3));
-
-%% Start measurement at t = 0.2 s
-d = d(t > 0.2);
-F = F(t > 0.2);
-t = t(t > 0.2); t = t - t(1);
-
-
-
-% The obtained time domain measurements are shown in Figure ref:fig:test_joints_meas_bend_time.
-
-
-%% Time Domain plots
-figure;
-tiledlayout(2, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
-
-ax1 = nexttile;
-plot(t, F);
-ylabel('Force [N]'); set(gca, 'XTickLabel',[]);
-
-ax2 = nexttile;
-plot(t, 1e3*d);
-hold off;
-xlabel('Time [s]'); ylabel('Displacement [mm]');
-
-linkaxes([ax1,ax2],'x');
-xlim([0,5]);
-
-
-
-% #+name: fig:test_joints_meas_bend_time
-% #+caption: Typical time domain measurements
-% #+RESULTS:
-% [[file:figs/test_joints_meas_bend_time.png]]
-
-% The displacement as a function of the force is then shown in Figure ref:fig:test_joints_meas_F_d.
-
-
-figure;
-plot(F, 1e3*d);
-xlabel('Force [N]'); ylabel('Displacement [mm]');
-xlim([0,6]); ylim([0,1]);
-
-
-
-% #+name: fig:test_joints_meas_F_d
-% #+caption: Typical measurement of the diplacement as a function of the applied force
-% #+RESULTS:
-% [[file:figs/test_joints_meas_F_d.png]]
-
-% The bending stiffness can be estimated by computing the slope of the curve in Figure ref:fig:test_joints_meas_F_d.
-% The bending stroke and the stiffness when touching the mechanical stop can also be estimated from the same figure.
-
-
-%% Determine the linear region and region when touching the mechanical stop
-% Find when the force sensor touches the flexible joint
-i_l_start = find(F > 0.3, 1, 'first');
-% Reset the measured diplacement at that point
-d = d - d(i_l_start);
-% Find then the maximum force is applied
-[~, i_s_stop] = max(F);
-% Linear region stops ~ when 90% of the stroke is reached
-i_l_stop = find(d > 0.9*d(i_s_stop), 1, 'first');
-% "Stop" region start ~1N before maximum force is applied
-i_s_start = find(F > max(F)-1, 1, 'first');
-
-%% Define variables for the two regions
-F_l = F(i_l_start:i_l_stop);
-d_l = d(i_l_start:i_l_stop);
-
-F_s = F(i_s_start:i_s_stop);
-d_s = d(i_s_start:i_s_stop);
-
-%% Fit the best straight line for the two regions
-fit_l = polyfit(F_l, d_l, 1);
-fit_s = polyfit(F_s, d_s, 1);
-
-%% Reset displacement based on fit
-d = d - fit_l(2);
-fit_s(2) = fit_s(2) - fit_l(2);
-fit_l(2) = 0;
-
-
-
-% The raw data as well as the fit corresponding to the two stiffnesses are shown in Figure ref:fig:test_joints_meas_F_d_lin_fit.
-
-
-figure;
-hold on;
-plot(F(1:i_s_stop), 1e3*d(1:i_s_stop), '.k')
-plot(F_l, 1e3*(F_l*fit_l(1) + fit_l(2)))
-plot(F_s, 1e3*(F_s*fit_s(1) + fit_s(2)))
-hold off;
-xlabel('Force [N]'); ylabel('Displacement [mm]');
-xlim([0,6]);
-
-
-
-% #+name: fig:test_joints_meas_F_d_lin_fit
-% #+caption: Typical measurement of the diplacement as a function of the applied force with estimated linear fits
-% #+RESULTS:
-% [[file:figs/test_joints_meas_F_d_lin_fit.png]]
-
-% Then, the bending stroke is estimated as crossing point between the two fitted lines:
-
-d_max = fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1));
-
-% Bending stiffness and bending stroke of all the flexible joints
-
-% Now, let's estimate the bending stiffness and stroke for all the flexible joints.
-
-
-%% Initialize variables
-kRx = zeros(1,16);
-kSx = zeros(1,16);
-Rmx = zeros(1,16);
-
-for i = 1:16
-    %% Load the data
-    load(['meas_stiff_flex_' num2str(i) '_x.mat'], 't', 'F', 'd');
-
-    %% Automatic Zero of the force
-    F = F - mean(F(t > 0.1 & t < 0.3));
-
-    %% Start measurement at t = 0.2 s
-    d = d(t > 0.2);
-    F = F(t > 0.2);
-    t = t(t > 0.2); t = t - t(1);
-
-    %% Estimate linear region and "stop" region
-    i_l_start = find(F > 0.3, 1, 'first');
-    d = d - d(i_l_start);
-    [~, i_s_stop] = max(F);
-    i_l_stop = find(d > 0.9*d(i_s_stop), 1, 'first');
-    i_s_start = find(F > max(F)-1, 1, 'first');
-
-    F_l = F(i_l_start:i_l_stop);
-    d_l = d(i_l_start:i_l_stop);
-
-    F_s = F(i_s_start:i_s_stop);
-    d_s = d(i_s_start:i_s_stop);
-
-    %% Straight line fit
-    fit_l = polyfit(F_l, d_l, 1);
-    fit_s = polyfit(F_s, d_s, 1);
-
-    %% Reset displacement based on fit
-    d = d - fit_l(2);
-    fit_s(2) = fit_s(2) - fit_l(2);
-    fit_l(2) = 0;
-
-    %% Estimated Stroke
-    d_max = fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1));
-
-    %% Save stiffnesses and stroke
-    kRx(i) = (h)^2/fit_l(1);
-    kSx(i) = (h)^2/fit_s(1);
-    Rmx(i) = atan2(d_max,h);
-
-end
-
-%% Initialize variables
-kRy = zeros(1,16);
-kSy = zeros(1,16);
-Rmy = zeros(1,16);
-
-for i = 1:16
-    %% Load the data
-    load(['meas_stiff_flex_' num2str(i) '_y.mat'], 't', 'F', 'd');
-
-    %% Automatic Zero of the force
-    F = F - mean(F(t > 0.1 & t < 0.3));
-
-    %% Start measurement at t = 0.2 s
-    d = d(t > 0.2);
-    F = F(t > 0.2);
-    t = t(t > 0.2); t = t - t(1);
-
-    %% Estimate linear region and "stop" region
-    i_l_start = find(F > 0.3, 1, 'first');
-    d = d - d(i_l_start);
-    [~, i_s_stop] = max(F);
-    i_l_stop = find(d > 0.9*d(i_s_stop), 1, 'first');
-    i_s_start = find(F > max(F)-1, 1, 'first');
-
-    F_l = F(i_l_start:i_l_stop);
-    d_l = d(i_l_start:i_l_stop);
-
-    F_s = F(i_s_start:i_s_stop);
-    d_s = d(i_s_start:i_s_stop);
-
-    %% Straight line fit
-    fit_l = polyfit(F_l, d_l, 1);
-    fit_s = polyfit(F_s, d_s, 1);
-
-    %% Reset displacement based on fit
-    d = d - fit_l(2);
-    fit_s(2) = fit_s(2) - fit_l(2);
-    fit_l(2) = 0;
-
-    %% Estimated Stroke
-    d_max = fit_l(1)*fit_s(2)/(fit_l(1) - fit_s(1));
-
-    %% Save stiffnesses and stroke
-    kRy(i) = (h)^2/fit_l(1);
-    kSy(i) = (h)^2/fit_s(1);
-    Rmy(i) = atan2(d_max,h);
-
-end
-
-% Analysis
-% The dispersion of the measured bending stiffness is shown in Figure ref:fig:test_joints_bend_stiff_hist and of the bending stroke in Figure ref:fig:test_joints_bend_stroke_hist.
-
-
-figure;
-hold on;
-histogram(kRx, 'DisplayName', '$k_{R_x}$')
-histogram(kRy, 'DisplayName', '$k_{R_y}$')
-hold off;
-xlabel('Bending Stiffness [Nm/rad]')
-legend();
-
-
-
-% #+name: fig:test_joints_bend_stiff_hist
-% #+caption: Histogram of the measured bending stiffness
-% #+RESULTS:
-% [[file:figs/test_joints_bend_stiff_hist.png]]
-
-
-figure;
-hold on;
-histogram(1e3*Rmx, 'DisplayName', '$k_{R_x}$')
-histogram(1e3*Rmy, 'DisplayName', '$k_{R_y}$')
-hold off;
-xlabel('Bending Stroke [mrad]')
-legend();
-
-
-
-% #+name: fig:test_joints_bend_stroke_hist
-% #+caption: Histogram of the measured bending stroke
-% #+RESULTS:
-% [[file:figs/test_joints_bend_stroke_hist.png]]
-
-% The relation between the measured beam thickness and the measured bending stiffness is shown in Figure ref:fig:test_joints_thickness_stiffness.
-
-
-load('flex_meas_dim.mat', 'meas_flex');
-figure;
-hold on;
-plot((meas_flex(:,1)+meas_flex(:,2))/2, kRx, 'o', 'DisplayName', '$x$')
-plot((meas_flex(:,3)+meas_flex(:,4))/2, kRy, 'o', 'DisplayName', '$y$')
-hold off;
-xlabel('Flexible Beam Thickness [$\mu m$]');
-ylabel('Bending Stiffness [Nm/rad]');
-legend('location', 'southeast');

test-bench-flexible-joints.org

@@ -623,7 +623,7 @@ With an estimated shear stiffness $k_s = 13\,N/\mu m$ from the finite element mo
 
 #+begin_src matlab
 %% Estimated error due to shear
-epsilon_s = 100*abs(1-1/(1 + kb/(ks*h^2))) % Error in %
+epsilon_s = 100*abs(1-1/(1 + kb/(ks*h^2))); % Error in %
 #+end_src
 
 **** Effect of load cell limited stiffness