Dear ITK users,I have another question. In order to visualize the deformation,
I want to replace the old node coordinates with the new one after generating
the solution, like in UpdateDisplacements() function in itkFemSolver.hxx 610-651
The code-block, which generates solution:m_FemSolver->Update();
//m_FemSolver->UpdateDisplacements();
const unsigned int invalidId= itk::fem::Element::InvalidDegreeOfFreedomID;
int numberOfNodes= m_FemSolver->GetInput()->GetNumberOfNodes();
for(int i= 0; i< numberOfNodes; i++)
{
itk::fem::Element::Node::Pointer node= m_FemSolver->GetInput()->GetNode(i);
std::cout<<"FEM Mesh Node: "<< node->GetGlobalNumber() << " Coordinates: "<<
node->GetCoordinates()<<std::endl;
for(unsigned int j= 0, dof; (dof= node->GetDegreeOfFreedom(j))!= invalidId;
j++)
{
cout <<"FEM Mesh Solution: "<< m_FemSolver->GetSolution(dof)<<endl;
}
}
gives me following output, that shows, I have some solutions for nodes 69 and
70:FEM Mesh Node: 68 Coordinates: -61.1098 -139.262 211.9
FEM Mesh Node: 69 Coordinates: -61.0742 -139.319 211.9
FEM Mesh Solution: 82.6806
FEM Mesh Solution: -32.5517
FEM Mesh Solution: 6.36021
FEM Mesh Node: 70 Coordinates: -61.3429 -138.832 211.9
FEM Mesh Solution: 0.798099
FEM Mesh Solution: 37.6971
FEM Mesh Solution: -4.36149
FEM Mesh Node: 71 Coordinates: -61.5039 -138.402 211.509
FEM Mesh Node: 72 Coordinates: -61.6232 -138.402 211.9
if I want to update the displacements://itk::fem::LinearSystemWrapper::Pointer
wrapper;
VectorType replacedNodeVector(3);
for(int i= 0; i< numberOfNodes; i++)
{
TetrahedronNodeType::Pointer node= m_FemSolver->GetInput()->GetNode(i);
//itk::fem::Element::Node::Pointer node=
m_Fem3dSolver->GetInput()->GetNode(i);
for(unsigned int j= 0; j< 3; j++)
{
std::cout<<"FEM Mesh Solution Vector:
"<<m_FemSolver->GetSolution(node->GetDegreeOfFreedom(j))<<std::endl;
replacedNodeVector[j] = node->GetCoordinates()[j] +
m_FemSolver->GetSolution(node->GetDegreeOfFreedom(j));
//replacedNodeVector[j] = node->GetCoordinates()[j] +
m_ls->GetSolutionValue(node->GetDegreeOfFreedom(j));
}
node->SetCoordinates(replacedNodeVector);
std::cout<<"FEM Mesh New Node: "<< node->GetGlobalNumber() << " Coordinates:
"<< node->GetCoordinates()<<std::endl;
}
I get as output:FEM Mesh New Node: 67 Coordinates: -61.0742 -139.262 211.797
FEM Mesh Solution Vector: 0
FEM Mesh Solution Vector: 0
FEM Mesh Solution Vector: 0
FEM Mesh New Node: 68 Coordinates: -61.1098 -139.262 211.9
FEM Mesh Solution Vector: 82.6806
FEM Mesh Solution Vector: -32.5517
FEM Mesh Solution Vector: 6.36021
FEM Mesh New Node: 69 Coordinates: 21.6064 -171.871 218.26
FEM Mesh Solution Vector: 0.798099
FEM Mesh Solution Vector: 37.6971
FEM Mesh Solution Vector: -4.36149
FEM Mesh New Node: 70 Coordinates: -60.5448 -101.135 207.539
FEM Mesh Solution Vector: 0
FEM Mesh Solution Vector: 0
FEM Mesh Solution Vector: 0
My questions are:
Do I do it correct for the replacements, if I sum the old node coordinate with
the generated solution vector for that spacial node?:replacedNodeVector[j] =
node->GetCoordinates()[j] +
m_FemSolver->GetSolution(node->GetDegreeOfFreedom(j));
In my test application, I apply some force at node 1( loadNode->SetNode(1); )
but in generated solution it shows that the force doesn’t act on node 1 and
nodes around it but for nodes 68, 69, 82, and 83. All the other solutions are
zero. My test mesh contains 86 nodes and 162 elements. If I want to apply the
force except node 1,eg at node 20, I get then segmentation fault.typedef
itk::fem::LoadNode LoadNodeType;
LoadNodeType::Pointer loadNode= LoadNodeType::New();
loadNode->SetElement(tetrahedronElement.GetPointer());
loadNode->SetGlobalNumber(3);
loadNode->SetNode(1);
vnl_vector<double> force(3);
force[0]= 110.0;
force[1]= 20.0;
force[2]= 31.0;
loadNode->SetForce(force);
m_FemObject->AddNextLoad(loadNode);
For detailed information please see runnable test code:vtkUGridToFemMesh.zip
|
|
|
| | |
|
|
|
| |
vtkUGridToFemMesh.zip
Shared with Dropbox | |
|
|
Thank you for any help!
The ITK community is transitioning from this mailing list to discourse.itk.org.
Please join us there!
______________________________________
Powered by www.kitware.com
Visit other Kitware open-source projects at
http://www.kitware.com/opensource/opensource.html
Kitware offers ITK Training Courses, for more information visit:
http://kitware.com/products/protraining.php
Please keep messages on-topic and check the ITK FAQ at:
http://www.itk.org/Wiki/ITK_FAQ
Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/insight-developers
