diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..e826a22 --- /dev/null +++ b/.gitignore @@ -0,0 +1,4 @@ +*.asv +Code/ICANSAT/V0/Result +Code/ICANSAT/V1/Result +Code/MAD/Result diff --git a/Code/AllClean.m b/Code/AllClean.m new file mode 100644 index 0000000..9c65a1b --- /dev/null +++ b/Code/AllClean.m @@ -0,0 +1,17 @@ +mfile_name = mfilename('fullpath'); +[pathstr,name,ext] = fileparts(mfile_name); +cd(pathstr); + +cd ICANSAT/V0/Result +delete *.txt +fprintf("ICANSAT V0 result files deleted.\n"); + +cd ../../V1/Result +delete *.txt +fprintf("ICANSAT V1 result files deleted.\n"); + +cd ../../../MAD/Result +delete *.txt +fprintf("MAD result files deleted.\n"); + +cd(pathstr); \ No newline at end of file diff --git a/Code/AllRun.m b/Code/AllRun.m new file mode 100644 index 0000000..1dc2882 --- /dev/null +++ b/Code/AllRun.m @@ -0,0 +1,14 @@ +mfile_name = mfilename('fullpath'); +[pathstr,name,ext] = fileparts(mfile_name); +cd(pathstr); + +cd ICANSAT/V0 +ICANSAT + +cd ../V1 +ICANSAT_V1 + +cd ../../MAD +MAD_Validation + +cd(pathstr); \ No newline at end of file diff --git a/Code/CalcAtmosQuantities.m b/Code/CalcAtmosQuantities.m new file mode 100644 index 0000000..8837c08 --- /dev/null +++ b/Code/CalcAtmosQuantities.m @@ -0,0 +1,70 @@ +function [t,p,rho] = CalcAtmosQuantities(h) + hg = h; % geometric alt in meters + R = 287; % Gas Constant + re = 6.356766e6; % radius of earth + g0 = 9.8066; % gravity + h = re/(re+hg)*hg; + + % Altitude in meters + h0 = 0.0; + h1 = 11000.0; + h2 = 25000.0; + h3 = 47000; + h4 = 53000; + h5 = 79000; + h6 = 90000; + h7 = 105000; + + % Temperature in Kelvin. + t0 = 288.16; + t1 = 216.66; + t2 = 216.66; + t3 = 282.66; + t4 = 282.66; + t5 = 165.66; + t6 = 165.66; + t7 = 225.66; + + %%equations + a01 = (t1-t0)/(h1-h0); % slope of T/h in K/m for 0 to 11 km + a23 = (t3-t2)/(h3-h2); % slope of T/h in K/m for 25 to 47 km + a45 = (t5-t4)/(h5-h4); % slope of T/h in K/m for 53 to 79 km + a67 = (t7-t6)/(h7-h6); % slope of T/h in K/m for 90 to 105 km + + % Pressure in N/m^2 + p0 = 101325.0; + p1 = p0*(t1/t0)^(-g0/(a01*R)); + p2 = p1*exp(-(g0/(R*t1))*(h2-h1)); + p3 = p2*(t3/t2)^(-g0/(a23*R)); + p4 = p3*exp(-(g0/(R*t3))*(h4-h3)); + p5 = p4*(t5/t4)^(-g0/(a45*R)); + p6 = p5*exp(-(g0/(R*t5))*(h6-h5)); + p7 = p6*(t7/t6)^(-g0/(a67*R)); + + if h < h1 + t = t0+a01*h; + p = p0*(t/t0)^(-g0/(a01*R)); + elseif h < h2 + t = t1; + p= p1*exp(-(g0/(R*t))*(h-h1)); + elseif h < h3 + t = t2+a23*(h-h2); + p = p2*(t/t2)^(-g0/(a23*R)); + elseif h < h4 + t = t3; + p = p3*exp(-(g0/(R*t))*(h-h3)); + elseif h < h5 + t = t4+a45*(h-h4); + p = p4*(t/t4)^(-g0/(a45*R)); + elseif h < h6 + t = t5; + p = p5*exp(-(g0/(R*t))*(h-h5)); + elseif h < h7 + t = t6+a67*(h-h6); + p = p6*(t/t6)^(-g0/(a67*R)); + else + t= t7; + p = p7*exp(-(g0/(R*t))*(h-h7)); + end + rho=p/(R*t); % Density in kg/m^3 +end \ No newline at end of file diff --git a/Code/DragInterpolations_ICANSAT.m b/Code/DragInterpolations_ICANSAT.m index 46dea55..b7ff755 100644 --- a/Code/DragInterpolations_ICANSAT.m +++ b/Code/DragInterpolations_ICANSAT.m @@ -1,7 +1,5 @@ function [coeffdrag] = DragInterpolations_ICANSAT(velocity, altitude, time, tb) -% global tb - drag_data = dlmread("ICANSAT/ICANSAT_DragData.txt", "\t",1,0); % Mach No @@ -13,7 +11,13 @@ cd_full = drag_data(:,3); % No base drag cd_zero = drag_data(:,2); -[temp, sound_speed, P, rho] = atmosisa(altitude); +R = 8.3145; +gamma = 1.4; +M_air = 0.0289645; + +[T,p,rho] = CalcAtmosQuantities(altitude); +sound_speed = sqrt(gamma*R*T/M_air); + mach_speed = velocity/sound_speed; if (time > tb) diff --git a/Code/DragInterpolations_MAD.m b/Code/DragInterpolations_MAD.m index 25221c0..0533f0f 100644 --- a/Code/DragInterpolations_MAD.m +++ b/Code/DragInterpolations_MAD.m @@ -27,7 +27,13 @@ cdzero_10000 = drag_data(:,3); % Interpolate for Cd values at two stated altitudes, then interpolate for % Cd value at a specific altitude -[temp, sound_speed, P, rho] = atmosisa(altitude); + +R = 8.3145; +gamma = 1.4; +M_air = 0.0289645; + +[T,p,rho] = CalcAtmosQuantities(altitude); +sound_speed = sqrt(gamma*R*T/M_air); mach_speed = velocity/sound_speed; % If propellant is still available diff --git a/Code/DragInterpolations_Spark.m b/Code/DragInterpolations_Spark.m index b4dbbe8..24218f8 100644 --- a/Code/DragInterpolations_Spark.m +++ b/Code/DragInterpolations_Spark.m @@ -2,7 +2,13 @@ function [coeffdrag] = DragInterpolations_Spark(velocity, altitude) spark_drag_data = dlmread("Spark/Spark_DragData.txt", "\t",1,0); -[temp, sound_speed, P, rho] = atmosisa(altitude); +R = 8.3145; +gamma = 1.4; +M_air = 0.0289645; + +[T,p,rho] = CalcAtmosQuantities(altitude); +sound_speed = sqrt(gamma*R*T/M_air); + mach_speed = velocity/sound_speed; mach = spark_drag_data(:,1); diff --git a/Code/ICANSAT/V1/ICANSAT_V1.m b/Code/ICANSAT/V1/ICANSAT_V1.m index bf97ba0..443b243 100644 --- a/Code/ICANSAT/V1/ICANSAT_V1.m +++ b/Code/ICANSAT/V1/ICANSAT_V1.m @@ -5,15 +5,15 @@ addpath('../..') % Df, Lr, Ls, ms, mp, tb icansat = SingleStageRocket(0.1567, 2.391, 4, 14, 1.160534, 4.0); -% RunTrajectorySim(icansat, 88); -% RunTrajectorySim(icansat, 86); -% RunTrajectorySim(icansat, 84); -% RunTrajectorySim(icansat, 82); -% RunTrajectorySim(icansat, 80); -% RunTrajectorySim(icansat, 78); -% RunTrajectorySim(icansat, 76); -% RunTrajectorySim(icansat, 74); -% RunTrajectorySim(icansat, 72); +RunTrajectorySim(icansat, 88); +RunTrajectorySim(icansat, 86); +RunTrajectorySim(icansat, 84); +RunTrajectorySim(icansat, 82); +RunTrajectorySim(icansat, 80); +RunTrajectorySim(icansat, 78); +RunTrajectorySim(icansat, 76); +RunTrajectorySim(icansat, 74); +RunTrajectorySim(icansat, 72); RunTrajectorySim(icansat, 70); function RunTrajectorySim(rocket, launch_angle) diff --git a/Code/SingleStageRocket.m b/Code/SingleStageRocket.m index 1e2d54b..eee3321 100644 --- a/Code/SingleStageRocket.m +++ b/Code/SingleStageRocket.m @@ -165,81 +165,9 @@ classdef SingleStageRocket < handle obj.CalculateMassFlow(length(resampled_time)); end - function [t,p,rho] = CalcAtmosQuantities(obj, h) - hg = h; % geometric alt in meters - R = 287; % Gas Constant - re = 6.356766e6; % radius of earth - g0 = 9.8066; % gravity - h = re/(re+hg)*hg; - - % Altitude in meters - h0 = 0.0; - h1 = 11000.0; - h2 = 25000.0; - h3 = 47000; - h4 = 53000; - h5 = 79000; - h6 = 90000; - h7 = 105000; - - % Temperature in Kelvin. - t0 = 288.16; - t1 = 216.66; - t2 = 216.66; - t3 = 282.66; - t4 = 282.66; - t5 = 165.66; - t6 = 165.66; - t7 = 225.66; - - %%equations - - a01 = (t1-t0)/(h1-h0); % slope of T/h in K/m for 0 to 11 km - a23 = (t3-t2)/(h3-h2); % slope of T/h in K/m for 25 to 47 km - a45 = (t5-t4)/(h5-h4); % slope of T/h in K/m for 53 to 79 km - a67 = (t7-t6)/(h7-h6); % slope of T/h in K/m for 90 to 105 km - - % Pressure in N/m^2 - p0 = 101325.0; - p1 = p0*(t1/t0)^(-g0/(a01*R)); - p2 = p1*exp(-(g0/(R*t1))*(h2-h1)); - p3 = p2*(t3/t2)^(-g0/(a23*R)); - p4 = p3*exp(-(g0/(R*t3))*(h4-h3)); - p5 = p4*(t5/t4)^(-g0/(a45*R)); - p6 = p5*exp(-(g0/(R*t5))*(h6-h5)); - p7 = p6*(t7/t6)^(-g0/(a67*R)); - - if h < h1 - t = t0+a01*h; - p = p0*(t/t0)^(-g0/(a01*R)); - elseif h < h2 - t = t1; - p= p1*exp(-(g0/(R*t))*(h-h1)); - elseif h < h3 - t = t2+a23*(h-h2); - p = p2*(t/t2)^(-g0/(a23*R)); - elseif h < h4 - t = t3; - p = p3*exp(-(g0/(R*t))*(h-h3)); - elseif h < h5 - t = t4+a45*(h-h4); - p = p4*(t/t4)^(-g0/(a45*R)); - elseif h < h6 - t = t5; - p = p5*exp(-(g0/(R*t))*(h-h5)); - elseif h < h7 - t = t6+a67*(h-h6); - p = p6*(t/t6)^(-g0/(a67*R)); - else - t= t7; - p = p7*exp(-(g0/(R*t))*(h-h7)); - end - rho=p/(R*t); % Density in kg/m^3 - end - % Thrust correction based on ASTOS manual function T = ThrustCorrection(obj, h, thrust) - [t,p,rho] = obj.CalcAtmosQuantities(h); + [t,p,rho] = CalcAtmosQuantities(h); T = thrust*1.095 - 0.014 * p; end @@ -277,7 +205,7 @@ classdef SingleStageRocket < handle % Density calculation function [T, rho] = CalculateRhoTemp(obj, h) - [T, p, rho] = obj.CalcAtmosQuantities(h); + [T,p,rho] = CalcAtmosQuantities(h); return; end @@ -427,15 +355,18 @@ classdef SingleStageRocket < handle % Write data file function WriteDataFile(obj) if obj.mdot_constant == true - result_dat_file = obj.ID + "_Result_" + obj.launch_angle + "deg.txt"; + result_dat_file = 'Result/' + obj.ID + "_Result_" + obj.launch_angle + "deg.txt"; else text = "VariableMassFlow"; - result_dat_file = obj.ID + "_Result_" + obj.launch_angle + "_deg_" + text + ".txt"; + result_dat_file = 'Result/' + obj.ID + "_Result_" + obj.launch_angle + "_deg_" + text + ".txt"; end result_matrix = [obj.t obj.x/1000 obj.y/1000 obj.M obj.thrust obj.DvDt obj.m obj.mdot obj.isp]; result_data = reshape(result_matrix, [], 9); dlmwrite(result_dat_file, 'Time(s),x(km),y(km),Mach,Thrust(N),Acceleration(m/s2),Mass(kg),mdot(kg/s),Isp (s)','delimiter','') dlmwrite(result_dat_file, result_data, '-append', 'delimiter', ','); + notify_text = 'Results written for ' + obj.ID + ' at ' + obj.launch_angle + " degrees.\n "; + fprintf(notify_text); + fprintf("-----------------------------------------------------------------------------------\n"); end end end \ No newline at end of file