diff --git a/CMake/External/CMakeLists.txt b/CMake/External/CMakeLists.txt
index 9121b16fc89ffac4e2402b0489a522e9ae165f2a..df3f072e564744883accc5646cf0410b5cd051a6 100644
--- a/CMake/External/CMakeLists.txt
+++ b/CMake/External/CMakeLists.txt
@@ -94,6 +94,7 @@ ExternalProject_Add( ${PROJECT_NAME}
     -D${PROJECT_NAME}_SUPERBUILD:BOOL=OFF
     -D${PROJECT_NAME}_BUILD_EXAMPLES:BOOL=${${PROJECT_NAME}_BUILD_EXAMPLES}
     -D${PROJECT_NAME}_BUILD_TESTING:BOOL=${${PROJECT_NAME}_BUILD_TESTING}
+    -D${PROJECT_NAME}_WRAP_CSHARP:BOOL=${${PROJECT_NAME}_WRAP_CSHARP}
     -D${PROJECT_NAME}_USE_OpenHaptics:BOOL=${${PROJECT_NAME}_USE_OpenHaptics}
     -D${PROJECT_NAME}_USE_MODEL_REDUCTION:BOOL=${${PROJECT_NAME}_USE_MODEL_REDUCTION}    
     -D${PROJECT_NAME}_ENABLE_AUDIO:BOOL=${${PROJECT_NAME}_ENABLE_AUDIO}
diff --git a/CMake/External/External_iMSTKData.cmake b/CMake/External/External_iMSTKData.cmake
index d849218257774f21519b021cf13e54f8609a577f..df40d3d47823bcd013044bddb9b781e43c2f4014 100644
--- a/CMake/External/External_iMSTKData.cmake
+++ b/CMake/External/External_iMSTKData.cmake
@@ -16,13 +16,14 @@ set(copy_data_command
   )
 
 include(imstkAddExternalProject)
+set(GIT_SHA "14824e3d53328ed6be481981959780f32881030b")
+set(DATA_URL "https://gitlab.kitware.com/iMSTK/imstk-data/-/archive/${GIT_SHA}/imstk-data-${GIT_SHA}.zip")
+
 imstk_add_external_project( iMSTKData
-  GIT_REPOSITORY  "https://gitlab.kitware.com/iMSTK/imstk-data.git"
-  GIT_TAG  "14824e3d53328ed6be481981959780f32881030b"
+  URL ${DATA_URL}
   CONFIGURE_COMMAND ${SKIP_STEP_COMMAND}
   BUILD_COMMAND ${SKIP_STEP_COMMAND}
   INSTALL_COMMAND COMMAND ${copy_data_command}
   DEPENDENCIES ""
-  #VERBOSE
 )
 
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 2115d505c8e461f577e32ee248e3b4d7f4732679..1d119e5e9d044addd7f7e1d638afb48b40d8d4e7 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -64,24 +64,41 @@ endif()
 #-----------------------------------------------------------------------------
 # Options
 #-----------------------------------------------------------------------------
+
+# SWIG 
+option(${PROJECT_NAME}_WRAP_CSHARP "Build iMSTK-C# wrapper code and lib" OFF)
+# SWIG will generate managed array using pinning if this option is ON. Otherwise using P/Invoke default array marshalling.
+option(${PROJECT_NAME}_SWIG_PINNED_ARRAY "Managed C# arrays using pinning" ON)
+mark_as_advanced(${PROJECT_NAME}_SWIG_PINNED_ARRAY)
+if (${PROJECT_NAME}_WRAP_CSHARP)
+  find_package(SWIG REQUIRED)
+  set(CMAKE_POSITION_INDEPENDENT_CODE ON)
+
+  include_directories(${CMAKE_CURRENT_SOURCE_DIR})
+
+  set(CMAKE_SWIG_FLAGS "")
+endif()
+
+# CUDA
 option(${PROJECT_NAME}_ENABLE_CUDA_BACKEND "Enable iMSTK CUDA backend" OFF)
 if (${PROJECT_NAME}_ENABLE_CUDA_BACKEND)
-	include(CheckLanguage)
-	check_language(CUDA)
-	if (CMAKE_CUDA_COMPILER)
-		enable_language(CUDA)
-		set(CMAKE_CUDA_SEPARABLE_COMPILATION ON)
-		set(CUDA_PROPAGATE_HOST_FLAGS ON)
-		if(NOT DEFINED CMAKE_CUDA_STANDARD)
-			set(CMAKE_CUDA_STANDARD 11)
-			set(CMAKE_CUDA_STANDARD_REQUIRED ON)
-		endif()
-	else()
-		message(STATUS "WARNING: CUDA compiler NOT FOUND; CUDA backend not enabled!")
-		set(${PROJECT_NAME}_ENABLE_CUDA_BACKEND OFF)
+  include(CheckLanguage)
+  check_language(CUDA)
+  if (CMAKE_CUDA_COMPILER)
+    enable_language(CUDA)
+	set(CMAKE_CUDA_SEPARABLE_COMPILATION ON)
+	set(CUDA_PROPAGATE_HOST_FLAGS ON)
+	if(NOT DEFINED CMAKE_CUDA_STANDARD)
+	  set(CMAKE_CUDA_STANDARD 11)
+	  set(CMAKE_CUDA_STANDARD_REQUIRED ON)
 	endif()
+  else()
+	message(STATUS "WARNING: CUDA compiler NOT FOUND; CUDA backend not enabled!")
+    set(${PROJECT_NAME}_ENABLE_CUDA_BACKEND OFF)
+  endif()
 endif ()
 
+# General Options
 option(${PROJECT_NAME}_BUILD_EXAMPLES "Build iMSTK examples" ON)
 option(${PROJECT_NAME}_BUILD_TESTING "Build iMSTK tests" ON)
 set(BUILD_TESTING OFF)
@@ -389,6 +406,10 @@ add_subdirectory(Source/Testing)
 add_subdirectory(Source/Filtering)
 add_subdirectory(Source/FilteringCore)
 
+if (${PROJECT_NAME}_WRAP_CSHARP)
+  add_subdirectory(Source/Wrappers)
+endif()
+
 #--------------------------------------------------------------------------
 # Add Examples subdirectories
 #--------------------------------------------------------------------------
diff --git a/Source/CollisionDetection/CollisionData/imstkCollisionData.h b/Source/CollisionDetection/CollisionData/imstkCollisionData.h
index 569d5ad58e5f9236e605bae010ff0cdd1b0d706f..690d231b32485cec5ee0dc1d908ed9b4228059c6 100644
--- a/Source/CollisionDetection/CollisionData/imstkCollisionData.h
+++ b/Source/CollisionDetection/CollisionData/imstkCollisionData.h
@@ -168,6 +168,29 @@ struct CollisionElement
         }
     }
 
+    CollisionElement& operator=(const CollisionElement& other)
+    {
+        m_type = other.m_type;
+        switch (m_type)
+        {
+        case CollisionElementType::Empty:
+            break;
+        case CollisionElementType::CellVertex:
+            m_element.m_CellVertexElement = other.m_element.m_CellVertexElement;
+            break;
+        case CollisionElementType::CellIndex:
+            m_element.m_CellIndexElement = other.m_element.m_CellIndexElement;
+            break;
+        case CollisionElementType::PointDirection:
+            m_element.m_PointDirectionElement = other.m_element.m_PointDirectionElement;
+            break;
+        case CollisionElementType::PointIndexDirection:
+            m_element.m_PointIndexDirectionElement = other.m_element.m_PointIndexDirectionElement;
+            break;
+        }
+        return *this;
+    }
+
     union Element
     {
         EmptyElement m_EmptyElement;
@@ -199,4 +222,4 @@ public:
     std::shared_ptr<Geometry>     geomA;
     std::shared_ptr<Geometry>     geomB;
 };
-}
\ No newline at end of file
+}
diff --git a/Source/CollisionDetection/Testing/imstkCollisionDataTest.cpp b/Source/CollisionDetection/Testing/imstkCollisionDataTest.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..c8076ffe523c858273768d6d12531f17b31d3f64
--- /dev/null
+++ b/Source/CollisionDetection/Testing/imstkCollisionDataTest.cpp
@@ -0,0 +1,175 @@
+#include "gtest/gtest.h"
+
+#include "imstkCollisionData.h"
+
+using namespace imstk;
+
+namespace {
+}
+
+struct imstkCollisionElementTest : public ::testing::Test
+{
+    imstkCollisionElementTest()
+    {
+        cv.pts[0] = Vec3d{ 1,2,3 };
+        cv.pts[1] = Vec3d{ 2,3,4 };
+        cv.pts[2] = Vec3d{ 3,4,5 };
+        cv.pts[3] = Vec3d{ 4,5,6 };
+        cv.size = 10;
+
+        ci.ids[0] = 7;
+        ci.ids[1] = 8;
+        ci.ids[2] = 9;
+        ci.ids[3] = 10;
+        ci.idCount = 11;
+        ci.cellType = IMSTK_VERTEX;
+
+        pd.pt = Vec3d{ 2,3,4 };
+        pd.dir = Vec3d{ 4,5,6 };
+        pd.penetrationDepth = 12.0;
+
+        pi.ptIndex = 13;
+        pi.dir = Vec3d{ 3,4,5 };
+        pi.penetrationDepth = 2.0;
+    };
+
+    CellVertexElement cv;
+    CellIndexElement ci;
+    PointDirectionElement pd;
+    PointIndexDirectionElement pi;
+};
+
+TEST_F(imstkCollisionElementTest, Constructor)
+{
+    {
+        CollisionElement e;
+        EXPECT_EQ(CollisionElementType::Empty, e.m_type);
+    }
+    {
+        CollisionElement e(cv);
+        EXPECT_EQ(CollisionElementType::CellVertex, e.m_type);
+        EXPECT_EQ(cv.size, e.m_element.m_CellVertexElement.size);
+        EXPECT_EQ(cv.pts[1], e.m_element.m_CellVertexElement.pts[1]);
+    }
+    {
+        CollisionElement e(ci);
+        EXPECT_EQ(CollisionElementType::CellIndex, e.m_type);
+        EXPECT_EQ(ci.idCount, e.m_element.m_CellIndexElement.idCount);
+        EXPECT_EQ(ci.ids[1], e.m_element.m_CellIndexElement.ids[1]);
+    }
+    {
+        CollisionElement e(pd);
+        EXPECT_EQ(CollisionElementType::PointDirection, e.m_type);
+        EXPECT_EQ(pd.dir, e.m_element.m_PointDirectionElement.dir);
+        EXPECT_EQ(pd.penetrationDepth, e.m_element.m_PointDirectionElement.penetrationDepth);
+    } 
+    {
+        CollisionElement e(pi);
+        EXPECT_EQ(CollisionElementType::PointIndexDirection, e.m_type);
+        EXPECT_EQ(pi.dir, e.m_element.m_PointIndexDirectionElement.dir);
+        EXPECT_EQ(pi.ptIndex, e.m_element.m_PointIndexDirectionElement.ptIndex);
+    }
+}
+
+TEST_F(imstkCollisionElementTest, CopyConstructor)
+{
+    {
+        CollisionElement old;
+        CollisionElement e(old);
+        EXPECT_EQ(CollisionElementType::Empty, e.m_type);
+    }
+    {
+        CollisionElement old(cv);
+        CollisionElement e(old);
+        EXPECT_EQ(CollisionElementType::CellVertex, e.m_type);
+        EXPECT_EQ(cv.size, e.m_element.m_CellVertexElement.size);
+        EXPECT_EQ(cv.pts[1], e.m_element.m_CellVertexElement.pts[1]);
+    }
+    {
+        CollisionElement old(ci);
+        CollisionElement e(old);
+        EXPECT_EQ(CollisionElementType::CellIndex, e.m_type);
+        EXPECT_EQ(ci.idCount, e.m_element.m_CellIndexElement.idCount);
+        EXPECT_EQ(ci.ids[1], e.m_element.m_CellIndexElement.ids[1]);
+    }
+    {
+        CollisionElement old(pd);
+        CollisionElement e(old);
+        EXPECT_EQ(CollisionElementType::PointDirection, e.m_type);
+        EXPECT_EQ(pd.dir, e.m_element.m_PointDirectionElement.dir);
+        EXPECT_EQ(pd.penetrationDepth, e.m_element.m_PointDirectionElement.penetrationDepth);
+    }
+    {
+        CollisionElement old(pi);
+        CollisionElement e(old);
+        EXPECT_EQ(CollisionElementType::PointIndexDirection, e.m_type);
+        EXPECT_EQ(pi.dir, e.m_element.m_PointIndexDirectionElement.dir);
+        EXPECT_EQ(pi.ptIndex, e.m_element.m_PointIndexDirectionElement.ptIndex);
+    }
+}
+
+TEST_F(imstkCollisionElementTest, DataAssignment)
+{
+        CollisionElement e;
+        e = cv;
+        EXPECT_EQ(CollisionElementType::CellVertex, e.m_type);
+        EXPECT_EQ(cv.size, e.m_element.m_CellVertexElement.size);
+        EXPECT_EQ(cv.pts[1], e.m_element.m_CellVertexElement.pts[1]);
+
+        e = ci;
+        EXPECT_EQ(CollisionElementType::CellIndex, e.m_type);
+        EXPECT_EQ(ci.idCount, e.m_element.m_CellIndexElement.idCount);
+        EXPECT_EQ(ci.ids[1], e.m_element.m_CellIndexElement.ids[1]);
+
+        e = pd;
+        EXPECT_EQ(CollisionElementType::PointDirection, e.m_type);
+        EXPECT_EQ(pd.dir, e.m_element.m_PointDirectionElement.dir);
+        EXPECT_EQ(pd.penetrationDepth, e.m_element.m_PointDirectionElement.penetrationDepth);
+
+        e = pi;
+        EXPECT_EQ(CollisionElementType::PointIndexDirection, e.m_type);
+        EXPECT_EQ(pi.dir, e.m_element.m_PointIndexDirectionElement.dir);
+        EXPECT_EQ(pi.ptIndex, e.m_element.m_PointIndexDirectionElement.ptIndex);
+}
+
+TEST_F(imstkCollisionElementTest, Assignment)
+{
+    {
+        CollisionElement old;
+        CollisionElement e(pi);
+        e = old;
+        EXPECT_EQ(CollisionElementType::Empty, e.m_type);
+    }
+    {
+        CollisionElement old(cv);
+        CollisionElement e;
+        e = old;
+        EXPECT_EQ(CollisionElementType::CellVertex, e.m_type);
+        EXPECT_EQ(cv.size, e.m_element.m_CellVertexElement.size);
+        EXPECT_EQ(cv.pts[1], e.m_element.m_CellVertexElement.pts[1]);
+    }
+    {
+        CollisionElement old(ci);
+        CollisionElement e;
+        e = old;
+        EXPECT_EQ(CollisionElementType::CellIndex, e.m_type);
+        EXPECT_EQ(ci.idCount, e.m_element.m_CellIndexElement.idCount);
+        EXPECT_EQ(ci.ids[1], e.m_element.m_CellIndexElement.ids[1]);
+    }
+    {
+        CollisionElement old(pd);
+        CollisionElement e;
+        e = old;
+        EXPECT_EQ(CollisionElementType::PointDirection, e.m_type);
+        EXPECT_EQ(pd.dir, e.m_element.m_PointDirectionElement.dir);
+        EXPECT_EQ(pd.penetrationDepth, e.m_element.m_PointDirectionElement.penetrationDepth);
+    }
+    {
+        CollisionElement old(pi);
+        CollisionElement e;
+        e = old;
+        EXPECT_EQ(CollisionElementType::PointIndexDirection, e.m_type);
+        EXPECT_EQ(pi.dir, e.m_element.m_PointIndexDirectionElement.dir);
+        EXPECT_EQ(pi.ptIndex, e.m_element.m_PointIndexDirectionElement.ptIndex);
+    }
+}
\ No newline at end of file
diff --git a/Source/Common/imstkDataArray.h b/Source/Common/imstkDataArray.h
index 7171597221de02444772c26a0fb983a3f62492c9..33aa00a05943896950648367c3f4bcbc1ba870b2 100644
--- a/Source/Common/imstkDataArray.h
+++ b/Source/Common/imstkDataArray.h
@@ -84,7 +84,7 @@ public:
         using iterator_category = std::forward_iterator_tag;
         using difference_type   = std::ptrdiff_t;
         using pointer   = T*;
-        using reference = T&;
+        using reference = const T&;
 
     public:
         const_iterator(pointer ptr) : ptr_(ptr) { }
@@ -98,7 +98,7 @@ public:
 
         self_type operator++(int junk) { ptr_++; return *this; }
 
-        const reference operator*() { return *ptr_; }
+        reference operator*() { return *ptr_; }
 
         const pointer operator->() { return ptr_; }
 
@@ -416,4 +416,4 @@ protected:
     bool m_mapped;
     T*   m_data;
 };
-}
\ No newline at end of file
+}
diff --git a/Source/Common/imstkMath.h b/Source/Common/imstkMath.h
index 4e4e01671786be293046fbb84a76c458f61dbf32..87a8fa4c10c9aa64760aa468cd881942a0efd469 100644
--- a/Source/Common/imstkMath.h
+++ b/Source/Common/imstkMath.h
@@ -53,22 +53,28 @@ make_unique(Args&& ... args)
 namespace imstk
 {
 // 2D vector
-using Vec2f = Eigen::Vector2f;
-using Vec2d = Eigen::Vector2d;
+// using Vec2f = Eigen::Vector2f;
+// using Vec2d = Eigen::Vector2d;
+using Vec2f = Eigen::Matrix<float, 2, 1>;
+using Vec2d = Eigen::Matrix<double, 2, 1>;
 using Vec2i = Eigen::Matrix<int, 2, 1>;
 using StdVectorOfVec2f = std::vector<Vec2f, Eigen::aligned_allocator<Vec2f>>;
 using StdVectorOfVec2d = std::vector<Vec2d, Eigen::aligned_allocator<Vec2d>>;
 
 // 3D vector
-using Vec3f = Eigen::Vector3f;
-using Vec3d = Eigen::Vector3d;
+// using Vec3f = Eigen::Vector3f;
+// using Vec3d = Eigen::Vector3d;
+using Vec3f = Eigen::Matrix<float, 3, 1>;
+using Vec3d = Eigen::Matrix<double, 3, 1>;
 using Vec3i = Eigen::Matrix<int, 3, 1>;
 using StdVectorOfVec3f = std::vector<Vec3f, Eigen::aligned_allocator<Vec3f>>;
 using StdVectorOfVec3d = std::vector<Vec3d, Eigen::aligned_allocator<Vec3d>>;
 
 // 4D vector
-using Vec4f = Eigen::Vector4f;
-using Vec4d = Eigen::Vector4d;
+// using Vec4f = Eigen::Vector4f;
+// using Vec4d = Eigen::Vector4d;
+using Vec4f = Eigen::Matrix<float, 4, 1>;
+using Vec4d = Eigen::Matrix<double, 4, 1>;
 using Vec4i = Eigen::Matrix<int, 4, 1>;
 using StdVectorOfVec4f = std::vector<Vec4f, Eigen::aligned_allocator<Vec4f>>;
 using StdVectorOfVec4d = std::vector<Vec4d, Eigen::aligned_allocator<Vec4d>>;
@@ -98,13 +104,13 @@ using Rotf = Eigen::AngleAxisf;
 using Rotd = Eigen::AngleAxisd;
 
 // 3x3 Matrix
-using Mat3f = Eigen::Matrix3f;
-using Mat3d = Eigen::Matrix3d;
+using Mat3f = Eigen::Matrix<float, 3, 3>;
+using Mat3d = Eigen::Matrix<double, 3, 3>;
 using StdVectorOfMat3d = std::vector<Mat3d, Eigen::aligned_allocator<Mat3d>>;
 
 // 4x4 Matrix
-using Mat4f = Eigen::Matrix4f;
-using Mat4d = Eigen::Matrix4d;
+using Mat4f = Eigen::Matrix<float, 4, 4> ;
+using Mat4d = Eigen::Matrix<double, 4, 4>;
 
 /// A dynamic size matrix of floats
 using Matrixf = Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic>;
diff --git a/Source/Common/imstkVecDataArray.h b/Source/Common/imstkVecDataArray.h
index 873126e3a814986db4d092913239ca557df68545..ea666aa34a45736290c822fad86ae772944ebedc 100644
--- a/Source/Common/imstkVecDataArray.h
+++ b/Source/Common/imstkVecDataArray.h
@@ -23,6 +23,7 @@
 
 #include "imstkDataArray.h"
 #include "imstkMath.h"
+#include "imstkLogger.h"
 //#include "imstkParallelReduce.h"
 
 namespace imstk
@@ -85,7 +86,7 @@ public:
         using iterator_category = std::forward_iterator_tag;
         using difference_type   = std::ptrdiff_t;
         using pointer   = VecType*;
-        using reference = VecType&;
+        using reference = const VecType&;
 
     public:
         const_iterator(pointer ptr) : ptr_(ptr) { }
@@ -99,7 +100,7 @@ public:
 
         self_type operator++(int junk) { ptr_++; return *this; }
 
-        const reference operator*() { return *ptr_; }
+        reference operator*() { return *ptr_; }
 
         const pointer operator->() { return ptr_; }
 
@@ -423,4 +424,4 @@ private:
     int      m_vecCapacity;
     VecType* m_dataCast;
 };
-}
\ No newline at end of file
+}
diff --git a/Source/Devices/imstkHapticDeviceClient.h b/Source/Devices/imstkHapticDeviceClient.h
index 73230d1b4435baf5d2196c2a95c0ff48dfee2a7f..2be1abea809b4aa7a37a6a1317043addff00135b 100644
--- a/Source/Devices/imstkHapticDeviceClient.h
+++ b/Source/Devices/imstkHapticDeviceClient.h
@@ -78,9 +78,13 @@ private:
     ///
     /// \brief Phantom omni device api callback
     ///
+#ifndef HDCALLBACK
+#define HDCALLBACK    
+#endif
+    typedef unsigned int HDCallbackCode;
     static HDCallbackCode HDCALLBACK hapticCallback(void* pData);
 
     HHD     m_handle; ///< device handle
     HDstate m_state;  ///< device reading state
 };
-}
\ No newline at end of file
+}
diff --git a/Source/DynamicalModels/ObjectModels/imstkFEMDeformableBodyModel.h b/Source/DynamicalModels/ObjectModels/imstkFEMDeformableBodyModel.h
index 9e2daaa1f67ee7d57690853cd2235e0c95139a51..59b13df204fdfdd43abed29e4aaffa859b15494c 100644
--- a/Source/DynamicalModels/ObjectModels/imstkFEMDeformableBodyModel.h
+++ b/Source/DynamicalModels/ObjectModels/imstkFEMDeformableBodyModel.h
@@ -71,8 +71,9 @@ struct FEMModelConfig
 ///
 class FEMDeformableBodyModel : public DynamicalModel<FeDeformBodyState>
 {
-using kinematicState = FeDeformBodyState;
-using System = NonLinearSystem<SparseMatrixd>;
+public:
+    using kinematicState = FeDeformBodyState;
+    using System = NonLinearSystem<SparseMatrixd>;
 
 public:
     ///
diff --git a/Source/DynamicalModels/ObjectModels/imstkRigidBodyModel2.h b/Source/DynamicalModels/ObjectModels/imstkRigidBodyModel2.h
index 252347d54a15e715472f4fef0a833263d5781377..13cd35df851dd520d7cc8b60c81911a25d34661b 100644
--- a/Source/DynamicalModels/ObjectModels/imstkRigidBodyModel2.h
+++ b/Source/DynamicalModels/ObjectModels/imstkRigidBodyModel2.h
@@ -158,4 +158,4 @@ protected:
 
     Eigen::VectorXd F;               // Reaction forces
 };
-}
\ No newline at end of file
+}
diff --git a/Source/Geometry/Implicit/imstkImplicitFunctionFiniteDifferenceFunctor.h b/Source/Geometry/Implicit/imstkImplicitFunctionFiniteDifferenceFunctor.h
index 1545901a2e9c854cb48b780e29a659850f6ec04f..c0d9798c9e56c012c1bee32445da973edd30b1cd 100644
--- a/Source/Geometry/Implicit/imstkImplicitFunctionFiniteDifferenceFunctor.h
+++ b/Source/Geometry/Implicit/imstkImplicitFunctionFiniteDifferenceFunctor.h
@@ -32,6 +32,8 @@ namespace imstk
 struct ImplicitFunctionGradient
 {
     public:
+        virtual ~ImplicitFunctionGradient() = default;
+
         virtual Vec3d operator()(const Vec3d& pos) const = 0;
 
         virtual void setDx(const Vec3d& dx)
@@ -57,6 +59,8 @@ struct ImplicitFunctionGradient
 struct ImplicitFunctionCentralGradient : public ImplicitFunctionGradient
 {
     public:
+        virtual ~ImplicitFunctionCentralGradient() = default;
+
         Vec3d operator()(const Vec3d& pos) const override
         {
             const ImplicitGeometry& funcRef = *m_func;
@@ -122,6 +126,8 @@ struct ImplicitFunctionForwardGradient : public ImplicitFunctionGradient
 struct ImplicitFunctionBackwardGradient : public ImplicitFunctionGradient
 {
     public:
+        virtual ~ImplicitFunctionBackwardGradient() = default;
+
         Vec3d operator()(const Vec3d& pos) const override
         {
             const ImplicitGeometry& funcRef      = *m_func;
@@ -142,6 +148,8 @@ struct ImplicitFunctionBackwardGradient : public ImplicitFunctionGradient
 struct StructuredForwardGradient : public ImplicitFunctionGradient
 {
     public:
+        virtual ~StructuredForwardGradient() = default;
+
         inline Vec3d operator()(const Vec3d& pos) const override
         {
             const SignedDistanceField& funcRef      = *static_cast<SignedDistanceField*>(m_func.get());
@@ -172,6 +180,8 @@ struct StructuredForwardGradient : public ImplicitFunctionGradient
 struct StructuredBackwardGradient : public ImplicitFunctionGradient
 {
     public:
+        virtual ~StructuredBackwardGradient() = default;
+
         inline Vec3d operator()(const Vec3d& pos) const override
         {
             const SignedDistanceField& funcRef      = *static_cast<SignedDistanceField*>(m_func.get());
@@ -270,4 +280,4 @@ struct ImplicitStructuredCurvature
     protected:
         Vec3i m_dxi;
 };
-}
\ No newline at end of file
+}
diff --git a/Source/Scene/imstkScene.cpp b/Source/Scene/imstkScene.cpp
index 631e164bf8d33e46ce264e44235b8af06c3e2005..fae6231a7b956ae1af67ab0d384ae98617f3f5d0 100644
--- a/Source/Scene/imstkScene.cpp
+++ b/Source/Scene/imstkScene.cpp
@@ -362,8 +362,9 @@ Scene::removeLight(const std::string& lightName)
 std::string
 Scene::getCameraName(const std::shared_ptr<Camera> cam) const
 {
+    using MapType = std::unordered_map<std::string, std::shared_ptr<Camera>>;
     auto i = std::find_if(m_cameras.begin(), m_cameras.end(),
-        [&cam](const NamedMap<Camera>::value_type& j) { return j.second == cam; });
+        [&cam](const MapType::value_type& j) { return j.second == cam; });
     if (i != m_cameras.end())
     {
         return i->first;
diff --git a/Source/Scene/imstkScene.h b/Source/Scene/imstkScene.h
index 96d99cc031bc948c2287675e11a407fd49a4f78f..cd6377feedb25a4a93d0c7a4c3a5a6539f17953a 100644
--- a/Source/Scene/imstkScene.h
+++ b/Source/Scene/imstkScene.h
@@ -78,10 +78,10 @@ struct SceneConfig
 ///
 class Scene : public EventObject
 {
-template<class T>
-using NamedMap = std::unordered_map<std::string, std::shared_ptr<T>>;
-
 public:
+    template<class T>
+    using NamedMap = std::unordered_map<std::string, std::shared_ptr<T>>;
+
     Scene(const std::string& name, std::shared_ptr<SceneConfig> config = std::make_shared<SceneConfig>());
     virtual ~Scene() override = default;
 
@@ -185,7 +185,8 @@ public:
     ///
     /// \brief Get and unordered map of cameras with names
     ///
-    const NamedMap<Camera>& getCameras() const { return m_cameras; }
+    const std::unordered_map<std::string, std::shared_ptr<Camera>>& 
+    getCameras() const { return m_cameras; }
 
     ///
     /// \brief Add light from the scene
@@ -288,10 +289,10 @@ protected:
 
     std::string m_name; ///> Name of the scene
     std::unordered_set<std::shared_ptr<SceneObject>> m_sceneObjects;
-    NamedMap<Light> m_lightsMap;
+    std::unordered_map<std::string, std::shared_ptr<Light>> m_lightsMap;
     std::shared_ptr<IBLProbe> m_globalIBLProbe = nullptr;
 
-    NamedMap<Camera> m_cameras;
+    std::unordered_map<std::string, std::shared_ptr<Camera>> m_cameras;
     std::shared_ptr<Camera> m_activeCamera;
 
     std::shared_ptr<CollisionGraph> m_collisionGraph;
diff --git a/Source/SimulationManager/imstkSimulationManager.h b/Source/SimulationManager/imstkSimulationManager.h
index a9835c1d0db4bb51fdb5a82f9b461192b2d28fc5..c40ddb8212ad9dc2e5a2a677e6a1feb74974c01b 100644
--- a/Source/SimulationManager/imstkSimulationManager.h
+++ b/Source/SimulationManager/imstkSimulationManager.h
@@ -53,6 +53,8 @@ public:
     };
 
 public:
+
+    SimulationManager() = default;
     virtual ~SimulationManager() override = default;
 
 public:
@@ -127,4 +129,4 @@ protected:
     double m_dt       = 0.0;    // Actual timestep
     int    m_numSteps = 0;
 };
-};
\ No newline at end of file
+};
diff --git a/Source/Wrappers/CMakeLists.txt b/Source/Wrappers/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..00824982f669003b7d4a54331768c351ff1b9e7e
--- /dev/null
+++ b/Source/Wrappers/CMakeLists.txt
@@ -0,0 +1,2 @@
+add_subdirectory(csharp)
+
diff --git a/Source/Wrappers/csharp/CMakeLists.txt b/Source/Wrappers/csharp/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..4705718402bf82c7468406f19dfe3a0b4b716eda
--- /dev/null
+++ b/Source/Wrappers/csharp/CMakeLists.txt
@@ -0,0 +1,52 @@
+include_directories(${CMAKE_CURRENT_SOURCE_DIR})
+include(UseSWIG)
+
+set_source_files_properties(SwigInterface/imstkCWrapper.i PROPERTIES CPLUSPLUS ON)
+
+swig_add_library(iMSTKCWrapper
+    TYPE SHARED
+    OUTPUT_DIR "${CMAKE_BINARY_DIR}/Source/Wrappers/csharp"
+    LANGUAGE csharp
+    SOURCES SwigInterface/imstkCWrapper.i)
+
+set_target_properties(iMSTKCWrapper PROPERTIES SWIG_COMPILE_OPTIONS "-namespace;imstk")
+if(iMSTK_SWIG_PINNED_ARRAY)
+    set_target_properties(iMSTKCWrapper PROPERTIES SWIG_COMPILE_DEFINITIONS "SWIG_PINNED_ARRAY")
+endif()
+if(iMSTK_USE_OpenHaptics)
+    set_target_properties(iMSTKCWrapper PROPERTIES SWIG_COMPILE_DEFINITIONS "iMSTK_USE_OpenHaptics")
+	target_compile_definitions(iMSTKCWrapper PRIVATE iMSTK_USE_OpenHaptics)
+endif()
+
+target_link_libraries(iMSTKCWrapper PUBLIC SimulationManager Filtering)
+target_compile_options(iMSTKCWrapper PRIVATE
+    $<$<OR:$<CXX_COMPILER_ID:Clang>,$<CXX_COMPILER_ID:AppleClang>,$<CXX_COMPILER_ID:GNU>>:
+        -Wall -Wno-unused-function>
+    $<$<CXX_COMPILER_ID:MSVC>:
+        -W4 -MP -wd4505 /bigobj>)
+
+# Under MSVC we can generate a project for the examples
+if (MSVC)
+  add_subdirectory(iMSTKCSharp)
+  if (${IMSTK_WRAPPER_BUILT})
+    add_subdirectory(Examples)
+  endif()
+endif()
+
+#-----------------------------------------------------------------------------
+# Install library
+#-----------------------------------------------------------------------------
+install( TARGETS iMSTKCWrapper EXPORT ${PROJECT_NAME}Targets
+    RUNTIME DESTINATION bin COMPONENT RuntimeLibraries
+    LIBRARY DESTINATION lib COMPONENT RuntimeLibraries
+    ARCHIVE DESTINATION lib COMPONENT Development)
+
+if (LINUX)
+    install(CODE "execute_process(COMMAND ${CMAKE_COMMAND} -E create_symlink iMSTKCWrapper.so ${CMAKE_INSTALL_PREFIX}/lib/libiMSTKCWrapper.so)")
+endif()
+
+# Install all ".cs" files
+install(
+    DIRECTORY "${CMAKE_BINARY_DIR}/Source/Wrappers/csharp/"
+    DESTINATION "${CMAKE_INSTALL_PREFIX}/include/iMSTKSharp/"
+    FILES_MATCHING PATTERN "*.cs")
diff --git a/Source/Wrappers/csharp/Examples/CMakeLists.txt b/Source/Wrappers/csharp/Examples/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..a09149a8815838e84b993c0117859a2a43de01fb
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/CMakeLists.txt
@@ -0,0 +1,37 @@
+project(CSharpExamples LANGUAGES CSharp)
+
+# Check https://stackoverflow.com/questions/52556785/code-generator-generating-its-own-cmake-files-and-targets/52714922#52714922
+# https://stackoverflow.com/questions/55713475/correct-use-of-vs-debugger-working-directory-etc
+
+function(add_example example_name sources)
+  
+  set(target "Example-${example_name}-CSharp")
+  add_executable(${target} ${sources})
+  set_property(TARGET ${target} PROPERTY VS_DOTNET_TARGET_FRAMEWORK_VERSION "v4.6.1")
+  set_property(TARGET ${target} PROPERTY VS_DOTNET_REFERENCES "System;iMSTKCSharp")
+  set_property(TARGET ${target} PROPERTY FOLDER Examples/CSharp)
+  
+  # Neither of these actually work, has to be set in the project explicitely 
+  set_property(TARGET ${target} PROPERTY VS_DEBUGGER_WORKING_DIRECTORY "${CMAKE_INSTALL_PREFIX}/bin")
+  set_property(TARGET ${target} PROPERTY VS_DEBUGGER_ENVIRONMENT "PATH=${CMAKE_INSTALL_PREFIX}/bin;%PATH%")
+  
+  add_custom_command(TARGET ${target} POST_BUILD
+                     COMMAND ${CMAKE_COMMAND} -E copy $<TARGET_FILE:${target}> ${CMAKE_INSTALL_PREFIX}/bin)
+  
+  add_dependencies(${target} iMSTKCSharp)
+  
+  target_link_libraries(${target} iMSTKCWrapper)
+endfunction()
+
+message (WARNING "Adding CSharp Examples, please note to run the examples from visual studio you need to set the working directory in the debug project properties to ${CMAKE_INSTALL_PREFIX}/bin")
+
+add_example(FemDeformable femDeformable.cs)
+add_example(PBDCloth pbdCloth.cs)
+add_example(PBDClothCollision pbdClothCollision.cs)
+add_example(PBDCollisionOneObject pbdCollisionOneObject.cs)
+add_example(PBDCutting pbdCutting.cs)
+add_example(PBDVolume pbdVolume.cs)
+add_example(RigidBody2 rigidBody2.cs)
+add_example(SDFHaptics sdfHaptics.cs)
+add_example(SPHFluid sphFluid.cs)
+add_example(TestGeometry testGeometry.cs)
\ No newline at end of file
diff --git a/Source/Wrappers/csharp/Examples/femDeformable.cs b/Source/Wrappers/csharp/Examples/femDeformable.cs
new file mode 100644
index 0000000000000000000000000000000000000000..af61fbdba386e5351171b8ba3dc7de2d5c31551e
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/femDeformable.cs
@@ -0,0 +1,164 @@
+using System;
+using imstk;
+
+
+public enum Geom
+{
+    Dragon,
+    Heart
+}
+
+public class Input
+{
+    public string meshFileName;
+    public VectorSizet fixedNodeIds;
+}
+
+// public class FeMain
+// {
+//
+// }
+
+public class FeDeformable
+{
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+        runFEMDeformable();
+    }
+
+    private static void runFEMDeformable()
+    {
+
+        input = new Input();
+        if (geom == Geom.Dragon)
+        {
+            input.meshFileName = "../data/asianDragon/asianDragon.veg";
+            input.fixedNodeIds = new VectorSizet(3);
+            input.fixedNodeIds.Add(50);
+            input.fixedNodeIds.Add(126);
+            input.fixedNodeIds.Add(177);
+        }
+        else if (geom == Geom.Heart)
+        {
+            input.meshFileName = "../data/textured_organs/heart_volume.vtk";
+            input.fixedNodeIds = new VectorSizet(13);
+            input.fixedNodeIds.Add(75);
+            input.fixedNodeIds.Add(82);
+            input.fixedNodeIds.Add(84);
+            input.fixedNodeIds.Add(94);
+            input.fixedNodeIds.Add(95);
+            input.fixedNodeIds.Add(105);
+            input.fixedNodeIds.Add(110);
+            input.fixedNodeIds.Add(124);
+            input.fixedNodeIds.Add(139);
+            input.fixedNodeIds.Add(150);
+            input.fixedNodeIds.Add(161);
+            input.fixedNodeIds.Add(171);
+            input.fixedNodeIds.Add(350);
+        }
+
+        // Construct the scene
+        Scene scene = new Scene("DeformableBodyFEM");
+        {
+            Camera cam = scene.getActiveCamera();
+            cam.setPosition(0.0, 2.0, -25.0);
+            cam.setFocalPoint(0.0, 0.0, 0.0);
+
+            // Load a tetrahedral mesh
+            TetrahedralMesh tetMesh = MeshIO.readTetrahedralMesh(input.meshFileName);
+            // CHECK(tetMesh != nullptr) << "Could not read mesh from file.";
+
+            // Scene object 1: fe-FeDeformableObject
+            FeDeformableObject deformableObj = makeFEDeformableObject(tetMesh);
+            scene.addSceneObject(deformableObj);
+
+            // Scene object 2: Plane
+            Plane planeGeom = new Plane();
+            planeGeom.setWidth(40.0);
+            planeGeom.setPosition(0.0, -8.0, 0.0);
+            CollidingObject planeObj = new CollidingObject("Plane");
+            planeObj.setVisualGeometry(planeGeom);
+            planeObj.setCollidingGeometry(planeGeom);
+            scene.addSceneObject(planeObj);
+
+            // Light
+            DirectionalLight light = new DirectionalLight();
+            light.setFocalPoint(new Vec3d(5.0, -8.0, -5.0));
+            light.setIntensity(1);
+            scene.addLight("light", light);
+        }
+
+        // Run the simulation
+        {
+            // Setup a viewer to render in its own thread
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene in its own thread
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.pause(); // Start simulation paused
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            driver.start();
+        }
+    }
+
+    private static FeDeformableObject makeFEDeformableObject(TetrahedralMesh tetMesh)
+    {
+        SurfaceMesh surfMesh = tetMesh.extractSurfaceMesh();
+        // surfMesh.flipNormals();
+
+        // Configure dynamic model
+        FEMDeformableBodyModel dynaModel = new FEMDeformableBodyModel();
+        FEMModelConfig config = new FEMModelConfig();
+        config.m_fixedNodeIds = input.fixedNodeIds;
+        dynaModel.configure(config);
+        //dynaModel.configure(iMSTK_DATA_ROOT "/asianDragon/asianDragon.config");
+
+        dynaModel.setTimeStepSizeType(TimeSteppingType.Fixed);
+        dynaModel.setModelGeometry(tetMesh);
+        BackwardEuler timeIntegrator = new BackwardEuler(0.01); // Create and add Backward Euler time integrator
+        dynaModel.setTimeIntegrator(timeIntegrator);
+
+        RenderMaterial mat = new RenderMaterial();
+        mat.setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+        mat.setPointSize(10.0f);
+        mat.setLineWidth(2.0f);
+        mat.setEdgeColor(Color.Orange);
+        VisualModel surfMeshModel = new VisualModel(surfMesh);
+        surfMeshModel.setRenderMaterial(mat);
+
+        // Scene object 1: Dragon
+        FeDeformableObject deformableObj = new FeDeformableObject("Dragon");
+        deformableObj.addVisualModel(surfMeshModel);
+        deformableObj.setPhysicsGeometry(tetMesh);
+        // Map simulated geometry to visual
+        deformableObj.setPhysicsToVisualMap(new OneToOneMap(tetMesh, surfMesh));
+        deformableObj.setDynamicalModel(dynaModel);
+
+        return deformableObj;
+    }
+
+    private const Geom geom = Geom.Heart;
+    private static Input input;
+}
+
+// femDefomrbale = new FeDeformable();
diff --git a/Source/Wrappers/csharp/Examples/pbdCloth.cs b/Source/Wrappers/csharp/Examples/pbdCloth.cs
new file mode 100644
index 0000000000000000000000000000000000000000..fd788123ee5eca739483cbe9c5c72b3bb1efdd1a
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/pbdCloth.cs
@@ -0,0 +1,235 @@
+using System;
+using System.Runtime.InteropServices;
+using imstk;
+
+public class PbdCloth
+{
+     public class CSReceiverFunc : ReceiverFunc {
+         public CSReceiverFunc(Action<KeyEvent> action) {
+             action_ = action;
+         }
+         public override void call(KeyEvent e) {
+             action_(e);
+         }
+         private Action<KeyEvent> action_;
+     }
+
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+
+        RunPbdCloth();
+    }
+
+    public static PbdObject makeClothObj(string name, double width, double height, int rowCount, int colCount)
+    {
+        PbdObject clothObj = new PbdObject(name);
+
+        // Setup the Geometry
+        SurfaceMesh clothMesh = makeClothGeometry(10.0, 10.0, 16, 16, 2.0);
+
+        // Setup the Parameters
+        PBDModelConfig pbdParams = new PBDModelConfig();
+        pbdParams.enableConstraint(PbdConstraint.Type.Distance, 1.0e2);
+        pbdParams.enableConstraint(PbdConstraint.Type.Dihedral, 1.0e1);
+        pbdParams.m_fixedNodeIds = new VectorSizet(2);
+        pbdParams.m_fixedNodeIds.Add(0);
+        pbdParams.m_fixedNodeIds.Add((uint)colCount - 1);
+        pbdParams.m_uniformMassValue = width * height / (rowCount * colCount);
+        pbdParams.m_gravity    = new Vec3d(0.0, -9.8, 0.0);
+        pbdParams.m_dt  = 0.005;
+        pbdParams.m_iterations = 5;
+
+        // Setup the Model
+        PbdModel pbdModel = new PbdModel();
+        pbdModel.setModelGeometry(clothMesh);
+        pbdModel.configure(pbdParams);
+
+        // Setup the VisualModel
+        RenderMaterial material = new RenderMaterial();
+        material.setBackFaceCulling(false);
+        material.setDisplayMode(RenderMaterial.DisplayMode.Surface);
+        // material.setDisplayMode(RenderMaterial.DisplayMode.Wireframe);
+        material.setShadingModel(RenderMaterial.ShadingModel.PBR);
+
+        setFleshTextures(material);
+        // setFabricTextures(material);
+        VisualModel visualModel = new VisualModel(clothMesh);
+        visualModel.setRenderMaterial(material);
+
+        // Setup the Object
+        clothObj.addVisualModel(visualModel);
+        clothObj.setPhysicsGeometry(clothMesh);
+        clothObj.setDynamicalModel(pbdModel);
+
+        return clothObj;
+    }
+
+    static SurfaceMesh makeClothGeometry(double width,
+                                         double height,
+                                         int    nRows,
+                                         int    nCols,
+                                         double uvScale)
+    {
+        SurfaceMesh clothMesh = new SurfaceMesh();
+
+        VecDataArray3d vertices = new VecDataArray3d(nRows * nCols);
+        double dy = width / (nCols - 1);
+        double dx = height / (nRows - 1);
+        for (int i = 0; i < nRows; ++i)
+        {
+            Vec3d xyz = new Vec3d();
+            for (int j = 0; j < nCols; j++)
+            {
+                xyz[0] = dx * i;
+                xyz[1] = 1.0;
+                xyz[2] = dy * j;
+                vertices[(uint)(i * nCols + j)] = xyz;
+            }
+        }
+
+        // Add connectivity data
+        VecDataArray3i indices = new VecDataArray3i();
+        for (int i = 0; i < nRows - 1; i++)
+        {
+            for (int j = 0; j < nCols - 1; j++)
+            {
+                int index1 = i * nCols + j;
+                int index2 = index1 + nCols;
+                int index3 = index1 + 1;
+                int index4 = index2 + 1;
+
+                // Interleave [/][\]
+                if (i % 2 != 0 || j % 2 != 0)
+                {
+                    indices.push_back(new Vec3i(index1, index2, index3));
+                    indices.push_back(new Vec3i(index4, index3, index2));
+                }
+                else
+                {
+                    indices.push_back(new Vec3i(index2, index4, index1));
+                    indices.push_back(new Vec3i(index4, index3, index1));
+                }
+            }
+        }
+
+        VecDataArray2f uvCoords = new VecDataArray2f(nRows * nCols);
+        for (uint i = 0; i < nRows; ++i)
+        {
+            for (uint j = 0; j < nCols; j++)
+            {
+                uvCoords[(uint)(i * nCols + j)] = new Vec2f((float)i / nRows, (float)j / nCols) * (float)uvScale;
+            }
+        }
+
+        clothMesh.initialize(vertices, indices);
+        clothMesh.setVertexTCoords("uvs", uvCoords);
+
+        return clothMesh;
+    }
+
+    public static void RunPbdCloth()
+    {
+        // Setup a scene
+        Scene scene = new Scene("PBDCloth");
+        PbdObject clothObj = makeClothObj("Cloth", 10.0, 10.0, 16, 16);
+        scene.addSceneObject(clothObj);
+
+        // Adjust camera
+        scene.getActiveCamera().setFocalPoint(0.0, -5.0, 5.0);
+        scene.getActiveCamera().setPosition(-15.0, -5.0, 25.0);
+
+        // Run the simulation
+        {
+            // Setup a viewer to render
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setExecutionType(Module.ExecutionType.ADAPTIVE);
+            sceneManager.setActiveScene(scene);
+            sceneManager.pause(); // Start simulation paused
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+            driver.setDesiredDt(0.001);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            Action<KeyEvent> receiverAction = (KeyEvent e) => {
+                // Set new textures
+                if (e.m_key == '1')
+                {
+                    setFleshTextures(clothObj.getVisualModel(0).getRenderMaterial());
+                }
+                else if (e.m_key == '2')
+                {
+                    setFabricTextures(clothObj.getVisualModel(0).getRenderMaterial());
+                }
+                // Darken the texture pixel values
+                else if (e.m_key == 'h')
+                {
+                    ImageData imageData = clothObj.getVisualModel(0).getRenderMaterial().getTexture(Texture.Type.Diffuse).getImageData();
+                    // std.shared_ptr<VecDataArray<unsigned char, 3>> scalars = std.dynamic_pointer_cast<VecDataArray<unsigned char, 3>>(imageData.getScalars());
+                    // VecDataArray3uc scalars = Utils.castToVecDataArray3uc(imageData.getScalars());
+                    VecDataArray3uc scalars = Utils.CastTo<VecDataArray3uc>(imageData.getScalars());
+                    // Console.WriteLine("{0} = {1}", scalars.size(), scalars2.size());
+                    byte[] newScalars = new byte[3*scalars.size()];
+                    scalars.getValues(newScalars);
+
+                    // Vec3uc scalarPtr = scalars.getPointer();
+                    for (int i = 0; i < newScalars.Length; i++)
+                    {
+                        // scalarPtr[i] = (scalarPtr[i].cast<double>() * 0.8).cast<unsigned char>();
+                        newScalars[i] = (byte)((double)newScalars[i] *0.8);
+                    }
+                    scalars.setValues(newScalars);
+                    // Console.WriteLine("'h' key is pressed to change the render material...");
+                    clothObj.getVisualModel(0).getRenderMaterial().getTexture(Texture.Type.Diffuse).postModified();
+                }
+            };
+            CSReceiverFunc receiverFunc = new CSReceiverFunc(receiverAction);
+            Utils.queueConnectKeyEvent(viewer.getKeyboardDevice(), Utils.KeyboardDeviceClient_getKeyPress_cb, sceneManager, receiverFunc);
+
+            driver.start();
+        }
+
+    }
+
+    private static void setFabricTextures(RenderMaterial material)
+    {
+        // ImageData diffuseTex = MeshIO.readImageData("../../../install/data/textures/fabricDiffuse.jpg");
+        ImageData diffuseTex = MeshIO.readImageData(dataPath + "textures/fabricDiffuse.jpg");
+        material.addTexture(new Texture(diffuseTex, Texture.Type.Diffuse));
+        ImageData normalTex = MeshIO.readImageData(dataPath + "textures/fabricNormal.jpg");
+        material.addTexture(new Texture(normalTex, Texture.Type.Normal));
+        ImageData ormTex = MeshIO.readImageData(dataPath + "textures/fabricORM.jpg");
+        material.addTexture(new Texture(ormTex, Texture.Type.ORM));
+    }
+
+    private static void setFleshTextures(RenderMaterial material)
+    {
+        ImageData diffuseTex = MeshIO.readImageData(dataPath + "textures/fleshDiffuse.jpg");
+        material.addTexture(new Texture(diffuseTex, Texture.Type.Diffuse));
+        ImageData normalTex = MeshIO.readImageData(dataPath + "textures/fleshNormal.jpg");
+        material.addTexture(new Texture(normalTex, Texture.Type.Normal));
+        ImageData ormTex = MeshIO.readImageData(dataPath + "textures/fleshORM.jpg");
+        material.addTexture(new Texture(ormTex, Texture.Type.ORM));
+    }
+
+    private static string dataPath = "../data/";
+
+}
diff --git a/Source/Wrappers/csharp/Examples/pbdClothCollision.cs b/Source/Wrappers/csharp/Examples/pbdClothCollision.cs
new file mode 100644
index 0000000000000000000000000000000000000000..13654219ace4ac9b0356ed334b719b6de202a2dc
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/pbdClothCollision.cs
@@ -0,0 +1,324 @@
+using System;
+using System.Runtime.InteropServices;
+using imstk;
+
+public class PbdCloth
+{
+     public class CSReceiverFunc : ReceiverFunc {
+         public CSReceiverFunc(Action<KeyEvent> action) {
+             action_ = action;
+         }
+         public override void call(KeyEvent e) {
+             action_(e);
+         }
+         private Action<KeyEvent> action_;
+     }
+
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+
+        RunPbdCloth();
+    }
+
+    public static PbdObject makeClothObj(string name, double width, double height, int rowCount, int colCount)
+    {
+        PbdObject clothObj = new PbdObject(name);
+
+        // Setup the Geometry
+        SurfaceMesh clothMesh = makeClothGeometry(width, height, rowCount, colCount, 2.0);
+
+        // Setup the Parameters
+        PBDModelConfig pbdParams = new PBDModelConfig();
+        pbdParams.enableConstraint(PbdConstraint.Type.Distance, 1.0e2);
+        pbdParams.enableConstraint(PbdConstraint.Type.Dihedral, 1.0e1);
+        pbdParams.m_uniformMassValue = width * height / (rowCount * colCount);
+        pbdParams.m_gravity    = new Vec3d(0.0, -9.8, 0.0);
+        pbdParams.m_dt  = 0.005;
+        pbdParams.m_iterations = 5;
+
+        // Setup the Model
+        PbdModel pbdModel = new PbdModel();
+        pbdModel.setModelGeometry(clothMesh);
+        pbdModel.configure(pbdParams);
+
+        // Setup the VisualModel
+        RenderMaterial material = new RenderMaterial();
+        material.setBackFaceCulling(false);
+        material.setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+        material.setColor(Color.Blue);
+
+        // setFabricTextures(material);
+        VisualModel visualModel = new VisualModel(clothMesh);
+        visualModel.setRenderMaterial(material);
+
+        // Setup the Object
+        clothObj.addVisualModel(visualModel);
+        clothObj.setPhysicsGeometry(clothMesh);
+        clothObj.setCollidingGeometry(clothMesh);
+        clothObj.setDynamicalModel(pbdModel);
+
+        return clothObj;
+    }
+
+    static SurfaceMesh makeClothGeometry(double width,
+                                         double height,
+                                         int    nRows,
+                                         int    nCols,
+                                         double uvScale)
+    {
+        SurfaceMesh clothMesh = new SurfaceMesh();
+
+        VecDataArray3d vertices = new VecDataArray3d(nRows * nCols);
+        double dy = width / (nCols - 1);
+        double dx = height / (nRows - 1);
+        Vec3d halfSize = new Vec3d(height, 0.0, width) * 0.5;
+        for (int i = 0; i < nRows; ++i)
+        {
+            Vec3d xyz = new Vec3d();
+            for (int j = 0; j < nCols; j++)
+            {
+                xyz[0] = dx * i;
+                xyz[1] = 0.05;
+                xyz[2] = dy * j - 1.0;
+                vertices[(uint)(i * nCols + j)] = xyz - halfSize;
+            }
+        }
+
+        // Add connectivity data
+        VecDataArray3i indices = new VecDataArray3i();
+        for (int i = 0; i < nRows - 1; i++)
+        {
+            for (int j = 0; j < nCols - 1; j++)
+            {
+                int index1 = i * nCols + j;
+                int index2 = index1 + nCols;
+                int index3 = index1 + 1;
+                int index4 = index2 + 1;
+
+                // Interleave [/][\]
+                if (i % 2 != 0 || j % 2 != 0)
+                {
+                    indices.push_back(new Vec3i(index1, index2, index3));
+                    indices.push_back(new Vec3i(index4, index3, index2));
+                }
+                else
+                {
+                    indices.push_back(new Vec3i(index2, index4, index1));
+                    indices.push_back(new Vec3i(index4, index3, index1));
+                }
+            }
+        }
+
+        VecDataArray2f uvCoords = new VecDataArray2f(nRows * nCols);
+        for (uint i = 0; i < nRows; ++i)
+        {
+            for (uint j = 0; j < nCols; j++)
+            {
+                uvCoords[(uint)(i * nCols + j)] = new Vec2f((float)i / nRows, (float)j / nCols) * (float)uvScale;
+            }
+        }
+
+        clothMesh.initialize(vertices, indices);
+        clothMesh.setVertexTCoords("uvs", uvCoords);
+
+        return clothMesh;
+    }
+
+    public static void RunPbdCloth()
+    {
+        Capsule capsule = new Capsule(new Vec3d(0.0, -4.0, 0.0), 2.0, 5.0, new Quatd(new Rotd(1.5708, new Vec3d(0.0, 0.0, 1.0))));
+        Sphere sphere = new Sphere(new Vec3d(0.0, -2.0, 0.0), 2.0);
+        OrientedBox cube = new OrientedBox(new Vec3d(0.0, -4.0, 0.0), new Vec3d(2.5, 2.5, 2.5));
+        Plane plane = new Plane(new Vec3d(0.0, -2.0, 0.0), new Vec3d(0.0, 1.0, 0.0));
+        plane.setWidth(20.0);
+
+        Geometry[] geometries = new Geometry[] {capsule, sphere, cube, plane};
+
+        // Setup a scene
+        Scene scene = new Scene("PBDClothCollision");
+        PbdObject clothObj = makeClothObj("Cloth", 10.0, 10.0, 16, 16);
+        scene.addSceneObject(clothObj);
+        CollidingObject collisionObj = new CollidingObject("test");
+        collisionObj.setCollidingGeometry(capsule);
+
+        for (int i=0; i<4; ++i)
+        {
+            VisualModel visualModel = new VisualModel(geometries[i]);
+            visualModel.getRenderMaterial().setBackFaceCulling(false);
+            visualModel.getRenderMaterial().setOpacity(0.5);
+            visualModel.hide();
+            collisionObj.addVisualModel(visualModel);
+        }
+        collisionObj.getVisualModel(0).show();
+        scene.addSceneObject(collisionObj);
+
+        PbdObjectCollision pbdInteraction = new PbdObjectCollision(clothObj, collisionObj, "PointSetToCapsuleCD");
+        pbdInteraction.setFriction(0.4);
+        pbdInteraction.setRestitution(0.0);
+        scene.getCollisionGraph().addInteraction(pbdInteraction);
+
+
+        // Adjust camera
+        scene.getActiveCamera().setFocalPoint(0.0, -2.0, 0.0);
+        scene.getActiveCamera().setPosition(5.0, 4.0, 18.0);
+
+        // Run the simulation
+        {
+            // Setup a viewer to render
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setExecutionType(Module.ExecutionType.ADAPTIVE);
+            sceneManager.setActiveScene(scene);
+            sceneManager.pause(); // Start simulation paused
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+            driver.setDesiredDt(0.001);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            Action<KeyEvent> receiverAction = (KeyEvent e) => {
+                // Switch to sphere and reset
+                if (e.m_key == '1')
+                {
+                    for (uint i = 0; i < 4; i++)
+                    {
+                        collisionObj.getVisualModel(i).hide();
+                    }
+                    collisionObj.getVisualModel(0).show();
+                    collisionObj.setCollidingGeometry(capsule);
+
+                    PointSetToCapsuleCD capsuleCD = new PointSetToCapsuleCD();
+                    capsuleCD.setInputGeometryA(clothObj.getCollidingGeometry());
+                    capsuleCD.setInputGeometryB(capsule);
+                    pbdInteraction.setCollisionDetection(capsuleCD);
+                    pbdInteraction.getCollisionHandlingA().setInputCollisionData(capsuleCD.getCollisionData());
+
+                    scene.buildTaskGraph();
+                    scene.reset();
+                }
+                // Switch to capsule and reset
+                else if (e.m_key == '2')
+                {
+                    for (uint i = 0; i < 4; i++)
+                    {
+                        collisionObj.getVisualModel(i).hide();
+                    }
+                    collisionObj.getVisualModel(1).show();
+                    collisionObj.setCollidingGeometry(sphere);
+
+                    PointSetToSphereCD sphereCD = new PointSetToSphereCD();
+                    sphereCD.setInputGeometryA(clothObj.getCollidingGeometry());
+                    sphereCD.setInputGeometryB(sphere);
+                    pbdInteraction.setCollisionDetection(sphereCD);
+                    pbdInteraction.getCollisionHandlingA().setInputCollisionData(sphereCD.getCollisionData());
+
+                    scene.buildTaskGraph();
+                    scene.reset();
+                }
+                // Switch to cube and reset
+                else if (e.m_key == '3')
+                {
+                    for (uint i = 0; i < 4; i++)
+                    {
+                        collisionObj.getVisualModel(i).hide();
+                    }
+                    collisionObj.getVisualModel(2).show();
+                    collisionObj.setCollidingGeometry(cube);
+
+                    PointSetToOrientedBoxCD cubeCD = new PointSetToOrientedBoxCD();
+                    cubeCD.setInputGeometryA(clothObj.getCollidingGeometry());
+                    cubeCD.setInputGeometryB(cube);
+                    pbdInteraction.setCollisionDetection(cubeCD);
+                    pbdInteraction.getCollisionHandlingA().setInputCollisionData(cubeCD.getCollisionData());
+
+                    scene.buildTaskGraph();
+                    scene.reset();
+                }
+                // Switch to plane and reset
+                else if (e.m_key == '4')
+                {
+                    for (uint i = 0; i < 4; i++)
+                    {
+                        collisionObj.getVisualModel(i).hide();
+                    }
+                    collisionObj.getVisualModel(3).show();
+                    collisionObj.setCollidingGeometry(plane);
+
+                    PointSetToPlaneCD planeCD = new PointSetToPlaneCD();
+                    planeCD.setInputGeometryA(clothObj.getCollidingGeometry());
+                    planeCD.setInputGeometryB(plane);
+                    pbdInteraction.setCollisionDetection(planeCD);
+                    pbdInteraction.getCollisionHandlingA().setInputCollisionData(planeCD.getCollisionData());
+
+                    scene.buildTaskGraph();
+                    scene.reset();
+                }
+                // Switch to sphere vs surface and reset
+                else if (e.m_key == '5')
+                {
+                    for (uint i = 0; i < 4; i++)
+                    {
+                        collisionObj.getVisualModel(i).hide();
+                    }
+                    collisionObj.getVisualModel(1).show();
+                    collisionObj.setCollidingGeometry(sphere);
+
+                    SurfaceMeshToSphereCD sphereCD = new SurfaceMeshToSphereCD();
+                    sphereCD.setInputGeometryA(clothObj.getCollidingGeometry());
+                    sphereCD.setInputGeometryB(sphere);
+                    pbdInteraction.setCollisionDetection(sphereCD);
+                    pbdInteraction.getCollisionHandlingA().setInputCollisionData(sphereCD.getCollisionData());
+
+                    scene.buildTaskGraph();
+                    scene.reset();
+                }
+                // Switch to sphere vs surface and reset
+                else if (e.m_key == '6')
+                {
+                    for (uint i = 0; i < 4; i++)
+                    {
+                        collisionObj.getVisualModel(i).hide();
+                    }
+                    collisionObj.getVisualModel(0).show();
+                    collisionObj.setCollidingGeometry(capsule);
+
+                    SurfaceMeshToCapsuleCD capsuleCD = new SurfaceMeshToCapsuleCD();
+                    capsuleCD.setInputGeometryA(clothObj.getCollidingGeometry());
+                    capsuleCD.setInputGeometryB(capsule);
+                    pbdInteraction.setCollisionDetection(capsuleCD);
+                    pbdInteraction.getCollisionHandlingA().setInputCollisionData(capsuleCD.getCollisionData());
+
+                    scene.buildTaskGraph();
+                    scene.reset();
+                }
+            };
+
+            CSReceiverFunc receiverFunc = new CSReceiverFunc(receiverAction);
+            Utils.queueConnectKeyEvent(viewer.getKeyboardDevice(), Utils.KeyboardDeviceClient_getKeyPress_cb, sceneManager, receiverFunc);
+
+            driver.start();
+        }
+
+    }
+
+    private static string dataPath = "../data/";
+
+}
+
diff --git a/Source/Wrappers/csharp/Examples/pbdCollisionOneObject.cs b/Source/Wrappers/csharp/Examples/pbdCollisionOneObject.cs
new file mode 100644
index 0000000000000000000000000000000000000000..0b6e5e168e2e8a558d97b3c5ec73e5cb227b3d0a
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/pbdCollisionOneObject.cs
@@ -0,0 +1,194 @@
+using imstk;
+
+public class PbdCollisionOneObject
+{
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+        runPbdCollisionOneObject();
+    }
+
+
+    private static void runPbdCollisionOneObject()
+    {
+        // Setup the scene
+        Scene scene = new Scene("PbdCollisionOneDragon");
+        {
+            scene.getActiveCamera().setPosition(0, 3.0, 20.0);
+            scene.getActiveCamera().setFocalPoint(0.0, -10.0, 0.0);
+
+            // set up the meshes
+            string tetMeshFileName  = dataPath + "asianDragon/asianDragon.veg";
+            TetrahedralMesh coarseTetMesh   = MeshIO.readTetrahedralMesh(tetMeshFileName);
+            string surfMeshFileName = dataPath + "asianDragon/asianDragon.obj";
+            // string surfMeshFileName = "/home/jianfeng/Documents/imstk/build_csharp/install/data/asianDragon/asianDragon.obj";
+            SurfaceMesh highResSurfMesh = MeshIO.readSurfaceMesh(surfMeshFileName);
+            SurfaceMesh coarseSurfMesh = coarseTetMesh.extractSurfaceMesh();
+
+            // set up visual model based on high res mesh
+            RenderMaterial material = new RenderMaterial();
+            material.setDisplayMode(RenderMaterial.DisplayMode.Surface);
+            material.setLineWidth(0.5f);
+            material.setEdgeColor(Color.Blue);
+            material.setShadingModel(RenderMaterial.ShadingModel.Phong);
+            VisualModel surfMeshModel = new VisualModel(highResSurfMesh);
+            surfMeshModel.setRenderMaterial(material);
+
+            // configure the deformable object
+            PbdObject deformableObj = new PbdObject("DeformableObj");
+            deformableObj.addVisualModel(surfMeshModel);
+            deformableObj.setCollidingGeometry(coarseSurfMesh);
+            deformableObj.setPhysicsGeometry(coarseTetMesh);
+            deformableObj.setPhysicsToCollidingMap(new OneToOneMap(coarseTetMesh, coarseSurfMesh));
+            deformableObj.setPhysicsToVisualMap(new TetraTriangleMap(coarseTetMesh, highResSurfMesh));
+
+            // Create model and object
+            PbdModel pbdModel = new PbdModel();
+            pbdModel.setModelGeometry(coarseTetMesh);
+
+            // configure model
+            PBDModelConfig pbdParams = new PBDModelConfig();
+
+            // FEM constraint
+            pbdParams.m_femParams.m_YoungModulus = youngModulus;
+            pbdParams.m_femParams.m_PoissonRatio = poissonRatio;
+            pbdParams.enableFEMConstraint(PbdConstraint.Type.FEMTet, PbdFEMConstraint.MaterialType.Corotation);
+
+            // Other parameters
+            // \todo use lumped mass
+            pbdParams.m_uniformMassValue = 1.0;
+            pbdParams.m_gravity    = new Vec3d(0, -10.0, 0);
+            pbdParams.m_dt  = timeStep;
+            pbdParams.m_iterations = maxIter;
+            // pbdParams.collisionParams.m_proximity = 0.3;
+            // pbdParams.collisionParams.m_stiffness = 0.1;
+
+            pbdModel.configure(pbdParams);
+            deformableObj.setDynamicalModel(pbdModel);
+
+            scene.addSceneObject(deformableObj);
+
+            // Build floor geometry
+            SurfaceMesh floorMesh = createUniformSurfaceMesh(100.0, 100.0, 2, 2);
+
+            RenderMaterial floorMaterial = new RenderMaterial();
+            floorMaterial.setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+            VisualModel floorVisualModel = new VisualModel(floorMesh);
+            floorVisualModel.setRenderMaterial(floorMaterial);
+
+            PbdObject floorObj = new PbdObject("Floor");
+            floorObj.setCollidingGeometry(floorMesh);
+            floorObj.setPhysicsGeometry(floorMesh);
+            floorObj.addVisualModel(floorVisualModel);
+
+            PbdModel floorPbdModel = new PbdModel();
+            floorPbdModel.setModelGeometry(floorMesh);
+
+            // configure model
+            PBDModelConfig floorPbdParams = new PBDModelConfig();
+            floorPbdParams.m_uniformMassValue = 0.0;
+            floorPbdParams.m_iterations       = 0;
+            // floorPbdParams.collisionParams.m_proximity = -0.1;
+
+            // Set the parameters
+            floorPbdModel.configure(floorPbdParams);
+            floorObj.setDynamicalModel(floorPbdModel);
+
+            scene.addSceneObject(floorObj);
+
+            // Collision
+            // scene.getCollisionGraph().addInteraction(Utils.makeObjectInteractionPair(deformableObj, floorObj,
+                                                                                 // InteractionType.PbdObjToPbdObjCollision, CollisionDetection.Type.MeshToMeshBruteForce));
+            PbdObjectCollision pbdInteraction = new PbdObjectCollision(deformableObj, floorObj, "PointSetToPlaneCD");
+            MeshToMeshBruteForceCD cd = new MeshToMeshBruteForceCD();
+            cd.setInputGeometryA(deformableObj.getCollidingGeometry());
+            cd.setInputGeometryB(floorObj.getCollidingGeometry());
+            pbdInteraction.setCollisionDetection(cd);
+            pbdInteraction.getCollisionHandlingA().setInputCollisionData(cd.getCollisionData());
+
+            scene.getCollisionGraph().addInteraction(pbdInteraction);
+            scene.buildTaskGraph();
+            scene.reset();
+
+            // Light
+            DirectionalLight light = new DirectionalLight();
+            light.setFocalPoint(new Vec3d(5, -8, -5));
+            light.setIntensity(1);
+            scene.addLight("light", light);
+        }
+
+        // Run the simulation
+        {
+            // Setup a viewer to render
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.pause(); // Start simulation paused
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            driver.start();
+        }
+    }
+
+    private static SurfaceMesh createUniformSurfaceMesh(double width, double height, int nRows, int nCols)
+    {
+        double dy = width / (double)(nCols - 1);
+        double dx = height / (double)(nRows - 1);
+
+        VecDataArray3d vertices = new VecDataArray3d(nRows * nCols);
+
+        for (int i = 0; i < nRows; ++i)
+        {
+            for (int j = 0; j < nCols; j++)
+            {
+                double y = (double)(dy * j);
+                double x = (double)(dx * i);
+                vertices[(uint)(i * nCols + j)] = new Vec3d(x - height * 0.5, -10.0, y - width * 0.5);
+            }
+        }
+
+        // c. Add connectivity data
+        VecDataArray3i triangles = new VecDataArray3i();
+        for (int i = 0; i < nRows - 1; ++i)
+        {
+            for (int j = 0; j < nCols - 1; j++)
+            {
+                triangles.push_back(new Vec3i(i * nCols + j, i * nCols + j + 1, (i + 1) * nCols + j));
+                triangles.push_back(new Vec3i((i + 1) * nCols + j + 1, (i + 1) * nCols + j, i * nCols + j + 1));
+            }
+        }
+
+        SurfaceMesh surfMesh = new SurfaceMesh();
+        surfMesh.initialize(vertices, triangles);
+
+        return surfMesh;
+    }
+
+    private static string dataPath = "../data/";
+    // parameters to play with
+    private const double youngModulus     = 1000.0;
+    private const double poissonRatio     = 0.3;
+    private const double timeStep         = 0.01;
+    private const double contactStiffness = 0.1;
+    private const int    maxIter = 5;
+}
+
diff --git a/Source/Wrappers/csharp/Examples/pbdCutting.cs b/Source/Wrappers/csharp/Examples/pbdCutting.cs
new file mode 100644
index 0000000000000000000000000000000000000000..aa6cd307fe9b7e65370a9923eea6b4b81e4e4cd4
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/pbdCutting.cs
@@ -0,0 +1,203 @@
+using System;
+using imstk;
+
+public class PbdCutting
+{
+     public class CSReceiverFunc : ReceiverFunc {
+         public CSReceiverFunc(Action<KeyEvent> action) {
+             action_ = action;
+         }
+         public override void call(KeyEvent e) {
+             action_(e);
+         }
+         private Action<KeyEvent> action_;
+     }
+
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+
+        RunPbcCutting();
+    }
+
+    public static PbdObject makeClothObj(string name, double width, double height, int rowCount, int colCount)
+    {
+        PbdObject clothObj = new PbdObject(name);
+
+        // Setup the Geometry
+        SurfaceMesh clothMesh = makeClothGeometry(width, height, nRows, nCols);
+
+        // Setup the Parameters
+        PBDModelConfig pbdParams = new PBDModelConfig();
+        pbdParams.enableConstraint(PbdConstraint.Type.Distance, 1.0e3);
+        pbdParams.enableConstraint(PbdConstraint.Type.Dihedral, 1.0e3);
+        pbdParams.m_fixedNodeIds = new VectorSizet(2);
+        pbdParams.m_fixedNodeIds.Add(0);
+        pbdParams.m_fixedNodeIds.Add((uint)colCount - 1);
+        pbdParams.m_uniformMassValue = width * height / (rowCount * colCount);
+        pbdParams.m_gravity    = new Vec3d(0.0, -9.8, 0.0);
+        pbdParams.m_dt  = 0.005;
+        pbdParams.m_iterations = 5;
+
+        // Setup the Model
+        PbdModel pbdModel = new PbdModel();
+        pbdModel.setModelGeometry(clothMesh);
+        pbdModel.configure(pbdParams);
+
+        // Setup the VisualModel
+        RenderMaterial material = new RenderMaterial();
+        material.setBackFaceCulling(false);
+        material.setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+
+        VisualModel visualModel = new VisualModel(clothMesh);
+        visualModel.setRenderMaterial(material);
+
+        // Setup the Object
+        clothObj.addVisualModel(visualModel);
+        clothObj.setPhysicsGeometry(clothMesh);
+        clothObj.setCollidingGeometry(clothMesh);
+        clothObj.setDynamicalModel(pbdModel);
+
+        return clothObj;
+    }
+
+    static SurfaceMesh makeClothGeometry(double width,
+                                         double height,
+                                         int    nRows,
+                                         int    nCols)
+    {
+        SurfaceMesh clothMesh = new SurfaceMesh("Cloth_SurfaceMesh");
+
+        VecDataArray3d vertices = new VecDataArray3d(nRows * nCols);
+        double dy = width / (nCols - 1);
+        double dx = height / (nRows - 1);
+        for (int i = 0; i < nRows; ++i)
+        {
+            Vec3d xyz = new Vec3d();
+            for (int j = 0; j < nCols; j++)
+            {
+                xyz[0] = dx * i;
+                xyz[1] = 1.0;
+                xyz[2] = dy * j;
+                vertices[(uint)(i * nCols + j)] = xyz;
+            }
+        }
+
+        // Add connectivity data
+        VecDataArray3i indices = new VecDataArray3i();
+        for (int i = 0; i < nRows - 1; i++)
+        {
+            for (int j = 0; j < nCols - 1; j++)
+            {
+                int index1 = i * nCols + j;
+                int index2 = index1 + nCols;
+                int index3 = index1 + 1;
+                int index4 = index2 + 1;
+
+                // Interleave [/][\]
+                if (i % 2 != 0 || j % 2 != 0)
+                {
+                    indices.push_back(new Vec3i(index1, index2, index3));
+                    indices.push_back(new Vec3i(index4, index3, index2));
+                }
+                else
+                {
+                    indices.push_back(new Vec3i(index2, index4, index1));
+                    indices.push_back(new Vec3i(index4, index3, index1));
+                }
+            }
+        }
+
+        clothMesh.initialize(vertices, indices);
+
+        return clothMesh;
+    }
+
+    public static void RunPbcCutting()
+    {
+        // Setup a scene
+        Scene scene = new Scene("PBDCutting");
+        SurfaceMesh cutGeom = makeClothGeometry(40, 40, 2, 2);
+        cutGeom.setTranslation(new Vec3d(-10.0, -20.0, 0.0));
+        cutGeom.updatePostTransformData();
+        CollidingObject cutObj = new CollidingObject("CuttingObject");
+        cutObj.setVisualGeometry(cutGeom);
+        cutObj.setCollidingGeometry(cutGeom);
+        cutObj.getVisualModel(0).getRenderMaterial().setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+        scene.addSceneObject(cutObj);
+
+        PbdObject clothObj = makeClothObj("Cloth", width, height, nRows, nCols);
+        scene.addSceneObject(clothObj);
+
+        // Add interaction pair for pbd cutting
+        PbdObjectCuttingPair cuttingPair = new PbdObjectCuttingPair(clothObj, cutObj);
+
+        // Device Sever
+        HapticDeviceManager server = new HapticDeviceManager();
+        HapticDeviceClient client = server.makeDeviceClient();
+
+        SceneObjectController controller = new SceneObjectController(cutObj, client);
+        scene.addController(controller);
+
+        // Adjust camera
+        scene.getActiveCamera().setPosition(new Vec3d(0.0, -50.0, 0.0));
+        scene.getActiveCamera().setFocalPoint(new Vec3d(100.0, 100.0, 100.0));
+
+        // Light
+        DirectionalLight light = new DirectionalLight();
+        light.setFocalPoint(new Vec3d(5.0, -8.0, -5.0));
+        light.setIntensity(1.0);
+        scene.addLight("light", light);
+
+        // Run the simulation
+        {
+            // Setup a viewer to render
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.setExecutionType(Module.ExecutionType.ADAPTIVE);
+            sceneManager.pause(); // Start simulation paused
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(server);
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+            driver.setDesiredDt(0.001);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            Action<KeyEvent> receiverAction = (KeyEvent e) => {
+                const int KEY_PRESS = 1;
+                // Set new textures
+                if (e.m_key == 'i' && e.m_keyPressType == KEY_PRESS)
+                {
+                    cuttingPair.apply();
+                }
+            };
+            CSReceiverFunc receiverFunc = new CSReceiverFunc(receiverAction);
+            Utils.queueConnectKeyEvent(viewer.getKeyboardDevice(), Utils.KeyboardDeviceClient_getKeyPress_cb, sceneManager, receiverFunc);
+
+            driver.start();
+        }
+    }
+	
+	private static double width = 50.0;
+	private static double height = 50.0;
+	private static int nRows = 12;
+	private static int nCols = 12;
+}
+
diff --git a/Source/Wrappers/csharp/Examples/pbdVolume.cs b/Source/Wrappers/csharp/Examples/pbdVolume.cs
new file mode 100644
index 0000000000000000000000000000000000000000..742ae2f16fcac0129aa60381380d40a55eae1453
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/pbdVolume.cs
@@ -0,0 +1,123 @@
+using System;
+using imstk;
+
+public class PbdVolume
+{
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+        runPbdDeformable();
+    }
+
+    private static PbdObject createAndAddPbdObject(string tetMeshName)
+    {
+        TetrahedralMesh tetMesh = MeshIO.readTetrahedralMesh(tetMeshName);
+        tetMesh.rotate(new Vec3d(1.0, 0.0, 0.0), -1.3, Geometry.TransformType.ApplyToData);
+        SurfaceMesh surfMesh = tetMesh.extractSurfaceMesh();
+        // surfMesh.flipNormals();
+
+        RenderMaterial material = new RenderMaterial();
+        material.setDisplayMode(RenderMaterial.DisplayMode.Surface);
+        material.setColor(new Color(220.0 / 255.0, 100.0 / 255.0, 70.0 / 255.0));
+        material.setMetalness(100.9f);
+        material.setRoughness(0.5f);
+        material.setEdgeColor(Color.Teal);
+        material.setShadingModel(RenderMaterial.ShadingModel.Phong);
+        material.setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+        VisualModel visualModel = new VisualModel();
+        visualModel.setGeometry(surfMesh);
+        visualModel.setRenderMaterial(material);
+
+        PbdObject deformableObj = new PbdObject("DeformableObject");
+        PbdModel pbdModel = new PbdModel();
+        pbdModel.setModelGeometry(tetMesh);
+
+        // Configure model
+        PBDModelConfig pbdParams = new PBDModelConfig();
+
+        // FEM constraint
+        pbdParams.m_femParams.m_YoungModulus = 500.0;
+        pbdParams.m_femParams.m_PoissonRatio = 0.3;
+        pbdParams.m_fixedNodeIds = new VectorSizet(13);
+        pbdParams.m_fixedNodeIds.Add(75);
+        pbdParams.m_fixedNodeIds.Add(82);
+        pbdParams.m_fixedNodeIds.Add(84);
+        pbdParams.m_fixedNodeIds.Add(94);
+        pbdParams.m_fixedNodeIds.Add(95);
+        pbdParams.m_fixedNodeIds.Add(105);
+        pbdParams.m_fixedNodeIds.Add(110);
+        pbdParams.m_fixedNodeIds.Add(124);
+        pbdParams.m_fixedNodeIds.Add(139);
+        pbdParams.m_fixedNodeIds.Add(150);
+        pbdParams.m_fixedNodeIds.Add(161);
+        pbdParams.m_fixedNodeIds.Add(171);
+        pbdParams.m_fixedNodeIds.Add(350);
+        pbdParams.enableFEMConstraint(PbdConstraint.Type.FEMTet, PbdFEMConstraint.MaterialType.StVK);
+
+        // Other parameters
+        pbdParams.m_uniformMassValue = 1.0;
+        pbdParams.m_gravity    = new Vec3d(0, -9.8, 0);
+        pbdParams.m_iterations = 6;
+        pbdParams.m_dt = 0.02;
+
+        // Set the parameters
+        pbdModel.configure(pbdParams);
+        pbdModel.setTimeStepSizeType(TimeSteppingType.Fixed);
+
+        deformableObj.setDynamicalModel(pbdModel);
+        deformableObj.addVisualModel(visualModel);
+        deformableObj.setPhysicsGeometry(tetMesh);
+        deformableObj.setPhysicsToVisualMap(new TetraTriangleMap(tetMesh, surfMesh));
+
+        return deformableObj;
+    }
+
+    private static void runPbdDeformable()
+    {
+        Scene scene = new Scene("PBDVolume");
+        scene.getActiveCamera().setPosition(0, 2.0, 15.0);
+
+        string tetMeshFileName = dataPath + "textured_organs/heart_volume.vtk";
+        // create and add a PBD object
+        scene.addSceneObject(createAndAddPbdObject(tetMeshFileName));
+
+        // Light
+        DirectionalLight light = new DirectionalLight();
+        light.setFocalPoint(new Vec3d(5, -8, -5));
+        light.setIntensity(1.1);
+        scene.addLight("light", light);
+
+        // Run the simulation
+        {
+            // Setup a viewer to render
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+            viewer.setBackgroundColors(new Color(0.3285, 0.3285, 0.6525), new Color(0.13836, 0.13836, 0.2748), true);
+
+            // Setup a scene manager to advance the scene
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.pause(); // Start simulation paused
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            driver.start();
+        }
+    }
+    private static string dataPath = "../data/";
+}
diff --git a/Source/Wrappers/csharp/Examples/rigidBody2.cs b/Source/Wrappers/csharp/Examples/rigidBody2.cs
new file mode 100644
index 0000000000000000000000000000000000000000..29646b2ff5648bee2c5f07417819104944e926f4
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/rigidBody2.cs
@@ -0,0 +1,198 @@
+using System;
+using imstk;
+
+public class RigidBody2
+{
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+        runRigidBody2();
+    }
+
+    private static void runRigidBody2()
+    {
+        Scene scene = new Scene("Rigid Body Dynamics");
+        RigidObject2 cubeObj = new RigidObject2("Cube");
+        {
+            // This model is shared among interacting rigid bodies
+            RigidBodyModel2 rbdModel = new RigidBodyModel2();
+            rbdModel.getConfig().m_gravity = new Vec3d(0.0, -2500.0, 0.0);
+            rbdModel.getConfig().m_maxNumIterations = 10;
+
+            // Create the first rbd, plane floor
+            CollidingObject planeObj = new CollidingObject("Plane");
+            {
+                // Subtract the sphere from the plane to make a crater
+                Plane planeGeom = new Plane();
+                planeGeom.setWidth(40.0);
+                Sphere sphereGeom = new Sphere();
+                sphereGeom.setRadius(25.0);
+                sphereGeom.setPosition(0.0, 10.0, 0.0);
+                CompositeImplicitGeometry compGeom = new CompositeImplicitGeometry();
+                compGeom.addImplicitGeometry(planeGeom, CompositeImplicitGeometry.GeometryBoolType.Union);
+                compGeom.addImplicitGeometry(sphereGeom, CompositeImplicitGeometry.GeometryBoolType.Difference);
+
+                // Rasterize the SDF into an image
+                ImplicitGeometryToImageData toImage = new ImplicitGeometryToImageData();
+                toImage.setInputGeometry(compGeom);
+                Vec6d bounds = new Vec6d();
+                bounds[0] = -20.0;
+                bounds[1] = 20.0;
+                bounds[2] = -20.0;
+                bounds[3] = 20.0;
+                bounds[4] = -20.0;
+                bounds[5] = 20.0;
+                toImage.setBounds(bounds);
+                toImage.setDimensions(new Vec3i(80, 80, 80));
+                toImage.update();
+
+                // Extract surface
+                SurfaceMeshFlyingEdges toSurfMesh = new SurfaceMeshFlyingEdges();
+                toSurfMesh.setInputImage(toImage.getOutputImage());
+                toSurfMesh.update();
+                toSurfMesh.getOutputMesh().flipNormals();
+
+                // Create the visual model
+                VisualModel visualModel = new VisualModel(toSurfMesh.getOutputMesh());
+
+                // Create the object
+                planeObj.addVisualModel(visualModel);
+                planeObj.setCollidingGeometry(compGeom);
+                //planeObj.getRigidBody().m_isStatic = true;
+                //planeObj.getRigidBody().m_mass     = 100.0;
+
+                scene.addSceneObject(planeObj);
+            }
+
+            // Create surface mesh cube (so we can use pointset for point.implicit collision)
+            {
+                OrientedBox cubeGeom = new OrientedBox(new Vec3d(0.0, 0.0, 0.0), new Vec3d(1.5, 3.0, 1.0));
+                SurfaceMesh surfMesh = Utils.toSurfaceMesh(cubeGeom);
+
+                SurfaceMeshSubdivide subdivide = new SurfaceMeshSubdivide();
+                subdivide.setInputMesh(surfMesh);
+                subdivide.setNumberOfSubdivisions(1);
+                subdivide.update();
+
+                // Create the visual model
+                VisualModel visualModel = new VisualModel(subdivide.getOutputMesh());
+                RenderMaterial mat = new RenderMaterial();
+                mat.setDisplayMode(RenderMaterial.DisplayMode.WireframeSurface);
+                mat.setLineWidth(2.0f);
+                mat.setColor(Color.Orange);
+                visualModel.setRenderMaterial(mat);
+
+                // Create the cube rigid object
+                cubeObj.setDynamicalModel(rbdModel);
+                cubeObj.setPhysicsGeometry(subdivide.getOutputMesh());
+                cubeObj.setCollidingGeometry(subdivide.getOutputMesh());
+                cubeObj.addVisualModel(visualModel);
+                cubeObj.getRigidBody().m_mass    = 100.0;
+                cubeObj.getRigidBody().m_initPos = new Vec3d(0.0, 8.0, 0.0);
+                Rotd rotd = new Rotd(0.4, new Vec3d(1.0, 0.0 ,0.0));
+                cubeObj.getRigidBody().m_initOrientation = new Quatd(new Rotd(0.4, new Vec3d(1.0, 0.0, 0.0)));
+                cubeObj.getRigidBody().m_intertiaTensor  = Mat3d.Identity();
+
+                scene.addSceneObject(cubeObj);
+            }
+
+            RigidObjectCollision rbdInteraction = new RigidObjectCollision(cubeObj, planeObj, "ImplicitGeometryToPointSetCD");
+            rbdInteraction.setFriction(0.0);
+            rbdInteraction.setStiffness(0.05);
+            scene.getCollisionGraph().addInteraction(rbdInteraction);
+            scene.getActiveCamera().setPosition(0.0, 40.0, 40.0);
+
+            // Light
+            DirectionalLight light = new DirectionalLight();
+            light.setIntensity(1.0);
+            scene.addLight("light", light);
+        }
+
+        // Run the simulation
+        {
+            // Setup a viewer to render in its own thread
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene in its own thread
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.setExecutionType(Module.ExecutionType.ADAPTIVE);
+            sceneManager.pause();
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+            driver.setDesiredDt(0.001);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            // LOG(INFO) << "Cube Controls:";
+            // LOG(INFO) << "----------------------------------------------------------------------";
+            // LOG(INFO) << " | i - forward movement";
+            // LOG(INFO) << " | j - left movement";
+            // LOG(INFO) << " | l - right movement";
+            // LOG(INFO) << " | k - backwards movement";
+            // LOG(INFO) << " | u - rotate left";
+            // LOG(INFO) << " | o - rotate right";
+
+            // Not perfectly thread safe movement lambda, ijkl movement instead of wasd because d is already used
+            // KeyboardDeviceClient keyDevice = viewer.getKeyboardDevice();
+            // Vec3d dx = new Vec3d();
+            // dx = scene.getActiveCamera().getPosition() - scene.getActiveCamera().getFocalPoint();
+            // connect<Event>(sceneManager, &SceneManager.postUpdate, [&](Event*)
+            //                {
+            //                new Vec3d extForce  = new Vec3d(0.0, 0.0, 0.0);
+            //                new Vec3d extTorque = new Vec3d(0.0, 0.0, 0.0);
+            //                // If w down, move forward
+            //                if (keyDevice.getButton('i') == KEY_PRESS)
+            //                {
+            //                extForce += new Vec3d(0.0, 0.0, -900.0);
+            //                }
+            //                if (keyDevice.getButton('k') == KEY_PRESS)
+            //                {
+            //                extForce += new Vec3d(0.0, 0.0, 900.0);
+            //                }
+            //                if (keyDevice.getButton('j') == KEY_PRESS)
+            //                {
+            //                extForce += new Vec3d(-900.0, 0.0, 0.0);
+            //                }
+            //                if (keyDevice.getButton('l') == KEY_PRESS)
+            //                {
+            //                extForce += new Vec3d(900.0, 0.0, 0.0);
+            //                }
+            //                if (keyDevice.getButton('u') == KEY_PRESS)
+            //                {
+            //                    extTorque += new Vec3d(0.0, 1.5, 0.0);
+            //                }
+            //                if (keyDevice.getButton('o') == KEY_PRESS)
+            //                {
+            //                    extTorque += new Vec3d(0.0, -1.5, 0.0);
+            //                }
+            //                *cubeObj.getRigidBody().m_force  = extForce;
+            //                *cubeObj.getRigidBody().m_torque = extTorque;
+            //                scene.getActiveCamera().setFocalPoint(cubeObj.getRigidBody().getPosition());
+            //                scene.getActiveCamera().setPosition(cubeObj.getRigidBody().getPosition() + dx);
+            //                });
+            // connect<Event>(sceneManager, &SceneManager.postUpdate, [&](Event*)
+            //                {
+            //                cubeObj.getRigidBodyModel2().getConfig().m_dt = sceneManager.getDt();
+            //                });
+
+            driver.start();
+        }
+    }
+}
+
+
diff --git a/Source/Wrappers/csharp/Examples/sdfHaptics.cs b/Source/Wrappers/csharp/Examples/sdfHaptics.cs
new file mode 100644
index 0000000000000000000000000000000000000000..04691a5a4047d524c4f207a2e82cc855da2feb81
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/sdfHaptics.cs
@@ -0,0 +1,116 @@
+using System;
+using imstk;
+
+public class SdfHaptics
+{
+     public class CSReceiverFunc : EventFunc {
+         public CSReceiverFunc(Action<Event> action) {
+             action_ = action;
+         }
+         public override void call(Event e) {
+             action_(e);
+         }
+         private Action<Event> action_;
+     }
+
+    public static void Main(string[] args)
+    {
+        Logger.startLogger();
+
+        Scene scene = new Scene("SDFHaptics");
+        SurfaceMesh  axesMesh = MeshIO.readSurfaceMesh("../data/axesPoly.vtk");
+        ImageData    sdfImage = MeshIO.readImageData("../data/stanfordBunny/stanfordBunny_SDF.nii");
+        SignedDistanceField sdf = new SignedDistanceField(sdfImage.cast((byte)Utils.IMSTK_DOUBLE));
+        {
+            scene.getActiveCamera().setPosition(-2.3, 23.81, 45.65);
+            scene.getActiveCamera().setFocalPoint(9.41, 8.45, 5.76);
+
+            CollidingObject bunnyObj = new CollidingObject("Bunny");
+            {
+                bunnyObj.setCollidingGeometry(sdf);
+
+                SurfaceMeshFlyingEdges isoExtract = new SurfaceMeshFlyingEdges();
+                isoExtract.setInputImage(sdfImage);
+                isoExtract.update();
+
+                isoExtract.getOutputMesh().flipNormals();
+                bunnyObj.setVisualGeometry(isoExtract.getOutputMesh());
+
+                scene.addSceneObject(bunnyObj);
+            }
+
+            SceneObject axesObj = new SceneObject("Axes");
+            {
+                axesObj.setVisualGeometry(axesMesh);
+                scene.addSceneObject(axesObj);
+            }
+
+            // Light (white)
+            DirectionalLight whiteLight = new DirectionalLight();
+            {
+                whiteLight.setDirection(new Vec3d(5.0, -8.0, -5.0));
+                whiteLight.setIntensity(1.0);
+                scene.addLight("whiteLight", whiteLight);
+            }
+        }
+
+        HapticDeviceManager hapticManager = new HapticDeviceManager();
+        HapticDeviceClient client = hapticManager.makeDeviceClient();
+
+        // Run the simulation
+        {
+            // Setup a viewer to render in its own thread
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+
+            // Setup a scene manager to advance the scene in its own thread
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.setExecutionType(Module.ExecutionType.ADAPTIVE);
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+            driver.addModule(hapticManager);
+
+            ImplicitFunctionCentralGradient centralGrad = new ImplicitFunctionCentralGradient();
+            centralGrad.setFunction(sdf);
+            centralGrad.setDx(sdf.getImage().getSpacing());
+
+            Action<Event> receiverAction = (Event e) => {
+				Vec3d tmpPos = client.getPosition();
+                Vec3d pos = tmpPos * 0.1 + new Vec3d(0.0, 0.1, 10.0);
+
+                client.update();
+                axesMesh.setTranslation(pos);
+                axesMesh.setRotation(client.getOrientation());
+                axesMesh.postModified();
+
+                double dx = sdf.getFunctionValue(pos);
+                if (dx < 0.0)
+                {
+                    Vec3d g = centralGrad.compute(pos);
+					Vec3d nrm_g = new Vec3d(-g[0]*dx*4.0, -g[1]*dx*4.0, -g[2]*dx*4.0); 
+                    client.setForce(nrm_g);
+                }
+            };
+
+            CSReceiverFunc receiverFunc = new CSReceiverFunc(receiverAction);
+            Utils.connectEvent(sceneManager, Utils.SceneManager_getPostUpdate_cb, receiverFunc);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            driver.start();
+        }
+    }
+}
diff --git a/Source/Wrappers/csharp/Examples/sphFluid.cs b/Source/Wrappers/csharp/Examples/sphFluid.cs
new file mode 100644
index 0000000000000000000000000000000000000000..1ccefd29ca93fa10cb3b30e00af5e0696a74e826
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/sphFluid.cs
@@ -0,0 +1,312 @@
+using System;
+using imstk;
+
+public class PbdCloth
+{
+    private static int SCENE_ID = 1;
+
+    public static void Main(string[] args)
+    {
+
+        // Setup logger (write to file and stdout)
+        Logger.startLogger();
+
+        double particleRadius = 0.1;
+
+        // Override particle radius for scene3 because particles in this scene were pre-generated using particle radius 0.08
+        if (SCENE_ID == 3)
+        {
+            particleRadius = 0.08;
+        }
+
+        Scene scene = new Scene("SPH Fluid");
+
+        SPHObject fluidObj = generateFluid(particleRadius);
+        CollidingObject[] solids = generateSolids(scene);
+        scene.addSceneObject(fluidObj);
+        for (int i = 0; i < solids.Length; i++)
+        {
+            scene.addSceneObject(solids[i]);
+        }
+
+        // Collision between fluid and solid objects
+        CollisionGraph collisionGraph = scene.getCollisionGraph();
+
+        collisionGraph.addInteraction(new SphObjectCollision(fluidObj, solids[0]));
+        collisionGraph.addInteraction(new SphObjectCollision(fluidObj, solids[1]));
+        collisionGraph.addInteraction(new SphObjectCollision(fluidObj, solids[2]));
+
+
+        // configure camera
+        scene.getActiveCamera().setPosition(-0.475, 8.116, -6.728);
+
+        // configure light (white)
+        DirectionalLight whiteLight = new DirectionalLight();
+        whiteLight.setFocalPoint(new Vec3d(5.0, -8.0, -5.0));
+        whiteLight.setIntensity(1.5);
+        scene.addLight("whiteLight", whiteLight);
+
+        // Run the simulation
+        {
+            // Setup a viewer to render
+            VTKViewer viewer = new VTKViewer("Viewer");
+            viewer.setActiveScene(scene);
+            viewer.setWindowTitle("SPH Fluid");
+            viewer.setSize(1920, 1080);
+            VTKTextStatusManager statusManager = viewer.getTextStatusManager();
+            statusManager.setStatusFontSize(VTKTextStatusManager.StatusType.Custom, 30);
+            statusManager.setStatusFontColor(VTKTextStatusManager.StatusType.Custom, Color.Red);
+
+            // Action<KeyEvent> receiverAction = (KeyEvent e) => {
+            //     statusManager.setCustomStatus("Number of particles: " +
+            //                                    std.to_string(fluidObj.getSPHModel().getCurrentState().getNumParticles()) +
+            //                                    "\nNumber of solids: " + std.to_string(solids.size()));
+            //
+            // };
+
+            // CSReceiverFunc eventFunc = new CSReceiverFunc(receiverAction);
+            // Utils.connectEvent(viewer, Utils.VTKViewer_getPostUpdate_cb, eventFunc);
+
+
+            // Setup a scene manager to advance the scene
+            SceneManager sceneManager = new SceneManager("Scene Manager");
+            sceneManager.setActiveScene(scene);
+            sceneManager.pause();
+
+            SimulationManager driver = new SimulationManager();
+            driver.addModule(viewer);
+            driver.addModule(sceneManager);
+
+            // Add mouse and keyboard controls to the viewer
+            {
+                MouseSceneControl mouseControl = new MouseSceneControl(viewer.getMouseDevice());
+                mouseControl.setSceneManager(sceneManager);
+                viewer.addControl(mouseControl);
+
+                KeyboardSceneControl keyControl = new KeyboardSceneControl(viewer.getKeyboardDevice());
+                keyControl.setSceneManager(new SceneManagerWeakPtr(sceneManager));
+                keyControl.setModuleDriver(new ModuleDriverWeakPtr(driver));
+                viewer.addControl(keyControl);
+            }
+
+            driver.start();
+        }
+    }
+
+    public static SPHObject generateFluid(double particleRadius)
+    {
+        VecDataArray3d particles = new VecDataArray3d();
+        switch (SCENE_ID)
+        {
+        case 1:
+            particles = generateSphereShapeFluid(particleRadius);
+            break;
+        case 2:
+            particles = generateBoxShapeFluid(particleRadius);
+            break;
+        case 3:
+            // particles = generateBunnyShapeFluid(particleRadius);
+            return null;
+            break;
+        default:
+            return null;
+        }
+
+
+        // Create a geometry object
+        PointSet geometry = new PointSet();
+        geometry.initialize(particles);
+
+        // Create a fluids object
+        SPHObject fluidObj = new SPHObject("Sphere");
+
+        // Create a visual model
+        VisualModel visualModel = new VisualModel(geometry);
+        RenderMaterial material = new RenderMaterial();
+        material.setDisplayMode(RenderMaterial.DisplayMode.Fluid);
+        //material.setDisplayMode(RenderMaterial.DisplayMode.Points);
+        if (material.getDisplayMode() == RenderMaterial.DisplayMode.Fluid)
+        {
+            material.setPointSize(0.1f);
+        }
+        else
+        {
+            material.setPointSize(20.0f);
+            material.setRenderPointsAsSpheres(true);
+            material.setColor(Color.Orange);
+        }
+        visualModel.setRenderMaterial(material);
+
+        // Create a physics model
+        SPHModel sphModel = new SPHModel();
+        sphModel.setModelGeometry(geometry);
+
+        // Configure model
+        SPHModelConfig sphParams = new SPHModelConfig(particleRadius);
+        sphParams.m_bNormalizeDensity = true;
+        if (SCENE_ID == 2)   // highly viscous fluid
+        {
+            sphParams.m_kernelOverParticleRadiusRatio = 6.0;
+            sphParams.m_surfaceTensionStiffness = 5.0;
+        }
+
+        if (SCENE_ID == 3)   // bunny-shaped fluid
+        {
+            sphParams.m_frictionBoundary = 0.3;
+        }
+
+        sphModel.configure(sphParams);
+        sphModel.setTimeStepSizeType(TimeSteppingType.RealTime);
+
+        // Add the component models
+        fluidObj.addVisualModel(visualModel);
+        fluidObj.setCollidingGeometry(geometry);
+        fluidObj.setDynamicalModel(sphModel);
+        fluidObj.setPhysicsGeometry(geometry);
+
+        return fluidObj;
+    }
+
+    ///
+    /// \brief Generate a sphere-shape fluid object
+    ///
+    public static VecDataArray3d generateSphereShapeFluid(double particleRadius)
+    {
+        double sphereRadius = 2.0;
+        Vec3d  sphereCenter = new Vec3d(0, 1, 0);
+        double  sphereRadiusSqr = sphereRadius * sphereRadius;
+        double  spacing = 2.0 * particleRadius;
+        int     N = (int)(2.0 * sphereRadius / spacing);              // Maximum number of particles in each dimension
+        // Vec3d lcorner = sphereCenter - new Vec3d(sphereRadius, sphereRadius, sphereRadius); // Cannot use auto here, due to Eigen bug
+        Vec3d lcorner = new Vec3d(sphereCenter[0] - sphereRadius, sphereCenter[1] - sphereRadius, sphereCenter[2] - sphereRadius); // Cannot use auto here, due to Eigen bug
+
+        VecDataArray3d particles = new VecDataArray3d();
+        particles.reserve(N * N * N);
+
+        for (int i = 0; i < N; ++i)
+        {
+            for (int j = 0; j < N; ++j)
+            {
+                for (int k = 0; k < N; ++k)
+                {
+                    Vec3d ppos = lcorner + new Vec3d(spacing * (double)(i), spacing * (double)(j), spacing * (double)(k));
+                    Vec3d cx = ppos - sphereCenter;
+                    double nrm = cx[0] * cx[0] + cx[1] * cx[1] + cx[2] * cx[2];
+                    if (nrm < sphereRadiusSqr)
+                    {
+                        particles.push_back(ppos);
+                    }
+                }
+            }
+        }
+
+        return particles;
+    }
+
+    ///
+    /// \brief Generate a box-shape fluid object
+    ///
+    public static VecDataArray3d generateBoxShapeFluid(double particleRadius)
+    {
+        double boxWidth = 4.0;
+        Vec3d  boxLowerCorner = new Vec3d(-2, -3, -2);
+
+        double spacing = 2.0 * particleRadius;
+        int N = (int)(boxWidth / spacing);
+
+        VecDataArray3d particles = new VecDataArray3d();
+        particles.reserve(N * N * N);
+
+        for (int i = 0; i < N; ++i)
+        {
+            for (int j = 0; j < N; ++j)
+            {
+                for (int k = 0; k < N; ++k)
+                {
+                    Vec3d ppos = boxLowerCorner + new Vec3d(spacing * i, spacing * j, spacing * k);
+                    particles.push_back(ppos);
+                }
+            }
+        }
+
+        return particles;
+    }
+
+    public static CollidingObject[] generateSolids(Scene scene)
+    {
+        switch (SCENE_ID)
+        {
+        case 1:
+            return generateSolidsScene1();
+        case 2:
+            return null; // To avoid warning
+            // return generateSolidsScene2();
+        case 3:
+            return null; // To avoid warning
+            // return generateSolidsScene3();
+        case 4:
+            return null; // To avoid warning
+            // return generateSolidsScene4(scene);
+        default:
+            return null; // To avoid warning
+        }
+    }
+
+    ///
+    /// \brief Generate two planes and a solid sphere
+    ///
+    public static CollidingObject[] generateSolidsScene1()
+    {
+        CollidingObject[] solids = new CollidingObject[3];
+
+        {
+            Plane geometry = new Plane();
+            geometry.setWidth(40.0);
+            geometry.setPosition(0.0, -6.0, 0.0);
+            geometry.setNormal(new Vec3d(0.0, 1.0, -0.5));
+
+            VisualModel visualModel = new VisualModel(geometry);
+            RenderMaterial material = new RenderMaterial();
+            material.setColor(Color.DarkGray);
+            visualModel.setRenderMaterial(material);
+
+            CollidingObject obj = new CollidingObject("Floor");
+            obj.addVisualModel(visualModel);
+            obj.setCollidingGeometry(geometry);
+            solids[0] = obj;
+        }
+        {
+            Plane geometry = new Plane();
+            geometry.setWidth(40.0);
+            geometry.setPosition(0.0, -6.0, 0.0);
+            geometry.setNormal(new Vec3d(0.0, 1.0, 1.0));
+
+            VisualModel visualModel = new VisualModel(geometry);
+            RenderMaterial material = new RenderMaterial();
+            material.setColor(Color.LightGray);
+            visualModel.setRenderMaterial(material);
+
+            CollidingObject obj = new CollidingObject("Back Plane");
+            obj.addVisualModel(visualModel);
+            obj.setCollidingGeometry(geometry);
+            solids[1] = obj;
+        }
+        {
+            Sphere geometry = new Sphere();
+            geometry.setRadius(2.0);
+            geometry.setPosition(0.0, -6.0, 0.0);
+
+            VisualModel visualModel = new VisualModel(geometry);
+            RenderMaterial material = new RenderMaterial();
+            material.setColor(Color.Red);
+            visualModel.setRenderMaterial(material);
+
+            CollidingObject obj = new CollidingObject("Sphere on Floor");
+            obj.addVisualModel(visualModel);
+            obj.setCollidingGeometry(geometry);
+            solids[2] = obj;
+        }
+
+        return solids;
+    }
+}
diff --git a/Source/Wrappers/csharp/Examples/testGeometry.cs b/Source/Wrappers/csharp/Examples/testGeometry.cs
new file mode 100644
index 0000000000000000000000000000000000000000..15adecdcb658f16cbbebf1cf1536c41c0d1e4e8a
--- /dev/null
+++ b/Source/Wrappers/csharp/Examples/testGeometry.cs
@@ -0,0 +1,78 @@
+using System;
+using imstk;
+
+
+public class PbdCloth
+{
+    public static void Main(string[] args)
+    {
+        // Write log to stdout and file
+        Logger.startLogger();
+        GeometryTest();
+    }
+
+    public static void GeometryTest()
+    {
+        VecDataArray4i conn = new VecDataArray4i(1);
+        Vec4i vec4i = new Vec4i();
+        vec4i[0] = 0;
+        vec4i[1] = 1;
+        vec4i[2] = 2;
+        vec4i[3] = 3;
+        conn[0] = vec4i;
+
+        VecDataArray3d coords = new VecDataArray3d(4);
+        Vec3d xyz = new Vec3d();
+        xyz[0] = 0.0;
+        xyz[1] = 0.0;
+        xyz[2] = 0.0;
+        coords[0] = xyz;
+        uint id = 0;
+        Console.WriteLine ("coords[{0}] = [{1}, {2}, {3}]", id, coords[id][0], coords[id][1], coords[id][2]);
+
+        xyz[0] = 1.0;
+        xyz[1] = 0.0;
+        xyz[2] = 0.0;
+        coords[1] = xyz;
+        id = 1;
+        Console.WriteLine ("coords[{0}] = [{1}, {2}, {3}]", id, coords[id][0], coords[id][1], coords[id][2]);
+
+        xyz[0] = 0.0;
+        xyz[1] = 1.0;
+        xyz[2] = 0.0;
+        coords[2] = xyz;
+        id = 2;
+        id = 2;
+        Console.WriteLine ("coords[{0}] = [{1}, {2}, {3}]", id, coords[id][0], coords[id][1], coords[id][2]);
+
+        xyz[0] = 0.0;
+        xyz[1] = 0.0;
+        xyz[2] = 1.0;
+        coords[3] = xyz;
+        id = 3;
+        Console.WriteLine ("coords[{0}] = [{1}, {2}, {3}]", id, coords[id][0], coords[id][1], coords[id][2]);
+
+        TetrahedralMesh tetMesh = new TetrahedralMesh();
+        tetMesh.initialize(coords, conn);
+        
+        Console.WriteLine ("vol = {0}", tetMesh.getVolume());
+        Console.WriteLine ("numberOfTets = {0}", tetMesh.getNumTetrahedra());
+        Vec4i conn_tet = tetMesh.getTetrahedronIndices(0);
+        Console.WriteLine ("conn[0] = [{0}, {1}, {2}, {3}]", conn_tet[0], conn_tet[1], conn_tet[2], conn_tet[3]);
+
+        {
+            VecDataArray4i conn2 = tetMesh.getTetrahedraIndices();
+            conn2[0] = new Vec4i(1, 2, 3, 0);
+            Vec4i new_conn_tet = tetMesh.getTetrahedronIndices(0);
+            Console.WriteLine ("conn[0] = [{0}, {1}, {2}, {3}]", new_conn_tet[0], new_conn_tet[1], new_conn_tet[2], new_conn_tet[3]);
+        }
+
+        tetMesh.print();
+
+        PointSet asianDragonMesh = MeshIO.read("/home/jianfeng/Documents/imstk/build_csharp_shared/install/data/asianDragon/asianDragon.veg");
+        // asianDragonMesh.print();
+        Console.WriteLine ("asianDragon.getVolume() = {0}", asianDragonMesh.getVolume());
+        Console.WriteLine ("asianDragon.getNumVertices() = {0}", asianDragonMesh.getNumVertices());
+    }
+}
+
diff --git a/Source/Wrappers/csharp/README.md b/Source/Wrappers/csharp/README.md
new file mode 100644
index 0000000000000000000000000000000000000000..8ef0d081eba80e615c719f67e2646d1c510724e3
--- /dev/null
+++ b/Source/Wrappers/csharp/README.md
@@ -0,0 +1,54 @@
+# General 
+
+Please note that the C# support is still a work in progress, while we will do our best to maintain the C/C++ layer of this part of iMSTK it might happen that some of the code in the wrapper comes out of sync with the rest of the system, we appreciate any reports or fixes that you want to contribute to help us out. The current version of the wrapper will be fully tested against the `5.0` release tag of imstk.
+
+## Components
+
+iMSTK uses the [simplified wrapper and interface generator (SWIG)](http://www.swig.org/) to provide an interface between C++ and C#. This enables us to use most of the C++ interface almost "as is" in C#. There are a few components to this wrapping proccess. The sources to generate the wrapper can be found under `Source/Wrappers/csharp` in the iMSTK source directory. Even on windows not all of the processes are automated yet.
+
+The directory contains the following items 
+    - SwigInterface: The swig `.i` files that tell swig how to build the interfaces for iMSTK, this is the C/C++ part of the wrapper
+    - iMSTKCSharp: A CMake project (Visual Studio only) that can build a combined dll with all the swig generated files, the C# part of the wrapper
+    - Examples: A set of iMSTK examples written in C#
+    - Testing: Unit tests for some critical parts of the wrapper (can't be built automatically atm)
+
+## Building
+
+Support for building differs between Linux and Visual Studio, as Visual Studio supports C# projects. Both require SWIG to be installed. On windows two different paths are available, one using visual studio and `.NET` libraries, or a more manual approach by building via `Mono`. In both cases the C/C++ layer has to be built first.
+
+# How to build iMSTK C/C++ wrappers?
+- Windows:
+    - download and install [swig](http://www.swig.org/download.html)
+    - Edit the system environment variables:
+        - create 2 "System Variables": `SWIG_DIR` points to `/path/to/swig/dir/`, and `SWIG_EXECUTABLE` points to `/path/to/swig/dir/swig.exe`
+        - add one entry to PATH variable, `/path/to/swig/dir`
+        - may have to restart the computer to apply all environment variables.
+    - Build imstk with `iMSTK_WRAP_CSHARP=ON`
+- Ubuntu
+    - sudo apt install swig
+    - Build imstk with `iMSTK_WRAP_CSHARP=ON`.
+
+The output of this step is a shared library `iMSTKCWrapper.[dll|sh]` (located in `/bin`), and a series of `.cs` files that together with the `iMSTKCWrapper` library now provide access to iMSTK C++ functionalies. Note that iMSTK depends on a variety of other libraries, these are not statically linked in the `iMSTKCWrapper` library, whenever writing C# code, the `iMSTKCWrapper` library and all pertinent iMSTK dependencies need to be accessible (e.g. by being on the library path).
+
+Some people like to compile all the `.cs` files into a separate library, the Visual Studio build does that in a separate target.
+
+# How ro run C# examples on Windows
+
+## Via Visual Studio .NET
+- The superbuild needs to be run twice for the C# Examples to be run from the install folder
+- If you are looking at the inner build, after building `iMSTKCWrapper`, you will need to run CMake "Configure" and "Generate" again to load the example projects.
+- Due to the dependency relation ship noted above you can't immediately run the example out of visual studio at the moment but you need to make sure that the "Working Directory" under each examples' project settings is set to the iMSTK install folder
+
+## Via Mono 
+- Install mono. See [https://www.mono-project.com/docs/getting-started/install/windows/](https://www.mono-project.com/docs/getting-started/install/windows/)
+- The C# examples run without building iMSTK csharp wrappers, but needs iMSTK built with OpenHaptics.
+    - build iMSTK with `iMSTK_WRAP_CSHARP=ON`
+    - There are several examples in `Source/Wrappers/csharp/Examples`. Taking pbdCloth.cs for an example, open a terminal, go to the imstk installation directory and execute 
+        - /path/to/mono/bin/csc '<imstkSource>\Source\Wrappers\csharp\Examples\pbdCloth.cs include\iMSTKCSharp\*.cs -out:bin\pbdCloth.exe
+        - bin\pbdCloth.exe can run directly.
+		
+# How to compile SWIG generated C# code into a dll, and run examples with it?
+- open a terminal, cd to the iMSTK installation directory
+- /path/to/mono/bin/csc -target:library include\iMSTKCSharp\*.cs\*.cs -out:bin\iMSTKCS.dll (This compiles the C# code into a dynamic lib, iMSTKCS.dll)
+- /path/to/mono/bin/csc -r:bin\iMSTKCS.dll <imstkSource>\Source\Wrappers\csharp\Examples\pbdCloth.cs -out:bin\pbdCloth.exe
+- bin\pbdCloth.exe can run directly.
diff --git a/Source/Wrappers/csharp/SwigInterface/callback.i b/Source/Wrappers/csharp/SwigInterface/callback.i
new file mode 100644
index 0000000000000000000000000000000000000000..1c24c7d1a5a8da0c1b08980838e462c23e95c8a6
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/callback.i
@@ -0,0 +1,28 @@
+%callback("%s_cb");
+/* std::string imstk::KeyboardDeviceClient::keyPress(); */
+std::string KeyboardDeviceClient_getKeyPress();
+std::string Module_getPostUpdate();
+std::string Module_getPreUpdate();
+std::string SceneManager_getPreUpdate();
+std::string SceneManager_getPostUpdate();
+%nocallback;
+
+%{
+std::string KeyboardDeviceClient_getKeyPress() {
+  return imstk::KeyboardDeviceClient::keyPress();
+}
+
+std::string Module_getPostUpdate() {
+  return imstk::Module::postUpdate();
+}
+std::string Module_getPreUpdate() {
+  return imstk::Module::preUpdate();
+}
+std::string SceneManager_getPostUpdate() {
+  return imstk::SceneManager::postUpdate();
+}
+std::string SceneManager_getPreUpdate() {
+  return imstk::SceneManager::preUpdate();
+}
+%}
+
diff --git a/Source/Wrappers/csharp/SwigInterface/common.i b/Source/Wrappers/csharp/SwigInterface/common.i
new file mode 100644
index 0000000000000000000000000000000000000000..40826d36a8a2bc805b608e43b51913b71c8ccec2
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/common.i
@@ -0,0 +1,587 @@
+%include <cpointer.i>
+%include <arrays_csharp.i>
+
+#ifdef SWIG_PINNED_ARRAY
+  %csmethodmodifiers imstk::VecDataArray::setValues "public unsafe";
+  %csmethodmodifiers imstk::VecDataArray::getValues "public unsafe";
+#endif
+
+#ifdef SWIG_PINNED_ARRAY
+  %apply unsigned char FIXED[] {unsigned char* val}
+  %apply unsigned char FIXED[] {const unsigned char* val}
+  %apply int FIXED[] {int* val}
+  %apply int FIXED[] {const int* val}
+  %apply float FIXED[] {float * val}
+  %apply float FIXED[] {const float * val}
+  %apply double FIXED[] {double* val}
+  %apply double FIXED[] {const double* val}
+#else
+  %apply unsigned char INPUT[] {const unsigned char* val}
+  %apply unsigned char OUTPUT[] {unsigned char* val}
+  %apply int INPUT[] {const int* val}
+  %apply int OUTPUT[] {int* val}
+  %apply float INPUT[] {const float * val}
+  %apply float OUTPUT[] {float * val}
+  %apply double INPUT[] {const double* val}
+  %apply double OUTPUT[] {double* val}
+#endif
+
+/*
+ * TODO: use %define to simplify the process of wrapping Eigen types.
+ */
+
+%pointer_functions(int, intPtr)
+%pointer_functions(float, floatPtr)
+%pointer_functions(double, doublePtr)
+%pointer_functions(unsigned char, ucharPtr)
+
+%rename(getValue) imstk::VecDataArray::operator[] (const size_t pos) const;
+%rename(setValue) imstk::VecDataArray::operator[] (const size_t pos);
+%rename(getValue) imstk::Vec::operator[] (const int pos) const;
+%rename(setValue) imstk::Vec::operator[] (const int pos);
+
+%rename(Rotf) imstk::imstkRotf;
+%rename(Rotd) imstk::imstkRotd;
+%rename(Quatf) imstk::imstkQuatf;
+%rename(Quatd) imstk::imstkQuatd;
+
+%define %extend_VecDataArray(T, N)
+%extend imstk::VecDataArray<T, N> {
+    void setValues(const T* val)
+    {
+        /* std::copy(val, val+$self->m_vecSize * N, $self->Base::m_data); */
+        T* data = $self->DataArray<T>::getPointer();
+        std::copy(val, val + $self->size() * N, data);
+    }
+
+    void setValues(const T* val, const int n)
+    {
+        T* data = $self->DataArray<T>::getPointer();
+        CHECK($self->size()*N >= n) << "number of values are larger than the array size";
+        std::copy(val, val + n, data);
+    }
+
+    void getValues(T* val)
+    {
+        T* data = $self->DataArray<T>::getPointer();
+        std::copy(data, data+$self->size()*N, val);
+    }
+};
+%enddef
+
+%extend_VecDataArray(float, 2)
+%extend_VecDataArray(double, 2)
+%extend_VecDataArray(int, 3)
+%extend_VecDataArray(double, 3)
+%extend_VecDataArray(unsigned char, 3)
+%extend_VecDataArray(int, 4)
+
+%typemap(cscode) imstk::imstkQuatf %{
+    public static implicit operator SWIGTYPE_p_Eigen__QuaternionT_float_t (Quatf cs_data) {
+        return cs_data.get();
+    }
+
+    public static implicit operator Quatf (SWIGTYPE_p_Eigen__QuaternionT_float_t eigen_data) {
+        return Utils.eigen2cs_Quatf(eigen_data);
+    }
+%}
+
+%typemap(cscode) imstk::imstkQuatd %{
+    public static implicit operator SWIGTYPE_p_Eigen__QuaternionT_double_t (Quatd cs_data) {
+        return cs_data.get();
+    }
+
+    public static implicit operator Quatd (SWIGTYPE_p_Eigen__QuaternionT_double_t eigen_data) {
+        return Utils.eigen2cs_Quatd(eigen_data);
+    }
+%}
+
+#ifdef SWIGCSHARP
+%inline %{
+namespace imstk {
+template<typename T, int N>
+class Vec {
+  public:
+    using EigenData = Eigen::Matrix<T, N, 1>;
+    Vec(const Vec& other) = default;
+    Vec() = default;
+    Vec(const T t0, const T t1) {
+        m_data[0] = t0;
+        m_data[1] = t1;
+    }
+    Vec(const T t0, const T t1, const T t2) {
+        m_data[0] = t0;
+        m_data[1] = t1;
+        m_data[2] = t2;
+    }
+    Vec(const T t0, const T t1, const T t2, const T t3) {
+        m_data[0] = t0;
+        m_data[1] = t1;
+        m_data[2] = t2;
+        m_data[3] = t3;
+    }
+    Vec(const T t0, const T t1, const T t2, const T t3, const T t4, const T t5) {
+        m_data[0] = t0;
+        m_data[1] = t1;
+        m_data[2] = t2;
+        m_data[3] = t3;
+        m_data[4] = t4;
+        m_data[5] = t5;
+    }
+    // Vec(int i0, int i1, int i2, int i3)
+    Vec(const EigenData& other) { m_data = other; }
+    inline T& operator[](const int pos) { return m_data[pos]; }
+    inline const T& operator[](const int pos) const { return m_data[pos]; }
+    inline const EigenData& get() const {return m_data;}
+  private:
+    EigenData m_data;
+};
+
+template <typename T, int N>
+class Mat {
+  public:
+    using EigenData = Eigen::Matrix<T, N, N>;
+    Mat (const EigenData& data) : m_data(data) {}
+    static Mat Identity () { return Mat (EigenData::Identity()); }
+    inline const EigenData& get() const {return m_data;}
+  private:
+    EigenData m_data;
+};
+
+class imstkRotf {
+public:
+    using EigenData = Eigen::AngleAxisf;
+    imstkRotf(float s, const imstk::Vec3f& vec): m_data(s, vec) { }
+    const EigenData& get() const { return m_data; }
+private:
+    EigenData m_data;
+};
+
+class imstkRotd {
+public:
+    using EigenData = Eigen::AngleAxisd;
+    imstkRotd(double s, const imstk::Vec3d& vec): m_data(s, vec) { }
+    inline const EigenData& get() const { return m_data; }
+private:
+    EigenData m_data;
+};
+
+class imstkQuatf {
+public:
+    using EigenData = Eigen::Quaternion<float>;
+    imstkQuatf(const EigenData& data) : m_data(data) { }
+    imstkQuatf(const imstkRotf& rot) : m_data(rot.get()) { }
+    inline const EigenData& get() const { return m_data; }
+private:
+    EigenData m_data;
+};
+
+class imstkQuatd {
+public:
+    using EigenData = Eigen::Quaternion<double>;
+    imstkQuatd(const EigenData& data) : m_data(data) { }
+    imstkQuatd(const imstkRotd& rot) : m_data(rot.get()) { }
+    inline const EigenData& get() const { return m_data; }
+private:
+    EigenData m_data;
+};
+
+template<typename T, int N>
+Vec<T, N> vec_add(const Vec<T, N>& u, const Vec<T, N>& v) 
+{
+    Vec<T, N> ret(u);
+    for (int i=0; i<N; ++i)
+    {
+        ret[i] += v[i];
+    }
+    return ret;
+}
+
+template<typename T, int N>
+Vec<T, N> vec_subtract(const Vec<T, N>& u, const Vec<T, N>& v) 
+{
+    Vec<T, N> ret(u);
+    for (int i=0; i<N; ++i)
+    {
+        ret[i] -= v[i];
+    }
+    return ret;
+}
+
+template<typename T, int N>
+Vec<T, N> vec_scale(const Vec<T, N>& v, const T c) 
+{
+    Vec<T, N> ret(v);
+    for (int i=0; i<N; ++i)
+    {
+        ret[i] *= c;
+    }
+    return ret;
+}
+
+template <typename T, int N>
+Vec<T, N> vec_to_eigen(typename VecDataArray<T, N>::VecType* v2, const Vec<T, N>& v1)
+{
+    /* *v2 = v1; */
+    *v2 = v1.get();
+    return v1;
+}
+
+Vec<double, 6> vec_to_eigen_6d(imstk::Vec6d* v2, const Vec<double, 6>& v1)
+{
+    /* *v2 = v1; */
+    *v2 = v1.get();
+    return v1;
+}
+
+template <typename T, int N>
+Vec<T, N> vec_from_eigen(const typename VecDataArray<T, N>::VecType& v1)
+{
+    return Vec<T, N>(v1);
+}
+
+Vec<double, 6> vec_from_eigen_6d(const imstk::Vec6d& v1)
+{
+    return Vec<double, 6>(v1);
+}
+
+template <typename T, int N>
+void vecdataarray_push_back(VecDataArray<T, N>& vec, const Vec<T, N>& vi)
+{
+    vec.push_back(vi.get());
+}
+
+template <typename csType, typename EigenType>
+csType eigen2cs(const EigenType& eigen_data) {
+    return csType(eigen_data);
+}
+
+template <typename csType, typename EigenType>
+EigenType cs2Eigen(const csType& cs_data, EigenType* eigen_data) {
+    return *eigen_data = cs_data.get();
+}
+
+} /* end of namespace imstk */
+
+%} // end of inline
+
+
+%typemap(cscode) imstk::Vec<int, 3> %{
+    public int this[int i] {
+        get => getValue(i);
+        set => Utils.intPtr_assign(setValue(i), value);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_int_3_1_t (Vec3i v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec3i(SWIGTYPE_p_Eigen__MatrixT_int_3_1_t eigen_v) {
+        return Utils.vec_from_eigen_3i(eigen_v);
+    }
+%}
+
+%typemap(cscode) imstk::Vec<int, 4> %{
+    public int this[int i] {
+        get => getValue(i);
+        set => Utils.intPtr_assign(setValue(i), value);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_int_4_1_t (Vec4i v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec4i(SWIGTYPE_p_Eigen__MatrixT_int_4_1_t eigen_v) {
+        return Utils.vec_from_eigen_4i(eigen_v);
+    }
+%}
+
+%typemap(cscode) imstk::Vec<float, 2> %{
+    public float this[int i] {
+        get => getValue(i);
+        set => Utils.floatPtr_assign(setValue(i), value);
+    }
+    public static Vec2f operator * (Vec2f v, float c) {
+        return Utils.vec_scale_2f(v, c);
+    }
+
+    public static Vec2f operator * (float c, Vec2f v) {
+        return Utils.vec_scale_2f(v, c);
+    }
+    public static Vec2f operator / (Vec2f v, float c) {
+        return Utils.vec_scale_2f(v, (float)1.0/c);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_float_2_1_t (Vec2f v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec2f(SWIGTYPE_p_Eigen__MatrixT_float_2_1_t eigen_v) {
+        return Utils.vec_from_eigen_2f(eigen_v);
+    }
+%}
+
+%typemap(cscode) imstk::Vec<double, 2> %{
+    public double this[int i] {
+        get => getValue(i);
+        set => Utils.doublePtr_assign(setValue(i), value);
+    }
+    public static Vec2d operator * (Vec2d v, double c) {
+        return Utils.vec_scale_2d(v, c);
+    }
+
+    public static Vec2d operator * (double c, Vec2d v) {
+        return Utils.vec_scale_2d(v, c);
+    }
+
+    public static Vec2d operator / (Vec2d v, double c) {
+        return Utils.vec_scale_2d(v, 1.0/c);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_double_2_1_t (Vec2d v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec2d(SWIGTYPE_p_Eigen__MatrixT_double_2_1_t eigen_v) {
+        return Utils.vec_from_eigen_2d(eigen_v);
+    }
+%}
+
+%typemap(cscode) imstk::Vec<double, 3> %{
+    public double this[int i] {
+        get => getValue(i);
+        set => Utils.doublePtr_assign(setValue(i), value);
+    }
+    public static Vec3d operator + (Vec3d u, Vec3d v) {
+        return Utils.vec_add_3d(u, v);
+    }
+    public static Vec3d operator - (Vec3d u, Vec3d v) {
+        return Utils.vec_subtract_3d(u, v);
+    }
+    public static Vec3d operator * (Vec3d v, double c) {
+        return Utils.vec_scale_3d(v, c);
+    }
+
+    public static Vec3d operator * (double c, Vec3d v) {
+        return Utils.vec_scale_3d(v, c);
+    }
+    public static Vec3d operator / (Vec3d v, double c) {
+        return Utils.vec_scale_3d(v, 1.0/c);
+    }
+
+    public static Vec3d operator / (double c, Vec3d v) {
+        return Utils.vec_scale_3d(v, 1.0/c);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_double_3_1_t (Vec3d v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec3d(SWIGTYPE_p_Eigen__MatrixT_double_3_1_t eigen_v) {
+        return Utils.vec_from_eigen_3d(eigen_v);
+    }
+    public Vec3d normalized() {
+        double nrm= getValue(0) * getValue(0) + getValue(1) * getValue(1) + getValue(2) * getValue(2);
+        return this / System.Math.Sqrt(nrm);
+    }
+%}
+
+%typemap(cscode) imstk::Vec<double, 6> %{
+    public double this[int i] {
+        get => getValue(i);
+        set => Utils.doublePtr_assign(setValue(i), value);
+    }
+    public static Vec6d operator * (Vec6d v, double c) {
+        return Utils.vec_scale_6d(v, c);
+    }
+
+    public static Vec6d operator * (double c, Vec6d v) {
+        return Utils.vec_scale_6d(v, c);
+    }
+    public static Vec6d operator / (Vec6d v, double c) {
+        return Utils.vec_scale_6d(v, 1.0/c);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_double_6_1_t (Vec6d v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec6d(SWIGTYPE_p_Eigen__MatrixT_double_6_1_t eigen_v) {
+        return Utils.vec_from_eigen_6d(eigen_v);
+    }
+%}
+
+%typemap(cscode) imstk::Vec<unsigned char, 3> %{
+    /* public unsigned char this[int i] { */
+    public byte this[int i] {
+        get => getValue(i);
+        set => Utils.ucharPtr_assign(setValue(i), value);
+    }
+
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_unsigned_char_3_1_t (Vec3uc v) {
+        return v.get(); 
+    }
+
+    public static implicit operator Vec3uc(SWIGTYPE_p_Eigen__MatrixT_unsigned_char_3_1_t eigen_v) {
+        return Utils.vec_from_eigen_3uc(eigen_v);
+    }
+%}
+
+%typemap(cscode) imstk::VecDataArray<int, 3> %{
+    public Vec3i this[uint i] {
+        get => Utils.vec_from_eigen_3i(getValue(i));
+        set => Utils.vec_to_eigen_3i(setValue(i), value);
+    }
+    public void push_back(Vec3i v) {
+        Utils.vecdataarray_push_back_3i(this, v);
+    }
+%}
+
+%typemap(cscode) imstk::VecDataArray<int, 4> %{
+    public Vec4i this[uint i] {
+        get => Utils.vec_from_eigen_4i(getValue(i));
+        set => Utils.vec_to_eigen_4i(setValue(i), value);
+    }
+    public void push_back(Vec4i v) {
+        Utils.vecdataarray_push_back_4i(this, v);
+    }
+
+%}
+
+%typemap(cscode) imstk::VecDataArray<float, 2> %{
+    public Vec2f this[uint i] {
+        get => Utils.vec_from_eigen_2f(getValue(i));
+        set => Utils.vec_to_eigen_2f(setValue(i), value);
+    }
+    public void push_back(Vec2f v) {
+        Utils.vecdataarray_push_back_2f(this, v);
+    }
+%}
+
+%typemap(cscode) imstk::VecDataArray<double, 2> %{
+    public Vec2d this[uint i] {
+        get => Utils.vec_from_eigen_2d(getValue(i));
+        set => Utils.vec_to_eigen_2d(setValue(i), value);
+    }
+    public void push_back(Vec2d v) {
+        Utils.vecdataarray_push_back_2d(this, v);
+    }
+%}
+
+%typemap(cscode) imstk::VecDataArray<double, 3> %{
+    public Vec3d this[uint i] {
+        get => Utils.vec_from_eigen_3d(getValue(i));
+        set => Utils.vec_to_eigen_3d(setValue(i), value);
+    }
+    public void push_back(Vec3d v) {
+        Utils.vecdataarray_push_back_3d(this, v);
+    }
+%}
+
+%typemap(cscode) imstk::VecDataArray<unsigned char, 3> %{
+    public Vec3uc this[uint i] {
+        get => Utils.vec_from_eigen_3uc(getValue(i));
+        set => Utils.vec_to_eigen_3uc(setValue(i), value);
+    }
+    public void push_back(Vec3uc v) {
+        Utils.vecdataarray_push_back_3uc(this, v);
+    }
+%}
+
+
+%typemap(cscode) imstk::Mat<double, 3> %{
+    public static implicit operator SWIGTYPE_p_Eigen__MatrixT_double_3_3_t (Mat3d cs_data) {
+        return cs_data.get();
+    }
+
+    public static implicit operator Mat3d (SWIGTYPE_p_Eigen__MatrixT_double_3_3_t eigen_data) {
+        return Utils.eigen2cs_Mat3d(eigen_data);
+    }
+%}
+
+/* Make the renamed function private */
+%csmethodmodifiers operator[] "private";
+
+#endif /* end of #ifdef SWIGCSHARP */
+
+%include "../../Common/imstkMacros.h"
+%include "../../Common/imstkTypes.h";
+%include "../../Common/imstkMath.h";
+%include "../../Common/imstkEventObject.h";
+%include "../../Common/imstkAbstractDataArray.h";
+%include "../../Common/imstkDataArray.h";
+%template(DataArrayi) imstk::DataArray<int>;
+%template(DataArrayf) imstk::DataArray<float>;
+%template(DataArrayd) imstk::DataArray<double>;
+%template(DataArrayuc) imstk::DataArray<unsigned char>;
+%include "../../Common/imstkLogger.h";
+%include "../../Common/imstkModule.h";
+%include "../../Common/imstkModuleDriver.h";
+%include "../../Common/imstkColor.h";
+%include "../../Common/imstkVecDataArray.h";
+%template(VecDataArray3i) imstk::VecDataArray<int, 3>;
+%template(VecDataArray4i) imstk::VecDataArray<int, 4>;
+%template(VecDataArray2f) imstk::VecDataArray<float, 2>;
+%template(VecDataArray2d) imstk::VecDataArray<double, 2>;
+%template(VecDataArray3d) imstk::VecDataArray<double, 3>;
+%template(VecDataArray3uc) imstk::VecDataArray<unsigned char, 3>;
+
+/*
+ * Instantiation of Vec and VecDataArray
+ */
+%template(Vec3i) imstk::Vec<int, 3>;
+%template(Vec4i) imstk::Vec<int, 4>;
+%template(Vec2f) imstk::Vec<float, 2>;
+%template(Vec2d) imstk::Vec<double, 2>;
+%template(Vec3d) imstk::Vec<double, 3>;
+%template(Vec6d) imstk::Vec<double, 6>;
+%template(Vec3uc) imstk::Vec<unsigned char, 3>;
+
+%template(vec_scale_2i) imstk::vec_scale<int, 2>;
+%template(vec_scale_3i) imstk::vec_scale<int, 3>;
+%template(vec_scale_2f) imstk::vec_scale<float, 2>;
+%template(vec_scale_2d) imstk::vec_scale<double, 2>;
+%template(vec_scale_3d) imstk::vec_scale<double, 3>;
+%template(vec_scale_6d) imstk::vec_scale<double, 6>;
+%template(vec_add_2i) imstk::vec_add<int, 2>;
+%template(vec_add_3i) imstk::vec_add<int, 3>;
+%template(vec_add_2f) imstk::vec_add<float, 2>;
+%template(vec_add_2d) imstk::vec_add<double, 2>;
+%template(vec_add_3d) imstk::vec_add<double, 3>;
+%template(vec_add_6d) imstk::vec_add<double, 6>;
+%template(vec_subtract_2i) imstk::vec_subtract<int, 2>;
+%template(vec_subtract_3i) imstk::vec_subtract<int, 3>;
+%template(vec_subtract_2f) imstk::vec_subtract<float, 2>;
+%template(vec_subtract_2d) imstk::vec_subtract<double, 2>;
+%template(vec_subtract_3d) imstk::vec_subtract<double, 3>;
+%template(vec_subtract_6d) imstk::vec_subtract<double, 6>;
+%template(vec_to_eigen_3i) imstk::vec_to_eigen<int, 3>;
+%template(vec_to_eigen_4i) imstk::vec_to_eigen<int, 4>;
+%template(vec_to_eigen_2f) imstk::vec_to_eigen<float, 2>;
+%template(vec_to_eigen_2d) imstk::vec_to_eigen<double, 2>;
+%template(vec_to_eigen_3d) imstk::vec_to_eigen<double, 3>;
+%template(vec_to_eigen_3uc) imstk::vec_to_eigen<unsigned char, 3>;
+%template(vec_from_eigen_3i) imstk::vec_from_eigen<int, 3>;
+%template(vec_from_eigen_4i) imstk::vec_from_eigen<int, 4>;
+%template(vec_from_eigen_2f) imstk::vec_from_eigen<float, 2>;
+%template(vec_from_eigen_2d) imstk::vec_from_eigen<double, 2>;
+%template(vec_from_eigen_3d) imstk::vec_from_eigen<double, 3>;
+%template(vec_from_eigen_3uc) imstk::vec_from_eigen<unsigned char, 3>;
+%template(vecdataarray_push_back_3i) imstk::vecdataarray_push_back<int, 3>;
+%template(vecdataarray_push_back_4i) imstk::vecdataarray_push_back<int, 4>;
+%template(vecdataarray_push_back_2f) imstk::vecdataarray_push_back<float, 2>;
+%template(vecdataarray_push_back_2d) imstk::vecdataarray_push_back<double, 2>;
+%template(vecdataarray_push_back_3d) imstk::vecdataarray_push_back<double, 3>;
+%template(vecdataarray_push_back_3uc) imstk::vecdataarray_push_back<unsigned char, 3>;
+%template(cs2Eigen_Quatf) imstk::cs2Eigen<imstk::imstkQuatf, Eigen::Quaternion<float> >;
+%template(cs2Eigen_Quatd) imstk::cs2Eigen<imstk::imstkQuatd, Eigen::Quaternion<double> >;
+%template(eigen2cs_Quatf) imstk::eigen2cs<imstk::imstkQuatf, Eigen::Quaternion<float> >;
+%template(eigen2cs_Quatd) imstk::eigen2cs<imstk::imstkQuatd, Eigen::Quaternion<double> >;
+%template(Mat3d) imstk::Mat<double, 3>;
+%template(cs2Eigen_Mat3f) imstk::cs2Eigen<imstk::Mat<float, 3>, Eigen::Matrix<float, 3, 3> >;
+%template(cs2Eigen_Mat3d) imstk::cs2Eigen<imstk::Mat<double, 3>, Eigen::Matrix<double, 3, 3> >;
+%template(eigen2cs_Mat3f) imstk::eigen2cs<imstk::Mat<float, 3>, Eigen::Matrix<float, 3, 3> >;
+%template(eigen2cs_Mat3d) imstk::eigen2cs<imstk::Mat<double, 3>, Eigen::Matrix<double, 3, 3> >;
+
+%template(queueConnectKeyEvent) imstk::queueConnect<imstk::KeyEvent>;
+%template(connectEvent) imstk::connect<imstk::Event>;
diff --git a/Source/Wrappers/csharp/SwigInterface/debug.i b/Source/Wrappers/csharp/SwigInterface/debug.i
new file mode 100644
index 0000000000000000000000000000000000000000..bc3e7692e0a4163faeff2fe4b803fd892aa8dfec
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/debug.i
@@ -0,0 +1,67 @@
+%inline %{
+#include <atomic>
+class DebugAndTestOnly {
+public:
+    static constexpr int NULL_PTR = -2;
+    static constexpr int CAST_FAILURE = -1;
+
+    static int use_count(void* objPtr, const std::string& typeName)
+    {
+        if (objPtr == nullptr) return NULL_PTR;
+
+        if ("Geometry" == typeName)
+        {
+            std::shared_ptr<imstk::Geometry>* obj = (std::shared_ptr<imstk::Geometry>*) objPtr; 
+            return obj->use_count();
+        }
+        else if ("PointSet" == typeName)
+        {
+            std::shared_ptr<imstk::PointSet>* obj = (std::shared_ptr<imstk::PointSet>*) objPtr; 
+            return obj->use_count();
+        }
+        else if ("SurfaceMesh" == typeName)
+        {
+            std::shared_ptr<imstk::SurfaceMesh>* obj = (std::shared_ptr<imstk::SurfaceMesh>*) objPtr; 
+            return obj->use_count();
+        }
+        else if ("TetrahedralMesh" == typeName)
+        {
+            std::shared_ptr<imstk::TetrahedralMesh>* obj = (std::shared_ptr<imstk::TetrahedralMesh>*) objPtr; 
+            return obj->use_count();
+        }
+
+        return CAST_FAILURE;
+    }
+
+    DebugAndTestOnly()
+    {
+        printf("DebugAndTestOnly: constructed\n");
+        m_numCreated += 1; 
+    }
+
+    ~DebugAndTestOnly()
+    {
+        printf("DebugAndTestOnly: destructed\n");
+        m_numDestroyed += 1;
+    }
+
+    static int getNumCreated()
+    {
+        return m_numCreated;
+    }
+
+    static int getNumDestroyed()
+    {
+        return m_numDestroyed;
+    }
+private:
+    static std::atomic<int> m_numCreated;
+    static std::atomic<int> m_numDestroyed;
+};
+
+std::atomic<int> DebugAndTestOnly::m_numCreated{0};
+std::atomic<int> DebugAndTestOnly::m_numDestroyed{0};
+
+%}
+
+%shared_ptr(DebugAndTestOnly)
diff --git a/Source/Wrappers/csharp/SwigInterface/ignored.i b/Source/Wrappers/csharp/SwigInterface/ignored.i
new file mode 100644
index 0000000000000000000000000000000000000000..4b01332360ebc3e59fd2d9812bf3f4d3555f64e1
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/ignored.i
@@ -0,0 +1,67 @@
+%ignore imstk::PbdCollisionConstraint;
+/* %ignore imstk::PbdFEMConstraint; */
+%ignore imstk::PbdModel::getIntegratePositionNode();
+%ignore imstk::PbdModel::getUpdateCollisionGeometryNode();
+%ignore imstk::PbdModel::getSolveNode();
+%ignore imstk::PbdModel::getUpdateVelocityNode();
+
+%ignore imstk::DataArray::iterator; /* fix the multiple-definition problem. */
+%ignore imstk::DataArray::const_iterator; /* fix the multiple-definition problem. */
+%ignore imstk::DataArray::begin(); /* fix the multiple-definition problem. */
+%ignore imstk::DataArray::cbegin() const; /* fix the multiple-definition problem. */
+%ignore imstk::DataArray::end(); /* fix the multiple-definition problem. */
+%ignore imstk::DataArray::cend() const; /* fix the multiple-definition problem. */
+%ignore imstk::VecDataArray::iterator; /* fix the multiple-definition problem. */
+%ignore imstk::VecDataArray::const_iterator;
+%ignore imstk::VecDataArray::begin(); /* fix the multiple-definition problem. */
+%ignore imstk::VecDataArray::cbegin() const; /* fix the multiple-definition problem. */
+%ignore imstk::VecDataArray::end(); /* fix the multiple-definition problem. */
+%ignore imstk::VecDataArray::cend() const; /* fix the multiple-definition problem. */
+%ignore imstk::VecDataArray::setData();
+%ignore imstk::stdSink;
+%ignore imstk::LogManager;
+%ignore imstk::Logger::Logger();
+%ignore imstk::Logger::getInstance();
+%ignore imstk::Logger::addStdoutSink();
+%ignore imstk::Logger::addFileSink;
+%ignore imstk::Logger::addSink;
+%ignore imstk::Logger::initialize();
+%ignore imstk::Logger::destroy();
+%ignore imstk::Log;
+%rename("%s") imstk::Logger::startLogger();
+
+%ignore imstk::InteractionPair::getTaskNodeInputs();
+%ignore imstk::InteractionPair::getTaskNodeOutputs();
+%ignore imstk::CollisionGraph::getInteractionPairMap();
+
+%ignore imstk::FeDeformableObject::getFEMModel();
+%ignore imstk::FEMDeformableBodyModel::initializeEigenMatrixFromVegaMatrix;
+
+%ignore imstk::RbdConstraint;
+%ignore imstk::CollisionHandling::getTaskNode();
+
+%ignore imstk::VTKTextStatusManager::getTextActor();
+%ignore imstk::AbstractVTKViewer::getVtkRenderWindow() const;
+
+%ignore imstk::GeometryUtils::coupleVtkDataArray;
+%ignore imstk::GeometryUtils::coupleVtkImageData;
+%ignore imstk::GeometryUtils::copyToVtkDataArray;
+%ignore imstk::GeometryUtils::copyToVtkImageData;
+%ignore imstk::GeometryUtils::copyToPointSet;
+%ignore imstk::GeometryUtils::copyToSurfaceMesh;
+%ignore imstk::GeometryUtils::copyToLineMesh;
+%ignore imstk::GeometryUtils::copyToVolumetricMesh;
+%ignore imstk::GeometryUtils::copyToVtkPointSet;
+%ignore imstk::GeometryUtils::copyToVtkPolyData;
+%ignore imstk::GeometryUtils::copyToVtkUnstructuredGrid;
+%ignore imstk::GeometryUtils::copyToVtkVecDataArray;
+%ignore imstk::GeometryUtils::copyToVtkPoints;
+%ignore imstk::GeometryUtils::copyToVtkCellArray;
+%ignore imstk::GeometryUtils::copyToDataArray;
+%ignore imstk::GeometryUtils::copyToVecDataArray;
+%ignore imstk::GeometryUtils::copyToImageData;
+
+%ignore imstk::CollisionElement::m_element;
+%ignore imstk::CollisionElement::m_type;
+%ignore imstk::CollisionData::elementsA;
+%ignore imstk::CollisionData::elementsB;
diff --git a/Source/Wrappers/csharp/SwigInterface/imstkCWrapper.i b/Source/Wrappers/csharp/SwigInterface/imstkCWrapper.i
new file mode 100644
index 0000000000000000000000000000000000000000..421f58e3b7a6e53e1587ede2abf3c5a49e7c9481
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/imstkCWrapper.i
@@ -0,0 +1,407 @@
+%module(directors="1") Utils
+#pragma SWIG nowarn=302,314,317,401,476,501,503,505,516,844,
+%{
+/* Common */
+#include "imstkMacros.h"
+#include "imstkMath.h"
+#include "imstkAbstractDataArray.h"
+#include "imstkDataArray.h"
+#include "imstkVecDataArray.h"
+#include "imstkLogger.h"
+#include "imstkModule.h"
+#include "imstkModuleDriver.h"
+#include "imstkColor.h"
+#include "imstkEventObject.h"
+#include "imstkTypes.h"
+
+/*
+ * DataStructures
+ */
+#include "imstkNeighborSearch.h"
+
+/* 
+ * Geometry 
+ */
+#include "imstkGeometry.h"
+#include "imstkGeometryUtilities.h"
+#include "imstkPointSet.h"
+#include "imstkSurfaceMesh.h"
+#include "imstkLineMesh.h"
+#include "imstkImageData.h"
+#include "imstkVolumetricMesh.h"
+#include "imstkTetrahedralMesh.h"
+#include "imstkHexahedralMesh.h"
+#include "imstkImplicitGeometry.h"
+#include "imstkAnalyticalGeometry.h"
+#include "imstkCompositeImplicitGeometry.h"
+#include "imstkPlane.h"
+#include "imstkSphere.h"
+#include "imstkOrientedBox.h"
+#include "imstkCapsule.h"
+#include "imstkGeometryUtilities.h"
+#include "imstkSignedDistanceField.h"
+#include "imstkImplicitFunctionFiniteDifferenceFunctor.h"
+
+/*
+ * GeometryMappers
+ */
+#include "imstkGeometryMap.h"
+#include "imstkOneToOneMap.h"
+#include "imstkTetraTriangleMap.h"
+
+/*
+ * Filter
+ */
+#include "imstkGeometryAlgorithm.h"
+#include "imstkSurfaceMeshSubdivide.h"
+#include "imstkImplicitGeometryToImageData.h"
+#include "imstkSurfaceMeshFlyingEdges.h"
+
+/* 
+ * MeshIO 
+ */
+#include "imstkMeshIO.h"
+
+/* 
+ * DynamicalModel 
+ */
+#include "imstkVectorizedState.h"
+#include "imstkPbdState.h"
+#include "imstkAbstractDynamicalModel.h"
+#include "imstkDynamicalModel.h"
+#include "imstkPbdModel.h"
+#include "imstkTimeIntegrator.h"
+#include "imstkBackwardEuler.h"
+#include "imstkPbdFEMConstraint.h"
+#include "imstkPbdCollisionConstraint.h"
+#include "imstkSPHBoundaryConditions.h"
+/* #include "imstkSPHHemorrhage.h" */
+#include "imstkInternalForceModelTypes.h"
+#include "imstkFEMDeformableBodyModel.h"
+#include "imstkRigidBodyState2.h"
+#include "imstkRigidBodyModel2.h"
+#include "imstkSPHState.h"
+#include "imstkSPHModel.h"
+
+/* 
+ * Rendering
+ */
+#include "imstkRenderMaterial.h"
+#include "imstkTexture.h"
+
+/*
+ * Constraints
+ */
+#include "imstkPbdConstraint.h"
+#include "imstkRbdConstraint.h"
+
+/* 
+ * SceneEntities
+ */
+#include "imstkSceneEntity.h"
+#include "imstkSceneObject.h"
+#include "imstkCollidingObject.h"
+#include "imstkDynamicObject.h"
+#include "imstkPbdObject.h"
+#include "imstkVisualModel.h"
+#include "imstkCamera.h"
+#include "imstkLight.h"
+#include "imstkDirectionalLight.h"
+#include "imstkFeDeformableObject.h"
+#include "imstkRigidObject2.h"
+#include "imstkSPHObject.h"
+
+/*
+ * CollisionDetection
+ */
+/* #include "imstkCollisionDetection.h" */
+#include "imstkCollisionData.h"
+#include "imstkCollisionDetectionAlgorithm.h"
+#include "imstkBidirectionalPlaneToSphereCD.h"
+#include "imstkCollisionDetectionAlgorithm.h"
+#include "imstkCollisionUtils.h"
+#include "imstkImplicitGeometryToPointSetCCD.h"
+#include "imstkImplicitGeometryToPointSetCD.h"
+#include "imstkMeshToMeshBruteForceCD.h"
+#include "imstkPointSetToCapsuleCD.h"
+#include "imstkPointSetToOrientedBoxCD.h"
+#include "imstkPointSetToPlaneCD.h"
+#include "imstkPointSetToSphereCD.h"
+#include "imstkSphereToCylinderCD.h"
+#include "imstkSphereToSphereCD.h"
+#include "imstkSurfaceMeshToCapsuleCD.h"
+#include "imstkSurfaceMeshToSphereCD.h"
+#include "imstkSurfaceMeshToSurfaceMeshCD.h"
+#include "imstkTetraToLineMeshCD.h"
+#include "imstkTetraToPointSetCD.h"
+#include "imstkUnidirectionalPlaneToSphereCD.h"
+
+
+/*
+ * CollisionHandling
+ */
+#include "imstkCollisionHandling.h"
+#include "imstkRigidBodyCH.h"
+/*
+ * Controller
+ */
+#include "imstkDeviceControl.h"
+#include "imstkMouseControl.h"
+#include "imstkKeyboardControl.h"
+#include "imstkTrackingDeviceControl.h"
+#include "imstkSceneObjectController.h"
+
+/*
+ * Scene
+ */
+#include "imstkScene.h"
+#include "imstkCollisionGraph.h"
+#include "imstkObjectInteractionPair.h"
+#include "imstkObjectInteractionFactory.h"
+#include "imstkCollisionPair.h"
+#include "imstkRigidObjectCollision.h"
+#include "imstkInteractionPair.h"
+#include "imstkObjectInteractionPair.h"
+#include "imstkPbdObjectCuttingPair.h"
+#include "imstkPbdObjectCollision.h"
+#include "imstkSphObjectCollision.h"
+
+/*
+ * SimulationManager
+ */
+#include "imstkModule.h"
+#include "imstkViewer.h"
+#include "imstkAbstractVTKViewer.h"
+#include "imstkVTKViewer.h"
+#include "imstkVTKTextStatusManager.h"
+#include "imstkSceneManager.h"
+#include "imstkSimulationManager.h"
+#include "imstkMouseSceneControl.h"
+#include "imstkKeyboardSceneControl.h"
+
+/*
+ * Devices
+ */
+#include "imstkDeviceClient.h"
+#include "imstkKeyboardDeviceClient.h"
+#ifdef iMSTK_USE_OpenHaptics
+#include "imstkHapticDeviceManager.h"
+#include "imstkHapticDeviceClient.h"
+#endif
+
+%} /* end of module */
+
+/*
+ * stl
+ */
+%include <stdint.i>
+%include <std_string.i>
+%include <std_vector.i>
+namespace std {
+  %template(VectorInt) vector<int>; 
+  %template(VectorSizet) vector<std::size_t>;
+  %template(VectorCollisionData) vector<imstk::CollisionElement>;
+}
+
+%include "shared_ptr_instantiation.i"
+%include "weak_ptr.i"
+%include "ignored.i"
+%include "modifiers.i"
+%include "type_cast.i"
+%include "std_function.i"
+%include "callback.i"
+
+%rename(compute) imstk::ImplicitFunctionGradient::operator();
+%rename(compute) imstk::ImplicitFunctionCentralGradient::operator();
+
+/*
+ * * * * * * * * * * * * * * * *
+ * list C/C++ declarations
+ * * * * * * * * * * * * * * * * *
+ */
+/*
+ * Common
+ */
+%include "common.i"
+
+/*
+ * DataStructures
+ */
+%include "../../../DataStructures/imstkNeighborSearch.h"
+
+/*
+ * Geometry
+ */
+%include "../../../Geometry/imstkGeometry.h";
+%include "../../../Geometry/Mesh/imstkPointSet.h"
+%include "../../../Geometry/Mesh/imstkSurfaceMesh.h"
+%include "../../../Geometry/Mesh/imstkLineMesh.h"
+%include "../../../Geometry/Mesh/imstkImageData.h"
+%include "../../../Geometry/Mesh/imstkVolumetricMesh.h"
+%include "../../../Geometry/Mesh/imstkTetrahedralMesh.h"
+%include "../../../Geometry/Mesh/imstkHexahedralMesh.h"
+%include "../../../Geometry/Implicit/imstkImplicitGeometry.h"
+%include "../../../Geometry/Implicit/imstkCompositeImplicitGeometry.h"
+%include "../../../Geometry/Analytic/imstkAnalyticalGeometry.h"
+%include "../../../Geometry/Analytic/imstkPlane.h"
+%include "../../../Geometry/Analytic/imstkSphere.h"
+%include "../../../Geometry/Analytic/imstkOrientedBox.h"
+%include "../../../Geometry/Analytic/imstkCapsule.h"
+%include "../../../Geometry/imstkGeometryUtilities.h"
+%include "../../../Geometry/Implicit/imstkSignedDistanceField.h"
+%include "../../../Geometry/Implicit/imstkImplicitFunctionFiniteDifferenceFunctor.h"
+
+/*
+ * GeometryMap
+ */
+%include "../../../GeometryMappers/imstkGeometryMap.h"
+%include "../../../GeometryMappers/imstkOneToOneMap.h"
+%include "../../../GeometryMappers/imstkTetraTriangleMap.h"
+
+/*
+ * Filtering
+ */
+%include "../../../FilteringCore/imstkGeometryAlgorithm.h"
+%include "../../../Filtering/imstkSurfaceMeshSubdivide.h"
+%include "../../../Filtering/imstkImplicitGeometryToImageData.h"
+%include "../../../Filtering/imstkSurfaceMeshFlyingEdges.h"
+
+/*
+ * MeshIO
+ */
+%include "../../../MeshIO/imstkMeshIO.h";
+%template(readImageData) imstk::MeshIO::read<imstk::ImageData>;
+%template(readSurfaceMesh) imstk::MeshIO::read<imstk::SurfaceMesh>;
+%template(readTetrahedralMesh) imstk::MeshIO::read<imstk::TetrahedralMesh>;
+
+/*
+ * Constraint
+ */
+%include "../../../Constraint/PbdConstraints/imstkPbdConstraint.h"
+%include "../../../Constraint/PbdConstraints/imstkPbdCollisionConstraint.h"
+%include "../../../Constraint/PbdConstraints/imstkPbdFEMConstraint.h"
+%include "../../../Constraint/RigidBodyConstraints/imstkRbdConstraint.h"
+
+/*
+ * DynamicalModel
+ */
+%include "../../../DynamicalModels/ObjectStates/imstkVectorizedState.h"
+%include "../../../DynamicalModels/ObjectStates/imstkPbdState.h"
+%include "../../../DynamicalModels/ObjectModels/imstkAbstractDynamicalModel.h"
+%include "../../../DynamicalModels/ObjectModels/imstkDynamicalModel.h"
+/* Instantiation of base class should be put before the derived class */
+%template(DynamicalModelPbdState) imstk::DynamicalModel<imstk::PbdState>;
+%include "../../../DynamicalModels/ObjectModels/imstkPbdModel.h"
+%template(DynamicalModelFeDeformBodyState) imstk::DynamicalModel<imstk::FeDeformBodyState>;
+%include "../../../DynamicalModels/InternalForceModel/imstkInternalForceModelTypes.h"
+%include "../../../DynamicalModels/ObjectModels/imstkFEMDeformableBodyModel.h"
+%include "../../../DynamicalModels/ObjectModels/imstkSPHBoundaryConditions.h"
+/* %include "../../../DynamicalModels/ObjectModels/imstkSPHHemorrhage.h" */
+%include "../../../DynamicalModels/TimeIntegrators/imstkTimeIntegrator.h"
+%include "../../../DynamicalModels/TimeIntegrators/imstkBackwardEuler.h"
+%include "../../../DynamicalModels/ObjectStates/imstkRigidBodyState2.h"
+%template(DynamicalModelRigidBodyState2) imstk::DynamicalModel<imstk::RigidBodyState2>;
+%include "../../../DynamicalModels/ObjectModels/imstkRigidBodyModel2.h"
+%include "../../../DynamicalModels/ObjectStates/imstkSPHState.h"
+%template(DynamicalModelSPHState) imstk::DynamicalModel<imstk::SPHState>;
+%include "../../../DynamicalModels/ObjectModels/imstkSPHModel.h"
+
+/* 
+ * Rendering 
+ */
+%include "../../../Rendering/Materials/imstkRenderMaterial.h";
+%include "../../../Rendering/Materials/imstkTexture.h";
+
+/*
+ * SceneEntities
+ */
+%include "../../../SceneEntities/imstkSceneEntity.h"
+%include "../../../SceneEntities/Objects/imstkSceneObject.h";
+%include "../../../SceneEntities/Objects/imstkCollidingObject.h";
+%include "../../../SceneEntities/Objects/imstkDynamicObject.h";
+%include "../../../SceneEntities/Objects/imstkPbdObject.h";
+%include "../../../SceneEntities/Objects/imstkVisualModel.h";
+%include "../../../SceneEntities/Objects/imstkFeDeformableObject.h";
+%include "../../../SceneEntities/Objects/imstkRigidObject2.h";
+%include "../../../SceneEntities/Objects/imstkSPHObject.h";
+%include "../../../SceneEntities/Camera/imstkCamera.h";
+%include "../../../SceneEntities/Lights/imstkLight.h";
+%include "../../../SceneEntities/Lights/imstkDirectionalLight.h";
+
+/*
+ * CollisionDetection
+ */
+/* %include "../../../CollisionDetection/CollisionDetection/imstkCollisionDetection.h"; */
+%include "../../../CollisionDetection/CollisionData/imstkCollisionData.h";
+%include "../../../CollisionDetection/CollisionDetection/imstkCollisionDetectionAlgorithm.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkBidirectionalPlaneToSphereCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkCollisionUtils.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkImplicitGeometryToPointSetCCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkImplicitGeometryToPointSetCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkMeshToMeshBruteForceCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkPointSetToCapsuleCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkPointSetToOrientedBoxCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkPointSetToPlaneCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkPointSetToSphereCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkSphereToCylinderCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkSphereToSphereCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkSurfaceMeshToCapsuleCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkSurfaceMeshToSphereCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkSurfaceMeshToSurfaceMeshCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkTetraToLineMeshCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkTetraToPointSetCD.h"
+%include "../../../CollisionDetection/CollisionDetection/imstkUnidirectionalPlaneToSphereCD.h"
+
+
+
+/*
+ * CollisionHandling
+ */ 
+%include "../../../CollisionHandling/imstkCollisionHandling.h";
+%include "../../../CollisionHandling/imstkRigidBodyCH.h";
+
+/* 
+ * Controllers
+ */
+%include "../../../Controllers/imstkDeviceControl.h"
+%include "../../../Controllers/imstkMouseControl.h"
+%include "../../../Controllers/imstkKeyboardControl.h"
+%include "../../../Controllers/imstkTrackingDeviceControl.h"
+%include "../../../Controllers/imstkSceneObjectController.h"
+
+/* 
+ * Scene
+ */
+%include "../../../Scene/imstkScene.h";
+%include "../../../Scene/imstkCollisionGraph.h";
+%include "../../../Scene/imstkInteractionPair.h"
+%include "../../../Scene/imstkObjectInteractionPair.h";
+%include "../../../Scene/imstkObjectInteractionFactory.h";
+%include "../../../Scene/imstkCollisionPair.h";
+%include "../../../Scene/imstkRigidObjectCollision.h";
+%include "../../../Scene/imstkPbdObjectCuttingPair.h"
+%include "../../../Scene/imstkPbdObjectCollision.h"
+%include "../../../Scene/imstkSphObjectCollision.h"
+
+/*
+ * SimulationManager
+ */
+%include "../../../SimulationManager/imstkViewer.h";
+%include "../../../SimulationManager/VTKRenderer/imstkAbstractVTKViewer.h";
+%include "../../../SimulationManager/VTKRenderer/imstkVTKViewer.h";
+%include "../../../SimulationManager/VTKRenderer/imstkVTKTextStatusManager.h";
+%include "../../../SimulationManager/imstkSceneManager.h"
+%include "../../../SimulationManager/imstkSimulationManager.h"
+%include "../../../SimulationManager/imstkMouseSceneControl.h"
+%include "../../../SimulationManager/imstkKeyboardSceneControl.h"
+
+/*
+ * Devices
+ */
+%include "../../../Devices/imstkDeviceClient.h"
+%include "../../../Devices/imstkKeyboardDeviceClient.h"
+#ifdef iMSTK_USE_OpenHaptics
+%include "../../../Devices/imstkHapticDeviceManager.h"
+%include "../../../Devices/imstkHapticDeviceClient.h"
+#endif
+
diff --git a/Source/Wrappers/csharp/SwigInterface/modifiers.i b/Source/Wrappers/csharp/SwigInterface/modifiers.i
new file mode 100644
index 0000000000000000000000000000000000000000..7a17babba7018e205ed96bc675312f89cf41386d
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/modifiers.i
@@ -0,0 +1,5 @@
+%csmethodmodifiers imstk::VecDataArray::getPointer "new public";
+%csmethodmodifiers imstk::AbstractVTKViewer::setSize "override public"
+%csmethodmodifiers imstk::VTKViewer::getScreenCaptureUtility "new public"
+
+
diff --git a/Source/Wrappers/csharp/SwigInterface/shared_ptr_instantiation.i b/Source/Wrappers/csharp/SwigInterface/shared_ptr_instantiation.i
new file mode 100644
index 0000000000000000000000000000000000000000..896a0d74f8ff0a327ba3d6eb6a732c70fd34cab0
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/shared_ptr_instantiation.i
@@ -0,0 +1,188 @@
+/*
+ * Instantiation of shared_ptr
+ */
+%include <std_shared_ptr.i>
+
+/*
+ * Common
+ */
+%shared_ptr(imstk::AbstractDataArray)
+%shared_ptr(imstk::DataArray<int>)
+%shared_ptr(imstk::DataArray<float>)
+%shared_ptr(imstk::DataArray<double>)
+%shared_ptr(imstk::DataArray<unsigned char>)
+%shared_ptr(imstk::VecDataArray<int, 3>)
+%shared_ptr(imstk::VecDataArray<int, 4>)
+%shared_ptr(imstk::VecDataArray<float, 2>)
+%shared_ptr(imstk::VecDataArray<double, 2>)
+%shared_ptr(imstk::VecDataArray<double, 3>)
+%shared_ptr(imstk::VecDataArray<unsigned char, 3>)
+%shared_ptr(imstk::ModuleDriver)
+%shared_ptr(imstk::Event)
+%shared_ptr(imstk::EventObject)
+%shared_ptr(imstk::ButtonEvent)
+%shared_ptr(imstk::KeyEvent)
+
+/* 
+ * Geometry 
+ */
+%shared_ptr(imstk::Geometry)
+%shared_ptr(imstk::PointSet)
+%shared_ptr(imstk::ImageData)
+%shared_ptr(imstk::LineMesh)
+%shared_ptr(imstk::PointSet)
+%shared_ptr(imstk::SurfaceMesh)
+%shared_ptr(imstk::VolumetricMesh)
+%shared_ptr(imstk::TetrahedralMesh)
+%shared_ptr(imstk::HexahedralMesh)
+%shared_ptr(imstk::ImplicitGeometry)
+%shared_ptr(imstk::AnalyticalGeometry)
+%shared_ptr(imstk::CompositeImplicitGeometry)
+%shared_ptr(imstk::Plane)
+%shared_ptr(imstk::Sphere)
+%shared_ptr(imstk::OrientedBox)
+%shared_ptr(imstk::Capsule)
+%shared_ptr(imstk::SignedDistanceField)
+
+/* 
+ * GeometryMap
+ */
+%shared_ptr(imstk::GeometryMap)
+%shared_ptr(imstk::OneToOneMap)
+%shared_ptr(imstk::TetraTriangleMap)
+
+/*
+ * FilteringCore
+ */
+%shared_ptr(imstk::GeometryAlgorithm)
+%shared_ptr(imstk::SurfaceMeshSubdivide)
+%shared_ptr(imstk::ImplicitGeometryToImageData)
+%shared_ptr(imstk::SurfaceMeshFlyingEdges)
+
+/* 
+ * DynamicalModel 
+ */
+%shared_ptr(imstk::PbdObject)
+%shared_ptr(imstk::PBDModelConfig)
+%shared_ptr(imstk::PbdCollisionConstraintConfig)
+%shared_ptr(imstk::PbdFEMConstraintConfig)
+%shared_ptr(imstk::PbdSolver)
+%shared_ptr(imstk::FeDeformBodyState)
+%shared_ptr(imstk::FEMModelConfig)
+%shared_ptr(imstk::FEMDeformableBodyModel)
+%shared_ptr(imstk::AbstractDynamicalModel)
+%shared_ptr(imstk::DynamicalModel<imstk::PbdState>)
+%shared_ptr(imstk::DynamicalModel<imstk::FeDeformBodyState>)
+%shared_ptr(imstk::DynamicalModel<imstk::SPHState>)
+%shared_ptr(imstk::RigidBodyState2)
+%shared_ptr(imstk::DynamicalModel<imstk::RigidBodyState2>)
+%shared_ptr(imstk::PbdModel)
+%shared_ptr(imstk::RigidBodyModel2Config)
+%shared_ptr(imstk::RigidBodyModel2)
+%shared_ptr(imstk::SPHModelConfig)
+%shared_ptr(imstk::SPHModel)
+%shared_ptr(imstk::TimeIntegrator)
+%shared_ptr(imstk::BackwardEuler)
+
+/* 
+ * Rendering
+ */
+%shared_ptr(imstk::RenderMaterial)
+%shared_ptr(imstk::Texture)
+
+/*
+ * Constraint
+ */
+%shared_ptr(imstk::RigidBody)
+
+/* 
+ * SceneEntities
+ */
+%shared_ptr(imstk::SceneEntity)
+%shared_ptr(imstk::SceneObject)
+%shared_ptr(imstk::CollidingObject)
+%shared_ptr(imstk::DynamicObject)
+%shared_ptr(imstk::PbdObject)
+%shared_ptr(imstk::FeDeformableObject)
+%shared_ptr(imstk::SPHObject)
+%shared_ptr(imstk::RigidObject2)
+%shared_ptr(imstk::VisualModel)
+%shared_ptr(imstk::Camera)
+%shared_ptr(imstk::Light)
+%shared_ptr(imstk::DirectionalLight)
+%shared_ptr(imstk::PointLight)
+%shared_ptr(imstk::SpotLight)
+
+/*
+ * CollisionDetection
+ */
+%shared_ptr(imstk::CollisionDetectionAlgorithm)
+%shared_ptr(imstk::BidirectionalPlaneToSphereCD)
+%shared_ptr(imstk::ImplicitGeometryToPointSetCCD)
+%shared_ptr(imstk::ImplicitGeometryToPointSetCD)
+%shared_ptr(imstk::MeshToMeshBruteForceCD)
+%shared_ptr(imstk::PointSetToCapsuleCD)
+%shared_ptr(imstk::PointSetToOrientedBoxCD)
+%shared_ptr(imstk::PointSetToSphereCD)
+%shared_ptr(imstk::PointSetToPlaneCD)
+%shared_ptr(imstk::SphereToCylinderCD)
+%shared_ptr(imstk::SphereToSphereCD)
+%shared_ptr(imstk::SurfaceMeshToCapsuleCD)
+%shared_ptr(imstk::SurfaceMeshToSphereCD)
+%shared_ptr(imstk::SurfaceMeshToSurfaceMeshCD)
+%shared_ptr(imstk::TetraToLineMeshCD)
+%shared_ptr(imstk::TetraToPointSetCD)
+%shared_ptr(imstk::UnidirectionalPlaneToSphereCD)
+
+/*
+ * CollisionHandling
+ */
+%shared_ptr(imstk::CollisionPair)
+%shared_ptr(imstk::CollisionHandling)
+%shared_ptr(imstk::RigidObjectCollidingCollisionPair)
+%shared_ptr(imstk::RigidBodyCH)
+
+/*
+ * Controller
+ */
+%shared_ptr(imstk::DeviceControl)
+%shared_ptr(imstk::MouseControl)
+%shared_ptr(imstk::KeyboardControl)
+%shared_ptr(imstk::TrackingDeviceControl)
+%shared_ptr(imstk::SceneObjectController)
+
+/*
+ * Scene
+ */
+%shared_ptr(imstk::SceneConfig)
+%shared_ptr(imstk::Scene)
+%shared_ptr(imstk::CollisionGraph)
+%shared_ptr(imstk::RigidObjectCollision)
+%shared_ptr(imstk::InteractionPair)
+%shared_ptr(imstk::ObjectInteractionPair)
+%shared_ptr(imstk::PbdObjectCuttingPair)
+%shared_ptr(imstk::PbdObjectCollision)
+%shared_ptr(imstk::SphObjectCollision)
+
+/*
+ * SimulationManager
+ */
+%shared_ptr(imstk::Module)
+%shared_ptr(imstk::Viewer)
+%shared_ptr(imstk::AbstractVTKViewer)
+%shared_ptr(imstk::VTKViewer)
+%shared_ptr(imstk::VTKTextStatusManager)
+%shared_ptr(imstk::SceneManager)
+%shared_ptr(imstk::SimulationManager)
+%shared_ptr(imstk::MouseSceneControl)
+%shared_ptr(imstk::KeyboardSceneControl)
+
+/*
+ * Devices
+ */
+%shared_ptr(imstk::DeviceClient)
+%shared_ptr(imstk::KeyboardDeviceClient)
+%shared_ptr(imstk::HapticDeviceClient)
+%shared_ptr(imstk::HapticDeviceManager)
+
+
diff --git a/Source/Wrappers/csharp/SwigInterface/std_function.i b/Source/Wrappers/csharp/SwigInterface/std_function.i
new file mode 100644
index 0000000000000000000000000000000000000000..e2b8c8c0db5d8ba22a477f4e8bca4cf16701c45a
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/std_function.i
@@ -0,0 +1,110 @@
+%{
+  #include <functional>
+%}
+
+// These are the things we actually use
+#define param(num, type) $typemap(cstype, type) arg##num
+#define unpack(num, type) arg##num
+#define lvalref(num, type) type&& arg##num
+#define rvalref(num, mytype) std::add_rvalue_reference<mytype>::type arg##num
+#define forward(num, type) std::forward<type>(arg##num)
+
+// This is the mechanics
+#define FE_0(...)
+#define FE_1(action, a1) action(0, a1)
+#define FE_2(action, a1, a2) action(0, a1), action(1, a2)
+#define FE_3(action, a1, a2, a3) action(0, a1), action(1, a2), action(2, a3)
+#define FE_4(action, a1, a2, a3, a4) action(0, a1), action(1, a2), action(2, a3), action(3, a4)
+#define FE_5(action, a1, a2, a3, a4, a5) action(0, a1), action(1, a2), action(2, a3), action(3, a4), action(4, a5)
+
+#define GET_MACRO(_1, _2, _3, _4, _5, NAME, ...) NAME
+%define FOR_EACH(action, ...) 
+  GET_MACRO(__VA_ARGS__, FE_5, FE_4, FE_3, FE_2, FE_1, FE_0)(action, __VA_ARGS__)
+%enddef
+
+%define %std_function(Name, Ret, ...)
+
+%feature("director") Name##Impl;
+%typemap(csclassmodifiers) Name##Impl "public abstract class";
+
+%{
+  struct Name##Impl {
+    virtual ~Name##Impl() {}
+    void do_nothing() {}
+    virtual Ret call(__VA_ARGS__) = 0;
+  };
+%}
+
+%csmethodmodifiers Name##Impl::call "public abstract";
+%typemap(csout) Ret Name##Impl::call ";" // Suppress the body of the abstract method
+
+struct Name##Impl {
+  virtual ~Name##Impl();
+protected:
+  virtual Ret call(__VA_ARGS__) = 0;
+};
+
+%typemap(maybereturn) SWIGTYPE "return ";
+%typemap(maybereturn) void "";
+%typemap(Functor) SWIGTYPE "Func"
+%typemap(Functor) void "Action"
+
+%typemap(csin) std::function<Ret(__VA_ARGS__)> "$csclassname.getCPtr($csclassname.makeNative($csinput))"
+%typemap(cscode) std::function<Ret(__VA_ARGS__)> %{
+
+  private class Name##ImplCB : Name##Impl {
+    public Name##ImplCB($csclassname func) {
+      func_ = func;
+    }
+
+    public override Ret call (FOR_EACH(param, __VA_ARGS__)) {
+      $typemap(maybereturn, Ret) func_.call(FOR_EACH(unpack, __VA_ARGS__));
+    }
+
+    protected $csclassname func_;
+  }
+
+  public Name() {
+    wrapper = new Name##ImplCB(this);
+    proxy = new $csclassname(wrapper);
+  }
+
+
+  public static $csclassname makeNative($csclassname func_in) {
+    if (null == func_in.wrapper) return func_in;
+    return func_in.proxy;
+  }
+
+  private Name##Impl wrapper;
+  private $csclassname proxy;
+%}
+
+%rename(Name) std::function<Ret(__VA_ARGS__)>;
+%rename(call) std::function<Ret(__VA_ARGS__)>::operator();
+
+namespace std {
+  struct function<Ret(__VA_ARGS__)> {
+    // Copy constructor
+    function<Ret(__VA_ARGS__)>(const std::function<Ret(__VA_ARGS__)>&);
+
+    // Call operator
+    virtual Ret operator()(__VA_ARGS__) const;
+
+    // Conversion constructor from function pointer
+    function<Ret(__VA_ARGS__)>(Ret(*const)(__VA_ARGS__));
+
+    %extend {
+      function<Ret(__VA_ARGS__)>(Name##Impl *in) {
+        return new std::function<Ret(__VA_ARGS__)>(
+          [=](FOR_EACH(lvalref, __VA_ARGS__)) {
+            return in->call(FOR_EACH(forward, __VA_ARGS__));
+          });
+      }
+    }
+  };
+}
+
+%enddef
+
+%std_function(EventFunc, void, imstk::Event*)
+%std_function(ReceiverFunc, void, imstk::KeyEvent*)
diff --git a/Source/Wrappers/csharp/SwigInterface/type_cast.i b/Source/Wrappers/csharp/SwigInterface/type_cast.i
new file mode 100644
index 0000000000000000000000000000000000000000..105e198b05858c207ef25f21bbfd17f7453a13cd
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/type_cast.i
@@ -0,0 +1,40 @@
+%inline %{
+/* template <typename Base, typename Derived> */
+/* Derived* type_cast(Base* base) { */
+/*     return dynamic_cast<Derived*>(base); */
+/* } */
+template <typename Base, typename Derived>
+std::shared_ptr<Derived> type_cast(std::shared_ptr<Base> base) {
+    return std::dynamic_pointer_cast<Derived>(base);
+}
+%}
+
+%template(castToRigidBodyCH) type_cast<imstk::CollisionHandling, imstk::RigidBodyCH>;
+%template(castToVecDataArray3uc) type_cast<imstk::AbstractDataArray, imstk::VecDataArray<unsigned char, 3>>;
+
+%template(castToPointSet) type_cast<imstk::Geometry, imstk::PointSet>;
+%template(castToLineMesh) type_cast<imstk::Geometry, imstk::LineMesh>;
+%template(castToLineMesh) type_cast<imstk::PointSet, imstk::LineMesh>;
+%template(castToSurfaceMesh) type_cast<imstk::Geometry, imstk::SurfaceMesh>;
+%template(castToSurfaceMesh) type_cast<imstk::PointSet, imstk::SurfaceMesh>;
+%template(castToTetrahedralMesh) type_cast<imstk::Geometry, imstk::TetrahedralMesh>;
+%template(castToTetrahedralMesh) type_cast<imstk::PointSet, imstk::TetrahedralMesh>;
+%template(castToHexahedral) type_cast<imstk::Geometry, imstk::HexahedralMesh>;
+%template(castToHexahedral) type_cast<imstk::PointSet, imstk::HexahedralMesh>;
+%template(castToImageData) type_cast<imstk::Geometry, imstk::ImageData>;
+%template(castToImageData) type_cast<imstk::PointSet, imstk::ImageData>;
+
+%pragma(csharp) modulecode=%{
+  public static T CastTo<T>(object from, bool cMemoryOwn=false)
+  {
+      System.Reflection.MethodInfo CPtrGetter = from.GetType().GetMethod("getCPtr", System.Reflection.BindingFlags.NonPublic | System.Reflection.BindingFlags.Static);
+      return CPtrGetter == null ? default(T) : (T) System.Activator.CreateInstance
+      (
+          typeof(T),
+          System.Reflection.BindingFlags.NonPublic | System.Reflection.BindingFlags.Instance,
+          null,
+          new object[] { ((System.Runtime.InteropServices.HandleRef) CPtrGetter.Invoke(null, new object[] { from })).Handle, cMemoryOwn },
+          null
+      );
+  }
+%}
diff --git a/Source/Wrappers/csharp/SwigInterface/weak_ptr.i b/Source/Wrappers/csharp/SwigInterface/weak_ptr.i
new file mode 100644
index 0000000000000000000000000000000000000000..ff767bf2da9addf615b648756cd64a1b3dfd028b
--- /dev/null
+++ b/Source/Wrappers/csharp/SwigInterface/weak_ptr.i
@@ -0,0 +1,28 @@
+namespace std {
+template<class Ty> class weak_ptr {
+public:
+    typedef Ty element_type;
+
+    weak_ptr();
+    weak_ptr(const weak_ptr&);
+    template<class Other>
+        weak_ptr(const weak_ptr<Other>&);
+    template<class Other>
+        weak_ptr(const shared_ptr<Other>&);
+
+    weak_ptr(const shared_ptr<Ty>&);
+
+
+    void swap(weak_ptr&);
+    void reset();
+
+    long use_count() const;
+    bool expired() const;
+    shared_ptr<Ty> lock() const;
+};
+} /* end of namespace std */
+
+
+%template(SceneManagerWeakPtr) std::weak_ptr<imstk::SceneManager>;
+%template(ModuleDriverWeakPtr) std::weak_ptr<imstk::ModuleDriver>;
+%template(SimulationManagerWeakPtr) std::weak_ptr<imstk::SimulationManager>;
diff --git a/Source/Wrappers/csharp/Testing/SharedPtrTest.cs b/Source/Wrappers/csharp/Testing/SharedPtrTest.cs
new file mode 100644
index 0000000000000000000000000000000000000000..88985446796b01c35370dcbb1850b5bfee119e21
--- /dev/null
+++ b/Source/Wrappers/csharp/Testing/SharedPtrTest.cs
@@ -0,0 +1,187 @@
+using System;
+using imstk;
+using NUnit.Framework;
+using System.Runtime.InteropServices;
+
+namespace imstkCSUnitTest
+{
+
+[TestFixture]
+public class SharedPtrTest
+{
+    // [Test]
+    // public void TestActor()
+    // {
+    //     Assert.Multiple(() =>
+    //     {
+    //         HandleRef surfMeshHandleRef;
+    //         {
+    //             SurfaceMesh surfMesh = new SurfaceMesh();
+    //             surfMeshHandleRef = SurfaceMesh.getCPtr(surfMesh);
+    //             int use_count = UtilsPINVOKE.DebugAndTestOnly_use_count(SurfaceMesh.getCPtr(surfMesh), "SurfaceMesh");
+    //             Assert.GreaterOrEqual(use_count, 1);
+    //         }
+    //
+    //         GC.Collect();
+    //         GC.WaitForPendingFinalizers();
+    //
+    //         // {
+    //         //     DebugAndTestOnly u = new DebugAndTestOnly();
+    //         //     Assert.AreEqual(1, DebugAndTestOnly.getNumCreated());
+    //         //     Assert.AreEqual(0, DebugAndTestOnly.getNumDestroyed());
+    //         // }
+    //
+    //         DebugAndTestOnly debugObj = createAndDestroy<DebugAndTestOnly>();
+    //         debugObj = null;
+    //         GC.Collect();
+    //         GC.WaitForPendingFinalizers();
+    //         Assert.AreEqual(1, DebugAndTestOnly.getNumCreated());
+    //         // Assert.AreEqual(1, DebugAndTestOnly.getNumDestroyed());
+    //     });
+    //
+    // }
+
+    [Test]
+    public void TestCast()
+    {
+        PbdModel pbdModel = createPbdModel();
+        Assert.Multiple(() =>
+        {
+            PointSet pointSet = Utils.castToPointSet(pbdModel.getModelGeometry());
+            Assert.That(pointSet, Is.Not.Null);
+            Assert.AreEqual(4, pointSet.getNumVertices());
+
+            SurfaceMesh surfMesh = Utils.castToSurfaceMesh(pbdModel.getModelGeometry());
+            Assert.That(surfMesh, Is.Not.Null);
+            Assert.AreEqual(2, surfMesh.getNumTriangles());
+
+            TetrahedralMesh tetMesh = Utils.castToTetrahedralMesh(Utils.castToPointSet(pbdModel.getModelGeometry()));
+            Assert.That(tetMesh, Is.Null);
+        });
+    }
+
+    [Test]
+    public void TestVertices()
+    {
+        Assert.Multiple(() =>
+        {
+            PbdModel pbdModel = createPbdModel();
+            PointSet pointSet = Utils.castToPointSet(pbdModel.getModelGeometry());
+            Assert.AreEqual(4, pointSet.getNumVertices());
+
+            VecDataArray3d vertices = pointSet.getVertexPositions();
+            Assert.AreEqual(4, vertices.size());
+            Assert.AreEqual(0.0, vertices[0][0], tol);
+            Assert.AreEqual(0.0, vertices[0][1], tol);
+            Assert.AreEqual(0.0, vertices[0][2], tol);
+
+            Assert.AreEqual(1.0, vertices[1][0], tol);
+            Assert.AreEqual(0.0, vertices[1][1], tol);
+            Assert.AreEqual(0.0, vertices[1][2], tol);
+
+            Assert.AreEqual(0.0, vertices[2][0], tol);
+            Assert.AreEqual(1.0, vertices[2][1], tol);
+            Assert.AreEqual(0.0, vertices[2][2], tol);
+
+            Assert.AreEqual(1.0, vertices[3][0], tol);
+            Assert.AreEqual(1.0, vertices[3][1], tol);
+            Assert.AreEqual(0.0, vertices[3][2], tol);
+
+            // make changes to vertices
+            vertices[0] = new Vec3d(0.1, 0.2, 0.3);
+
+            PointSet pointSet2 = Utils.castToPointSet(pbdModel.getModelGeometry());
+            VecDataArray3d changedVertices = pointSet.getVertexPositions();
+            Assert.AreEqual(0.1, changedVertices[0][0], tol);
+            Assert.AreEqual(0.2, changedVertices[0][1], tol);
+            Assert.AreEqual(0.3, changedVertices[0][2], tol);
+
+            Assert.AreEqual(1.0, changedVertices[1][0], tol);
+            Assert.AreEqual(0.0, changedVertices[1][1], tol);
+            Assert.AreEqual(0.0, changedVertices[1][2], tol);
+
+            Assert.AreEqual(0.0, changedVertices[2][0], tol);
+            Assert.AreEqual(1.0, changedVertices[2][1], tol);
+            Assert.AreEqual(0.0, changedVertices[2][2], tol);
+
+            Assert.AreEqual(1.0, changedVertices[3][0], tol);
+            Assert.AreEqual(1.0, changedVertices[3][1], tol);
+            Assert.AreEqual(0.0, changedVertices[3][2], tol);
+        });
+    }
+        
+    [Test]
+    public void TestIndices()
+    {
+        Assert.Multiple(() =>
+        {
+            PbdModel pbdModel = createPbdModel();
+            SurfaceMesh surfMesh = Utils.castToSurfaceMesh(pbdModel.getModelGeometry());
+            Assert.That(surfMesh, Is.Not.Null);
+            Assert.AreEqual(2, surfMesh.getNumTriangles());
+
+            VecDataArray3i indices = surfMesh.getTriangleIndices();
+            Assert.AreEqual(0, indices[0][0]);
+            Assert.AreEqual(1, indices[0][1]);
+            Assert.AreEqual(2, indices[0][2]);
+
+            Assert.AreEqual(1, indices[1][0]);
+            Assert.AreEqual(3, indices[1][1]);
+            Assert.AreEqual(2, indices[1][2]);
+
+            // make changes to indices
+            indices[0] = new Vec3i(1, 2, 0);
+
+            SurfaceMesh surfMesh2 = Utils.castToSurfaceMesh(pbdModel.getModelGeometry());
+            Assert.That(surfMesh2, Is.Not.Null);
+            Assert.AreEqual(2, surfMesh2.getNumTriangles());
+
+            VecDataArray3i changedIndices = surfMesh2.getTriangleIndices();
+            Assert.AreEqual(2, changedIndices.size());
+            Assert.AreEqual(1, changedIndices[0][0]);
+            Assert.AreEqual(2, changedIndices[0][1]);
+            Assert.AreEqual(0, changedIndices[0][2]);
+
+            Assert.AreEqual(1, changedIndices[1][0]);
+            Assert.AreEqual(3, changedIndices[1][1]);
+            Assert.AreEqual(2, changedIndices[1][2]);
+        });
+    }
+
+    public static PbdModel createPbdModel()
+    {
+        PbdModel pbdModel = new PbdModel();
+        SurfaceMesh surfMesh = createSurfaceMesh();
+        pbdModel.setModelGeometry(surfMesh);
+        return pbdModel;
+    }
+
+    public static SurfaceMesh createSurfaceMesh()
+    {
+
+        VecDataArray3i conn = new VecDataArray3i(2);
+        conn[0] = new Vec3i(0, 1, 2);
+        conn[1] = new Vec3i(1, 3, 2);
+
+        VecDataArray3d coords = new VecDataArray3d(4);
+        coords[0] = new Vec3d(0.0, 0.0, 0.0);
+        coords[1] = new Vec3d(1.0, 0.0, 0.0);
+        coords[2] = new Vec3d(0.0, 1.0, 0.0);
+        coords[3] = new Vec3d(1.0, 1.0, 0.0);
+
+        SurfaceMesh surfMesh = new SurfaceMesh();
+        surfMesh.initialize(coords, conn);
+
+        return surfMesh;
+    }
+
+    private static T createAndDestroy<T>() where T : new()
+    {
+        T u = new T();
+        return u;
+    }
+
+
+    private static double tol = 1e-14;
+}
+}
diff --git a/Source/Wrappers/csharp/Testing/VecDataArrayTest.cs b/Source/Wrappers/csharp/Testing/VecDataArrayTest.cs
new file mode 100644
index 0000000000000000000000000000000000000000..febe480c8655aef99414d0a803563cb267caa7bc
--- /dev/null
+++ b/Source/Wrappers/csharp/Testing/VecDataArrayTest.cs
@@ -0,0 +1,147 @@
+using System;
+using imstk;
+using NUnit.Framework;
+
+
+namespace imstkCSUnitTest
+{
+[TestFixture]
+public class VecDataArrayTest
+{
+    [Test]
+    public void TestConstructor()
+    {
+        // default constructors
+        {
+            VecDataArray3d u = new VecDataArray3d();
+            Assert.AreEqual(0, u.size());
+            Assert.AreEqual(3, u.getNumberOfComponents());
+            Assert.GreaterOrEqual(u.getCapacity(), 3*0);
+            Assert.AreEqual(Utils.IMSTK_DOUBLE, u.getScalarType());
+        }
+
+        {
+            VecDataArray3d u = new VecDataArray3d(4);
+            Assert.AreEqual(4, u.size());
+            Assert.AreEqual(3, u.getNumberOfComponents());
+            Assert.GreaterOrEqual(u.getCapacity(), 3*4);
+        }
+
+    }
+
+    [Test]
+    public void TestResize()
+    {
+        Assert.Multiple(() =>
+        {
+            VecDataArray3d u = new VecDataArray3d(4);
+            u.resize(5);
+            Assert.AreEqual(5, u.size());
+            Assert.AreEqual(3, u.getNumberOfComponents());
+            Assert.GreaterOrEqual(u.getCapacity(), 3*u.size());
+        
+        });
+    }
+
+    [Test]
+    public void TestFill()
+    {
+        Assert.Multiple(() =>
+        {
+            VecDataArray3d u = new VecDataArray3d(2);
+            Vec3d val = new Vec3d(0.0, 1.0, 2.0);
+            u.fill(val);
+            Assert.AreEqual(u[0][0], 0.0, tol);
+            Assert.AreEqual(u[0][1], 1.0, tol);
+            Assert.AreEqual(u[0][2], 2.0, tol);
+
+            Assert.AreEqual(u[1][0], 0.0, tol);
+            Assert.AreEqual(u[1][1], 1.0, tol);
+            Assert.AreEqual(u[1][2], 2.0, tol);
+        });
+    }
+
+    [Test]
+    public void TestEraseAndPushBack()
+    {
+        Assert.Multiple(() =>
+        {
+            VecDataArray3d u = new VecDataArray3d(4);
+            u[0] = new Vec3d(0.0, 1.0, 2.0);
+            u[1] = new Vec3d(3.0, 4.0, 5.0);
+            u[2] = new Vec3d(6.0, 7.0, 8.0);
+            u[3] = new Vec3d(9.0, 10.0, 11.0);
+
+            u.erase(1);
+            Assert.AreEqual(3, u.size());
+            Assert.GreaterOrEqual(u.getCapacity(), u.size()*3);
+
+            Assert.AreEqual(0.0, u[0][0], tol);
+            Assert.AreEqual(1.0, u[0][1], tol);
+            Assert.AreEqual(2.0, u[0][2], tol);
+            Assert.AreEqual(6.0, u[1][0], tol);
+            Assert.AreEqual(7.0, u[1][1], tol);
+            Assert.AreEqual(8.0, u[1][2], tol);
+            Assert.AreEqual(9.0, u[2][0], tol);
+            Assert.AreEqual(10.0, u[2][1], tol);
+            Assert.AreEqual(11.0, u[2][2], tol);
+
+            u.push_back(new Vec3d(12.0, 13.0, 14.0));
+            Assert.AreEqual(4, u.size());
+            Assert.GreaterOrEqual(u.getCapacity(), u.size()*3);
+            Assert.AreEqual(0.0, u[0][0], tol);
+            Assert.AreEqual(1.0, u[0][1], tol);
+            Assert.AreEqual(2.0, u[0][2], tol);
+
+            Assert.AreEqual(6.0, u[1][0], tol);
+            Assert.AreEqual(7.0, u[1][1], tol);
+            Assert.AreEqual(8.0, u[1][2], tol);
+
+            Assert.AreEqual(9.0, u[2][0], tol);
+            Assert.AreEqual(10.0, u[2][1], tol);
+            Assert.AreEqual(11.0, u[2][2], tol);
+
+            Assert.AreEqual(12.0, u[3][0], tol);
+            Assert.AreEqual(13.0, u[3][1], tol);
+            Assert.AreEqual(14.0, u[3][2], tol);
+        });
+    }
+
+    [Test]
+    public void TestSetValuesAndGetValues()
+    {
+        Assert.Multiple(() =>
+        {
+            const int sz = 99;
+            const int numData = (int)sz * 3;
+            const double eps = 0.01;
+            double[] values = new double[numData];
+            for (uint i=0; i<values.Length; ++i)
+            {
+                values[i] = (double)i + eps;
+            }
+
+            VecDataArray3d u = new VecDataArray3d(sz);
+            u.setValues(values);
+
+            for (uint i=0; i<sz; ++i)
+            {
+                Assert.AreEqual(u[i][0], i*3+0 + eps, tol);
+                Assert.AreEqual(u[i][1], i*3+1 + eps, tol);
+                Assert.AreEqual(u[i][2], i*3+2 + eps, tol);
+            }
+
+            double[] toValues = new double[numData];
+
+            u.getValues(toValues);
+
+            for (int i=0; i<numData; ++i)
+            {
+                Assert.AreEqual(toValues[i], values[i], tol);
+            }
+        });
+    }
+
+    private static double tol = 1e-14;
+}
+}
diff --git a/Source/Wrappers/csharp/Testing/VecTest.cs b/Source/Wrappers/csharp/Testing/VecTest.cs
new file mode 100644
index 0000000000000000000000000000000000000000..309bafcb58b51a72eabc467895a1b04fdebd2131
--- /dev/null
+++ b/Source/Wrappers/csharp/Testing/VecTest.cs
@@ -0,0 +1,228 @@
+
+using System;
+using imstk;
+using NUnit.Framework;
+using System.Runtime.InteropServices;
+
+namespace imstkCSUnitTest
+{
+
+[TestFixture]
+public class Vec3dTest
+{
+    [Test]
+    public void TestConstructorAndIndexing()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec3d u = new Vec3d(0.1, 0.2);
+            Assert.AreEqual(0.1, u[0], tol);
+            Assert.AreEqual(0.2, u[1], tol);
+
+            u[2] = 100.0;
+            Assert.AreEqual(100.0, u[2], tol);
+
+            Vec3d v = new Vec3d(0.1, 0.2, 0.3);
+            Assert.AreEqual(0.1, v[0], tol);
+            Assert.AreEqual(0.2, v[1], tol);
+            Assert.AreEqual(0.3, v[2], tol);
+            v[2] = 100.0 + 1e-12;
+            Assert.AreEqual(100.0+1e-12, v[2], tol);
+        });
+    }
+
+    [Test]
+    public void TestOperators()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec3d v = new Vec3d(0.1, 0.2, 0.3);
+            // *
+            const double factor = 2.0 + 1e-6;
+            v = v*factor;
+            Assert.AreEqual(0.1*factor, v[0], tol);
+            Assert.AreEqual(0.2*factor, v[1], tol);
+            Assert.AreEqual(0.3*factor, v[2], tol);
+
+            v = factor*v;
+            Assert.AreEqual(0.1*factor*factor, v[0], tol);
+            Assert.AreEqual(0.2*factor*factor, v[1], tol);
+            Assert.AreEqual(0.3*factor*factor, v[2], tol);
+
+            // /
+            v = v / factor;
+            Assert.AreEqual(0.1*factor, v[0], tol);
+            Assert.AreEqual(0.2*factor, v[1], tol);
+            Assert.AreEqual(0.3*factor, v[2], tol);
+
+            // +
+            Vec3d u = new Vec3d(0.1, 0.2, 0.3);
+            v = u + new Vec3d(0.9, 0.8, 0.7);
+            Assert.AreEqual(1.0, v[0], tol);
+            Assert.AreEqual(1.0, v[1], tol);
+            Assert.AreEqual(1.0, v[2], tol);
+
+            // -
+            v = u - new Vec3d(0.1-1e-6, 0.2-1e-7, 0.3-1e-8);
+            Assert.AreEqual(1e-6, v[0], tol);
+            Assert.AreEqual(1e-7, v[1], tol);
+            Assert.AreEqual(1e-8, v[2], tol);
+
+            // normalized
+            u = new Vec3d(3.0, 4.0, 0.0);
+            Vec3d w = u.normalized();
+            Assert.AreEqual(0.6, w[0], tol);
+            Assert.AreEqual(0.8, w[1], tol);
+            Assert.AreEqual(0.0, w[2], tol);
+
+        });
+    }
+
+    private static double tol = 1e-14;
+}
+
+[TestFixture]
+public class Vec2fTest
+{
+    [Test]
+    public void TestConstructorAndIndexing()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec2f u = new Vec2f(0.1f, 0.2f);
+            Assert.AreEqual(0.1f, u[0], tol);
+            Assert.AreEqual(0.2f, u[1], tol);
+
+            u[1] = 100.0f + 1e-4f;
+            Assert.AreEqual(100.0f+1e-4f, u[1], tol);
+        });
+    }
+
+    [Test]
+    public void TestOperators()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec2f v = new Vec2f(0.1f, 0.2f);
+            // *
+            const float factor = 2.0f + 1e-4f;
+            v = factor*v;
+            Assert.AreEqual(0.1f*factor, v[0], tol);
+            Assert.AreEqual(0.2f*factor, v[1], tol);
+
+            // /
+            v = v / factor;
+            Assert.AreEqual(0.1f, v[0], tol);
+            Assert.AreEqual(0.2f, v[1], tol);
+        });
+    }
+
+    private static double tol = 1e-6;
+}
+
+[TestFixture]
+public class Vec3iTest
+{
+    [Test]
+    public void TestConstructorAndIndexing()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec3i u = new Vec3i(1, 2);
+            Assert.AreEqual(1, u[0]);
+            Assert.AreEqual(2, u[1]);
+
+            u[2] = 100;
+            Assert.AreEqual(100, u[2]);
+        });
+    }
+}
+
+[TestFixture]
+public class Vec4iTest
+{
+    [Test]
+    public void TestConstructorAndIndexing()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec4i u = new Vec4i(1, 2, 3, 4);
+            Assert.AreEqual(1, u[0]);
+            Assert.AreEqual(2, u[1]);
+            Assert.AreEqual(3, u[2]);
+            Assert.AreEqual(4, u[3]);
+
+            u[1] = 100;
+            Assert.AreEqual(100, u[1]);
+        });
+    }
+}
+
+
+[TestFixture]
+public class Vec6dTest
+{
+    [Test]
+    public void TestConstructorAndIndexing()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec6d u = new Vec6d(0.1, 0.2, 0.3, 0.4, 0.5, 0.6);
+            Assert.AreEqual(0.1, u[0], tol);
+            Assert.AreEqual(0.2, u[1], tol);
+            Assert.AreEqual(0.3, u[2], tol);
+            Assert.AreEqual(0.4, u[3], tol);
+            Assert.AreEqual(0.5, u[4], tol);
+            Assert.AreEqual(0.6, u[5], tol);
+
+            u[2] = 100.0;
+            Assert.AreEqual(100.0, u[2], tol);
+
+            Vec6d v = new Vec6d(0.1, 0.2, 0.3);
+            Assert.AreEqual(0.1, v[0], tol);
+            Assert.AreEqual(0.2, v[1], tol);
+            Assert.AreEqual(0.3, v[2], tol);
+            v[5] = 100.0 + 1e-12;
+            Assert.AreEqual(100.0+1e-12, v[5], tol);
+        });
+    }
+
+    [Test]
+    public void TestOperators()
+    {
+        Assert.Multiple(() =>
+        {
+            Vec6d v = new Vec6d(0.1, 0.2, 0.3, 0.4, 0.5, 0.6);
+            // *
+            const double factor = 2.0 + 1e-6;
+            v = v*factor;
+            Assert.AreEqual(0.1*factor, v[0], tol);
+            Assert.AreEqual(0.2*factor, v[1], tol);
+            Assert.AreEqual(0.3*factor, v[2], tol);
+            Assert.AreEqual(0.4*factor, v[3], tol);
+            Assert.AreEqual(0.5*factor, v[4], tol);
+            Assert.AreEqual(0.6*factor, v[5], tol);
+
+            v = factor*v;
+            Assert.AreEqual(0.1*factor*factor, v[0], tol);
+            Assert.AreEqual(0.2*factor*factor, v[1], tol);
+            Assert.AreEqual(0.3*factor*factor, v[2], tol);
+            Assert.AreEqual(0.4*factor*factor, v[3], tol);
+            Assert.AreEqual(0.5*factor*factor, v[4], tol);
+            Assert.AreEqual(0.6*factor*factor, v[5], tol);
+
+            // /
+            v = v / factor;
+            Assert.AreEqual(0.1*factor, v[0], tol);
+            Assert.AreEqual(0.2*factor, v[1], tol);
+            Assert.AreEqual(0.3*factor, v[2], tol);
+            Assert.AreEqual(0.4*factor, v[3], tol);
+            Assert.AreEqual(0.5*factor, v[4], tol);
+            Assert.AreEqual(0.6*factor, v[5], tol);
+        });
+    }
+
+    private static double tol = 1e-14;
+}
+}
+
diff --git a/Source/Wrappers/csharp/Testing/run_test.sh b/Source/Wrappers/csharp/Testing/run_test.sh
new file mode 100644
index 0000000000000000000000000000000000000000..3f25bc418dc537aedeb293e5ebf16316ad5da3cd
--- /dev/null
+++ b/Source/Wrappers/csharp/Testing/run_test.sh
@@ -0,0 +1,11 @@
+# How to run NUnit on Ubuntu: https://stackoverflow.com/questions/31038629/run-nunit-test-on-ubuntu-from-command-line.
+
+export MONO_PATH=$(pwd)/NUnit.ConsoleRunner.3.12.0/tools:$(pwd)/NUnit.3.13.2/lib/net40
+export LD_LIBRARY_PATH=/home/jianfeng/Documents/imstk/build_csharp/install/lib:$LD_LIBRARY_PATH
+
+testName=$1
+rm -f ${testName}.dll
+iMSTKSharp=/home/jianfeng/Documents/imstk/build_csharp/install/include/iMSTKSharp
+mcs ${testName}.cs ${iMSTKSharp}/*.cs -unsafe -debug -target:library -r:NUnit.3.13.2/lib/net40/nunit.framework.dll -out:${testName}.dll
+mono ./NUnit.ConsoleRunner.3.12.0/tools/nunit3-console.exe ${testName}.dll -noresult
+
diff --git a/Source/Wrappers/csharp/iMSTKCSharp/CMakeLists.txt b/Source/Wrappers/csharp/iMSTKCSharp/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..7693c7ab62e762b3e278d31fb9380b153e564f4b
--- /dev/null
+++ b/Source/Wrappers/csharp/iMSTKCSharp/CMakeLists.txt
@@ -0,0 +1,17 @@
+project(iMSTKCSharp VERSION 0.1.0 LANGUAGES CSharp)
+
+# Check https://stackoverflow.com/questions/52556785/code-generator-generating-its-own-cmake-files-and-targets/52714922#52714922
+
+file(GLOB csharp_files "${CMAKE_INSTALL_PREFIX}/include/iMSTKSharp/*.cs")
+
+if (csharp_files)
+	set(IMSTK_WRAPPER_BUILT ON CACHE INTERNAL "")
+	add_library(${PROJECT_NAME} SHARED ${csharp_files})
+	set_property(TARGET ${PROJECT_NAME} PROPERTY VS_DOTNET_TARGET_FRAMEWORK_VERSION "v4.6.1")
+	set_property(TARGET ${PROJECT_NAME} PROPERTY VS_DOTNET_REFERENCES "System")
+	add_custom_command(TARGET ${PROJECT_NAME} POST_BUILD
+                   COMMAND ${CMAKE_COMMAND} -E copy $<TARGET_FILE:${PROJECT_NAME}> ${CMAKE_INSTALL_PREFIX}/bin)
+else ()
+	set(IMSTK_WRAPPER_BUILT OFF CACHE INTERNAL "")
+	message(WARNING "C# files not found, rerun configure after building all of immstk or just the swig iMSTKCWrapper")
+endif()
\ No newline at end of file