From 2b195ee93ec34b2fa56549759275404a5111268e Mon Sep 17 00:00:00 2001 From: tianchiTJ <149181511@qq.com> Date: Wed, 4 Nov 2020 09:48:04 -0800 Subject: [PATCH] Correct the names of nodal phases --- include/node.tcc | 2 +- include/node_base.h | 6 +- include/node_twophase.tcc | 92 +- include/solvers/mpm_explicit_twophase.tcc | 32 +- tests/nodes/node_test.cc | 18 +- tests/nodes/node_twophase_test.cc | 996 +++++++++++----------- tests/particles/particle_twophase_test.cc | 70 +- 7 files changed, 608 insertions(+), 608 deletions(-) diff --git a/include/node.tcc b/include/node.tcc index e4840f1a0..5583ea056 100644 --- a/include/node.tcc +++ b/include/node.tcc @@ -626,7 +626,7 @@ void mpm::Node momentum = property_handle_->property("momenta", prop_id_, *mitr, Tdim); const Eigen::Matrix change_in_momenta = - velocity_.col(mpm::NodePhase::nSolid) * mass - momentum; + velocity_.col(mpm::NodePhase::NSolid) * mass - momentum; property_handle_->update_property("change_in_momenta", prop_id_, *mitr, change_in_momenta, Tdim); } diff --git a/include/node_base.h b/include/node_base.h index baae25e1a..e5cf021c8 100644 --- a/include/node_base.h +++ b/include/node_base.h @@ -19,10 +19,10 @@ namespace mpm { //! Particle phases enum NodePhase : unsigned int { - nSolid = 0, - nLiquid = 1, + NSolid = 0, + NLiquid = 1, nGas = 2, - nMixture = 0 + NMixture = 0 }; //! NodeBase base class for nodes diff --git a/include/node_twophase.tcc b/include/node_twophase.tcc index 3d153a0eb..d1f6f033e 100644 --- a/include/node_twophase.tcc +++ b/include/node_twophase.tcc @@ -19,26 +19,26 @@ bool mpm::Node:: compute_acceleration_velocity_twophase_explicit(double dt) noexcept { bool status = false; const double tolerance = 1.0E-15; - if (this->mass(mpm::NodePhase::nSolid) > tolerance && - this->mass(mpm::NodePhase::nLiquid) > tolerance) { + if (this->mass(mpm::NodePhase::NSolid) > tolerance && + this->mass(mpm::NodePhase::NLiquid) > tolerance) { // Compute drag force VectorDim drag_force = drag_force_coefficient_.cwiseProduct( - velocity_.col(mpm::NodePhase::nLiquid) - - velocity_.col(mpm::NodePhase::nSolid)); + velocity_.col(mpm::NodePhase::NLiquid) - + velocity_.col(mpm::NodePhase::NSolid)); // Acceleration of pore fluid (momentume balance of fluid phase) - this->acceleration_.col(mpm::NodePhase::nLiquid) = - (this->external_force_.col(mpm::NodePhase::nLiquid) + - this->internal_force_.col(mpm::NodePhase::nLiquid) - drag_force) / - this->mass_(mpm::NodePhase::nLiquid); + this->acceleration_.col(mpm::NodePhase::NLiquid) = + (this->external_force_.col(mpm::NodePhase::NLiquid) + + this->internal_force_.col(mpm::NodePhase::NLiquid) - drag_force) / + this->mass_(mpm::NodePhase::NLiquid); // Acceleration of solid skeleton (momentume balance of mixture) - this->acceleration_.col(mpm::NodePhase::nSolid) = - (this->external_force_.col(mpm::NodePhase::nMixture) + - this->internal_force_.col(mpm::NodePhase::nMixture) - - this->mass_(mpm::NodePhase::nLiquid) * - this->acceleration_.col(mpm::NodePhase::nLiquid)) / - this->mass_(mpm::NodePhase::nSolid); + this->acceleration_.col(mpm::NodePhase::NSolid) = + (this->external_force_.col(mpm::NodePhase::NMixture) + + this->internal_force_.col(mpm::NodePhase::NMixture) - + this->mass_(mpm::NodePhase::NLiquid) * + this->acceleration_.col(mpm::NodePhase::NLiquid)) / + this->mass_(mpm::NodePhase::NSolid); // Apply friction constraints this->apply_friction_constraints(dt); @@ -52,14 +52,14 @@ bool mpm::Node:: // Set a threshold for (unsigned i = 0; i < Tdim; ++i) { - if (std::abs(velocity_.col(mpm::NodePhase::nSolid)(i)) < tolerance) - velocity_.col(mpm::NodePhase::nSolid)(i) = 0.; - if (std::abs(acceleration_.col(mpm::NodePhase::nSolid)(i)) < tolerance) - acceleration_.col(mpm::NodePhase::nSolid)(i) = 0.; - if (std::abs(velocity_.col(mpm::NodePhase::nLiquid)(i)) < tolerance) - velocity_.col(mpm::NodePhase::nLiquid)(i) = 0.; - if (std::abs(acceleration_.col(mpm::NodePhase::nLiquid)(i)) < tolerance) - acceleration_.col(mpm::NodePhase::nLiquid)(i) = 0.; + if (std::abs(velocity_.col(mpm::NodePhase::NSolid)(i)) < tolerance) + velocity_.col(mpm::NodePhase::NSolid)(i) = 0.; + if (std::abs(acceleration_.col(mpm::NodePhase::NSolid)(i)) < tolerance) + acceleration_.col(mpm::NodePhase::NSolid)(i) = 0.; + if (std::abs(velocity_.col(mpm::NodePhase::NLiquid)(i)) < tolerance) + velocity_.col(mpm::NodePhase::NLiquid)(i) = 0.; + if (std::abs(acceleration_.col(mpm::NodePhase::NLiquid)(i)) < tolerance) + acceleration_.col(mpm::NodePhase::NLiquid)(i) = 0.; } status = true; } @@ -74,36 +74,36 @@ bool mpm::Node:: bool status = false; const double tolerance = 1.0E-15; - if (this->mass(mpm::NodePhase::nSolid) > tolerance && - this->mass(mpm::NodePhase::nLiquid) > tolerance) { + if (this->mass(mpm::NodePhase::NSolid) > tolerance && + this->mass(mpm::NodePhase::NLiquid) > tolerance) { // Compute drag force VectorDim drag_force = drag_force_coefficient_.cwiseProduct( - velocity_.col(mpm::NodePhase::nLiquid) - - velocity_.col(mpm::NodePhase::nSolid)); + velocity_.col(mpm::NodePhase::NLiquid) - + velocity_.col(mpm::NodePhase::NSolid)); // Unbalanced force of liquid phase auto unbalanced_force_liquid = - this->external_force_.col(mpm::NodePhase::nLiquid) + - this->internal_force_.col(mpm::NodePhase::nLiquid) - drag_force; + this->external_force_.col(mpm::NodePhase::NLiquid) + + this->internal_force_.col(mpm::NodePhase::NLiquid) - drag_force; // Acceleration of liquid phase (momentume balance of fluid phase) - this->acceleration_.col(mpm::NodePhase::nLiquid) = + this->acceleration_.col(mpm::NodePhase::NLiquid) = (unbalanced_force_liquid - damping_factor * unbalanced_force_liquid.norm() * - this->velocity_.col(mpm::NodePhase::nLiquid).cwiseSign()) / - this->mass(mpm::NodePhase::nLiquid); + this->velocity_.col(mpm::NodePhase::NLiquid).cwiseSign()) / + this->mass(mpm::NodePhase::NLiquid); // Unbalanced force of solid phase auto unbalanced_force_solid = - this->external_force_.col(mpm::NodePhase::nMixture) + - this->internal_force_.col(mpm::NodePhase::nMixture) - - this->mass_(mpm::NodePhase::nLiquid) * - this->acceleration_.col(mpm::NodePhase::nLiquid); + this->external_force_.col(mpm::NodePhase::NMixture) + + this->internal_force_.col(mpm::NodePhase::NMixture) - + this->mass_(mpm::NodePhase::NLiquid) * + this->acceleration_.col(mpm::NodePhase::NLiquid); // Acceleration of solid phase (momentume balance of mixture) - this->acceleration_.col(mpm::NodePhase::nSolid) = + this->acceleration_.col(mpm::NodePhase::NSolid) = (unbalanced_force_solid - damping_factor * unbalanced_force_solid.norm() * - this->velocity_.col(mpm::NodePhase::nSolid).cwiseSign()) / - this->mass(mpm::NodePhase::nSolid); + this->velocity_.col(mpm::NodePhase::NSolid).cwiseSign()) / + this->mass(mpm::NodePhase::NSolid); // Apply friction constraints this->apply_friction_constraints(dt); @@ -117,14 +117,14 @@ bool mpm::Node:: // Set a threshold for (unsigned i = 0; i < Tdim; ++i) { - if (std::abs(velocity_.col(mpm::NodePhase::nSolid)(i)) < tolerance) - velocity_.col(mpm::NodePhase::nSolid)(i) = 0.; - if (std::abs(acceleration_.col(mpm::NodePhase::nSolid)(i)) < tolerance) - acceleration_.col(mpm::NodePhase::nSolid)(i) = 0.; - if (std::abs(velocity_.col(mpm::NodePhase::nLiquid)(i)) < tolerance) - velocity_.col(mpm::NodePhase::nLiquid)(i) = 0.; - if (std::abs(acceleration_.col(mpm::NodePhase::nLiquid)(i)) < tolerance) - acceleration_.col(mpm::NodePhase::nLiquid)(i) = 0.; + if (std::abs(velocity_.col(mpm::NodePhase::NSolid)(i)) < tolerance) + velocity_.col(mpm::NodePhase::NSolid)(i) = 0.; + if (std::abs(acceleration_.col(mpm::NodePhase::NSolid)(i)) < tolerance) + acceleration_.col(mpm::NodePhase::NSolid)(i) = 0.; + if (std::abs(velocity_.col(mpm::NodePhase::NLiquid)(i)) < tolerance) + velocity_.col(mpm::NodePhase::NLiquid)(i) = 0.; + if (std::abs(acceleration_.col(mpm::NodePhase::NLiquid)(i)) < tolerance) + acceleration_.col(mpm::NodePhase::NLiquid)(i) = 0.; } status = true; } diff --git a/include/solvers/mpm_explicit_twophase.tcc b/include/solvers/mpm_explicit_twophase.tcc index debb780db..822e8394e 100644 --- a/include/solvers/mpm_explicit_twophase.tcc +++ b/include/solvers/mpm_explicit_twophase.tcc @@ -161,29 +161,29 @@ bool mpm::MPMExplicitTwoPhase::solve() { // MPI all reduce nodal mass for solid phase mesh_->template nodal_halo_exchange( std::bind(&mpm::NodeBase::mass, std::placeholders::_1, - mpm::NodePhase::nSolid), + mpm::NodePhase::NSolid), std::bind(&mpm::NodeBase::update_mass, std::placeholders::_1, - false, mpm::NodePhase::nSolid, std::placeholders::_2)); + false, mpm::NodePhase::NSolid, std::placeholders::_2)); // MPI all reduce nodal momentum for solid phase mesh_->template nodal_halo_exchange, Tdim>( std::bind(&mpm::NodeBase::momentum, std::placeholders::_1, - mpm::NodePhase::nSolid), + mpm::NodePhase::NSolid), std::bind(&mpm::NodeBase::update_momentum, - std::placeholders::_1, false, mpm::NodePhase::nSolid, + std::placeholders::_1, false, mpm::NodePhase::NSolid, std::placeholders::_2)); // MPI all reduce nodal mass for liquid phase mesh_->template nodal_halo_exchange( std::bind(&mpm::NodeBase::mass, std::placeholders::_1, - mpm::NodePhase::nLiquid), + mpm::NodePhase::NLiquid), std::bind(&mpm::NodeBase::update_mass, std::placeholders::_1, - false, mpm::NodePhase::nLiquid, std::placeholders::_2)); + false, mpm::NodePhase::NLiquid, std::placeholders::_2)); // MPI all reduce nodal momentum for liquid phase mesh_->template nodal_halo_exchange, Tdim>( std::bind(&mpm::NodeBase::momentum, std::placeholders::_1, - mpm::NodePhase::nLiquid), + mpm::NodePhase::NLiquid), std::bind(&mpm::NodeBase::update_momentum, - std::placeholders::_1, false, mpm::NodePhase::nLiquid, + std::placeholders::_1, false, mpm::NodePhase::NLiquid, std::placeholders::_2)); } #endif @@ -263,31 +263,31 @@ bool mpm::MPMExplicitTwoPhase::solve() { // MPI all reduce external force of mixture mesh_->template nodal_halo_exchange, Tdim>( std::bind(&mpm::NodeBase::external_force, std::placeholders::_1, - mpm::NodePhase::nMixture), + mpm::NodePhase::NMixture), std::bind(&mpm::NodeBase::update_external_force, - std::placeholders::_1, false, mpm::NodePhase::nMixture, + std::placeholders::_1, false, mpm::NodePhase::NMixture, std::placeholders::_2)); // MPI all reduce external force of pore fluid mesh_->template nodal_halo_exchange, Tdim>( std::bind(&mpm::NodeBase::external_force, std::placeholders::_1, - mpm::NodePhase::nLiquid), + mpm::NodePhase::NLiquid), std::bind(&mpm::NodeBase::update_external_force, - std::placeholders::_1, false, mpm::NodePhase::nLiquid, + std::placeholders::_1, false, mpm::NodePhase::NLiquid, std::placeholders::_2)); // MPI all reduce internal force of mixture mesh_->template nodal_halo_exchange, Tdim>( std::bind(&mpm::NodeBase::internal_force, std::placeholders::_1, - mpm::NodePhase::nMixture), + mpm::NodePhase::NMixture), std::bind(&mpm::NodeBase::update_internal_force, - std::placeholders::_1, false, mpm::NodePhase::nMixture, + std::placeholders::_1, false, mpm::NodePhase::NMixture, std::placeholders::_2)); // MPI all reduce internal force of pore liquid mesh_->template nodal_halo_exchange, Tdim>( std::bind(&mpm::NodeBase::internal_force, std::placeholders::_1, - mpm::NodePhase::nLiquid), + mpm::NodePhase::NLiquid), std::bind(&mpm::NodeBase::update_internal_force, - std::placeholders::_1, false, mpm::NodePhase::nLiquid, + std::placeholders::_1, false, mpm::NodePhase::NLiquid, std::placeholders::_2)); // MPI all reduce drag force diff --git a/tests/nodes/node_test.cc b/tests/nodes/node_test.cc index 7552300eb..8b2d4968d 100644 --- a/tests/nodes/node_test.cc +++ b/tests/nodes/node_test.cc @@ -159,13 +159,13 @@ TEST_CASE("Node is checked for 1D case", "[node][1D]") { // Check pressure constraints SECTION("Check nodal pressure constraints") { // Check assign pressure constraint - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nSolid, 8000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NSolid, 8000, nullptr) == true); // Check apply pressure constraint REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nSolid)); + node->apply_pressure_constraint(mpm::NodePhase::NSolid)); // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(8000).epsilon(Tolerance)); } } @@ -636,13 +636,13 @@ TEST_CASE("Node is checked for 2D case", "[node][2D]") { // Check pressure constraints SECTION("Check nodal pressure constraints") { // Check assign pressure constraint - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nSolid, 8000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NSolid, 8000, nullptr) == true); // Check apply pressure constraint REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nSolid)); + node->apply_pressure_constraint(mpm::NodePhase::NSolid)); // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(8000).epsilon(Tolerance)); } } @@ -1250,13 +1250,13 @@ TEST_CASE("Node is checked for 3D case", "[node][3D]") { // Check pressure constraints SECTION("Check nodal pressure constraints") { // Check assign pressure constraint - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nSolid, 8000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NSolid, 8000, nullptr) == true); // Check apply pressure constraint REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nSolid)); + node->apply_pressure_constraint(mpm::NodePhase::NSolid)); // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(8000).epsilon(Tolerance)); } } diff --git a/tests/nodes/node_twophase_test.cc b/tests/nodes/node_twophase_test.cc index eebb0182c..d5fa22f1f 100644 --- a/tests/nodes/node_twophase_test.cc +++ b/tests/nodes/node_twophase_test.cc @@ -121,111 +121,111 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { std::make_shared>(id, coords); // Check mass - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); double solid_mass = 100.5; double liquid_mass = 200.5; // Update mass to 100.5 and 200.5 REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(true, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(true, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.5).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(200.5).epsilon(Tolerance)); // Update mass to 201 and 401 REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(true, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(true, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(201.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(401.0).epsilon(Tolerance)); // Assign mass to 100 and 200 solid_mass = 100.; liquid_mass = 200.; REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(200.0).epsilon(Tolerance)); SECTION("Check nodal pressure") { // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); double pressure = 1000.7; double pore_pressure = 2000.7; // Update pressure to 1000.7 and 2000.7 - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NSolid, solid_mass * pressure)); - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NLiquid, liquid_mass * pore_pressure)); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.7).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.7).epsilon(Tolerance)); // Update pressure to 2001.4 and 4001.4 - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NSolid, solid_mass * pressure)); - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NLiquid, liquid_mass * pore_pressure)); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(2001.4).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(4001.4).epsilon(Tolerance)); // Assign pressure to 1000 and 2000 pressure = 1000.; pore_pressure = 2000.; - node->assign_pressure(mpm::NodePhase::nSolid, pressure); - node->assign_pressure(mpm::NodePhase::nLiquid, pore_pressure); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + node->assign_pressure(mpm::NodePhase::NSolid, pressure); + node->assign_pressure(mpm::NodePhase::NLiquid, pore_pressure); + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.0).epsilon(Tolerance)); // Assign mass to 0 solid_mass = 0.; liquid_mass = 0.; REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); // Try to update pressure to 2000, should throw and keep to 1000. - node->assign_pressure(mpm::NodePhase::nSolid, pressure); - node->assign_pressure(mpm::NodePhase::nLiquid, pore_pressure); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + node->assign_pressure(mpm::NodePhase::NSolid, pressure); + node->assign_pressure(mpm::NodePhase::NLiquid, pore_pressure); + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.0).epsilon(Tolerance)); // Check pressure constraints SECTION("Check nodal pressure constraints") { // Check assign pressure constraint - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nSolid, 8000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NSolid, 8000, nullptr) == true); - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nLiquid, 7000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NLiquid, 7000, nullptr) == true); // Check apply pressure constraint REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nSolid)); + node->apply_pressure_constraint(mpm::NodePhase::NSolid)); REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nLiquid)); + node->apply_pressure_constraint(mpm::NodePhase::NLiquid)); // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(8000).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(7000).epsilon(Tolerance)); } } @@ -237,45 +237,45 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { // Check current external force is zero for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Update force to 10.0 REQUIRE_NOTHROW( - node->update_external_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::nLiquid, + node->update_external_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } // Update force to 20.0 REQUIRE_NOTHROW( - node->update_external_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::nLiquid, + node->update_external_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Assign force as 10.0 REQUIRE_NOTHROW( - node->update_external_force(false, mpm::NodePhase::nMixture, force)); + node->update_external_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } @@ -283,7 +283,7 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { // Set external force to zero force.setZero(); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nMixture, force)); + false, mpm::NodePhase::NMixture, force)); // concentrated force std::shared_ptr ffunction = nullptr; @@ -291,38 +291,38 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { const unsigned Direction = 0; // Check external force for (unsigned i = 0; i < Dim; ++i) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->assign_concentrated_force(mpm::NodePhase::nMixture, + REQUIRE(node->assign_concentrated_force(mpm::NodePhase::NMixture, Direction, concentrated_force, ffunction) == true); double current_time = 0.0; - node->apply_concentrated_force(mpm::NodePhase::nMixture, current_time); + node->apply_concentrated_force(mpm::NodePhase::NMixture, current_time); for (unsigned i = 0; i < Dim; ++i) { if (i == Direction) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(concentrated_force).epsilon(Tolerance)); else - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); } // Check for incorrect direction / phase const unsigned wrong_dir = 4; - REQUIRE(node->assign_concentrated_force(mpm::NodePhase::nMixture, + REQUIRE(node->assign_concentrated_force(mpm::NodePhase::NMixture, wrong_dir, concentrated_force, ffunction) == false); // Check again to ensure value hasn't been updated for (unsigned i = 0; i < Dim; ++i) { if (i == Direction) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(concentrated_force).epsilon(Tolerance)); else - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); } } @@ -335,45 +335,45 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { // Check current internal force is zero for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Update force to 10.0 REQUIRE_NOTHROW( - node->update_internal_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::nLiquid, + node->update_internal_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } // Update force to 20.0 REQUIRE_NOTHROW( - node->update_internal_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::nLiquid, + node->update_internal_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Assign force as 10.0 REQUIRE_NOTHROW( - node->update_internal_force(false, mpm::NodePhase::nMixture, force)); + node->update_internal_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_internal_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } } @@ -422,12 +422,12 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { double liquid_mass = 100.; // Update mass to 100.5 REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(solid_mass).epsilon(Tolerance)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(liquid_mass).epsilon(Tolerance)); // Check internal force @@ -436,14 +436,14 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { for (unsigned i = 0; i < force.size(); ++i) force(i) = 10. * i; // Update force to 10.0 REQUIRE_NOTHROW( - node->update_internal_force(false, mpm::NodePhase::nMixture, force)); + node->update_internal_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_internal_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); // Internal force for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(force(i)).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * force(i)).epsilon(Tolerance)); } @@ -451,13 +451,13 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { for (unsigned i = 0; i < force.size(); ++i) force(i) = 5. * i; // Update force to 10.0 REQUIRE_NOTHROW( - node->update_external_force(false, mpm::NodePhase::nMixture, force)); + node->update_external_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(force(i)).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * force(i)).epsilon(Tolerance)); } @@ -482,9 +482,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { solid_acceleration << 0.; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -492,9 +492,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { Eigen::Matrix solid_velocity = solid_acceleration * dt; Eigen::Matrix liquid_velocity = liquid_acceleration * dt; for (unsigned i = 0; i < solid_velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(solid_velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(liquid_velocity(i)).epsilon(Tolerance)); } @@ -508,9 +508,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { solid_acceleration[0] = 0.5 * solid_acceleration[0]; liquid_acceleration[0] = 0.5 * liquid_acceleration[0]; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -522,9 +522,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { solid_acceleration[0] = 0.; liquid_acceleration[0] = 0.; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -538,9 +538,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { solid_velocity[0] = 10.5; liquid_velocity[0] = 10.5; for (unsigned i = 0; i < solid_velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(solid_velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(liquid_velocity(i)).epsilon(Tolerance)); } @@ -548,19 +548,19 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { solid_acceleration[0] = 0.; liquid_acceleration[0] = 0.; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } // Exception check when mass is zero // Update mass to 0. - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nSolid, 0.)); - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nLiquid, 0.)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NSolid, 0.)); + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NLiquid, 0.)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.).epsilon(Tolerance)); REQUIRE(node->compute_acceleration_velocity_twophase_explicit(dt) == false); @@ -573,82 +573,82 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { // Check initial momentum for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Check update momentum to 10 REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(true, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(true, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Check update momentum to 20 REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(true, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(true, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(20.).epsilon(Tolerance)); } // Check assign momentum to 10 REQUIRE_NOTHROW( - node->update_momentum(false, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(false, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(false, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(false, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Check zero velocity for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Check mass double mass = 0.; - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nSolid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NSolid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nLiquid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NLiquid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); // Compute and check velocity this should throw zero mass node->compute_velocity(); mass = 100.; // Update mass to 100. - REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::nSolid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::NSolid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.).epsilon(Tolerance)); - REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::nLiquid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::NLiquid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(100.).epsilon(Tolerance)); // Compute and check velocity node->compute_velocity(); for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.1).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.1).epsilon(Tolerance)); } @@ -661,9 +661,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { Eigen::Matrix velocity; velocity << 0.1; for (unsigned i = 0; i < velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } @@ -673,9 +673,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { // Check apply constraints velocity << 10.5; for (unsigned i = 0; i < velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } } @@ -686,21 +686,21 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { for (unsigned i = 0; i < acceleration.size(); ++i) acceleration(i) = 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } - REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::NSolid, acceleration)); - REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::NLiquid, 0.5 * acceleration)); for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(5.).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(2.5).epsilon(Tolerance)); } @@ -713,9 +713,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { acceleration.resize(Dim); acceleration << 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * acceleration(i)).epsilon(Tolerance)); } @@ -725,9 +725,9 @@ TEST_CASE("Twophase Node is checked for 1D case", "[node][1D][2Phase]") { // Check apply constraints acceleration << 0.0; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } } @@ -866,135 +866,135 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { std::make_shared>(id, coords); // Check mass - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); double solid_mass = 100.5; double liquid_mass = 200.5; // Update mass to 100.5 and 200.5 REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(true, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(true, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.5).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(200.5).epsilon(Tolerance)); // Update mass to 201 and 401 REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(true, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(true, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(201.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(401.0).epsilon(Tolerance)); // Assign mass to 100 and 200 solid_mass = 100.; liquid_mass = 200.; REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(200.0).epsilon(Tolerance)); SECTION("Check nodal pressure") { // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); double pressure = 1000.7; double pore_pressure = 2000.7; // Update pressure to 1000.7 and 2000.7 - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NSolid, solid_mass * pressure)); - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NLiquid, liquid_mass * pore_pressure)); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.7).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.7).epsilon(Tolerance)); // Update pressure to 2001.4 and 4001.4 - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NSolid, solid_mass * pressure)); - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NLiquid, liquid_mass * pore_pressure)); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(2001.4).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(4001.4).epsilon(Tolerance)); // Assign pressure to 1000 and 2000 pressure = 1000.; pore_pressure = 2000.; - node->assign_pressure(mpm::NodePhase::nSolid, pressure); - node->assign_pressure(mpm::NodePhase::nLiquid, pore_pressure); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + node->assign_pressure(mpm::NodePhase::NSolid, pressure); + node->assign_pressure(mpm::NodePhase::NLiquid, pore_pressure); + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.0).epsilon(Tolerance)); // Assign mass to 0 solid_mass = 0.; liquid_mass = 0.; REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); // Try to update pressure to 2000, should throw and keep to 1000. - node->assign_pressure(mpm::NodePhase::nSolid, pressure); - node->assign_pressure(mpm::NodePhase::nLiquid, pore_pressure); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + node->assign_pressure(mpm::NodePhase::NSolid, pressure); + node->assign_pressure(mpm::NodePhase::NLiquid, pore_pressure); + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.0).epsilon(Tolerance)); // Check pressure constraints SECTION("Check nodal pressure constraints") { // Check assign pressure constraint - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nSolid, 8000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NSolid, 8000, nullptr) == true); - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nLiquid, 7000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NLiquid, 7000, nullptr) == true); // Check apply pressure constraint REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nSolid)); + node->apply_pressure_constraint(mpm::NodePhase::NSolid)); REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nLiquid)); + node->apply_pressure_constraint(mpm::NodePhase::NLiquid)); // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(8000).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(7000).epsilon(Tolerance)); } } SECTION("Check volume") { // Check volume - REQUIRE(node->volume(mpm::NodePhase::nMixture) == + REQUIRE(node->volume(mpm::NodePhase::NMixture) == Approx(0.0).epsilon(Tolerance)); double volume = 100.5; // Update volume to 100.5 REQUIRE_NOTHROW( - node->update_volume(true, mpm::NodePhase::nMixture, volume)); - REQUIRE(node->volume(mpm::NodePhase::nMixture) == + node->update_volume(true, mpm::NodePhase::NMixture, volume)); + REQUIRE(node->volume(mpm::NodePhase::NMixture) == Approx(100.5).epsilon(Tolerance)); // Update volume to 201 REQUIRE_NOTHROW( - node->update_volume(true, mpm::NodePhase::nMixture, volume)); - REQUIRE(node->volume(mpm::NodePhase::nMixture) == + node->update_volume(true, mpm::NodePhase::NMixture, volume)); + REQUIRE(node->volume(mpm::NodePhase::NMixture) == Approx(201.0).epsilon(Tolerance)); // Assign volume to 100 volume = 100.; REQUIRE_NOTHROW( - node->update_volume(false, mpm::NodePhase::nMixture, volume)); - REQUIRE(node->volume(mpm::NodePhase::nMixture) == + node->update_volume(false, mpm::NodePhase::NMixture, volume)); + REQUIRE(node->volume(mpm::NodePhase::NMixture) == Approx(100.0).epsilon(Tolerance)); } @@ -1005,45 +1005,45 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check current external force is zero for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Update force to 10.0 REQUIRE_NOTHROW( - node->update_external_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::nLiquid, + node->update_external_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } // Update force to 20.0 REQUIRE_NOTHROW( - node->update_external_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::nLiquid, + node->update_external_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Assign force as 10.0 REQUIRE_NOTHROW( - node->update_external_force(false, mpm::NodePhase::nMixture, force)); + node->update_external_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } @@ -1051,7 +1051,7 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Set external force to zero force.setZero(); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nMixture, force)); + false, mpm::NodePhase::NMixture, force)); // Concentrated force std::shared_ptr ffunction = nullptr; @@ -1059,38 +1059,38 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { const unsigned Direction = 0; // Check traction for (unsigned i = 0; i < Dim; ++i) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->assign_concentrated_force(mpm::NodePhase::nMixture, + REQUIRE(node->assign_concentrated_force(mpm::NodePhase::NMixture, Direction, concentrated_force, ffunction) == true); double current_time = 0.0; - node->apply_concentrated_force(mpm::NodePhase::nMixture, current_time); + node->apply_concentrated_force(mpm::NodePhase::NMixture, current_time); for (unsigned i = 0; i < Dim; ++i) { if (i == Direction) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(concentrated_force).epsilon(Tolerance)); else - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); } // Check for incorrect direction / phase const unsigned wrong_dir = 4; - REQUIRE(node->assign_concentrated_force(mpm::NodePhase::nMixture, + REQUIRE(node->assign_concentrated_force(mpm::NodePhase::NMixture, wrong_dir, concentrated_force, ffunction) == false); // Check again to ensure value hasn't been updated for (unsigned i = 0; i < Dim; ++i) { if (i == Direction) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(concentrated_force).epsilon(Tolerance)); else - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); } } @@ -1103,45 +1103,45 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check current internal force is zero for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Update force to 10.0 REQUIRE_NOTHROW( - node->update_internal_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::nLiquid, + node->update_internal_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } // Update force to 20.0 REQUIRE_NOTHROW( - node->update_internal_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::nLiquid, + node->update_internal_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Assign force as 10.0 REQUIRE_NOTHROW( - node->update_internal_force(false, mpm::NodePhase::nMixture, force)); + node->update_internal_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_internal_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } } @@ -1190,12 +1190,12 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { double liquid_mass = 100.; // Update mass to 100. REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(solid_mass).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(liquid_mass).epsilon(Tolerance)); // Check internal force @@ -1204,14 +1204,14 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { for (unsigned i = 0; i < force.size(); ++i) force(i) = 10. * i; // Update force to 10.0 REQUIRE_NOTHROW( - node->update_internal_force(false, mpm::NodePhase::nMixture, force)); + node->update_internal_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_internal_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); // Internal force for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(force(i)).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * force(i)).epsilon(Tolerance)); } @@ -1219,13 +1219,13 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { for (unsigned i = 0; i < force.size(); ++i) force(i) = 5. * i; // Update force to 10.0 REQUIRE_NOTHROW( - node->update_external_force(false, mpm::NodePhase::nMixture, force)); + node->update_external_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(force(i)).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * force(i)).epsilon(Tolerance)); } @@ -1250,9 +1250,9 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { solid_acceleration << 0., 0.075; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -1260,9 +1260,9 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { Eigen::Matrix solid_velocity = solid_acceleration * dt; Eigen::Matrix liquid_velocity = liquid_acceleration * dt; for (unsigned i = 0; i < solid_velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(solid_velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(liquid_velocity(i)).epsilon(Tolerance)); } @@ -1278,9 +1278,9 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { solid_velocity << 10.5, 0.03; liquid_velocity << 20.5, 1.03; for (unsigned i = 0; i < solid_velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(solid_velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(liquid_velocity(i)).epsilon(Tolerance)); } @@ -1288,9 +1288,9 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { solid_acceleration.setZero(); liquid_acceleration.setZero(); for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -1302,18 +1302,18 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { solid_acceleration << 0., 0.; liquid_acceleration << 0., 0.; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } // Exception check when mass is zero - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nSolid, 0.)); - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nLiquid, 0.)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NSolid, 0.)); + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NLiquid, 0.)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.).epsilon(Tolerance)); REQUIRE(node->compute_acceleration_velocity_twophase_explicit(dt) == false); @@ -1329,82 +1329,82 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check initial momentum for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Check update momentum to 10 REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(true, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(true, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Check update momentum to 20 REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(true, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(true, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(20.).epsilon(Tolerance)); } // Check assign momentum to 10 REQUIRE_NOTHROW( - node->update_momentum(false, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(false, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(false, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(false, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Check zero velocity for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Check mass double mass = 0.; - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nSolid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NSolid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nLiquid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NLiquid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); // Compute and check velocity this should throw zero mass node->compute_velocity(); mass = 100.; // Update mass to 100.5 - REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::nSolid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::NSolid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.).epsilon(Tolerance)); - REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::nLiquid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::NLiquid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(100.).epsilon(Tolerance)); // Compute and check velocity node->compute_velocity(); for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.1).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.1).epsilon(Tolerance)); } @@ -1413,20 +1413,20 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { for (unsigned i = 0; i < acceleration.size(); ++i) acceleration(i) = 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } - REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::NSolid, acceleration)); - REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::NLiquid, acceleration)); for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(5.).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } @@ -1439,18 +1439,18 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { Eigen::Matrix velocity; velocity << 0.1, 0.1; for (unsigned i = 0; i < velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } // Check acceleration before constraints acceleration << 5., 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } @@ -1466,17 +1466,17 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check apply constraints velocity << -12.5, 0.1; for (unsigned i = 0; i < velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } acceleration << 0., 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } } @@ -1501,35 +1501,35 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check applied velocity constraints in the global coordinates velocity << -9.583478335521184, -8.025403099849004; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } // Check that the velocity is as specified in local coordinate REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(0) == + node->velocity(mpm::NodePhase::NSolid))(0) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(0) == + node->velocity(mpm::NodePhase::NLiquid))(0) == Approx(-12.5).epsilon(Tolerance)); // Check applied constraints on acceleration in the global coordinates acceleration << -0.396139826697847, 0.472101061636807; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } // Check that the acceleration is 0 in local coordinate REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(0) == + node->acceleration(mpm::NodePhase::NSolid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(0) == + node->acceleration(mpm::NodePhase::NLiquid))(0) == Approx(0).epsilon(Tolerance)); } @@ -1555,47 +1555,47 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check applied velocity constraints in the global coordinates velocity << -14.311308834766370, 2.772442864323454; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } // Check that the velocity is as specified in local coordinate REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(0) == + node->velocity(mpm::NodePhase::NSolid))(0) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(1) == + node->velocity(mpm::NodePhase::NSolid))(1) == Approx(7.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(0) == + node->velocity(mpm::NodePhase::NLiquid))(0) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(1) == + node->velocity(mpm::NodePhase::NLiquid))(1) == Approx(7.5).epsilon(Tolerance)); // Check applied constraints on acceleration in the global coordinates acceleration << 0, 0; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } // Check that the acceleration is 0 in local coordinate REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(0) == + node->acceleration(mpm::NodePhase::NSolid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(1) == + node->acceleration(mpm::NodePhase::NSolid))(1) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(0) == + node->acceleration(mpm::NodePhase::NLiquid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(1) == + node->acceleration(mpm::NodePhase::NLiquid))(1) == Approx(0).epsilon(Tolerance)); } @@ -1611,7 +1611,7 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check apply constraints acceleration << 4., 5.; for (unsigned i = 0; i < acceleration.size(); ++i) - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } @@ -1633,14 +1633,14 @@ TEST_CASE("Twophase Node is checked for 2D case", "[node][2D][2Phase]") { // Check applied constraints on acceleration in the global coordinates acceleration << 4.905579787672637, 4.920772034660430; for (unsigned i = 0; i < Dim; ++i) - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); // Check the acceleration in local coordinate acceleration << 6.920903430595146, 0.616284167162194; for (unsigned i = 0; i < Dim; ++i) REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(i) == + node->acceleration(mpm::NodePhase::NSolid))(i) == Approx(acceleration(i)).epsilon(Tolerance)); } } @@ -1779,111 +1779,111 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { std::make_shared>(id, coords); // Check mass - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); double solid_mass = 100.5; double liquid_mass = 200.5; // Update mass to 100.5 and 200.5 REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(true, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(true, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.5).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(200.5).epsilon(Tolerance)); // Update mass to 201 and 401 REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(true, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(true, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(true, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(201.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(401.0).epsilon(Tolerance)); // Assign mass to 100 and 200 solid_mass = 100.; liquid_mass = 200.; REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(200.0).epsilon(Tolerance)); SECTION("Check nodal pressure") { // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); double pressure = 1000.7; double pore_pressure = 2000.7; // Update pressure to 1000.7 and 2000.7 - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NSolid, solid_mass * pressure)); - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NLiquid, liquid_mass * pore_pressure)); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.7).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.7).epsilon(Tolerance)); // Update pressure to 2001.4 and 4001.4 - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NSolid, solid_mass * pressure)); - REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_mass_pressure(mpm::NodePhase::NLiquid, liquid_mass * pore_pressure)); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(2001.4).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(4001.4).epsilon(Tolerance)); // Assign pressure to 1000 and 2000 pressure = 1000.; pore_pressure = 2000.; - node->assign_pressure(mpm::NodePhase::nSolid, pressure); - node->assign_pressure(mpm::NodePhase::nLiquid, pore_pressure); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + node->assign_pressure(mpm::NodePhase::NSolid, pressure); + node->assign_pressure(mpm::NodePhase::NLiquid, pore_pressure); + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.0).epsilon(Tolerance)); // Assign mass to 0 solid_mass = 0.; liquid_mass = 0.; REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); // Try to update pressure to 2000, should throw and keep to 1000. - node->assign_pressure(mpm::NodePhase::nSolid, pressure); - node->assign_pressure(mpm::NodePhase::nLiquid, pore_pressure); - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + node->assign_pressure(mpm::NodePhase::NSolid, pressure); + node->assign_pressure(mpm::NodePhase::NLiquid, pore_pressure); + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(1000.0).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(2000.0).epsilon(Tolerance)); // Check pressure constraints SECTION("Check nodal pressure constraints") { // Check assign pressure constraint - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nSolid, 8000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NSolid, 8000, nullptr) == true); - REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::nLiquid, 7000, + REQUIRE(node->assign_pressure_constraint(mpm::NodePhase::NLiquid, 7000, nullptr) == true); // Check apply pressure constraint REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nSolid)); + node->apply_pressure_constraint(mpm::NodePhase::NSolid)); REQUIRE_NOTHROW( - node->apply_pressure_constraint(mpm::NodePhase::nLiquid)); + node->apply_pressure_constraint(mpm::NodePhase::NLiquid)); // Check pressure - REQUIRE(node->pressure(mpm::NodePhase::nSolid) == + REQUIRE(node->pressure(mpm::NodePhase::NSolid) == Approx(8000).epsilon(Tolerance)); - REQUIRE(node->pressure(mpm::NodePhase::nLiquid) == + REQUIRE(node->pressure(mpm::NodePhase::NLiquid) == Approx(7000).epsilon(Tolerance)); } } @@ -1895,45 +1895,45 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check current external force is zero for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Update force to 10.0 REQUIRE_NOTHROW( - node->update_external_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::nLiquid, + node->update_external_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } // Update force to 20.0 REQUIRE_NOTHROW( - node->update_external_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::nLiquid, + node->update_external_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_external_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Assign force as 10.0 REQUIRE_NOTHROW( - node->update_external_force(false, mpm::NodePhase::nMixture, force)); + node->update_external_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } @@ -1941,7 +1941,7 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Set external force to zero force.setZero(); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nMixture, force)); + false, mpm::NodePhase::NMixture, force)); // Concentrated force std::shared_ptr ffunction = nullptr; @@ -1949,38 +1949,38 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { const unsigned Direction = 0; // Check traction for (unsigned i = 0; i < Dim; ++i) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->assign_concentrated_force(mpm::NodePhase::nMixture, + REQUIRE(node->assign_concentrated_force(mpm::NodePhase::NMixture, Direction, concentrated_force, ffunction) == true); double current_time = 0.0; - node->apply_concentrated_force(mpm::NodePhase::nMixture, current_time); + node->apply_concentrated_force(mpm::NodePhase::NMixture, current_time); for (unsigned i = 0; i < Dim; ++i) { if (i == Direction) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(concentrated_force).epsilon(Tolerance)); else - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); } // Check for incorrect direction / phase const unsigned wrong_dir = 4; - REQUIRE(node->assign_concentrated_force(mpm::NodePhase::nMixture, + REQUIRE(node->assign_concentrated_force(mpm::NodePhase::NMixture, wrong_dir, concentrated_force, ffunction) == false); // Check again to ensure value hasn't been updated for (unsigned i = 0; i < Dim; ++i) { if (i == Direction) - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(concentrated_force).epsilon(Tolerance)); else - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); } } @@ -1993,45 +1993,45 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check current internal force is zero for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Update force to 10.0 REQUIRE_NOTHROW( - node->update_internal_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::nLiquid, + node->update_internal_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } // Update force to 20.0 REQUIRE_NOTHROW( - node->update_internal_force(true, mpm::NodePhase::nMixture, force)); - REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::nLiquid, + node->update_internal_force(true, mpm::NodePhase::NMixture, force)); + REQUIRE_NOTHROW(node->update_internal_force(true, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Assign force as 10.0 REQUIRE_NOTHROW( - node->update_internal_force(false, mpm::NodePhase::nMixture, force)); + node->update_internal_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_internal_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } } @@ -2080,12 +2080,12 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { double liquid_mass = 100.; // Update mass to 100. REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nSolid, solid_mass)); + node->update_mass(false, mpm::NodePhase::NSolid, solid_mass)); REQUIRE_NOTHROW( - node->update_mass(false, mpm::NodePhase::nLiquid, liquid_mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + node->update_mass(false, mpm::NodePhase::NLiquid, liquid_mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(solid_mass).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(liquid_mass).epsilon(Tolerance)); // Check internal force @@ -2094,14 +2094,14 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { for (unsigned i = 0; i < force.size(); ++i) force(i) = 10. * i; // Update force to 10.0 REQUIRE_NOTHROW( - node->update_internal_force(false, mpm::NodePhase::nMixture, force)); + node->update_internal_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_internal_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); // Internal force for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->internal_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NMixture)(i) == Approx(force(i)).epsilon(Tolerance)); - REQUIRE(node->internal_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->internal_force(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * force(i)).epsilon(Tolerance)); } @@ -2109,13 +2109,13 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { for (unsigned i = 0; i < force.size(); ++i) force(i) = 5. * i; // Update force to 10.0 REQUIRE_NOTHROW( - node->update_external_force(false, mpm::NodePhase::nMixture, force)); + node->update_external_force(false, mpm::NodePhase::NMixture, force)); REQUIRE_NOTHROW(node->update_external_force( - false, mpm::NodePhase::nLiquid, 0.5 * force)); + false, mpm::NodePhase::NLiquid, 0.5 * force)); for (unsigned i = 0; i < force.size(); ++i) { - REQUIRE(node->external_force(mpm::NodePhase::nMixture)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NMixture)(i) == Approx(force(i)).epsilon(Tolerance)); - REQUIRE(node->external_force(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->external_force(mpm::NodePhase::NLiquid)(i) == Approx(0.5 * force(i)).epsilon(Tolerance)); } @@ -2141,9 +2141,9 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check acceleration for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -2151,9 +2151,9 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { Eigen::Matrix solid_velocity = solid_acceleration * dt; Eigen::Matrix liquid_velocity = liquid_acceleration * dt; for (unsigned i = 0; i < solid_velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(solid_velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(liquid_velocity(i)).epsilon(Tolerance)); } @@ -2171,9 +2171,9 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { solid_velocity << 10.5, 0.03, 5.13; liquid_velocity << 20.5, 1.03, 7.13; for (unsigned i = 0; i < solid_velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(solid_velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(liquid_velocity(i)).epsilon(Tolerance)); } @@ -2181,9 +2181,9 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { solid_acceleration.setZero(); liquid_acceleration.setZero(); for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -2195,9 +2195,9 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { solid_acceleration << 0., 0., 0.; liquid_acceleration << 0., 0., 0.; for (unsigned i = 0; i < solid_acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(solid_acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(liquid_acceleration(i)).epsilon(Tolerance)); } @@ -2208,11 +2208,11 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { true); // Exception check when mass is zero - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nSolid, 0.)); - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nLiquid, 0.)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NSolid, 0.)); + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NLiquid, 0.)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.).epsilon(Tolerance)); REQUIRE(node->compute_acceleration_velocity_twophase_explicit(dt) == false); @@ -2228,90 +2228,90 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check initial momentum for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Check update momentum to 10 REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(true, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(true, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Check update momentum to 20 REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(true, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(true, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(true, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(20.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(20.).epsilon(Tolerance)); } // Check assign momentum to 10 REQUIRE_NOTHROW( - node->update_momentum(false, mpm::NodePhase::nSolid, momentum)); + node->update_momentum(false, mpm::NodePhase::NSolid, momentum)); REQUIRE_NOTHROW( - node->update_momentum(false, mpm::NodePhase::nLiquid, momentum)); + node->update_momentum(false, mpm::NodePhase::NLiquid, momentum)); for (unsigned i = 0; i < momentum.size(); ++i) { - REQUIRE(node->momentum(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NSolid)(i) == Approx(10.).epsilon(Tolerance)); - REQUIRE(node->momentum(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->momentum(mpm::NodePhase::NLiquid)(i) == Approx(10.).epsilon(Tolerance)); } // Check zero velocity for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Check mass double mass = 0.; - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nSolid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NSolid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(0.0).epsilon(Tolerance)); - REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::nLiquid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE_NOTHROW(node->update_mass(false, mpm::NodePhase::NLiquid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(0.0).epsilon(Tolerance)); // Compute and check velocity this should throw zero mass node->compute_velocity(); mass = 100.; // Update mass to 100.5 - REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::nSolid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nSolid) == + REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::NSolid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NSolid) == Approx(100.).epsilon(Tolerance)); - REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::nLiquid, mass)); - REQUIRE(node->mass(mpm::NodePhase::nLiquid) == + REQUIRE_NOTHROW(node->update_mass(true, mpm::NodePhase::NLiquid, mass)); + REQUIRE(node->mass(mpm::NodePhase::NLiquid) == Approx(100.).epsilon(Tolerance)); // Check zero velocity for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } // Compute and check velocity node->compute_velocity(); for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(0.1).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(0.1).epsilon(Tolerance)); } @@ -2320,20 +2320,20 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { for (unsigned i = 0; i < acceleration.size(); ++i) acceleration(i) = 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(0.).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(0.).epsilon(Tolerance)); } - REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::nSolid, + REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::NSolid, acceleration)); - REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::nLiquid, + REQUIRE_NOTHROW(node->update_acceleration(true, mpm::NodePhase::NLiquid, acceleration)); for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(5.).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(5.).epsilon(Tolerance)); } @@ -2341,18 +2341,18 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { Eigen::Matrix velocity; velocity << 0.1, 0.1, 0.1; for (unsigned i = 0; i < velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } // Check acceleration before constraints acceleration << 5., 5., 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } SECTION("Check Cartesian velocity constraints") { @@ -2367,17 +2367,17 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check apply constraints velocity << -12.5, 0.1, 0.1; for (unsigned i = 0; i < velocity.size(); ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } acceleration << 0., 5., 5.; for (unsigned i = 0; i < acceleration.size(); ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } } @@ -2403,47 +2403,47 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check apply constraints velocity << -14.5068204271, -0.1432759442, 1.4260971922; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } // Check that the velocity is as specified in local coordinate REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(0) == + node->velocity(mpm::NodePhase::NSolid))(0) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(1) == + node->velocity(mpm::NodePhase::NSolid))(1) == Approx(7.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(0) == + node->velocity(mpm::NodePhase::NLiquid))(0) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(1) == + node->velocity(mpm::NodePhase::NLiquid))(1) == Approx(7.5).epsilon(Tolerance)); // Check applied constraints on acceleration in the global coordinates acceleration << 0.1998888554, -1.1336260315, 1.9937880031; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } // Check that the acceleration is 0 in local coordinate REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(0) == + node->acceleration(mpm::NodePhase::NSolid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(1) == + node->acceleration(mpm::NodePhase::NSolid))(1) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(0) == + node->acceleration(mpm::NodePhase::NLiquid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(1) == + node->acceleration(mpm::NodePhase::NLiquid))(1) == Approx(0).epsilon(Tolerance)); } @@ -2471,53 +2471,53 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check apply constraints velocity << 13.351984588153375, -5.717804716697730, 10.572663655835457; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->velocity(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NSolid)(i) == Approx(velocity(i)).epsilon(Tolerance)); - REQUIRE(node->velocity(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->velocity(mpm::NodePhase::NLiquid)(i) == Approx(velocity(i)).epsilon(Tolerance)); } // Check that the velocity is as specified in local coordinate REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(0) == + node->velocity(mpm::NodePhase::NSolid))(0) == Approx(10.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(1) == + node->velocity(mpm::NodePhase::NSolid))(1) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nSolid))(2) == + node->velocity(mpm::NodePhase::NSolid))(2) == Approx(7.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(0) == + node->velocity(mpm::NodePhase::NLiquid))(0) == Approx(10.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(1) == + node->velocity(mpm::NodePhase::NLiquid))(1) == Approx(-12.5).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->velocity(mpm::NodePhase::nLiquid))(2) == + node->velocity(mpm::NodePhase::NLiquid))(2) == Approx(7.5).epsilon(Tolerance)); // Check apply constraints acceleration << 0, 0, 0; for (unsigned i = 0; i < Dim; ++i) { - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); - REQUIRE(node->acceleration(mpm::NodePhase::nLiquid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NLiquid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } // Check that the acceleration is 0 in local coordinate REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(0) == + node->acceleration(mpm::NodePhase::NSolid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(1) == + node->acceleration(mpm::NodePhase::NSolid))(1) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(0) == + node->acceleration(mpm::NodePhase::NLiquid))(0) == Approx(0).epsilon(Tolerance)); REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nLiquid))(1) == + node->acceleration(mpm::NodePhase::NLiquid))(1) == Approx(0).epsilon(Tolerance)); } @@ -2533,7 +2533,7 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check apply constraints acceleration << 3.939339828220179, 3.939339828220179, 5.; for (unsigned i = 0; i < acceleration.size(); ++i) - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); } @@ -2556,14 +2556,14 @@ TEST_CASE("Twophase Node is checked for 3D case", "[node][3D][2Phase]") { // Check applied constraints on acceleration in the global coordinates acceleration << 4.602895052828914, 4.492575657560740, 4.751301246937935; for (unsigned i = 0; i < Dim; ++i) - REQUIRE(node->acceleration(mpm::NodePhase::nSolid)(i) == + REQUIRE(node->acceleration(mpm::NodePhase::NSolid)(i) == Approx(acceleration(i)).epsilon(Tolerance)); // Check the acceleration in local coordinate acceleration << 6.878925666702865, 3.365244416454818, 2.302228080558999; for (unsigned i = 0; i < Dim; ++i) REQUIRE((inverse_rotation_matrix * - node->acceleration(mpm::NodePhase::nSolid))(i) == + node->acceleration(mpm::NodePhase::NSolid))(i) == Approx(acceleration(i)).epsilon(Tolerance)); } } diff --git a/tests/particles/particle_twophase_test.cc b/tests/particles/particle_twophase_test.cc index 58d71ee69..024efbe5a 100644 --- a/tests/particles/particle_twophase_test.cc +++ b/tests/particles/particle_twophase_test.cc @@ -913,12 +913,12 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", // Check nodal mass for (unsigned i = 0; i < nodes.size(); ++i) { // Solid phase - REQUIRE(nodes.at(i)->mass(mpm::NodePhase::nSolid) == + REQUIRE(nodes.at(i)->mass(mpm::NodePhase::NSolid) == Approx(nodal_mass.at(i) * (1 - particle->porosity())) .epsilon(Tolerance)); // Liquid phase REQUIRE( - nodes.at(i)->mass(mpm::NodePhase::nLiquid) == + nodes.at(i)->mass(mpm::NodePhase::NLiquid) == Approx(nodal_mass.at(i) * particle->porosity()).epsilon(Tolerance)); } @@ -937,11 +937,11 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", for (unsigned i = 0; i < nodal_momentum.rows(); ++i) for (unsigned j = 0; j < nodal_momentum.cols(); ++j) { // Solid phase - REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::nSolid)(j) == + REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::NSolid)(j) == Approx(nodal_momentum(i, j) * (1 - particle->porosity())) .epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::nLiquid)(j) == + REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::NLiquid)(j) == Approx(nodal_momentum(i, j) * particle->porosity()) .epsilon(Tolerance)); } @@ -957,10 +957,10 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", for (unsigned i = 0; i < nodal_velocity.rows(); ++i) for (unsigned j = 0; j < nodal_velocity.cols(); ++j) { // Solid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nSolid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NSolid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nLiquid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NLiquid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); } @@ -974,11 +974,11 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", for (unsigned i = 0; i < nodes.size(); ++i) { // Solid phase REQUIRE_NOTHROW(nodes.at(i)->update_momentum( - false, mpm::NodePhase::nSolid, + false, mpm::NodePhase::NSolid, nodal_momentum.row(i) * (1 - particle->porosity()))); // Liquid phase REQUIRE_NOTHROW(nodes.at(i)->update_momentum( - false, mpm::NodePhase::nLiquid, + false, mpm::NodePhase::NLiquid, nodal_momentum.row(i) * particle->porosity())); } @@ -995,10 +995,10 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", for (unsigned i = 0; i < nodal_velocity.rows(); ++i) for (unsigned j = 0; j < nodal_velocity.cols(); ++j) { // Solid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nSolid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NSolid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nLiquid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NLiquid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); } @@ -1069,7 +1069,7 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", // Check nodal body force for (unsigned i = 0; i < body_force.rows(); ++i) for (unsigned j = 0; j < body_force.cols(); ++j) - REQUIRE(nodes[i]->external_force(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->external_force(mpm::NodePhase::NSolid)[j] == Approx(body_force(i, j)).epsilon(Tolerance)); // Check traction force @@ -1101,7 +1101,7 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", // Check nodal traction force for (unsigned i = 0; i < traction_force.rows(); ++i) for (unsigned j = 0; j < traction_force.cols(); ++j) - REQUIRE(nodes[i]->external_force(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->external_force(mpm::NodePhase::NSolid)[j] == Approx(traction_force(i, j)).epsilon(Tolerance)); // Reset traction particle->assign_traction(direction, @@ -1111,7 +1111,7 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", // Check nodal external force for (unsigned i = 0; i < traction_force.rows(); ++i) for (unsigned j = 0; j < traction_force.cols(); ++j) - REQUIRE(nodes[i]->external_force(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->external_force(mpm::NodePhase::NSolid)[j] == Approx(body_force(i, j)).epsilon(Tolerance)); // Internal force @@ -1130,7 +1130,7 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", // Check nodal internal force for (unsigned i = 0; i < internal_force.rows(); ++i) for (unsigned j = 0; j < internal_force.cols(); ++j) - REQUIRE(nodes[i]->internal_force(mpm::NodePhase::nMixture)[j] == + REQUIRE(nodes[i]->internal_force(mpm::NodePhase::NMixture)[j] == Approx(internal_force(i, j)).epsilon(Tolerance)); // Internal force @@ -1170,10 +1170,10 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", for (unsigned i = 0; i < nodal_velocity.rows(); ++i) for (unsigned j = 0; j < nodal_velocity.cols(); ++j) { // Solid phase - REQUIRE(nodes[i]->velocity(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->velocity(mpm::NodePhase::NSolid)[j] == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes[i]->velocity(mpm::NodePhase::nLiquid)[j] == + REQUIRE(nodes[i]->velocity(mpm::NodePhase::NLiquid)[j] == Approx(nodal_liquid_velocity(i, j)).epsilon(Tolerance)); } @@ -1194,10 +1194,10 @@ TEST_CASE("TwoPhase Particle is checked for 2D case", for (unsigned i = 0; i < nodal_acceleration.rows(); ++i) for (unsigned j = 0; j < nodal_acceleration.cols(); ++j) { // Solid phase - REQUIRE(nodes[i]->acceleration(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->acceleration(mpm::NodePhase::NSolid)[j] == Approx(nodal_acceleration(i, j)).epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes[i]->acceleration(mpm::NodePhase::nLiquid)[j] == + REQUIRE(nodes[i]->acceleration(mpm::NodePhase::NLiquid)[j] == Approx(nodal_liquid_acceleration(i, j)).epsilon(Tolerance)); } // Approx(nodal_velocity(i, j) / dt).epsilon(Tolerance)); @@ -2350,11 +2350,11 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal mass for (unsigned i = 0; i < nodes.size(); ++i) { // Solid phase - REQUIRE(nodes.at(i)->mass(mpm::NodePhase::nSolid) == + REQUIRE(nodes.at(i)->mass(mpm::NodePhase::NSolid) == Approx(nodal_mass.at(i) * (1 - particle->porosity()))); // Liquid phase REQUIRE( - nodes.at(i)->mass(mpm::NodePhase::nLiquid) == + nodes.at(i)->mass(mpm::NodePhase::NLiquid) == Approx(nodal_mass.at(i) * particle->porosity()).epsilon(Tolerance)); } @@ -2378,11 +2378,11 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", for (unsigned i = 0; i < nodal_momentum.rows(); ++i) for (unsigned j = 0; j < nodal_momentum.cols(); ++j) { // Solid phase - REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::nSolid)(j) == + REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::NSolid)(j) == Approx(nodal_momentum(i, j) * (1 - particle->porosity())) .epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::nLiquid)(j) == + REQUIRE(nodes.at(i)->momentum(mpm::NodePhase::NLiquid)(j) == Approx(nodal_momentum(i, j) * particle->porosity()) .epsilon(Tolerance)); } @@ -2403,10 +2403,10 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", for (unsigned i = 0; i < nodal_velocity.rows(); ++i) for (unsigned j = 0; j < nodal_velocity.cols(); ++j) { // Solid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nSolid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NSolid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nLiquid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NLiquid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); } @@ -2424,11 +2424,11 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", for (unsigned i = 0; i < nodes.size(); ++i) { // Solid phase REQUIRE_NOTHROW(nodes.at(i)->update_momentum( - false, mpm::NodePhase::nSolid, + false, mpm::NodePhase::NSolid, nodal_momentum.row(i) * (1 - particle->porosity()))); // Liquid phase REQUIRE_NOTHROW(nodes.at(i)->update_momentum( - false, mpm::NodePhase::nLiquid, + false, mpm::NodePhase::NLiquid, nodal_momentum.row(i) * particle->porosity())); } @@ -2449,10 +2449,10 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", for (unsigned i = 0; i < nodal_velocity.rows(); ++i) for (unsigned j = 0; j < nodal_velocity.cols(); ++j) { // Solid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nSolid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NSolid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); // Liquid phase - REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::nLiquid)(j) == + REQUIRE(nodes.at(i)->velocity(mpm::NodePhase::NLiquid)(j) == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); } @@ -2527,7 +2527,7 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal body force for (unsigned i = 0; i < body_force.rows(); ++i) for (unsigned j = 0; j < body_force.cols(); ++j) - REQUIRE(nodes[i]->external_force(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->external_force(mpm::NodePhase::NSolid)[j] == Approx(body_force(i, j)).epsilon(Tolerance)); // Check traction force @@ -2562,7 +2562,7 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal traction force for (unsigned i = 0; i < traction_force.rows(); ++i) for (unsigned j = 0; j < traction_force.cols(); ++j) - REQUIRE(nodes[i]->external_force(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->external_force(mpm::NodePhase::NSolid)[j] == Approx(traction_force(i, j)).epsilon(Tolerance)); // Reset traction particle->assign_traction(direction, @@ -2572,7 +2572,7 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal external force for (unsigned i = 0; i < traction_force.rows(); ++i) for (unsigned j = 0; j < traction_force.cols(); ++j) - REQUIRE(nodes[i]->external_force(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->external_force(mpm::NodePhase::NSolid)[j] == Approx(body_force(i, j)).epsilon(Tolerance)); // Internal force @@ -2595,12 +2595,12 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal internal force for (unsigned i = 0; i < internal_force.rows(); ++i) for (unsigned j = 0; j < internal_force.cols(); ++j) - REQUIRE(nodes[i]->internal_force(mpm::NodePhase::nMixture)[j] == + REQUIRE(nodes[i]->internal_force(mpm::NodePhase::NMixture)[j] == Approx(internal_force(i, j)).epsilon(Tolerance)); // Calculate nodal acceleration and velocity for (const auto& node : nodes) - node->compute_acceleration_velocity(mpm::NodePhase::nSolid, dt); + node->compute_acceleration_velocity(mpm::NodePhase::NSolid, dt); // Check nodal velocity // clang-format off @@ -2616,7 +2616,7 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal velocity for (unsigned i = 0; i < nodal_velocity.rows(); ++i) for (unsigned j = 0; j < nodal_velocity.cols(); ++j) - REQUIRE(nodes[i]->velocity(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->velocity(mpm::NodePhase::NSolid)[j] == Approx(nodal_velocity(i, j)).epsilon(Tolerance)); // Check nodal acceleration @@ -2634,7 +2634,7 @@ TEST_CASE("TwoPhase Particle is checked for 3D case", // Check nodal acceleration for (unsigned i = 0; i < nodal_acceleration.rows(); ++i) for (unsigned j = 0; j < nodal_acceleration.cols(); ++j) - REQUIRE(nodes[i]->acceleration(mpm::NodePhase::nSolid)[j] == + REQUIRE(nodes[i]->acceleration(mpm::NodePhase::NSolid)[j] == Approx(nodal_acceleration(i, j)).epsilon(Tolerance)); // Approx(nodal_velocity(i, j) / dt).epsilon(Tolerance));