From 3e2ecd16fd4f18b62d209f179425a60f7d0d0e2f Mon Sep 17 00:00:00 2001 From: Richard Christie Date: Mon, 4 Jul 2022 14:14:53 +1200 Subject: [PATCH] Add FieldMeshIntegral get set mesh API For #230. Note FieldMeshIntegralSquares has FieldMeshIntegral as a base class so these new API also work for it. --- CHANGELOG.txt | 3 + CMakeLists.txt | 2 +- src/api/opencmiss/zinc/fieldmeshoperators.h | 18 +++++ src/api/opencmiss/zinc/fieldmeshoperators.hpp | 10 +++ .../computed_field_mesh_operators.cpp | 45 ++++++++++- tests/fieldmodule/fieldmeshoperators.cpp | 74 +++++++++++++++++++ 6 files changed, 149 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.txt b/CHANGELOG.txt index b95e5c7ca..cf4fe28ee 100644 --- a/CHANGELOG.txt +++ b/CHANGELOG.txt @@ -1,5 +1,8 @@ CHANGE LOG: OpenCMISS-Zinc Library +v3.8.0 +Add missing type-specific field API. + v3.7.0 Support optimisation at time for time-varying fields. Add Region change notifier, for tree structure changes. diff --git a/CMakeLists.txt b/CMakeLists.txt index a58545e07..4ed2c49db 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -17,7 +17,7 @@ if(POLICY CMP0086) endif() # This is the project name and shows up in IDEs -project(Zinc VERSION 3.7.0 LANGUAGES C CXX) +project(Zinc VERSION 3.8.0 LANGUAGES C CXX) # Set this to the development release or release candidate version. set(Zinc_DEVELOPER_VERSION ".dev0") diff --git a/src/api/opencmiss/zinc/fieldmeshoperators.h b/src/api/opencmiss/zinc/fieldmeshoperators.h index 9a0820c52..1dd6b8401 100644 --- a/src/api/opencmiss/zinc/fieldmeshoperators.h +++ b/src/api/opencmiss/zinc/fieldmeshoperators.h @@ -87,6 +87,24 @@ ZINC_C_INLINE cmzn_field_id cmzn_field_mesh_integral_base_cast( ZINC_API int cmzn_field_mesh_integral_destroy( cmzn_field_mesh_integral_id *mesh_integral_field_address); +/** + * Get the mesh integrated over. + * @param mesh_integral_field Handle to the mesh integral field to query. + * @return Handle to the mesh, or NULL/invalid handle if invalid field. + */ +ZINC_API cmzn_mesh_id cmzn_field_mesh_integral_get_mesh( + cmzn_field_mesh_integral_id mesh_integral_field); + +/** + * Set the mesh integrated over. + * @param mesh_integral_field Handle to the mesh integral field to modify. + * @param mesh The mesh to integrate over. + * @return Result OK on success, or ERROR_ARGUMENT if invalid mesh or mesh is + * not from this region. + */ +ZINC_API int cmzn_field_mesh_integral_set_mesh( + cmzn_field_mesh_integral_id mesh_integral_field, cmzn_mesh_id mesh); + /** * Get the numbers of quadrature points in each element axis. * @see cmzn_field_mesh_integral_set_numbers_of_points diff --git a/src/api/opencmiss/zinc/fieldmeshoperators.hpp b/src/api/opencmiss/zinc/fieldmeshoperators.hpp index a4b3cf77d..14853ab6a 100644 --- a/src/api/opencmiss/zinc/fieldmeshoperators.hpp +++ b/src/api/opencmiss/zinc/fieldmeshoperators.hpp @@ -36,6 +36,16 @@ class FieldMeshIntegral : public Field Field(reinterpret_cast(field_mesh_integral_id)) { } + Mesh getMesh() const + { + return Mesh(cmzn_field_mesh_integral_get_mesh(this->getDerivedId())); + } + + int setMesh(const Mesh& mesh) const + { + return cmzn_field_mesh_integral_set_mesh(this->getDerivedId(), mesh.getId()); + } + int getNumbersOfPoints(int valuesCount, int *valuesOut) const { return cmzn_field_mesh_integral_get_numbers_of_points(getDerivedId(), diff --git a/src/computed_field/computed_field_mesh_operators.cpp b/src/computed_field/computed_field_mesh_operators.cpp index 1ce5bf7b8..5fe381a6b 100644 --- a/src/computed_field/computed_field_mesh_operators.cpp +++ b/src/computed_field/computed_field_mesh_operators.cpp @@ -116,12 +116,31 @@ class Computed_field_mesh_integral : public Computed_field_core char* get_command_string(); - cmzn_mesh_id getMesh() + /** @return Non-accessed mesh */ + cmzn_mesh_id getMesh() const { return mesh; } - int getNumbersOfPoints(int valuesCount, int *values) + /** @return Result OK on success, ERROR_ARGUMENT if mesh missing or from wrong region */ + int setMesh(cmzn_mesh *meshIn) + { + if ((!meshIn) || (cmzn_mesh_get_region_internal(meshIn) != this->field->getRegion())) + { + display_message(ERROR_MESSAGE, "FieldMeshIntegral setMesh: Invalid mesh"); + return CMZN_ERROR_ARGUMENT; + } + if (!cmzn_mesh_match(this->mesh, meshIn)) + { + cmzn_mesh *oldMesh = this->mesh; + this->mesh = cmzn_mesh_access(meshIn); + cmzn_mesh_destroy(&oldMesh); + this->field->setChanged(); + } + return CMZN_OK; + } + + int getNumbersOfPoints(int valuesCount, int *values) const { if ((0 == valuesCount) || ((0 < valuesCount) && values)) { @@ -727,6 +746,28 @@ inline Computed_field_mesh_integral *Computed_field_mesh_integral_core_cast( reinterpret_cast(mesh_integral_field)->core)); } +cmzn_mesh_id cmzn_field_mesh_integral_get_mesh( + cmzn_field_mesh_integral_id mesh_integral_field) +{ + if (mesh_integral_field) + { + Computed_field_mesh_integral *mesh_integral_core = Computed_field_mesh_integral_core_cast(mesh_integral_field); + return cmzn_mesh_access(mesh_integral_core->getMesh()); + } + return nullptr; +} + +int cmzn_field_mesh_integral_set_mesh( + cmzn_field_mesh_integral_id mesh_integral_field, cmzn_mesh_id mesh) +{ + if (mesh_integral_field) + { + Computed_field_mesh_integral *mesh_integral_core = Computed_field_mesh_integral_core_cast(mesh_integral_field); + return mesh_integral_core->setMesh(mesh); + } + return CMZN_ERROR_ARGUMENT; +} + int cmzn_field_mesh_integral_get_numbers_of_points( cmzn_field_mesh_integral_id mesh_integral_field, int values_count, int *values) diff --git a/tests/fieldmodule/fieldmeshoperators.cpp b/tests/fieldmodule/fieldmeshoperators.cpp index 5a26c4319..f56c60d42 100644 --- a/tests/fieldmodule/fieldmeshoperators.cpp +++ b/tests/fieldmodule/fieldmeshoperators.cpp @@ -64,12 +64,15 @@ TEST(ZincFieldMeshIntegral, quadrature) FieldMeshIntegral volumeField = zinc.fm.createFieldMeshIntegral(integrandField, coordinateField, mesh3d); EXPECT_TRUE(volumeField.isValid()); + EXPECT_EQ(mesh3d, volumeField.getMesh()); EXPECT_EQ(Element::QUADRATURE_RULE_GAUSSIAN, volumeField.getElementQuadratureRule()); FieldMeshIntegral surfaceAreaField = zinc.fm.createFieldMeshIntegral(integrandField, coordinateField, exteriorMesh2d); EXPECT_TRUE(surfaceAreaField.isValid()); + EXPECT_EQ(exteriorMesh2d, surfaceAreaField.getMesh()); EXPECT_EQ(Element::QUADRATURE_RULE_GAUSSIAN, surfaceAreaField.getElementQuadratureRule()); FieldMeshIntegral lengthField = zinc.fm.createFieldMeshIntegral(integrandField, coordinateField, mesh1d); EXPECT_TRUE(lengthField.isValid()); + EXPECT_EQ(mesh1d, lengthField.getMesh()); EXPECT_EQ(Element::QUADRATURE_RULE_GAUSSIAN, lengthField.getElementQuadratureRule()); const double scaleFactors[3] = { 1.5, 0.75, 2.0 }; @@ -265,6 +268,71 @@ TEST(ZincFieldMeshIntegral, quadrature) } } +TEST(ZincFieldMeshIntegral, cube_setMesh) +{ + ZincTestSetupCpp zinc; + int result; + + EXPECT_EQ(RESULT_OK, result = zinc.root_region.readFile( + TestResources::getLocation(TestResources::FIELDMODULE_CUBE_RESOURCE))); + + Mesh mesh3d = zinc.fm.findMeshByDimension(3); + EXPECT_EQ(1, mesh3d.getSize()); + Mesh mesh2d = zinc.fm.findMeshByDimension(2); + EXPECT_EQ(6, mesh2d.getSize()); + Mesh mesh1d = zinc.fm.findMeshByDimension(1); + EXPECT_EQ(12, mesh1d.getSize()); + + Field coordinates = zinc.fm.findFieldByName("coordinates"); + EXPECT_TRUE(coordinates.isValid()); + + const double oneValue = 1.0; + Field one = zinc.fm.createFieldConstant(1, &oneValue); + EXPECT_TRUE(one.isValid()); + + FieldMeshIntegral meshIntegral = zinc.fm.createFieldMeshIntegral(one, coordinates, mesh3d); + EXPECT_TRUE(meshIntegral.isValid()); + EXPECT_EQ(mesh3d, meshIntegral.getMesh()); + + Fieldcache fieldcache = zinc.fm.createFieldcache(); + const double TOL = 1.0E-12; + double volume; + EXPECT_EQ(RESULT_OK, meshIntegral.evaluateReal(fieldcache, 1, &volume)); + EXPECT_NEAR(1.0, volume, TOL); + + EXPECT_EQ(RESULT_OK, meshIntegral.setMesh(mesh2d)); + EXPECT_EQ(mesh2d, meshIntegral.getMesh()); + double area; + EXPECT_EQ(RESULT_OK, meshIntegral.evaluateReal(fieldcache, 1, &area)); + EXPECT_NEAR(6.0, area, TOL); + + EXPECT_EQ(RESULT_OK, meshIntegral.setMesh(mesh1d)); + EXPECT_EQ(mesh1d, meshIntegral.getMesh()); + double length; + EXPECT_EQ(RESULT_OK, meshIntegral.evaluateReal(fieldcache, 1, &length)); + EXPECT_NEAR(12.0, length, TOL); + + FieldElementGroup elementGroup2 = zinc.fm.createFieldElementGroup(mesh2d); + EXPECT_TRUE(elementGroup2.isValid()); + MeshGroup meshGroup2 = elementGroup2.getMeshGroup(); + EXPECT_TRUE(meshGroup2.isValid()); + EXPECT_EQ(0, meshGroup2.getSize()); + EXPECT_EQ(RESULT_OK, meshIntegral.setMesh(meshGroup2)); + EXPECT_EQ(meshGroup2, meshIntegral.getMesh()); + EXPECT_EQ(RESULT_OK, meshIntegral.evaluateReal(fieldcache, 1, &area)); + EXPECT_NEAR(0.0, area, TOL); + + Element face2 = mesh2d.findElementByIdentifier(2); + EXPECT_EQ(RESULT_OK, meshGroup2.addElement(face2)); + EXPECT_EQ(RESULT_OK, meshIntegral.evaluateReal(fieldcache, 1, &area)); + EXPECT_NEAR(1.0, area, TOL); + + Element face5 = mesh2d.findElementByIdentifier(5); + EXPECT_EQ(RESULT_OK, meshGroup2.addElement(face5)); + EXPECT_EQ(RESULT_OK, meshIntegral.evaluateReal(fieldcache, 1, &area)); + EXPECT_NEAR(2.0, area, TOL); +} + TEST(ZincFieldMeshIntegral, midpoint_quadrature_large_numbers) { ZincTestSetupCpp zinc; @@ -349,8 +417,14 @@ TEST(ZincFieldMeshIntegralSquares, quadrature) FieldMeshIntegralSquares integralSquaresField = zinc.fm.createFieldMeshIntegralSquares(integrandField, coordinateField, mesh3d); EXPECT_TRUE(integralSquaresField.isValid()); + EXPECT_EQ(mesh3d, integralSquaresField.getMesh()); EXPECT_EQ(Element::QUADRATURE_RULE_GAUSSIAN, integralSquaresField.getElementQuadratureRule()); + // test can cast to MeshIntegral field and query it, as it's a base class of this field + FieldMeshIntegral integralField = integralSquaresField.castMeshIntegral(); + EXPECT_TRUE(integralField.isValid()); + EXPECT_EQ(mesh3d, integralField.getMesh()); + const double expectedValue[] = { 3.0 + 1.0/6.0, 4.0*(3.0 + 1.0/6.0) }; double value[2]; const double tolerance = 1.0E-12;