BUG: Failure to created eigen vectors

Clears runtime error that shows failure to create eigen vectors in problemState.
Small formatting changws.

Base/ForceModel/imstkForceModelConfig.cpp

@@ -28,7 +28,7 @@ ForceModelConfig::ForceModelConfig(const std::string &configFileName) : m_loadSu
 {
     if (configFileName.empty())
     {
-        LOG(INFO) << "WARNING: Empty configuration filename." << std::endl;
+        LOG(INFO) << "WARNING: Empty configuration filename.";
         return;
     }
     else
@@ -56,9 +56,13 @@ ForceModelConfig::parseConfig(const std::string &configFileName)
     vegaConfigFileOptions.addOptionOptional(optNameList.inversionThresholdName.c_str(), &optList.inversionThreshold, optList.inversionThreshold);
     vegaConfigFileOptions.addOptionOptional(optNameList.numberOfThreadsName.c_str(), &optList.numberOfThreads, optList.numberOfThreads);
 
-
     // Parse the configuration file
     if (vegaConfigFileOptions.parseOptions(configFileName.data()) != 0)
+    {
+        LOG(WARNING) << "ForceModelConfig::parseConfig - Unable to load the configuration file";
+        return false;
+    }
+    else
     {
         this->m_vegaConfigFileName = configFileName;
         m_loadSuccessful = true;
@@ -66,10 +70,6 @@ ForceModelConfig::parseConfig(const std::string &configFileName)
         // Print option variables
         vegaConfigFileOptions.printOptions();
     }
-    else
-    {
-        return false;
-    }
 
     // Store parsed string values
     this->m_stringsOptionMap.emplace(optNameList.femMethodName, optList.femMethod);
@@ -107,7 +107,7 @@ ForceModelConfig::getForceModelType()
     {
         return ForceModelType::Linear;
     }
-    else if (this->m_stringsOptionMap["femMethod"] == "Invertible")
+    else if (this->m_stringsOptionMap["femMethod"] == "InvertibleFEM")
     {
         return ForceModelType::Invertible;
     }
@@ -141,4 +141,14 @@ ForceModelConfig::getHyperelasticMaterialType()
     }
 }
 
+void
+ForceModelConfig::print()
+{
+    LOG(INFO) << "Floating point type options:\n";
+    for (auto const &fo : m_floatsOptionMap)
+    {
+        LOG(INFO) << fo.first << ": " << fo.second;
+    }
+}
+
 } // imstk

Base/ForceModel/imstkForceModelConfig.h

@@ -138,7 +138,7 @@ public:
     ~ForceModelConfig(){};
 
     ///
-    /// \brief Parse the deformable objeyc config file
+    /// \brief Parse the deformable object config file
     ///
     bool parseConfig(const std::string &configFileName);
 
@@ -176,6 +176,11 @@ public:
     ///
     HyperElasticMaterialType getHyperelasticMaterialType();
 
+    ///
+    /// \brief Print the parsed options
+    ///
+    void print();
+
 protected:
     std::string m_vegaConfigFileName;                      ///> Store configuration file.
     std::map<std::string, double> m_floatsOptionMap;       ///> Map for floating point

Base/SceneElements/Controllers/imstkTrackingController.h

@@ -27,7 +27,8 @@
 
 #include <memory>
 
-namespace imstk {
+namespace imstk
+{
 
 ///
 /// \class TrackingController
@@ -103,6 +104,7 @@ protected:
     Quatd m_rotationOffset = Quatd::Identity();   //!< Rotation concatenated to the device rotation
     unsigned char m_invertFlags = 0x00;           //!< Invert flags to be masked with TrackingController::InvertFlag
 };
+
 }
 
 #endif // ifndef imstkTrackingController_h

Base/SceneElements/Light/imstkLight.h

@@ -95,6 +95,7 @@ protected:
     std::string m_name;
 
 };
+
 }
 
 #endif // ifndef imstkLight_h

Base/SceneElements/Objects/imstkCollidingObject.h

@@ -34,17 +34,28 @@ class GeometryMap;
 class CollidingObject : public SceneObject
 {
 public:
-
+    ///
+    /// \brief
+    ///
     CollidingObject(std::string name) : SceneObject(name)
     {
         m_type = SceneObject::Type::Static;
     }
 
+    ///
+    /// \brief
+    ///
     ~CollidingObject() = default;
 
+    ///
+    /// \brief
+    ///
     std::shared_ptr<Geometry> getCollidingGeometry() const;
     void setCollidingGeometry(std::shared_ptr<Geometry> geometry);
 
+    ///
+    /// \brief
+    ///
     std::shared_ptr<GeometryMap> getCollidingToVisualMap() const;
     void setCollidingToVisualMap(std::shared_ptr<GeometryMap> map);
 

Base/SceneElements/Objects/imstkDeformableBodyModel.cpp

@@ -88,8 +88,10 @@ DeformableBodyModel::configure(const std::string& configFileName)
 }
 
 void
-DeformableBodyModel::initialize()
+DeformableBodyModel::initialize(std::shared_ptr<Geometry> geometry)
 {
+    setModelGeometry(geometry);
+
     // prerequisite of for successfully initializing
     if (!m_forceModelGeometry || !m_forceModelConfiguration)
     {
@@ -97,7 +99,7 @@ DeformableBodyModel::initialize()
         return;
     }
 
-    m_vegaPhysicsMesh = VegaMeshReader::readVegaMesh(m_forceModelConfiguration->getStringOptionsMap().at("asianDragon/asianDragon.veg"));// remove hard coded name after test
+    m_vegaPhysicsMesh = VegaMeshReader::readVegaMesh("asianDragon/asianDragon.veg");// remove hard coded name after test
 
     initializeForceModel();
     initializeMassMatrix();
@@ -116,10 +118,10 @@ DeformableBodyModel::initialize()
 void
 DeformableBodyModel::loadInitialStates()
 {
-    // Initialize the states
-    m_initialState->initialize(m_numDOF);
-    m_previousState->initialize(m_numDOF);
-    m_currentState->initialize(m_numDOF);
+    // For now the initial states are set to zero
+    m_initialState = std::make_shared<kinematicState>(m_numDOF);
+    m_initialState = std::make_shared<kinematicState>(m_numDOF);
+    m_initialState = std::make_shared<kinematicState>(m_numDOF);
 }
 
 void
@@ -162,6 +164,8 @@ DeformableBodyModel::initializeForceModel()
     const float g = m_forceModelConfiguration->getFloatsOptionsMap().at("gravity");
     const bool isGravityPresent = (g > 0) ? true : false;
 
+    m_numDOF = m_vegaPhysicsMesh->getNumVertices() * 3;
+
     switch (m_forceModelConfiguration->getForceModelType())
     {
     case ForceModelType::StVK:

Base/SceneElements/Objects/imstkDeformableBodyModel.h

@@ -105,7 +105,7 @@ public:
     ///
     /// \brief Initialize the deformable body model
     ///
-    void initialize();
+    void initialize(std::shared_ptr<Geometry> geometry);
 
     ///
     /// \brief Load the initial conditions of the deformable object
@@ -145,14 +145,13 @@ public:
     ///
     /// \brief Compute the RHS of the resulting linear system
     ///
-    void computeImplicitSystemRHS(const kinematicState& prevState,
-                                  kinematicState& newState);
+    void computeImplicitSystemRHS(const kinematicState& prevState, kinematicState& newState);
 
     ///
     /// \brief Compute the LHS of the resulting linear system
     ///
-    void computeImplicitSystemLHS(const kinematicState& prevState,
-                                  kinematicState& newState);
+    void computeImplicitSystemLHS(const kinematicState& prevState, kinematicState& newState);
+
     ///
     /// \brief Initialize explicit external forces
     ///

Base/SceneElements/Objects/imstkProblemState.h

@@ -41,19 +41,27 @@ public:
     ///
     /// \brief Constructor
     ///
-    ProblemState() = default;
-    ProblemState(const Vectord& u, const Vectord& v, const Vectord& a);
+    ProblemState(){};
+    ProblemState(const size_t size){ initialize(size); };
+    //ProblemState(const Vectord& u, const Vectord& v, const Vectord& a);
 
     ///
     /// \brief Destructor
     ///
     ~ProblemState(){}
 
+    ///
+    /// \brief Initialize the problem state
+    ///
     void initialize(const size_t numDof)
     {
         m_q.resize(numDof);
         m_qDot.resize(numDof);
         m_qDotDot.resize(numDof);
+
+        m_q.setZero();
+        m_qDot.setZero();
+        m_qDotDot.setZero();
     };
 
     ///
@@ -66,6 +74,9 @@ public:
         m_qDotDot = a;
     }
 
+    ///
+    /// \brief Set the state to a given one
+    ///
     void setState(std::shared_ptr<ProblemState> rhs)
     {
         m_q = rhs->getQ();

Base/SceneElements/Objects/imstkSceneObject.h

@@ -28,7 +28,8 @@
 #include <memory>
 #include <string>
 
-namespace imstk {
+namespace imstk
+{
 
 //class Geometry;
 //class GeometryMap;

Base/SceneElements/Objects/imstkVirtualCouplingObject.h

@@ -33,7 +33,9 @@ class GeometryMap;
 class VirtualCouplingObject : public CollidingObject, public TrackingController
 {
 public:
-
+    ///
+    /// \brief Constructor
+    ///
     VirtualCouplingObject(std::string name,
                           std::shared_ptr<DeviceClient> deviceClient = nullptr,
                           double scaling = 1.0) :
@@ -43,10 +45,13 @@ public:
         m_type = Type::VirtualCoupling;
     }
 
+    ///
+    /// \brief Destructor
+    ///
     ~VirtualCouplingObject() = default;
 
     ///
-    /// \brief Initialise offset based on object geometry
+    /// \brief Initialize offset based on object geometry
     ///
     void initOffsets();
 

Examples/Sandbox/main.cpp

@@ -369,15 +369,15 @@ void testReadMesh()
     auto scene = sdk->createNewScene("SceneTestMesh");
 
     // Read surface mesh
-    auto objMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/asianDragon/asianDragon.obj");
+    /*auto objMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/asianDragon/asianDragon.obj");
     auto plyMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/Cube/models/cube.ply");
     auto stlMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/Cube/models/cube.stl");
     auto vtkMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/Cube/models/cube.vtk");
-    auto vtpMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/Cube/models/cube.vtp");
+    auto vtpMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/Cube/models/cube.vtp");*/
 
     // Read volumetricMesh
-    auto vtkMesh2 = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/AVM/nidus-model/nidus10KTet.vtk");
-    auto vegaMesh = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/asianDragon/asianDragon.veg");
+    //auto vtkMesh2 = imstk::MeshReader::read("/home/virtualfls/Projects/IMSTK/resources/AVM/nidus-model/nidus10KTet.vtk");
+    auto vegaMesh = imstk::MeshReader::read("asianDragon.veg");
 
     // Extract surface mesh
     auto volumeMesh = std::dynamic_pointer_cast<imstk::VolumetricMesh>(vegaMesh); // change to any volumetric mesh above
@@ -801,7 +801,7 @@ void testDeformableBody()
     auto scene = sdk->createNewScene("DeformableBodyTest");
 
     // b. Load a tetrahedral mesh
-    auto tetMesh = imstk::MeshReader::read("asianDragon.veg");
+    auto tetMesh = imstk::MeshReader::read("asianDragon/asianDragon.veg");
 
     // c. Extract the surface mesh
     auto surfMesh = std::make_shared<imstk::SurfaceMesh>();
@@ -826,9 +826,8 @@ void testDeformableBody()
 
     // Configure dynamic model
     auto dynaModel = std::make_shared<DeformableBodyModel>(DynamicalModel::Type::elastoDynamics);
-    dynaModel->configure("dragon.config");
-    dynaModel->setModelGeometry(volTetMesh);
-    dynaModel->initialize();
+    dynaModel->configure("asianDragon.config");
+    dynaModel->initialize(volTetMesh);
     auto timeIntegrator = std::make_shared<BackwardEuler>(0.01);// Create and add Backward Euler time integrator
     dynaModel->setTimeIntegrator(timeIntegrator);
 
@@ -856,17 +855,17 @@ void testDeformableBody()
 
 
     // create a nonlinear system
-    auto nlSystem = std::make_shared<NonLinearSystem>();
-    nlSystem->setFunction(dynaModel->getFunction());
-    nlSystem->setJacobian(dynaModel->getFunctionGradient());
+    //auto nlSystem = std::make_shared<NonLinearSystem>();
+    //nlSystem->setFunction(dynaModel->getFunction());
+    //nlSystem->setJacobian(dynaModel->getFunctionGradient());
 
-    // create a linear solver
-    auto cgLinSolver = std::make_shared<ConjugateGradient>();
+    //// create a linear solver
+    //auto cgLinSolver = std::make_shared<ConjugateGradient>();
 
-    // create a non-linear solver
-    auto nlSolver = std::make_shared<NewtonMethod>();
-    nlSolver->setLinearSolver(cgLinSolver);
-    nlSolver->setSystem(nlSystem);
+    //// create a non-linear solver
+    //auto nlSolver = std::make_shared<NewtonMethod>();
+    //nlSolver->setLinearSolver(cgLinSolver);
+    //nlSolver->setSystem(nlSystem);
 
     // Run the simulation
     sdk->setCurrentScene("DeformableBodyTest");