Skip to content

Commit

Permalink
Merge pull request #17 from alanw0/master
Browse files Browse the repository at this point in the history
Basic unit-test infrastructure.
  • Loading branch information
spdomin authored Oct 14, 2016
2 parents 636a674 + 2571a9c commit d54ac14
Show file tree
Hide file tree
Showing 6 changed files with 219 additions and 0 deletions.
11 changes: 11 additions & 0 deletions CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -98,6 +98,8 @@ add_library (nalu ${SOURCE} ${HEADER})
target_link_libraries(nalu ${Trilinos_LIBRARIES})
target_link_libraries(nalu ${YAML_LIBRARY})

file (GLOB UNIT_TESTS_SOURCES unit_tests/*.C)

set(nalu_ex_name "naluX")
message("CMAKE_BUILD_TYPE = ${CMAKE_BUILD_TYPE}")
if (CMAKE_BUILD_TYPE STREQUAL "DEBUG")
Expand All @@ -108,4 +110,13 @@ endif()

add_executable(${nalu_ex_name} nalu.C)
target_link_libraries(${nalu_ex_name} nalu)

set(utest_ex_name "unittestX")
if (CMAKE_BUILD_TYPE STREQUAL "DEBUG")
set(utest_ex_name "unittestXd")
endif()

add_executable(${utest_ex_name} unit_tests.C ${UNIT_TESTS_SOURCES})
target_link_libraries(${utest_ex_name} nalu)

MESSAGE("\nAnd CMake says...:")
31 changes: 31 additions & 0 deletions unit_tests.C
Original file line number Diff line number Diff line change
@@ -0,0 +1,31 @@
/*------------------------------------------------------------------------*/
/* Copyright 2014 Sandia Corporation. */
/* This software is released under the license detailed */
/* in the file, LICENSE, which is located in the top-level Nalu */
/* directory structure */
/*------------------------------------------------------------------------*/

#include <gtest/gtest.h> // for InitGoogleTest, etc
#include <mpi.h> // for MPI_Comm_rank, MPI_Finalize, etc

#include "include/NaluEnv.h"

int gl_argc = 0;
char** gl_argv = 0;

int main(int argc, char **argv)
{
MPI_Init(&argc, &argv);
//NaluEnv will call MPI_Finalize for us.
sierra::nalu::NaluEnv::self();

testing::InitGoogleTest(&argc, argv);

gl_argc = argc;
gl_argv = argv;

int returnVal = RUN_ALL_TESTS();

return returnVal;
}

71 changes: 71 additions & 0 deletions unit_tests/UnitTest1ElemCoordCheck.C
Original file line number Diff line number Diff line change
@@ -0,0 +1,71 @@
#include <gtest/gtest.h>
#include <limits>

#include <stk_util/parallel/Parallel.hpp>
#include <stk_mesh/base/MetaData.hpp>
#include <stk_mesh/base/BulkData.hpp>
#include <stk_mesh/base/Bucket.hpp>
#include <stk_mesh/base/CoordinateSystems.hpp>
#include <stk_mesh/base/FieldBase.hpp>
#include <stk_mesh/base/Field.hpp>
#include <stk_mesh/base/GetEntities.hpp>

#include "UnitTestUtils.h"

namespace {

unsigned count_locally_owned_elems(const stk::mesh::BulkData& bulk)
{
return stk::mesh::count_selected_entities(bulk.mesh_meta_data().locally_owned_part(),
bulk.buckets(stk::topology::ELEM_RANK));
}

void verify_elems_are_unit_cubes(const stk::mesh::BulkData& bulk)
{
typedef stk::mesh::Field<double,stk::mesh::Cartesian> CoordFieldType;
CoordFieldType* coordField = bulk.mesh_meta_data().get_field<CoordFieldType>(stk::topology::NODE_RANK, "coordinates");
EXPECT_TRUE(coordField != nullptr);

stk::mesh::EntityVector elems;
stk::mesh::get_entities(bulk, stk::topology::ELEM_RANK, elems);

const double tolerance = std::numeric_limits<double>::epsilon();

for(stk::mesh::Entity elem : elems) {
double minX = 999.99, minY = 999.99, minZ = 999.99;
double maxX = 0, maxY = 0, maxZ = 0;
const stk::mesh::Entity* nodes = bulk.begin_nodes(elem);
unsigned numNodes = bulk.num_nodes(elem);
for(unsigned i=0; i<numNodes; ++i) {
double* nodeCoords = stk::mesh::field_data(*coordField, nodes[i]);
minX = std::min(minX, nodeCoords[0]);
minY = std::min(minY, nodeCoords[1]);
minZ = std::min(minZ, nodeCoords[2]);
maxX = std::max(maxX, nodeCoords[0]);
maxY = std::max(maxY, nodeCoords[1]);
maxZ = std::max(maxZ, nodeCoords[2]);
}

EXPECT_NEAR(1.0, (maxX-minX), tolerance);
EXPECT_NEAR(1.0, (maxY-minY), tolerance);
EXPECT_NEAR(1.0, (maxZ-minZ), tolerance);
}
}

}//namespace

TEST(Basic, CheckCoords1Elem)
{
stk::ParallelMachine comm = MPI_COMM_WORLD;

unsigned spatialDimension = 3;
stk::mesh::MetaData meta(spatialDimension);
stk::mesh::BulkData bulk(meta, comm);

unit_test_utils::fill_mesh_1_elem_per_proc_hex8(bulk);

EXPECT_EQ(1u, count_locally_owned_elems(bulk));

verify_elems_are_unit_cubes(bulk);
}

76 changes: 76 additions & 0 deletions unit_tests/UnitTestHexSCVDeterminant.C
Original file line number Diff line number Diff line change
@@ -0,0 +1,76 @@
#include <gtest/gtest.h>
#include <limits>

#include <stk_util/parallel/Parallel.hpp>
#include <stk_mesh/base/MetaData.hpp>
#include <stk_mesh/base/BulkData.hpp>
#include <stk_mesh/base/Bucket.hpp>
#include <stk_mesh/base/CoordinateSystems.hpp>
#include <stk_mesh/base/FieldBase.hpp>
#include <stk_mesh/base/Field.hpp>
#include <stk_mesh/base/GetEntities.hpp>

#include <master_element/MasterElement.h>

#include "UnitTestUtils.h"

namespace {

void check_HexSCV_determinant(const stk::mesh::BulkData& bulk)
{
typedef stk::mesh::Field<double,stk::mesh::Cartesian> CoordFieldType;
CoordFieldType* coordField = bulk.mesh_meta_data().get_field<CoordFieldType>(stk::topology::NODE_RANK, "coordinates");
EXPECT_TRUE(coordField != nullptr);

stk::mesh::EntityVector elems;
stk::mesh::get_entities(bulk, stk::topology::ELEM_RANK, elems);

const double tolerance = std::numeric_limits<double>::epsilon();

stk::topology hex8 = stk::topology::HEX_8;

const unsigned spatialDim = bulk.mesh_meta_data().spatial_dimension();
std::vector<double> hex8_node_coords(hex8.num_nodes()*spatialDim, 0.0);
std::vector<double> hex8_scvolumes(hex8.num_nodes(), 0.0);

sierra::nalu::HexSCV hexSCV;
double error[1] = {0};
for(stk::mesh::Entity elem : elems) {
EXPECT_EQ(hex8, bulk.bucket(elem).topology());

const stk::mesh::Entity* nodes = bulk.begin_nodes(elem);
unsigned numNodes = bulk.num_nodes(elem);
unsigned counter = 0;
for(unsigned i=0; i<numNodes; ++i) {
double* nodeCoords = stk::mesh::field_data(*coordField, nodes[i]);

for(unsigned d=0; d<spatialDim; ++d) {
hex8_node_coords[counter++] = nodeCoords[d];
}
}

hexSCV.determinant(1, hex8_node_coords.data(), hex8_scvolumes.data(), error);

EXPECT_EQ(0, error[0]);

for(unsigned i=0; i<hex8.num_nodes(); ++i) {
EXPECT_NEAR(0.125, hex8_scvolumes[i], tolerance);
}
}
}

}//namespace

TEST(HexSCV, determinant)
{
stk::ParallelMachine comm = MPI_COMM_WORLD;

unsigned spatialDimension = 3;
stk::mesh::MetaData meta(spatialDimension);
stk::mesh::BulkData bulk(meta, comm);

unit_test_utils::fill_mesh_1_elem_per_proc_hex8(bulk);

check_HexSCV_determinant(bulk);
}

21 changes: 21 additions & 0 deletions unit_tests/UnitTestUtils.C
Original file line number Diff line number Diff line change
@@ -0,0 +1,21 @@
#include <string>

#include <stk_mesh/base/BulkData.hpp>
#include <stk_io/StkMeshIoBroker.hpp>

namespace unit_test_utils {

void fill_mesh_1_elem_per_proc_hex8(stk::mesh::BulkData& bulk)
{
int nprocs = bulk.parallel_size();
std::string meshSpec = "generated:1x1x"+std::to_string(nprocs);

stk::io::StkMeshIoBroker io(bulk.parallel());
io.set_bulk_data(bulk);
io.add_mesh_database(meshSpec, stk::io::READ_MESH);
io.create_input_mesh();
io.populate_bulk_data();
}

}

9 changes: 9 additions & 0 deletions unit_tests/UnitTestUtils.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,9 @@

#include <stk_mesh/base/BulkData.hpp>

namespace unit_test_utils {

void fill_mesh_1_elem_per_proc_hex8(stk::mesh::BulkData& bulk);

}

0 comments on commit d54ac14

Please sign in to comment.