From 9b807d53a4e074f0c8f7ea7d53e473479479a189 Mon Sep 17 00:00:00 2001 From: Travis Hunter Date: Thu, 19 Oct 2023 10:58:01 -0600 Subject: [PATCH] Fixed US Standard Atmosphere class and tests --- sim/USStandardAtmosphere.cpp | 59 ++++++++++++++++++++----- sim/USStandardAtmosphere.h | 1 + sim/tests/USStandardAtmosphereTests.cpp | 42 +++++++++--------- utils/BinMap.cpp | 2 +- 4 files changed, 72 insertions(+), 32 deletions(-) diff --git a/sim/USStandardAtmosphere.cpp b/sim/USStandardAtmosphere.cpp index 83e4753..54231ba 100644 --- a/sim/USStandardAtmosphere.cpp +++ b/sim/USStandardAtmosphere.cpp @@ -36,13 +36,13 @@ utils::BinMap initTemps() utils::BinMap initLapseRates() { utils::BinMap map; - map.insert(std::make_pair(0.0, -0.0065)); + map.insert(std::make_pair(0.0, 0.0065)); map.insert(std::make_pair(11000.0, 0.0)); - map.insert(std::make_pair(20000.0, 0.001)); - map.insert(std::make_pair(32000.0, 0.0028)); + map.insert(std::make_pair(20000.0, -0.001)); + map.insert(std::make_pair(32000.0, -0.0028)); map.insert(std::make_pair(47000.0, 0.0)); - map.insert(std::make_pair(51000.0, -0.0028)); - map.insert(std::make_pair(71000.0, -0.002)); + map.insert(std::make_pair(51000.0, 0.0028)); + map.insert(std::make_pair(71000.0, 0.002)); return map; } @@ -61,9 +61,24 @@ utils::BinMap initDensities() return map; } +utils::BinMap initPressures() +{ + utils::BinMap map; + map.insert(std::make_pair(0.0, 101325)); + map.insert(std::make_pair(11000.0, 22632.1)); + map.insert(std::make_pair(20000.0, 5474.89)); + map.insert(std::make_pair(32000.0, 868.02)); + map.insert(std::make_pair(47000.0, 110.91)); + map.insert(std::make_pair(51000.0, 66.94)); + map.insert(std::make_pair(71000.0, 3.96)); + + return map; +} + utils::BinMap USStandardAtmosphere::temperatureLapseRate(initLapseRates()); utils::BinMap USStandardAtmosphere::standardTemperature(initTemps()); utils::BinMap USStandardAtmosphere::standardDensity(initDensities()); +utils::BinMap USStandardAtmosphere::standardPressure(initPressures()); @@ -81,28 +96,50 @@ double USStandardAtmosphere::getDensity(double altitude) { if(utils::math::floatingPointEqual(temperatureLapseRate[altitude], 0.0)) { - return standardDensity[altitude] * std::exp((-Constants::g0 * Constants::airMolarMass * (altitude - standardDensity.getBinBase(altitude))) + return standardDensity[altitude] * std::exp( + (-Constants::g0 * Constants::airMolarMass * (altitude - standardDensity.getBinBase(altitude))) / (Constants::Rstar * standardTemperature[altitude])); } else { - double base = (standardTemperature[altitude] - temperatureLapseRate[altitude] * (altitude - standardDensity.getBinBase(altitude))) / standardTemperature[altitude]; + double base = (standardTemperature[altitude] - temperatureLapseRate[altitude] * + (altitude - standardDensity.getBinBase(altitude))) / standardTemperature[altitude]; double exponent = (Constants::g0 * Constants::airMolarMass) / (Constants::Rstar * temperatureLapseRate[altitude]) - 1.0; return standardDensity[altitude] * std::pow(base, exponent); - } } -double USStandardAtmosphere::getTemperature(double /*altitude*/) +double USStandardAtmosphere::getTemperature(double altitude) { + double baseTemp = standardTemperature[altitude]; + double baseAltitude = standardTemperature.getBinBase(altitude); + return baseTemp - (altitude - baseAltitude) * temperatureLapseRate[altitude]; + return 0.0; } -double USStandardAtmosphere::getPressure(double /* altitude */) +double USStandardAtmosphere::getPressure(double altitude) { - return 0.0; + if(utils::math::floatingPointEqual(temperatureLapseRate[altitude], 0.0)) + { + return standardPressure[altitude] * std::exp( + (-Constants::g0 * Constants::airMolarMass * (altitude - standardPressure.getBinBase(altitude))) + / (Constants::Rstar * standardTemperature[altitude])); + + } + else + { + double base = (standardTemperature[altitude] - temperatureLapseRate[altitude] * + (altitude - standardPressure.getBinBase(altitude))) / standardTemperature[altitude]; + + double exponent = (Constants::g0 * Constants::airMolarMass) / + (Constants::Rstar * temperatureLapseRate[altitude]); + + return standardPressure[altitude] * std::pow(base, exponent); + } + } } // namespace sim diff --git a/sim/USStandardAtmosphere.h b/sim/USStandardAtmosphere.h index 2835d81..4b18083 100644 --- a/sim/USStandardAtmosphere.h +++ b/sim/USStandardAtmosphere.h @@ -32,6 +32,7 @@ private: static utils::BinMap temperatureLapseRate; static utils::BinMap standardTemperature; static utils::BinMap standardDensity; + static utils::BinMap standardPressure; }; diff --git a/sim/tests/USStandardAtmosphereTests.cpp b/sim/tests/USStandardAtmosphereTests.cpp index 1424232..a1375a4 100644 --- a/sim/tests/USStandardAtmosphereTests.cpp +++ b/sim/tests/USStandardAtmosphereTests.cpp @@ -6,24 +6,26 @@ TEST(USStandardAtmosphereTests, DensityTests) { sim::USStandardAtmosphere atmosphere; - EXPECT_DOUBLE_EQ(atmosphere.getDensity(0.0), 1.225); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(1000.0), 1.112); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(2000.0), 1.007); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(3000.0), 0.9093); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(4000.0), 0.8194); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(5000.0), 0.7364); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(6000.0), 0.6601); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(7000.0), 0.5900); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(8000.0), 0.5258); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(9000.0), 0.4647); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(10000.0), 0.4135); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(15000.0), 0.1948); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(20000.0), 0.08891); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(25000.0), 0.04008); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(30000.0), 0.01841); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(40000.0), 0.003996); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(50000.0), 0.001027); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(60000.0), 0.0003097); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(70000.0), 0.00008283); - EXPECT_DOUBLE_EQ(atmosphere.getDensity(80000.0), 0.00001846); + // Test that the calucated values are with 1% of the published values in the NOAA report + EXPECT_NEAR(atmosphere.getDensity(0.0) / 1.225, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(1000.0) / 1.112, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(2000.0) / 1.007, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(3000.0) / 0.9093, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(4000.0) / 0.8194, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(5000.0) / 0.7364, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(6000.0) / 0.6601, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(7000.0) / 0.5900, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(8000.0) / 0.5258, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(9000.0) / 0.4647, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(10000.0) / 0.4135, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(15000.0) / 0.1948, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(20000.0) / 0.08891, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(25000.0) / 0.039466, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(30000.0) / 0.018012, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(40000.0) / 0.0038510, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(50000.0) / 0.00097752, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(60000.0) / 0.00028832, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(70000.0) / 0.000074243, 1.0, 0.01); + EXPECT_NEAR(atmosphere.getDensity(80000.0) / 0.000015701, 1.0, 0.01); + } diff --git a/utils/BinMap.cpp b/utils/BinMap.cpp index 95d17d6..5ebe68b 100644 --- a/utils/BinMap.cpp +++ b/utils/BinMap.cpp @@ -61,7 +61,7 @@ double BinMap::operator[](double key) // Increment it and start searching If we reach the end without finding an existing key // greater than our search term, then we've just hit the last bin and return that iter++; - double retVal = bins.end()->second; + double retVal = bins.back().second; while(iter != bins.end()) { if(key < iter->first)