...
| Drawio | ||||||||||||||||||||||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
|
Force and Stall Analysis [WIP]
For analyzing stall, simulation is performed at maximum aileron deflection and at expected cruising speed. Initially, a steady state simulation was run but the plots of residuals, lift & drag struggled to converge and had a lot of noiseoscillations. A Hence, a transient simulation is also performed, which led to a less oscillating values. Pictures are shown below:
...
Something to note is that the result shown in Fig. 5 was an assembly geometry where the aileron and wing were separate bodies. There was some gap between the aileron and the wing. To simplify the meshing/discretization process, a single body CAD was created (used for simulation in fig. 6 and 7)
...
There are no large vorticities above the aileron, hence no significant sign of stall. conducted..
Mesh details:
Face sizing: 0.002 m applied only on the wing geometry
Surface mesh: enable curvature only, don’t need proximity feature
Boundary layers: last_ratio method with 5 layers. First layer height is 0.0001
Volume mesh: poly-hexcore method, everything else default. Poly-hexcore to reduce total cell count.
The minimum ortho quality of the mesh generated is 0.05 which isn’t the best, will mostly affect cd accuracy. But the average quality is pretty high so it’s acceptable.
Aileron deflection | Pressure plot | Velocity pathline | Normal force (surface integrated) | Drag coefficient | Lift Coefficient |
|---|---|---|---|---|---|
25 degrees DOWN |  |  | Top_face = -5.86 N Bot_face = 3.29 N Total: 9.15 N | 0.01519 | 0.4368 |
25 degrees UP |  |
| Top = 2.175 N Bot = -5.04 Total: 7.22 N | 0.008457 | -0.01134 |
Using the force results obtained, torque about the hinge line is calculated. It is assumed that the force acts at the mid-point of the chord-line.
For down deflection, the torque is 0.0268 kg*cm. For up deflection, the torque is 0.0212 kg*cm. These values are pretty low compared to what online calculators give. So, an online calculator is used to double check…
Link: RC Airplane Calculator (radiocontrolinfo.com)
...
Some servo candidates:
Servo Model | Torque (kg*cm) | Weight (g) | Recommended RC Plane application | Cost $$$ |
|---|---|---|---|---|
HS-311 | 3.0 ~ 3.5 | 43.0 | 1.8 kg | |
HS 645 MG | 7.7 ~ 9.6 | 55.2 | 6.8 kg | |
DS 843MG | 4 - 4.8 | 11 g | $16.39 |
Roll rate and time calculation script:
| Code Block |
|---|
%SI (m, kg, m/s, N)
%General
MASS_TOTAL = 5.137;
rho = 1.225;
V = 25;
Ixx = 0.06; %moment of inertia about the roll axis (kg*m^2)
Sv = 0.0423;
Sh = 0.0666;
%Wing geometry
b = 1.5;
fuselage_width = 0.132;
C = 0.23;
Sw = 0.345;
taper = 1;
MAC = 0.23;
NACA = 4412;
Cr = C; %
%Control surface effectiveness param data points
CA_C_ratio = [0.004918032786885199,0.03442622950819667,0.06967213114754095, 0.09918032786885245, 0.12213114754098359, 0.19344262295081968, 0.28196721311475414, 0.4024590163934427, 0.5016393442622952, 0.5434426229508198, 0.6778688524590166
];
TAU = [-0.002272727272726982, 0.10000000000000009, 0.20000000000000007, 0.26363636363636367, 0.3045454545454547, 0.4, 0.5, 0.6068181818181819, 0.6886363636363637, 0.7181818181818183, 0.8136363636363637
];
%Aileron inputs
span_ratio = 0.4;
chord_ratio = 0.25;
bi = 1.5*0.5-span_ratio*1.5*0.5;
bo = 1.5*0.5;
Ca = chord_ratio*C;
tau = interpolate_tau(CA_C_ratio, TAU, chord_ratio); %Aileron efficiency paramter
delta_max = [25, 25]; % [UP, DOWN]
CLw_slope = 6.03;%cl per rad
C_DR = 0.9; %avg drag from wing, horizontal tail, vertical tail
yD = (b-fuselage_width)*0.25;
phi_req = 90;
t_req = 1.7;
%EQUATIONS
eq_Cl_deltaA = @(yi, yo) ((0.5*yo^2 + (2/3)*yo^3*(taper-1)/b)-(0.5*yi^2 + ...
(2/3)*yi^3*(taper-1)/b))*(2*CLw_slope*tau*Cr)/(Sw*b); %Cl derivative
eq_Cl = @(x) x*delta_max(1,1)*pi/180; %Cl - rolling moment coefficient
eq_LA = @(Cl) 0.5*rho*V^2*Sw*Cl*b; %Lift contribution from aileron
eq_Pss = @(LA) sqrt((2*LA)/(rho*(Sw+Sv+Sh)*C_DR*yD^3));
eq_phi1 = @(Pss) (Ixx)*(log(Pss^2))/(rho*yD^3 * (Sw+Sv+Sh)*C_DR);
eq_roll_rate_derivative = @(phi1, Pss) 0.5*Pss^2/phi1;
%Assume that phi1 is going to be larger than phi_req for GA
eq_t2 = @(phi, P_dot) sqrt(2*phi/P_dot);
%Testing the equation
Cl_deltaA = eq_Cl_deltaA(bi, bo);
Cl = eq_Cl(Cl_deltaA);
LA = eq_LA(Cl);
Pss = eq_Pss(LA);
phi1 = eq_phi1(Pss);
P_dot = eq_roll_rate_derivative(phi1, Pss);
tss = sqrt(2*phi1/P_dot);
if (phi1 > phi_req)
t2 = eq_t2(phi_req*pi/180, P_dot);
else
delta_t = (phi_req - phi1)/Pss;
t2 = tss + delta_t;
end
err = (t2-t_req)/t_req *100;
%Plotting some points
time_v_phi(0.01, 5, tss,P_dot,Pss,phi1);
function graph_t_v_phi = time_v_phi(time_step,max_time,tss,P_dot,Pss,phi1)
%Time vs phi
t_nonlinear = 0:time_step:tss;
phi = 0.5.*t_nonlinear.^2 .* P_dot;
plot(t_nonlinear, phi);
title("Time vs. Bank Position");
xlabel('Time (s)');
ylabel('Bank angle (deg)');
hold on
%Linear region
t_linear = tss:time_step:max_time;
b = phi1 - Pss*tss;
phi_lin = b + Pss.*t_linear;
plot(t_linear, phi_lin);
hold off
end
function tau = interpolate_tau(x_points, y_points, x)
data_x = x_points(:);
data_y = y_points(:);
if x < min(data_x) || x > max(data_x)
error('Outside the range!');
end
%Find indices of 2 points, 1st point less than x, 2nd point greater
indx1 = find(data_x <= x, 1, 'last');
indx2 = find(data_x >= x, 1, 'first');
%If x matches a data point, return corresponding y
if data_x(indx1) == x_points
tau = data_y(indx1);
end
%Linear interpolation
x1 = data_x(indx1);
x2 = data_x(indx2);
y1 = data_y(indx1);
y2 = data_y(indx2);
tau = y1 + (y2 - y1) * (x - x1) / (x2 - x1);
end
|