Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Linear Wave Convergence Test #305

Merged
merged 3 commits into from
Jul 31, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
152 changes: 150 additions & 2 deletions src/system_tests/mhd_system_tests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@
// Local includes
#include "../io/io.h"
#include "../system_tests/system_tester.h"
#include "../utils/testing_utilities.h"

// =============================================================================
// Test Suite: tMHDSYSTEMLinearWavesParameterizedAngle
Expand All @@ -32,15 +33,15 @@ class tMHDSYSTEMLinearWavesParameterizedAngle : public ::testing::TestWithParam<

protected:
systemTest::SystemTestRunner waveTest;
inline static std::unordered_map<std::string, double> high_res_l2norms;

void setLaunchParams(double const &waveSpeed, double const &rEigenVec_rho, double const &rEigenVec_MomentumX,
double const &rEigenVec_MomentumY, double const &rEigenVec_MomentumZ, double const &rEigenVec_E,
double const &rEigenVec_Bx, double const &rEigenVec_By, double const &rEigenVec_Bz,
double const &pitch, double const &yaw, double const &domain, int const &domain_direction,
double const &vx = 0.0)
double const &vx = 0.0, size_t const &N = 32)
{
// Constant for all tests
size_t const N = 32;
double const gamma = 5. / 3.;
double const tOut = 2 * domain / waveSpeed;

Expand Down Expand Up @@ -154,6 +155,8 @@ TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, FastMagnetosonicWaveRightMovingC
#elif defined(PPMC)
waveTest.runL1ErrorTest(6.11E-8, 5.5E-8);
#endif // PCM

high_res_l2norms["fast_" + std::to_string(domain_direction)] = waveTest.getL2Norm();
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, FastMagnetosonicWaveLeftMovingCorrectInputExpectCorrectOutput)
Expand Down Expand Up @@ -228,6 +231,8 @@ TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, SlowMagnetosonicWaveRightMovingC
#elif defined(PPMC)
waveTest.runL1ErrorTest(1.45E-9, 1.3E-9);
#endif // PCM

high_res_l2norms["slow_" + std::to_string(domain_direction)] = waveTest.getL2Norm();
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, SlowMagnetosonicWaveLeftMovingCorrectInputExpectCorrectOutput)
Expand Down Expand Up @@ -301,6 +306,8 @@ TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, AlfvenWaveRightMovingCorrectInpu
#elif defined(PPMC)
waveTest.runL1ErrorTest(1.95e-09, 2.16e-09);
#endif // PCM

high_res_l2norms["alfven_" + std::to_string(domain_direction)] = waveTest.getL2Norm();
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, AlfvenWaveLeftMovingCorrectInputExpectCorrectOutput)
Expand Down Expand Up @@ -375,6 +382,147 @@ TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, MHDContactWaveCorrectInputExpect
#elif defined(PPMC)
waveTest.runL1ErrorTest(1.41e-09, 1.5E-09);
#endif // PCM

high_res_l2norms["contact_" + std::to_string(domain_direction)] = waveTest.getL2Norm();
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, FastMagnetosonicWaveExpectSecondOrderConvergence)
{
// Get the test parameters
auto [pitch, yaw, domain, domain_direction] = GetParam();

// Specific to this test
double const waveSpeed = 2.;
std::vector<int> const numTimeSteps = {107, 102, 110};

double const prefix = 1. / (2 * std::sqrt(5));
double const rEigenVec_rho = prefix * 2;
double const rEigenVec_MomentumX = prefix * 4;
double const rEigenVec_MomentumY = prefix * -2;
double const rEigenVec_MomentumZ = prefix * 0;
double const rEigenVec_Bx = prefix * 0;
double const rEigenVec_By = prefix * 4;
double const rEigenVec_Bz = prefix * 0;
double const rEigenVec_E = prefix * 9;

// Set the launch parameters
setLaunchParams(waveSpeed, rEigenVec_rho, rEigenVec_MomentumX, rEigenVec_MomentumY, rEigenVec_MomentumZ, rEigenVec_E,
rEigenVec_Bx, rEigenVec_By, rEigenVec_Bz, pitch, yaw, domain, domain_direction, 0.0, 16);

// Set the number of timesteps
waveTest.setFiducialNumTimeSteps(numTimeSteps[domain_direction - 1]);

// Run the wave
waveTest.runL1ErrorTest(3.0E-8, 4.0E-8);

// Check the scaling
double const low_res_l2norm = waveTest.getL2Norm();
testingUtilities::checkResults(4.0, low_res_l2norm / high_res_l2norms["fast_" + std::to_string(domain_direction)], "",
0.17);
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, SlowMagnetosonicWaveExpectSecondOrderConvergence)
{
// Get the test parameters
auto [pitch, yaw, domain, domain_direction] = GetParam();

// Specific to this test
double const waveSpeed = 0.5;
std::vector<int> const numTimeSteps = {427, 407, 440};

double const prefix = 1. / (2 * std::sqrt(5));
double const rEigenVec_rho = prefix * 4;
double const rEigenVec_MomentumX = prefix * 2;
double const rEigenVec_MomentumY = prefix * 4;
double const rEigenVec_MomentumZ = prefix * 0;
double const rEigenVec_Bx = prefix * 0;
double const rEigenVec_By = prefix * -2;
double const rEigenVec_Bz = prefix * 0;
double const rEigenVec_E = prefix * 3;

// Set the launch parameters
setLaunchParams(waveSpeed, rEigenVec_rho, rEigenVec_MomentumX, rEigenVec_MomentumY, rEigenVec_MomentumZ, rEigenVec_E,
rEigenVec_Bx, rEigenVec_By, rEigenVec_Bz, pitch, yaw, domain, domain_direction, 0.0, 16);

// Set the number of timesteps
waveTest.setFiducialNumTimeSteps(numTimeSteps[domain_direction - 1]);

// Run the wave
waveTest.runL1ErrorTest(3.0E-8, 4.0E-8);

// Check the scaling
double const low_res_l2norm = waveTest.getL2Norm();
testingUtilities::checkResults(4.0, low_res_l2norm / high_res_l2norms["slow_" + std::to_string(domain_direction)], "",
0.17);
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, AlfvenWaveExpectSecondOrderConvergence)
{
// Get the test parameters
auto [pitch, yaw, domain, domain_direction] = GetParam();

// Specific to this test
double const waveSpeed = 1.0;
std::vector<int> const numTimeSteps = {214, 204, 220};

double const rEigenVec_rho = 0;
double const rEigenVec_MomentumX = 0;
double const rEigenVec_MomentumY = 0;
double const rEigenVec_MomentumZ = -1;
double const rEigenVec_Bx = 0;
double const rEigenVec_By = 0;
double const rEigenVec_Bz = 1;
double const rEigenVec_E = 0;

// Set the launch parameters
setLaunchParams(waveSpeed, rEigenVec_rho, rEigenVec_MomentumX, rEigenVec_MomentumY, rEigenVec_MomentumZ, rEigenVec_E,
rEigenVec_Bx, rEigenVec_By, rEigenVec_Bz, pitch, yaw, domain, domain_direction, 0.0, 16);

// Set the number of timesteps
waveTest.setFiducialNumTimeSteps(numTimeSteps[domain_direction - 1]);

// Run the wave
waveTest.runL1ErrorTest(3.0E-8, 4.0E-8);

// Check the scaling
double const low_res_l2norm = waveTest.getL2Norm();
testingUtilities::checkResults(4.0, low_res_l2norm / high_res_l2norms["alfven_" + std::to_string(domain_direction)],
"", 0.17);
}

TEST_P(tMHDSYSTEMLinearWavesParameterizedAngle, MHDContactWaveExpectSecondOrderConvergence)
{
// Get the test parameters
auto [pitch, yaw, domain, domain_direction] = GetParam();

// Specific to this test
double const waveSpeed = 1.0;
std::vector<int> const numTimeSteps = {321, 310, 327};

double const rEigenVec_rho = 1;
double const rEigenVec_MomentumX = 1;
double const rEigenVec_MomentumY = 0;
double const rEigenVec_MomentumZ = 0;
double const rEigenVec_Bx = 0;
double const rEigenVec_By = 0;
double const rEigenVec_Bz = 0;
double const rEigenVec_E = 0.5;
double const velocityX = waveSpeed;

// Set the launch parameters
setLaunchParams(waveSpeed, rEigenVec_rho, rEigenVec_MomentumX, rEigenVec_MomentumY, rEigenVec_MomentumZ, rEigenVec_E,
rEigenVec_Bx, rEigenVec_By, rEigenVec_Bz, pitch, yaw, domain, domain_direction, velocityX, 16);

// Set the number of timesteps
waveTest.setFiducialNumTimeSteps(numTimeSteps[domain_direction - 1]);

// Run the wave
waveTest.runL1ErrorTest(3.0E-8, 4.0E-8);

// Check the scaling
double const low_res_l2norm = waveTest.getL2Norm();
testingUtilities::checkResults(4.0, low_res_l2norm / high_res_l2norms["contact_" + std::to_string(domain_direction)],
"", 0.17);
}

INSTANTIATE_TEST_SUITE_P(, tMHDSYSTEMLinearWavesParameterizedAngle,
Expand Down
16 changes: 8 additions & 8 deletions src/system_tests/system_tester.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -108,7 +108,7 @@ void systemTest::SystemTestRunner::runTest(bool const &compute_L2_norm_only, dou
<< std::endl;

// Compute the L1 Error.
double L2Norm = 0;
_L2Norm = 0;
double maxError = 0;
// Loop over the datasets to be tested
for (auto const &dataSetName : _fiducialDataSetNames) {
Expand Down Expand Up @@ -187,14 +187,14 @@ void systemTest::SystemTestRunner::runTest(bool const &compute_L2_norm_only, dou

if (compute_L2_norm_only) {
L1_error /= static_cast<double>(testDims[0] * testDims[1] * testDims[2]);
L2Norm += L1_error * L1_error;
_L2Norm += L1_error * L1_error;
}
}

if (compute_L2_norm_only) {
// Check the L2 Norm
L2Norm = std::sqrt(L2Norm);
EXPECT_LT(L2Norm, maxAllowedL1Error) << "the norm of the L1 error vector has exceeded the allowed value";
_L2Norm = std::sqrt(_L2Norm);
EXPECT_LT(_L2Norm, maxAllowedL1Error) << "the norm of the L1 error vector has exceeded the allowed value";

// Check the Max Error
EXPECT_LT(maxError, maxAllowedError) << "The maximum error has exceeded the allowed value";
Expand Down Expand Up @@ -270,7 +270,7 @@ void systemTest::SystemTestRunner::runL1ErrorTest(double const &maxAllowedL1Erro
<< std::endl;

// Loop over the datasets to be tested
double L2Norm = 0;
_L2Norm = 0;
double maxError = 0;
for (auto const &dataSetName : _fiducialDataSetNames) {
if (dataSetName == "GasEnergy") {
Expand Down Expand Up @@ -312,16 +312,16 @@ void systemTest::SystemTestRunner::runL1ErrorTest(double const &maxAllowedL1Erro
}

L1_error /= static_cast<double>(initialDims[0] * initialDims[1] * initialDims[2]);
L2Norm += L1_error * L1_error;
_L2Norm += L1_error * L1_error;

// Perform the correctness check
EXPECT_LT(L1_error, maxAllowedL1Error)
<< "the L1 error for the " << dataSetName << " data has exceeded the allowed value";
}

// Check the L2 Norm
L2Norm = std::sqrt(L2Norm);
EXPECT_LT(L2Norm, maxAllowedL1Error) << "the norm of the L1 error vector has exceeded the allowed value";
_L2Norm = std::sqrt(_L2Norm);
EXPECT_LT(_L2Norm, maxAllowedL1Error) << "the norm of the L1 error vector has exceeded the allowed value";

// Check the Max Error
EXPECT_LT(maxError, maxAllowedError) << "The maximum error has exceeded the allowed value";
Expand Down
10 changes: 10 additions & 0 deletions src/system_tests/system_tester.h
Original file line number Diff line number Diff line change
Expand Up @@ -106,6 +106,13 @@ class systemTest::SystemTestRunner
*/
std::string getChollaSettingsFilePath() { return _chollaSettingsPath; };

/*!
* \brief Get the L2Norm
*
* \return double The L2Norm of the last run test
*/
double getL2Norm() { return _L2Norm; };

/*!
* \brief Get the Output Directory object
*
Expand Down Expand Up @@ -304,6 +311,9 @@ class systemTest::SystemTestRunner
/// appear to differ from NVIDIA/GCC/XL by roughly 1E-12
double _fixedEpsilon = 5.0E-12;

/// The L2 norm of the error vector
double _L2Norm;

/// Flag to indicate if a fiducial HDF5 data file is being used or a
/// programmatically generated H5File object. `true` = use a file, `false` =
/// use generated H5File object
Expand Down