Integration of OpenCASCADE for aortic valve particle generation

.github/workflows/ci.yml

@@ -37,7 +37,10 @@ jobs:
             pkg-config `# for installing libraries with vcpkg`\
             git \
             cmake \
-            ninja-build
+            ninja-build \
+            libfontconfig1-dev `# From here required for vcpkg opencascade`\
+            libx11-dev \
+            libgl-dev
 
       - uses: hendrikmuhs/ccache-action@v1.2
         with:
@@ -61,15 +64,17 @@ jobs:
             simbody \
             gtest \
             xsimd \
-            pybind11
+            pybind11 \
+            opencascade
 
       - name: Generate buildsystem using float (No test)
         run: |
           cmake -G Ninja \
             -D CMAKE_BUILD_TYPE=Release \
-            -D SPHINXSYS_USE_FLOAT=ON \
             -D CMAKE_TOOLCHAIN_FILE="${{github.workspace}}/vcpkg/scripts/buildsystems/vcpkg.cmake" \
             -D CMAKE_C_COMPILER_LAUNCHER=ccache -D CMAKE_CXX_COMPILER_LAUNCHER=ccache \
+            -D SPHINXSYS_USE_FLOAT=ON \
+            -D SPHINXSYS_MODULE_OPENCASCADE=ON \
             -D SPHINXSYS_CI=ON \
             -S ${{github.workspace}} \
             -B ${{github.workspace}}/build
@@ -97,7 +102,10 @@ jobs:
             pkg-config `# for installing libraries with vcpkg`\
             git \
             cmake \
-            ninja-build
+            ninja-build \
+            libfontconfig1-dev `# From here required for vcpkg opencascade`\
+            libx11-dev \
+            libgl-dev
 
       - uses: hendrikmuhs/ccache-action@v1.2
         with:
@@ -121,16 +129,18 @@ jobs:
             simbody \
             gtest \
             xsimd \
-            pybind11
+            pybind11 \
+            opencascade
 
       - name: Generate buildsystem using double
         run: |
           cmake -G Ninja \
             -D CMAKE_BUILD_TYPE=Release \
-            -D SPHINXSYS_USE_FLOAT=OFF \
             -D CMAKE_TOOLCHAIN_FILE="${{github.workspace}}/vcpkg/scripts/buildsystems/vcpkg.cmake" \
             -D CMAKE_C_COMPILER_LAUNCHER=ccache -D CMAKE_CXX_COMPILER_LAUNCHER=ccache \
             -D SPHINXSYS_CI=ON \
+            -D SPHINXSYS_USE_FLOAT=OFF \
+            -D SPHINXSYS_MODULE_OPENCASCADE=ON \
             -S ${{github.workspace}} \
             -B ${{github.workspace}}/build
 
@@ -209,7 +219,7 @@ jobs:
             simbody `
             gtest `
             xsimd `
-            pybind11
+            pybind11 
 
       - uses: ilammy/msvc-dev-cmd@v1
 
@@ -297,15 +307,16 @@ jobs:
 
       - name: Install dependencies
         run: |
-          ${{github.workspace}}/vcpkg/vcpkg install --clean-after-build \
+          ${{github.workspace}}/vcpkg/vcpkg install --clean-after-build --allow-unsupported \
             eigen3 \
             tbb \
             boost-program-options \
             boost-geometry \
             gtest \
             simbody \
             xsimd \
-            pybind11
+            pybind11 \
+            opencascade
 
       - name: Generate buildsystem
         run: |
@@ -356,4 +367,4 @@ jobs:
         if: ${{ steps.fourth-try.outcome == 'failure' }} 
         run: |
           cd build 
-          ctest --rerun-failed --output-on-failure
+          ctest --rerun-failed --output-on-failure

CMakeLists.txt

@@ -11,6 +11,7 @@ option(SPHINXSYS_3D "Build sphinxsys_3d library" ON)
 option(SPHINXSYS_BUILD_TESTS "Build tests" ON)
 option(SPHINXSYS_DEVELOPER_MODE "Developer mode has more flags active for code quality" ON)
 option(SPHINXSYS_USE_FLOAT "Build using float (single-precision floating-point format) as primary type" OFF)
+option(SPHINXSYS_MODULE_OPENCASCADE "Build extension relying on OpenCASCADE" OFF)
 
 # ------ Global properties (Some cannot be set on INTERFACE targets)
 set(CMAKE_VERBOSE_MAKEFILE OFF CACHE BOOL "Enable verbose compilation commands for Makefile and Ninja" FORCE) # Extra fluff needed for Ninja: https://github.com/ninja-build/ninja/issues/900

modules/CMakeLists.txt

@@ -1,3 +1,7 @@
-add_library(structural_simulation_module structural_simulation/structural_simulation_class.cpp structural_simulation/structural_simulation_class.h)
-target_include_directories(structural_simulation_module PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/structural_simulation> $<INSTALL_INTERFACE:include/structural_simulation>)
-target_link_libraries(structural_simulation_module PUBLIC sphinxsys_3d)
+SUBDIRLIST(SUBDIRS ${CMAKE_CURRENT_SOURCE_DIR})
+
+foreach(subdir ${SUBDIRS})
+    if(EXISTS ${CMAKE_CURRENT_SOURCE_DIR}/${subdir}/CMakeLists.txt)
+	    add_subdirectory(${subdir})
+    endif()
+endforeach()

modules/opencascade/CMakeLists.txt

@@ -0,0 +1,29 @@
+if(NOT SPHINXSYS_MODULE_OPENCASCADE OR ${CMAKE_SYSTEM_NAME} MATCHES "Windows")
+    return()
+endif()
+
+find_package(OpenCASCADE CONFIG REQUIRED)
+# Must be shared library if OpenCASCADE is built static otherwise coplyeft OpenCASCADE license propagates
+get_target_property(OCCT_LIBRARY_TYPE TKernel TYPE) # https://cmake.org/cmake/help/latest/prop_tgt/TYPE.html
+if(OCCT_LIBRARY_TYPE MATCHES "STATIC_LIBRARY")
+    set(SPHINXSYS_OPENCASCADE_LIBRARY_TYPE SHARED)
+endif()
+add_library(sphinxsys_opencascade
+    ${SPHINXSYS_OPENCASCADE_LIBRARY_TYPE}
+    opencascade/relax_dynamics_surface.cpp 
+    opencascade/relax_dynamics_surface.h 
+    opencascade/surface_shape.cpp 
+    opencascade/surface_shape.h 
+    opencascade/vector.cpp 
+    opencascade/vector.h
+)
+target_include_directories(sphinxsys_opencascade PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/opencascade>)
+target_link_libraries(sphinxsys_opencascade PUBLIC sphinxsys_3d)
+
+# Must list all dependent packages given here https://dev.opencascade.org/doc/refman/html/toolkit_tkstep.html because improper CMake target link visibility
+target_link_libraries(sphinxsys_opencascade PRIVATE TKSTEP TKSTEPAttr TKSTEP209 TKSTEPBase TKXSBase TKShHealing TKTopAlgo TKGeomAlgo TKBRep TKGeomBase) 
+target_link_libraries(sphinxsys_opencascade PUBLIC TKG3d TKG2d TKMath TKernel) 
+
+if(SPHINXSYS_BUILD_TESTS)
+    add_subdirectory(tests)
+endif()

modules/opencascade/opencascade/relax_dynamics_surface.cpp

@@ -0,0 +1,108 @@
+#include "relax_dynamics.h"
+#include "base_particles.hpp"
+#include "particle_generator_lattice.h"
+#include "level_set_shape.h"
+#include "surface_shape.h"
+#include "relax_dynamics_surface.h"
+
+#include <opencascade/GeomAPI_ProjectPointOnSurf.hxx>
+#include <opencascade/gp_Pnt.hxx>
+//========================================================================================================//
+namespace SPH
+{
+	//=====================================================================================================//
+	namespace relax_dynamics
+	{
+        //=================================================================================================//
+        ShapeSurfaceBounding2::ShapeSurfaceBounding2(RealBody &real_body_)
+            : LocalDynamics(real_body_), RelaxDataDelegateSimple(real_body_), pos_(particles_->pos_) /*, surface_shape_(DynamicCast<SurfaceShape>(this, sph_body_.body_shape_))*/
+        {
+            shape_ = real_body_.body_shape_;
+
+        }
+        //=================================================================================================//
+        void ShapeSurfaceBounding2::update(size_t index_i, Real dt)
+        {
+            pos_[index_i] = shape_->findClosestPoint(pos_[index_i]);
+
+        }
+        //=================================================================================================//
+        RelaxationStepInnerFirstHalf::
+            RelaxationStepInnerFirstHalf(BaseInnerRelation &inner_relation, bool level_set_correction)
+            : BaseDynamics<void>(inner_relation.getSPHBody()), real_body_(inner_relation.real_body_),
+              inner_relation_(inner_relation)
+        {
+            relaxation_acceleration_inner_ =
+                makeUnique<InteractionDynamics<RelaxationAccelerationInner>>(inner_relation);
+        }
+        //=================================================================================================//
+        void RelaxationStepInnerFirstHalf::exec(Real dt)
+        {
+            real_body_->updateCellLinkedList();
+            inner_relation_.updateConfiguration();
+            relaxation_acceleration_inner_->exec();
+        }
+
+        //=================================================================================================//
+        RelaxationStepInnerSecondHalf::
+            RelaxationStepInnerSecondHalf(BaseInnerRelation &inner_relation, bool level_set_correction)
+            : BaseDynamics<void>(inner_relation.getSPHBody()), real_body_(inner_relation.real_body_),
+              get_time_step_square_(*real_body_), update_particle_position_(*real_body_),
+              surface_bounding_(*real_body_)
+        {
+        }
+        //=================================================================================================//
+        void RelaxationStepInnerSecondHalf::exec(Real dt)
+        {
+            Real dt_square = get_time_step_square_.exec();
+            update_particle_position_.exec(dt_square);
+            surface_bounding_.exec();
+        }
+
+		//=================================================================================================//
+        SurfaceNormalDirection::SurfaceNormalDirection(SPHBody &sph_body)
+            : RelaxDataDelegateSimple(sph_body), LocalDynamics(sph_body), 
+            surface_shape_(DynamicCast<SurfaceShape>(this, sph_body.body_shape_)), 
+            pos_(particles_->pos_), n_(*particles_->getVariableByName<Vecd>("NormalDirection")) {}
+
+        //=================================================================================================//
+        void SurfaceNormalDirection::update(size_t index_i, Real dt)
+        {
+            Extrema_ExtAlgo Algo = Extrema_ExtAlgo_Tree;
+            gp_Vec tangent_u, tangent_v;
+            gp_Vec norm;
+            Vecd normal_direction_;
+
+
+            gp_Pnt point1 = EigenToOcct(pos_[index_i]);
+            GeomAPI_ProjectPointOnSurf pointonsurf(point1, surface_shape_->surface_, Algo);
+
+
+            Standard_Real u;
+            Standard_Real v;
+            gp_Pnt Point;
+
+            pointonsurf.Parameters(1, u, v);
+
+            surface_shape_->surface_->D1(u, v, Point, tangent_u, tangent_v);
+            norm = tangent_u.Crossed(tangent_v);
+            normal_direction_ = OcctVecToEigen(norm);
+            n_[index_i] = normal_direction_.normalized();
+        }
+        //=================================================================================================//
+        ConstrainSurfaceBodyRegion::
+            ConstrainSurfaceBodyRegion(BodyPartByParticle &body_part)
+            : BaseLocalDynamics<BodyPartByParticle>(body_part), RelaxDataDelegateSimple(sph_body_),
+              acc_(particles_->acc_)
+        {
+        }
+        //=================================================================================================//
+        void ConstrainSurfaceBodyRegion::update(size_t index_i, Real dt)
+        {
+            acc_[index_i] = Vecd::Zero();
+        }
+        //=================================================================================================//
+    }
+	//=================================================================================================//
+}
+//=================================================================================================//

modules/opencascade/opencascade/relax_dynamics_surface.h

@@ -0,0 +1,123 @@
+/* -------------------------------------------------------------------------*
+ *								SPHinXsys									*
+ * -------------------------------------------------------------------------*
+ * SPHinXsys (pronunciation: s'finksis) is an acronym from Smoothed Particle*
+ * Hydrodynamics for industrial compleX systems. It provides C++ APIs for	*
+ * physical accurate simulation and aims to model coupled industrial dynamic*
+ * systems including fluid, solid, multi-body dynamics and beyond with SPH	*
+ * (smoothed particle hydrodynamics), a meshless computational method using	*
+ * particle discretization.													*
+ *																			*
+ * SPHinXsys is partially funded by German Research Foundation				*
+ * (Deutsche Forschungsgemeinschaft) DFG HU1527/6-1, HU1527/10-1,			*
+ *  HU1527/12-1 and HU1527/12-4													*
+ *                                                                          *
+ * Portions copyright (c) 2017-2022 Technical University of Munich and		*
+ * the authors' affiliations.												*
+ *                                                                          *
+ * Licensed under the Apache License, Version 2.0 (the "License"); you may  *
+ * not use this file except in compliance with the License. You may obtain a*
+ * copy of the License at http://www.apache.org/licenses/LICENSE-2.0.       *
+ *                                                                          *
+ * ------------------------------------------------------------------------*/
+/**
+ * @file 	relax_dynamics.h
+ * @brief 	This is the classes of particle relaxation in order to produce body fitted
+ * 			initial particle distribution.
+ * @author	Chi Zhang and Xiangyu Hu
+ */
+
+#ifndef RELAX_DYNAMICS_SURFACE_H
+#define RELAX_DYNAMICS_SURFACE_H
+
+#include "sphinxsys.h"
+#include "vector.h"
+
+namespace SPH
+{
+
+    class SurfaceShape;
+
+	namespace relax_dynamics
+    {
+     class ShapeSurfaceBounding2 : public LocalDynamics,
+                                  public RelaxDataDelegateSimple
+    {
+      public:
+        ShapeSurfaceBounding2(RealBody &real_body_);
+        virtual ~ShapeSurfaceBounding2(){};
+        void update(size_t index_i, Real dt = 0.0);
+
+      protected:
+        StdLargeVec<Vecd> &pos_;
+        Shape *shape_;
+    };
+
+
+      class RelaxationStepInnerFirstHalf : public BaseDynamics<void>
+    {
+      public:
+        explicit RelaxationStepInnerFirstHalf(BaseInnerRelation &inner_relation,
+                                              bool level_set_correction = false);
+        virtual ~RelaxationStepInnerFirstHalf(){};
+        virtual void exec(Real dt = 0.0) override;
+
+
+      protected:
+        RealBody *real_body_;
+        BaseInnerRelation &inner_relation_;
+        UniquePtr<BaseDynamics<void>> relaxation_acceleration_inner_;
+    };
+
+    class RelaxationStepInnerSecondHalf : public BaseDynamics<void>
+    {
+      public:
+        explicit RelaxationStepInnerSecondHalf(BaseInnerRelation &inner_relation,
+                                               bool level_set_correction = false);
+        virtual ~RelaxationStepInnerSecondHalf(){};
+        SimpleDynamics<ShapeSurfaceBounding2> &SurfaceBounding() { return surface_bounding_; };
+        virtual void exec(Real dt = 0.0) override;
+
+
+      protected:
+        RealBody *real_body_;
+        ReduceDynamics<GetTimeStepSizeSquare> get_time_step_square_;
+        SimpleDynamics<UpdateParticlePosition> update_particle_position_;
+        SimpleDynamics<ShapeSurfaceBounding2> surface_bounding_;
+    };
+
+    /**
+     * @class SurfaceNormalDirection
+     * @brief get the normal direction of surface particles.
+     */
+    class SurfaceNormalDirection : public RelaxDataDelegateSimple, public LocalDynamics
+    {
+      public:
+        explicit SurfaceNormalDirection(SPHBody &sph_body);
+        virtual ~SurfaceNormalDirection(){};
+        void update(size_t index_i, Real dt = 0.0);
+
+      protected:
+        SurfaceShape *surface_shape_;
+        StdLargeVec<Vecd> &pos_, &n_;
+    };
+
+
+    /**@class ConstrainSuefaceBodyRegion
+     * @brief Fix the position surafce body part.
+     */
+    class ConstrainSurfaceBodyRegion : public BaseLocalDynamics<BodyPartByParticle>, public RelaxDataDelegateSimple
+    {
+      public:
+        ConstrainSurfaceBodyRegion(BodyPartByParticle &body_part);
+        virtual ~ConstrainSurfaceBodyRegion(){};
+        void update(size_t index_i, Real dt = 0.0);
+
+      protected:
+        StdLargeVec<Vecd> &acc_;
+    };
+
+
+    }
+  }
+#endif // RELAX_DYNAMICS_H