===== Example: Shape optimization of cross-hinge =====

This example provides a shape optimization of a single cross flexure by using matlab's ''fmincon'' function. For this optimization a crosshinge is considered with a mass and inertia connected to node 3. Furthermore, the optimization objective is to maximize the first parasitic frequency of the system (note: rotation of the cross-hinge is prescribed to de desired rotation angle. Therefore, the first parasitic frequency of the system is given by the first eigenfrequency). A rotation of 0.5rad, a maximum stress of 300MPa and and a maximum actuation force of 2.5Nm is considered.

The cost function of the optimization is defined by the inverse of the first eigenfrequency. Additionally, a limit on the maximum stress and actuation force is implemented by applying a penalty on the value of the cost function. The matlab function which calls spacarlight and computes the cost function is defined by:

<code matlab run_sim.m>
function [cost, out] = run_sim(P, Silent)% EXAMPLE SCRIPT FOR RUNNING SPACAR LIGHT
if nargin == 1
    Silent = true; %run in silent mode when called by optimizer
end 
    
% DESIGN PARAMETERS
H = P(1);
B = P(2);
w = P(3);
t = P(4);

% Optimization parameters
stressmax = 300e6; %Max stress [MPa]
Mmax = 2.5; %Max actuation moment [Nm]
Stroke = 0.5; %Max stroke [Rad]


%% NODE POSITIONS
%           x y z
nodes = [   0 0 0;    %node 1  
            0 H 0;    %node 2  
            B H 0;    %node 3
            B 0 0];   %node 4

        
%% ELEMENT CONNECTIVITY
%               p   q
elements = [    1   3;  %element 1
                2   3;  %element 2
                2   4]; %element 3
        
%% NODE PROPERTIES  
%node 1
nprops(1).fix               = true;         %Fix node 1

%node 3
nprops(3).rot_z           = Stroke;
nprops(3).mass            = 1;
nprops(3).mominertia      = [1e-3 0 0 1e-3 0 1e-3];

%node 4
nprops(4).fix               = true;         %Fix node 4


%% ELEMENT PROPERTIES
%Property set 1
eprops(1).elems    = [1 3];            %Add this set of properties to elements 1 and 3
eprops(1).emod     = 210e9;            %E-modulus [Pa]
eprops(1).smod     = 70e9;             %G-modulus [Pa]
eprops(1).dens     = 7800;             %Density [kg/m^3]
eprops(1).cshape   = 'rect';           %Rectangular cross-section
eprops(1).dim      = [w t];   %Width: 50 mm, thickness: 0.2 mm
eprops(1).orien    = [0 0 1];          %Orientation of the cross-section as a vector pointing along "width-direction"
eprops(1).nbeams   = 2;                %Number of beams used to model this element 
eprops(1).flex     = 1:6;              %Consider all deformation modes (1 to 6) as flexible
eprops(1).color    = 'grey';           %Color
eprops(1).opacity  = 0.7;              %Opacity
eprops(1).cw       = true;             %Constrain warping at nodes (typical condition for flexure clamped at both ends). For a "free" end, warping is typically not constrained

%Property set 2
eprops(2).elems    = 2;                %Add this set of properties to element 2
eprops(2).dens     = 3000;             %Density [kg/m^3]
eprops(2).cshape   = 'rect';           %Rectangular cross-section
eprops(2).dim      = [w 10e-3];    %Width: 50 mm, thickness: 10 mm
eprops(2).orien    = [0 0 1];          %Orientation of the cross-section as a vector pointing along "width-direction"
eprops(2).nbeams   = 1;                %1 beam for simulating this element (as it is rigid an no more elements are required)
eprops(2).color    = 'darkblue';


%% OPTIONAL ARGUMENTS
opt.filename    = 'crosshinge';     %Filename
opt.calccompl   = false;            %Disable calculation of compliance matrices (can reduce computation time for large simulations)
opt.silent      = Silent;            %Run in silent mode

%% CALL SPACAR_LIGHT AND COMPUTE COST FUNCTION
try
    out = spacarlight(nodes, elements, nprops, eprops, opt);
    cost = 1/min(out.freq(1,:)); %cost function (inverse of first frequency )
    stress = max(out.stressmax); %max stress
    Mact = max(out.node(3).Mreac(3,:));%max actuation moment
    
    %add penalty if limits are exceeded
    if stress > stressmax
        cost = cost * (1 + (stress-stressmax)/stressmax)^5;
    end
    if Mact > Mmax
        cost = cost * (1 + (Mact-Mmax)/Mmax)^5;
    end
catch %if sumulation fails, cost is inf
    cost = inf;
end
</code>

An example to conduct the optimization is provided below, taking into account additional boundary conditions on the design parameters of the cross-hinge. See matlab's documentation on ''fmincon'' for more information.

<code matlab optimization.m>
%filename
simulation = @run_sim;

%starting point
x0 = [0.1;  0.1;  50e-3; 0.2e-3];

%constraints
b1 = [0.05; 0.05; 25e-3; 0.1e-3];   % Lower bound (b1 < x)
b2 = [0.2;  0.2;  75e-3; 1e-3  ];       % Upper bound (x < b2)

A = [-eye(4); eye(4)];
b = [-b1;b2];

%optimization
x_opt = fmincon(simulation, x0, A, b);

%show results
[cost, out] = run_sim(x_opt,false);
</code>
//Note: prescribing input rotation is only supported for SPACAR light version 1.27 or higher.//