From 2defd94c5b208763d4f4cdedd6a07c1b49911f7c Mon Sep 17 00:00:00 2001 From: Aditeya Shukla Date: Wed, 7 Feb 2024 09:37:22 -0500 Subject: [PATCH] new PRD code (#13) --- docs/src/examples/payload_range.md | 136 ++++++++++------------ example/example_PRD.jl | 181 ++++++++++------------------- src/IO/outputs.jl | 73 ++++++++++++ 3 files changed, 194 insertions(+), 196 deletions(-) diff --git a/docs/src/examples/payload_range.md b/docs/src/examples/payload_range.md index fc8bb3ae..fdb721be 100644 --- a/docs/src/examples/payload_range.md +++ b/docs/src/examples/payload_range.md @@ -4,9 +4,9 @@ ## Choosing a design mission -To plot a payload-range diagram with a fleet of missions you must first load a model with >=1 non design mission +To plot a payload-range diagram with a fleet of missions you must first load any aircraft model -Start by choosing a design mission. Your design mission should be what you want the second corner point in your Payload Range plot to be. Once you have a chosen a specific design range and payload weight (For eg: 3500 nmi and 195 pax) you can add it to the input toml file: `default_input.toml` +Start by choosing a design mission. Your design mission should be what you want the second corner point in your Payload Range plot to be. Once you have a chosen a specific design range and payload weight (For eg: 3500 nmi and 195 pax) you can add it to the input toml file for eg: `default_input.toml` ```toml [Mission] @@ -19,7 +19,9 @@ Start by choosing a design mission. Your design mission should be what you want ``` ## Julia script for Payload Range Diagram + Start the script importing `TASOPT.jl`, `PyPlot` and `index.inc` and then loading the default `aircraft` model. + ```julia # Import modules using PyPlot @@ -27,7 +29,7 @@ using TASOPT # you can optionally define # const tas = TASOPT # to use as a shorthand -include("../src/misc/index.inc") +include(joinpath(TASOPT.__TASOPTroot__, "./src/misc/index.inc")) # import indices for calling parameters # Load aircraft using default module @@ -35,80 +37,66 @@ ac = TASOPT.read_aircraft_model(joinpath(TASOPT.__TASOPTroot__, "../example/PRD_ time_wsize = @elapsed size_aircraft!(ac) ``` -Initialize some variables for mission range and payloads +One way is to call the `PayloadRange` function: ```julia -nmis = size(ac.parm)[2] # Get number of missions -MAX_PAYLOAD = 230 * 215 * 4.448222 # Starting Payload, #Passengers * 215 lb * lb_to_N -# Ranges: -Range_arr_nmi = LinRange(100, 1000, nmis-1) # Starting range array (in nmi) -Range_arr = Range_arr_nmi .* 1852.0 # Convert to SI units -Max_range = Range_arr[nmis-1] -curr_range = Range_arr[1] -prev_range = 0.1 +TASOPT.PayloadRange(ac) ``` -## Initialize variables for mission iterations +If you want a more customizable diagram, first initialize some variables for mission range and payloads ```julia -NPSS = Base.Process -NPSS_PT = true -itermax = 15 -initeng = 0 -Ldebug = false -saveOD = false -SEC_B = false -x_range = [] -y_payload = [] -y_Wemtpy = [] +# Copy aircraft structure as we will be changing range and payloads +ac = deepcopy(ac_og) +# Make an array of ranges to plot +RangeArray = ac.parm[imRange] * LinRange(0.1,1.2,Rpts) +# Store constant values to compare later +Wmax = ac.parg[igWMTO] +Fuelmax = ac.parg[igWfmax] +Wempty = ac.parg[igWMTO] - ac.parg[igWfuel] - ac.parg[igWpay] +# Arrays for plotting +RangesToPlot = [] +PayloadToPlot = [] +maxPay = ac.parm[imWpay ] ``` ## Main iteration loop ```julia -while MAX_PAYLOAD >= 0 - println("-----------------------") - println("Starting new set of missions, \n Ranges: ", Range_arr, "\n Payload: ", MAX_PAYLOAD./(215*4.448222)) - try - # Set non design mission ranges and payloads - ac.parm[imRange, 2:5] .= Range_arr - ac.parm[imWpay, 2:5] .= MAX_PAYLOAD - # Analyze each mission - for mi in 2:size(ac.parm)[2] - prev_range = curr_range - curr_range = ac.parm[imRange, mi] - # Call TASOPT woper function - @views TASOPT.woper(ac.pari,ac.parg,ac.parm[:,mi:mi],ac.para[:,:,mi:mi],ac.pare[:,:,mi:mi], ac.para[:,:,1:1],ac.pare[:,:,1:1], itermax,initeng, NPSS_PT, NPSS) - # Check if mission fuel and MTO weight is greater than design fuel and MTO weight - if (ac.parm[imWfuel,mi] - ac.parg[igWfuel] > 1) || (ac.parm[imWTO,mi] - ac.parg[igWMTO] > 1 ) - printstyled("Mission designed beyond capacity!", "\n"; color=:red) - println([ac.parg[igWfuel], ac.parm[imWfuel,mi], ac.parg[igWMTO], ac.parm[imWTO,mi]]) - throw(UndefVarError([ac.parg[igWfuel], ac.parm[imWfuel,mi], ac.parg[igWMTO], ac.parm[imWTO,mi]])) - end - - # Add to dataframes - append!(x_range, curr_range) - append!(y_payload, MAX_PAYLOAD) - append!(y_WTO, ac.parm[imWTO,mi]) - - # If in section B of Payload Range plot, decrease payload by 1 passenger - if SEC_B - println("decreasing payload since Section B") - MAX_PAYLOAD = max(0, MAX_PAYLOAD-(1*215*4.448222)) - end - end - println("PRD converged for all: Increasing Ranges") - prev_range = curr_range - Max_range = prev_range+ (500* 1852.0) #Increase range array by 500 nmi +for Range = RangeArray + # Make an array of payloads to plot + Payloads = (maxPay) * LinRange(1, 0.1, Ppts) + ac.parm[imRange] = Range + for mWpay = Payloads + println("Checking for Range (nmi): ",Range/1852.0, " and Pax = ", mWpay/(215*4.44822)) + ac.parm[imWpay ] = mWpay + # Try woper after setting new range and payload + try + @views TASOPT.woper(ac.pari,ac.parg,ac.parm[:,1:1],ac.para[:,:,1:1],ac.pare[:,:,1:1], ac.para[:,:,1:1],ac.pare[:,:,1:1], itermax,0.0) + # woper success: store maxPay, break loop + WTO = Wempty + mWpay + ac.parm[imWfuel] + mWfuel = ac.parm[imWfuel] + + # Compare with previously stored constants + if WTO > Wmax || mWfuel > Fuelmax || WTO < 0.0 || mWfuel < 0.0 + WTO = 0.0 + mWfuel = 0.0 + println("Max out error!") + else + maxPay = mWpay + println("Converged - moving to next range...") + break + end catch - println("PRD failed: Decreasing Payload and range") - prev_range = curr_range - Max_range = prev_range+ ((Max_range-prev_range)*0.75) - MAX_PAYLOAD = MAX_PAYLOAD-pax2N(5) - SEC_B = true + println("Not Converged - moving to lower payload...") end - Range_arr = LinRange(prev_range, Max_range, nmis-1) - println("-----------------------") + end + append!(RangesToPlot, Range) + if OEW + append!(PayloadToPlot, maxPay+Wempty) + else + append!(PayloadToPlot, maxPay) + end end ``` @@ -116,15 +104,13 @@ end ```julia using PyPlot -figure() -y_OEW = y_Wemtpy .+ y_payload -plot(x_range ./ (1852.0*1000), y_OEW./ (4.448222* 1000), linestyle="-", color="b", label="OEW + Payload ") - -xlabel("Range (1000 nmi)") -ylabel("Weight (1000 lbs)") -title("Payload Range Plot") - -legend() -grid() -savefig("./PayloadRangeExample.png") +fig, ax = plt.subplots(figsize=(8,5), dpi = 300) +ax.plot(RangesToPlot ./ (1000*1852.0), PayloadToPlot./ (9.8*1000), linestyle="-", color="b", label="Payload ") +ax.set_xlabel("Range (1000 nmi)") +ax.set_ylabel("Weight (1000 kg)") +ax.legend() +ax.set_title("Payload Range Plot") +ax.grid() + +fig.savefig("./PayloadRangeExample.png") ``` diff --git a/example/example_PRD.jl b/example/example_PRD.jl index 64a45e0f..0b607b80 100644 --- a/example/example_PRD.jl +++ b/example/example_PRD.jl @@ -12,129 +12,68 @@ using TASOPT include("../src/misc/index.inc") # import indices for calling parameters -function km2nmi(km) - return km./1852.0 -end - -function nmi2km(nmi) - return nmi.*1852.0 -end - -function pax2N(pax) - pax.*215*4.448222 -end - -function N2pax(N) - N./(215*4.448222) -end - -# Alogrithm: - -# Input: Range (arbitrary), WMpay -# Size the aircraft -# For a certain range of Ranges - # Try woper for each mission other than design - # If converges then append to data struct - # Increase range of ranges - # Reduce Payload if in section B - # Else if error/ non converge - # Reduce payload - # Activate section B - # End loop if difference between min range and max range below threshold - # OR payload <= 0 - -function PayloadRange() - lbf_to_N = 4.448222 - MAX_PAYLOAD = pax2N(230) - - Range_arr_nmi = LinRange(100, 1000, 4) - Range_arr = Range_arr_nmi .* 1852.0 - Max_range = Range_arr[4] - - curr_range = Range_arr[1] - prev_range = 0.1 - - # 5. Initialze some variables - NPSS = Base.Process - NPSS_PT = true - itermax = 15 - initeng = 0 - - x_range = [] - y_payload = [] - y_WTO = [] - y_Wfuel = [] - y_Wemtpy = [] - - SEC_B = false - - # Load default model - ac = TASOPT.read_aircraft_model(joinpath(TASOPT.__TASOPTroot__, "../example/PRD_input.toml")) - - size_aircraft!(ac) - - while MAX_PAYLOAD >= 0 - println("-----------------------") - println("Starting new set of missions, \n Ranges: ", km2nmi(Range_arr), "\n Payload: ", N2pax(MAX_PAYLOAD)) - - try - ac.parm[imRange, 2:5] .= Range_arr - ac.parm[imWpay, 2:5] .= MAX_PAYLOAD - - for mi in 2:size(ac.parm)[2] - prev_range = curr_range - curr_range = ac.parm[imRange, mi] - println("Analysing Mission Number: ",mi) - println("Range: ",curr_range) - @views TASOPT.woper(ac.pari,ac.parg,ac.parm[:,mi:mi],ac.para[:,:,mi:mi],ac.pare[:,:,mi:mi], ac.para[:,:,1:1],ac.pare[:,:,1:1], itermax,initeng, NPSS_PT, NPSS) - if (ac.parm[imWfuel,mi] - ac.parg[igWfuel] > 1) || (ac.parm[imWTO,mi] - ac.parg[igWMTO] > 1 ) - printstyled("Mission designed beyond capacity!", "\n"; color=:red) - println([ac.parg[igWfuel], ac.parm[imWfuel,mi], ac.parg[igWMTO], ac.parm[imWTO,mi]]) - throw(UndefVarError([ac.parg[igWfuel], ac.parm[imWfuel,mi], ac.parg[igWMTO], ac.parm[imWTO,mi]])) - end - #---- max TO weight - WMTO = ac.parg[igWMTO] - # ---- zero-fuel weight for this mission - Wzero = WMTO-ac.parg[igWfuel]-ac.parg[igWpay] - - append!(y_Wemtpy, Wzero) - append!(x_range, curr_range) - append!(y_payload, MAX_PAYLOAD) - append!(y_Wfuel, ac.parm[imWfuel,mi]) - append!(y_WTO, ac.parm[imWTO,mi]) - if SEC_B - println("decreasing payload since Section B") - MAX_PAYLOAD = max(0, MAX_PAYLOAD-pax2N(1)) - end +function PayloadRange(ac_og, Rpts = 20, Ppts = 20, filename = "./example/PayloadRangeExample.png", OEW = false, itermax = 20.0) + # Copy aircraft structure as we will be changing range and payloads + ac = deepcopy(ac_og) + # Make an array of ranges to plot + RangeArray = ac.parm[imRange] * LinRange(0.1,1.2,Rpts) + # Store constant values to compare later + Wmax = ac.parg[igWMTO] + Fuelmax = ac.parg[igWfmax] + Wempty = ac.parg[igWMTO] - ac.parg[igWfuel] - ac.parg[igWpay] + # Arrays for plotting + RangesToPlot = [] + PayloadToPlot = [] + maxPay = ac.parm[imWpay ] + + for Range = RangeArray + # Make an array of payloads to plot + Payloads = (maxPay) * LinRange(1, 0.1, Ppts) + ac.parm[imRange] = Range + for mWpay = Payloads + println("Checking for Range (nmi): ",Range/1852.0, " and Pax = ", mWpay/(215*4.44822)) + ac.parm[imWpay ] = mWpay + # Try woper after setting new range and payload + try + @views TASOPT.woper(ac.pari,ac.parg,ac.parm[:,1:1],ac.para[:,:,1:1],ac.pare[:,:,1:1], ac.para[:,:,1:1],ac.pare[:,:,1:1], itermax,0.0) + # woper success: store maxPay, break loop + WTO = Wempty + mWpay + ac.parm[imWfuel] + mWfuel = ac.parm[imWfuel] + + # Compare with previously stored constants + if WTO > Wmax || mWfuel > Fuelmax || WTO < 0.0 || mWfuel < 0.0 + WTO = 0.0 + mWfuel = 0.0 + println("Max out error!") + else + maxPay = mWpay + println("Converged - moving to next range...") + break + end + catch + println("Not Converged - moving to lower payload...") end - println("PRD converged for all: Increasing Ranges") - prev_range = curr_range - Max_range = prev_range+(500* 1852.0) - catch - println("Decreasing Payload and range") - prev_range = curr_range - Max_range = prev_range+ ((Max_range-prev_range)*0.75) - MAX_PAYLOAD = MAX_PAYLOAD-pax2N(5) - SEC_B = true end - Range_arr = LinRange(prev_range, Max_range, 4) - println("-----------------------", abs(Max_range-prev_range)) + append!(RangesToPlot, Range) + if OEW + append!(PayloadToPlot, maxPay+Wempty) + else + append!(PayloadToPlot, maxPay) + end end - - - - figure() - y_OEW = y_Wemtpy .+ y_payload - plot(x_range ./ (1000*1852.0), y_OEW./ (9.8*1000), linestyle="-", color="b", label="OEW + Payload ") - - xlabel("Range (1000 nmi)") - ylabel("Weight (1000 kg)") - title("Payload Range Plot") - legend() - grid() - - savefig("./example/PayloadRangeExample.png") - + println(RangesToPlot) + println(PayloadToPlot) + fig, ax = plt.subplots(figsize=(8,5), dpi = 300) + ax.plot(RangesToPlot ./ (1000*1852.0), PayloadToPlot./ (9.8*1000), linestyle="-", color="b", label="Payload ") + ax.set_xlabel("Range (1000 nmi)") + ax.set_ylabel("Weight (1000 kg)") + ax.legend() + ax.set_title("Payload Range Plot") + ax.grid() + + fig.savefig(filename) end -PayloadRange() \ No newline at end of file +# Load default model +ac = TASOPT.read_aircraft_model(joinpath(TASOPT.__TASOPTroot__, "../example/PRD_input.toml")) +size_aircraft!(ac) \ No newline at end of file diff --git a/src/IO/outputs.jl b/src/IO/outputs.jl index b9815eb8..d63bbb04 100644 --- a/src/IO/outputs.jl +++ b/src/IO/outputs.jl @@ -1624,3 +1624,76 @@ function arrange_seats(seats_per_row, Rfuse, end return yseats, symmetric_seats end # function arrange_seats + +""" + PayloadRange(ac_og, Rpts, Ppts, filename, OEW, itermax) + +Function to plot a payload range diagram for an aircraft + +!!! details "🔃 Inputs and Outputs" + **Inputs:** + - `ac_og::aircraft`: Aircraft structure for payload range diagram. + - `Rpts::Int64`: Density of ranges to be plot (Optional). + - `Ppts::Int64`: Density of payloads to be plot (Optional). + - `filename::String`: filename string for the plot to be stored (Optional). + - `OEW::Boolean`: Whether to have OEW+Payload on the y-axis or just Payload (Optional). + - `itermax::Int64`: Max Iterations for woper loop (Optional). +""" + +function PayloadRange(ac_og, Rpts = 20, Ppts = 20, filename = "PayloadRangeDiagram.png", OEW = false, itermax = 20.0) + ac = deepcopy(ac_og) + RangeArray = ac.parm[imRange] * LinRange(0.1,1.2,Rpts) + maxPay = 0 + + Wmax = ac.parg[igWMTO] + Fuelmax = ac.parg[igWfmax] + Wempty = ac.parg[igWMTO] - ac.parg[igWfuel] - ac.parg[igWpay] + + RangesToPlot = [] + PayloadToPlot = [] + maxPay = ac.parm[imWpay ] + + for Range = RangeArray + Payloads = (maxPay) * LinRange(1, 0.1, Ppts) + ac.parm[imRange] = Range + for mWpay = Payloads + println("Checking for Range (nmi): ",Range/1852.0, " and Pax = ", mWpay/(215*4.44822)) + ac.parm[imWpay ] = mWpay + try + @views TASOPT.woper(ac.pari,ac.parg,ac.parm[:,1:1],ac.para[:,:,1:1],ac.pare[:,:,1:1], ac.para[:,:,1:1],ac.pare[:,:,1:1], itermax,0.0) + # woper success: store maxPay, break loop + WTO = Wempty + mWpay + ac.parm[imWfuel] + mWfuel = ac.parm[imWfuel] + + if WTO > Wmax || mWfuel > Fuelmax || WTO < 0.0 || mWfuel < 0.0 + WTO = 0.0 + mWfuel = 0.0 + println("Max out error!") + else + maxPay = mWpay + println("Converged - moving to next range...") + break + end + catch + println("Not Converged - moving to lower payload...") + end + end + append!(RangesToPlot, Range) + if OEW + append!(PayloadToPlot, maxPay+Wempty) + else + append!(PayloadToPlot, maxPay) + end + end + println(RangesToPlot) + println(PayloadToPlot) + fig, ax = plt.subplots(figsize=(8,5), dpi = 300) + ax.plot(RangesToPlot ./ (1000*1852.0), PayloadToPlot./ (9.8*1000), linestyle="-", color="b", label="Payload ") + ax.set_xlabel("Range (1000 nmi)") + ax.set_ylabel("Weight (1000 kg)") + ax.legend() + ax.set_title("Payload Range Plot") + ax.grid() + + fig.savefig(filename) +end \ No newline at end of file