Fixed US Standard Atmosphere class and tests

This commit is contained in:
Travis Hunter 2023-10-19 10:58:01 -06:00
parent 9391b85a07
commit 9b807d53a4
4 changed files with 72 additions and 32 deletions

View File

@ -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

View File

@ -32,6 +32,7 @@ private:
static utils::BinMap temperatureLapseRate;
static utils::BinMap standardTemperature;
static utils::BinMap standardDensity;
static utils::BinMap standardPressure;
};

View File

@ -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);
}

View File

@ -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)