clear
syms wb x real

% Parameters
lb = 5.5e-2; % beam length
lf = 70e-3; % flexure length
wf = 1e-3; % flexure width
theta = deg2rad(60); % angle between padding beams
r_egg = 21.5e-3; % egg radius along horizontal plane
dx_beam = 1e-3; % max deflection of one padding beam
dx_beam_eff = dx_beam * (1 + cos(theta / 2));
xrange = 7.088e-2; % distance between gripper peak and base padding
egg_contact = r_egg * (1 + 1 / sin(theta / 2)); % distance from base padding to egg contact
x_contact = xrange - egg_contact;
x_max = x_contact + dx_beam_eff;

% Define material and constants
t = 0.00635; % thickness of ABS in meters
E = 1.9e9; % Young's modulus for ABS
gamma = 0.85147;
ktheta = 2.65;

% Derived parameters
Ib = t * wb^3 / 12; % Moment of inertia for beam
If = t * wf^3 / 12; % Moment of inertia for flexure
Lcp = (1 - gamma) * lf / 2; % characteristic length of flexure
B = asin(x / lf); % flexure angle as a function of x

% Stiffness constants
kbeam = 192 * E * Ib / lb^3; % stiffness of base padding beam
kbeam2 = kbeam * sin(theta / 2) / 2; % stiffness of single tilted padding beam
kgripper = 1 / (1 / kbeam + 1 / (2 * kbeam2)); % total gripper stiffness

% Flexure stiffness
kt = gamma * ktheta * E * If / lf;

% Energy stored in the flexure
Ustored = 2 * kt * B^2;
Fstored = diff(Ustored, x);
Kstored = diff(Fstored, x);

% Energy in the gripper (exerted force)
x_eff = @(x_in) max(0, x_in - x_contact); % effective displacement
K_exert = @(wb_in) double(subs(kgripper, wb, wb_in)); % convert kgripper to numeric based on wb
F_exert = @(x_in, wb_in) K_exert(wb_in) * x_eff(x_in); % force exerted on x based on wb

% Setting up the grid for x and wb
n = 50;
x_vals = linspace(x_contact, x_max, n);
wb_vals = linspace(0, 5e-3, n);
[WB, X] = meshgrid( wb_vals,x_vals);

% Compute Z values using F_exert over the grid
Z = arrayfun(@(x, wb) F_exert(x, wb), X, WB);

% 3D surface plot
figure;
surf(X,WB, Z);
xlabel('x (m)');
ylabel('wb (m)');
zlabel('F_{exert}(x, wb) (N)');
title('3D Plot of F_{exert}(x, wb)');
shading interp;
colormap winter;
colorbar;
view(3);