add element to tetra mesh sim problem

i forgot to mention that i am using sofa as a library in c++ code.

the part of scn file that i posted is written out for the frame from the c++ code.

i’m not sure what version number it is. this was before sofa was moved to github and before some big fork that was mentioned.

just to make it cleaner i am using:

MechanicalObject
TetrahedronSetTopologyContainer
EulerImplicitSolver
CGLinearSolver
UniformMass
TetrahedralCorotationalFEMForceField
TetrahedronSetGeometryAlgorithms
TetrahedronSetTopologyAlgorithms
TetrahedronSetTopologyModifier
TetrahedronModel
FixedConstraint template=”Vec3d” indices=”20 11 12 13″
UniformMass template=”Vec3d” name=”mass2″ mass=”10000″ totalmass=”210000″ localRange=”11 13″
UniformMass template=”Vec3d” name=”mass3″ mass=”10000″ totalmass=”210000″ localRange=”20 20″

the component responsible for topological change should be TetrahedronSetTopologyModifier

i have 3 uniformmass just so i can set the nodes at the end to be heavy so the bar deforms.

it seems different without tet2tri node. it seems tetras are stretching. maybe it doesnt have the stiffness applied to them?

https://www.youtube.com/watch?v=N87LdmiUD9Y

** just a side note, the forum keeps preventing me from posting. it forbids me then i have to go to cleantalk and delete my IP, but it still doesnt work. do you know what is wrong with the forum? **

that scn file i posted was only to show what components i was using and i didn’t want to make a large post with so many points and tetra information.

i am using sofa as a c++ library.

so i created a small demo that adds a tetra to a tetra mesh. i modified sofa/applications/tutorials/oneTetrahedron.cpp

this creates a *.obj file for each frame. at frame 2 it creates a new tetra. you can see the new point of the tetra not behaving as expected.

you need to change YOURPATH to a path on your system.

this was compiled on ubuntu 16.04 LTS (virtual machine in this case) using sofa 16.12

/******************************************************************************
*       SOFA, Simulation Open-Framework Architecture, development version     *
*                (c) 2006-2017 INRIA, USTL, UJF, CNRS, MGH                    *
*                                                                             *
* This program is free software; you can redistribute it and/or modify it     *
* under the terms of the GNU General Public License as published by the Free  *
* Software Foundation; either version 2 of the License, or (at your option)   *
* any later version.                                                          *
*                                                                             *
* This program is distributed in the hope that it will be useful, but WITHOUT *
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or       *
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for    *
* more details.                                                               *
*                                                                             *
* You should have received a copy of the GNU General Public License along     *
* with this program. If not, see <http://www.gnu.org/licenses/>.              *
*******************************************************************************
* Authors: The SOFA Team and external contributors (see Authors.txt)          *
*                                                                             *
* Contact information: contact@sofa-framework.org                             *
******************************************************************************/

#include <SofaComponentCommon/initComponentCommon.h>
#include <SofaComponentBase/initComponentBase.h>
#include <SofaComponentGeneral/initComponentGeneral.h>
#include <SofaComponentAdvanced/initComponentAdvanced.h>
#include <SofaComponentMisc/initComponentMisc.h>

#include <SofaGraphComponent/Gravity.h>
#include <SofaImplicitOdeSolver/EulerImplicitSolver.h>
#include <SofaBaseLinearSolver/CGLinearSolver.h>
#include <sofa/component/typedef/Sofa_typedef.h>
//#include <SofaOpenglVisual/OglModel.h>
#include <SofaBaseVisual/VisualStyle.h>

#include <sofa/core/objectmodel/Context.h>
#include <sofa/core/VecId.h>

#include <sofa/gui/GUIManager.h>
#include <sofa/gui/Main.h>

#include <sofa/helper/ArgumentParser.h>
#include <sofa/helper/system/FileRepository.h>
#include <sofa/helper/system/glut.h>

#include <SofaSimulationTree/init.h>
#include <SofaSimulationTree/GNode.h>
#include <SofaSimulationTree/TreeSimulation.h>

#include <SofaBaseTopology/TetrahedronSetTopologyContainer.h>
#include <SofaBaseTopology/TetrahedronSetTopologyModifier.h>

#include <iostream>
#include <fstream>

using namespace sofa::simulation::tree;
using namespace sofa;
using sofa::simulation::Node;
using sofa::component::odesolver::EulerImplicitSolver;
typedef component::linearsolver::CGLinearSolver<component::linearsolver::GraphScatteredMatrix, component::linearsolver::GraphScatteredVector> CGLinearSolver;
using sofa::component::topology::MeshTopology;
//using sofa::component::visualmodel::OglModel;
using sofa::core::objectmodel::Data;
using sofa::helper::ReadAccessor;
using sofa::helper::WriteAccessor;
using sofa::core::VecId;
using sofa::core::objectmodel::New;

int main(int argc, char** argv)
{

	sofa::simulation::setSimulation(new sofa::simulation::tree::TreeSimulation());
    sofa::simulation::Node::SPtr root = sofa::simulation::getSimulation()->createNewGraph("root");
    root->setGravity( Coord3(0,-10,0) );

    simulation::Node::SPtr nodeTetra = root->createChild("tetramesh");

    MechanicalObject3::SPtr DOF = sofa::core::objectmodel::New<MechanicalObject3>();
    nodeTetra->addObject(DOF);
    DOF->resize(21);
    DOF->setName("DOF");

    sofa::component::topology::TetrahedronSetTopologyContainer::SPtr topology = sofa::core::objectmodel::New<sofa::component::topology::TetrahedronSetTopologyContainer>();

    WriteAccessor<Data<VecCoord3> > xx = *DOF->write(VecId::position());

    xx[0] = Coord3(0.000000,  0.000000,  0.000000);
    xx[1] = Coord3(0.200000,  0.000000,  0.000000);
    xx[2] = Coord3(0.000000,  0.200000,  0.000000);
    xx[3] = Coord3(0.200000,  0.200000,  0.000000);
    xx[4] = Coord3(0.000000,  0.000000,  1.200000);

    xx[5] = Coord3(0.200000,  0.000000,  1.200000);
    xx[6] = Coord3(0.000000,  0.200000,  1.200000);
    xx[7] = Coord3(0.200000,  0.200000,  1.200000);
    xx[8] = Coord3(0.200000,  0.200000,  0.965324);
    xx[9] = Coord3(0.200000,  0.200000,  0.595933);

    xx[10] = Coord3(0.200000,  0.200000,  0.246792);
    xx[11] = Coord3(0.000000,  0.200000,  0.958366);
    xx[12] = Coord3(0.000000,  0.200000,  0.611182);
    xx[13] = Coord3(0.000000,  0.200000,  0.247715);
    xx[14] = Coord3(0.200000,  0.000000,  0.958366);

    xx[15] = Coord3(0.200000,  0.000000,  0.611182);
    xx[16] = Coord3(0.200000,  0.000000,  0.247715);
    xx[17] = Coord3(0.000000,  0.000000,  1.006005);
    xx[18] = Coord3(0.000000,  0.000000,  0.746736);
    xx[19] = Coord3(0.000000,  0.000000,  0.444605);

    xx[20] = Coord3(0.000000,  0.000000,  0.206183);

    topology->addTetra(17,8,18,14);
    topology->addTetra(12,8,15,18);
    topology->addTetra(18,11,12,8);
    topology->addTetra(20,10,13,2);
    topology->addTetra(6,5,17,4);
    topology->addTetra(5,8,6,17);
    topology->addTetra(20,3,1,10);
    topology->addTetra(12,9,19,15);
    topology->addTetra(15,8,12,9);
    topology->addTetra(18,12,19,15);

    topology->addTetra(9,16,19,15);
    topology->addTetra(18,8,17,11);
    topology->addTetra(20,16,19,10);
    topology->addTetra(20,3,2,0);
    topology->addTetra(16,20,1,10);
    topology->addTetra(18,14,8,15);
    topology->addTetra(19,9,13,10);
    topology->addTetra(14,5,17,8);
    topology->addTetra(6,8,5,7);
    topology->addTetra(20,3,10,2);

    topology->addTetra(12,9,13,19);
    topology->addTetra(11,8,17,6);
    topology->addTetra(3,20,1,0);
    topology->addTetra(20,19,13,10);
    topology->addTetra(19,9,10,16);

	EulerImplicitSolver::SPtr solver = sofa::core::objectmodel::New<EulerImplicitSolver>();
    solver->setName("solver");
    nodeTetra->addObject(solver);

    CGLinearSolver::SPtr linearSolver = New<CGLinearSolver>();
    linearSolver->setName("linearSolver");
    linearSolver->f_maxIter.setValue(2000);
    linearSolver->f_tolerance.setValue(1e-15);
    linearSolver->f_smallDenominatorThreshold.setValue(1e-15);
    nodeTetra->addObject(linearSolver);

    topology->setName("Container");
    topology->setSrc("@"+DOF->getName(), DOF.get());
    nodeTetra->addObject(topology);

	sofa::component::topology::TetrahedronSetTopologyModifier::SPtr tettopomod = sofa::core::objectmodel::New<sofa::component::topology::TetrahedronSetTopologyModifier>();
	tettopomod->setName("tettopomod");
	nodeTetra->addObject(tettopomod);

	TetrahedronSetTopologyAlgorithms3::SPtr tettopoalgo = sofa::core::objectmodel::New<TetrahedronSetTopologyAlgorithms3>();
	tettopoalgo->setName("tettopoalgo");
	nodeTetra->addObject(tettopoalgo);

	TetrahedronSetGeometryAlgorithms3::SPtr tetgeoalgo = sofa::core::objectmodel::New<TetrahedronSetGeometryAlgorithms3>();
	tetgeoalgo->setName("tetgeoalgo");
	nodeTetra->addObject(tetgeoalgo);

    UniformMass3::SPtr mass = sofa::core::objectmodel::New<UniformMass3>();
    nodeTetra->addObject(mass);
    mass->setMass(1);
    mass->setName("mass");

	FixedConstraint3::SPtr fixed1 = sofa::core::objectmodel::New<FixedConstraint3>();
	fixed1->removeConstraint(0);

	fixed1->addConstraint(0);
	fixed1->addConstraint(1);
	fixed1->addConstraint(2);
	fixed1->addConstraint(3);

	nodeTetra->addObject(fixed1);

	TetrahedralCorotationalFEMForceField3::SPtr fem = sofa::core::objectmodel::New<TetrahedralCorotationalFEMForceField3>();
	fem->setName("FEM");
	fem->setSrc("@"+DOF->getName(), DOF.get());
	fem->setMethod(2);
	fem->setUpdateStiffnessMatrix(true);
	fem->setComputeGlobalMatrix(true);
	fem->setYoungModulus(5e+5);
	fem->setPoissonRatio(0.1);
	nodeTetra->addObject(fem);

	sofa::simulation::getSimulation()->init(root.get());

	std::ofstream myfile;
    std::string f1 = "YOURPATH";
    std::string f2 = ".obj";
    string filename = "";

	for(int i=0;i<10;++i)
	{
		sofa::simulation::getSimulation()->animate(root.get(),0.01);

		if(i==2)
		{

			sofa::core::topology::BaseMeshTopology::Tetra tetra = topology->getTetra(18);

            sofa::component::topology::TetrahedronSetTopologyModifier* topoMod;
            topology->getContext()->get(topoMod);

            //add point
            topoMod->addPoints(1,true);

            MechanicalObject3::WriteVecCoord x_pos = DOF->writePositions();
            MechanicalObject3::WriteVecCoord x_restpos = DOF->writeRestPositions();

            int idx = DOF->getSize();

            x_pos[idx-1] = Coord3( topology->getPX(tetra[0]),topology->getPY(tetra[0])+0.2,topology->getPZ(tetra[0]) );
            x_restpos[idx-1] = Coord3( 0.2,1.0,1.2 );

            //add tetra
			helper::vector<core::topology::Topology::Tetrahedron > vitems;
			vitems.resize(1);
			vitems[0][0] = 6;
			vitems[0][1] = 8;
			vitems[0][2] = 7;
			vitems[0][3] = idx-1;
			topoMod->addTetrahedra(vitems);

            topoMod->notifyEndingEvent();
            topoMod->propagateTopologicalChanges();

            sofa::simulation::tree::getSimulation()->init( root.get() );

		}

		std::ostringstream filestring;
		filestring << f1 << i+1 << f2;
		filename = filestring.str();
		myfile.open(filename.c_str());

		for(int a=0;a<DOF->getSize();++a)
		{
			myfile << "v " << topology->getPX(a) << " " << topology->getPY(a) << " " << topology->getPZ(a) << "\n";
		}

		myfile.close();

	}
	
    return 0;
}

i just load the obj files into a 3d application and you can see the point shoot off into space.

here is a different example:

https://www.youtube.com/watch?v=0MMVekCvaAg

i was just trying to think what could be wrong. is there some corresponding rest mesh that i have to update? do i need to update the stiffness? maybe the order of the vertices in the tetra?

do you know what could be the issue?