[prev in list] [next in list] [prev in thread] [next in thread] 

List:       insight-users
Subject:    [ITK-users] itk::fem updating displacements
From:       Yusuf OEZBEK via Insight-users <insight-users () itk ! org>
Date:       2017-11-09 19:59:10
Message-ID: 587898442.10039537.1510257550772 () mail ! yahoo ! com
[Download RAW message or body]

[Attachment #2 (multipart/alternative)]


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!


[Attachment #5 (text/html)]

<html xmlns="http://www.w3.org/1999/xhtml" xmlns:v="urn:schemas-microsoft-com:vml" \
xmlns:o="urn:schemas-microsoft-com:office:office"><head><!--[if gte mso \
9]><xml><o:OfficeDocumentSettings><o:AllowPNG/><o:PixelsPerInch>96</o:PixelsPerInch></o:OfficeDocumentSettings></xml><![endif]--></head><body><div \
style="color:#000; background-color:#fff; font-family:Helvetica Neue, Helvetica, \
Arial, Lucida Grande, sans-serif;font-size:10px"><div \
id="yui_3_16_0_ym19_1_1510251141194_6598"><font size="2">Dear ITK users,</font></div> \
<div id="yui_3_16_0_ym19_1_1510251141194_6599"><font size="2">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<br \
id="yui_3_16_0_ym19_1_1510251141194_6600"> The code-block, which generates \
solution:</font></div> <pre id="yui_3_16_0_ym19_1_1510251141194_6601"><font \
id="yui_3_16_0_ym19_1_1510251141194_6707" size="2"><code \
id="yui_3_16_0_ym19_1_1510251141194_6602">m_FemSolver-&gt;Update(); \
//m_FemSolver-&gt;UpdateDisplacements(); const unsigned int invalidId= \
itk::fem::Element::InvalidDegreeOfFreedomID; int numberOfNodes= \
m_FemSolver-&gt;GetInput()-&gt;GetNumberOfNodes();

for(int i= 0; i&lt; numberOfNodes; i++)
{
  itk::fem::Element::Node::Pointer node= m_FemSolver-&gt;GetInput()-&gt;GetNode(i);
  std::cout&lt;&lt;"FEM Mesh Node: "&lt;&lt; node-&gt;GetGlobalNumber() &lt;&lt; " \
Coordinates: "&lt;&lt; node-&gt;GetCoordinates()&lt;&lt;std::endl;  
  for(unsigned int j= 0, dof; (dof= node-&gt;GetDegreeOfFreedom(j))!= invalidId; j++)
  {
    cout &lt;&lt;"FEM Mesh Solution: "&lt;&lt; \
m_FemSolver-&gt;GetSolution(dof)&lt;&lt;endl;  }
}
</code></font></pre>
<div id="yui_3_16_0_ym19_1_1510251141194_6603"><font \
id="yui_3_16_0_ym19_1_1510251141194_6665" size="2">gives me following output, that \
shows, I have some solutions for nodes 69 and 70:</font></div> <pre \
id="yui_3_16_0_ym19_1_1510251141194_6604"><font \
id="yui_3_16_0_ym19_1_1510251141194_6664" size="2"><code \
id="yui_3_16_0_ym19_1_1510251141194_6605">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
</code></font></pre>
<div id="yui_3_16_0_ym19_1_1510251141194_6606"><font size="2">if I want to update the \
displacements:</font></div> <pre id="yui_3_16_0_ym19_1_1510251141194_6607"><font \
id="yui_3_16_0_ym19_1_1510251141194_6644" size="2"><code \
id="yui_3_16_0_ym19_1_1510251141194_6608">//itk::fem::LinearSystemWrapper::Pointer \
wrapper; VectorType replacedNodeVector(3);
for(int i= 0; i&lt; numberOfNodes; i++)
{
  TetrahedronNodeType::Pointer node= m_FemSolver-&gt;GetInput()-&gt;GetNode(i);
  //itk::fem::Element::Node::Pointer node= \
m_Fem3dSolver-&gt;GetInput()-&gt;GetNode(i);

  for(unsigned int j= 0; j&lt; 3; j++)
  {
    std::cout&lt;&lt;"FEM Mesh Solution Vector: \
"&lt;&lt;m_FemSolver-&gt;GetSolution(node-&gt;GetDegreeOfFreedom(j))&lt;&lt;std::endl;
  replacedNodeVector[j] = node-&gt;GetCoordinates()[j] +  \
m_FemSolver-&gt;GetSolution(node-&gt;GetDegreeOfFreedom(j));  //replacedNodeVector[j] \
= node-&gt;GetCoordinates()[j] + \
m_ls-&gt;GetSolutionValue(node-&gt;GetDegreeOfFreedom(j));  }
  node-&gt;SetCoordinates(replacedNodeVector);
  std::cout&lt;&lt;"FEM Mesh New Node: "&lt;&lt; node-&gt;GetGlobalNumber() &lt;&lt; \
" Coordinates: "&lt;&lt; node-&gt;GetCoordinates()&lt;&lt;std::endl; }
</code></font></pre>
<div id="yui_3_16_0_ym19_1_1510251141194_6609"><font size="2">I get as \
output:</font></div> <pre id="yui_3_16_0_ym19_1_1510251141194_6610"><font \
id="yui_3_16_0_ym19_1_1510251141194_6643" size="2"><code \
id="yui_3_16_0_ym19_1_1510251141194_6611">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
</code></font></pre>
<div id="yui_3_16_0_ym19_1_1510251141194_6612"><font \
id="yui_3_16_0_ym19_1_1510251141194_6720" size="2">My questions are:<br \
id="yui_3_16_0_ym19_1_1510251141194_6613"> Do I do it correct for the replacements, \
if I sum the old node  coordinate with the generated solution vector for that spacial \
node?:</font></div> <pre id="yui_3_16_0_ym19_1_1510251141194_6614"><font \
id="yui_3_16_0_ym19_1_1510251141194_6736" size="2"><code \
id="yui_3_16_0_ym19_1_1510251141194_6615">replacedNodeVector[j] = \
node-&gt;GetCoordinates()[j] + \
m_FemSolver-&gt;GetSolution(node-&gt;GetDegreeOfFreedom(j)); </code></font></pre>
<div id="yui_3_16_0_ym19_1_1510251141194_6616"><font \
id="yui_3_16_0_ym19_1_1510251141194_6662" size="2">In my test application, I apply \
some force at node 1(  loadNode-&gt;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.</font></div>
<pre id="yui_3_16_0_ym19_1_1510251141194_6617"><font \
id="yui_3_16_0_ym19_1_1510251141194_6659" size="2"><code \
id="yui_3_16_0_ym19_1_1510251141194_6618">typedef itk::fem::LoadNode LoadNodeType; \
LoadNodeType::Pointer loadNode= LoadNodeType::New(); \
loadNode-&gt;SetElement(tetrahedronElement.GetPointer()); \
loadNode-&gt;SetGlobalNumber(3); loadNode-&gt;SetNode(1);

vnl_vector&lt;double&gt; force(3);
force[0]= 110.0;
force[1]= 20.0;
force[2]= 31.0;
loadNode-&gt;SetForce(force);
m_FemObject-&gt;AddNextLoad(loadNode);
</code></font></pre><div id="yui_3_16_0_ym19_1_1510251141194_6718">
<font size="2">For detailed information please see runnable test \
code:</font></div><div dir="ltr" id="yui_3_16_0_ym19_1_1510251141194_6719"><font \
id="yui_3_16_0_ym19_1_1510251141194_6755" size="2"><a \
id="yui_3_16_0_ym19_1_1510251141194_6754" \
class="enhancr2_c997ba3c-d147-b7a8-faff-66dfef0cd3dd" \
href="https://www.dropbox.com/s/n551i11lxv22m3b/vtkUGridToFemMesh.zip?dl=0">vtkUGridToFemMesh.zip</a></font></div><div><br></div><div \
id="enhancr2_c997ba3c-d147-b7a8-faff-66dfef0cd3dd" class="yahoo-link-enhancr-card  \
ymail-preserve-class ymail-preserve-style" \
style="max-width:400px;font-family:'Helvetica Neue', Helvetica, Arial, sans-serif;" \
data-url="https://www.dropbox.com/s/n551i11lxv22m3b/vtkUGridToFemMesh.zip?dl=0" \
data-type="yenhancr" data-category="website" data-embed-url="" data-size="medium" \
dir="ltr" contenteditable="false"> <a id="yui_3_16_0_ym19_1_1510251141194_6782" \
href="https://www.dropbox.com/s/n551i11lxv22m3b/vtkUGridToFemMesh.zip?dl=0" \
style="text-decoration:none !important; color: #000 !important;" \
class="yahoo-enhancr-cardlink" target="_blank" rel="noreferrer"> <table \
id="yui_3_16_0_ym19_1_1510251141194_6781" class="card-wrapper yahoo-ignore-table" \
style="max-width:400px;" border="0" cellpadding="0" cellspacing="0"> <tbody \
id="yui_3_16_0_ym19_1_1510251141194_6780"><tr \
id="yui_3_16_0_ym19_1_1510251141194_6779"> <td \
id="yui_3_16_0_ym19_1_1510251141194_6778" width="400"> <table \
id="yui_3_16_0_ym19_1_1510251141194_6777" class="card yahoo-ignore-table" \
style="max-width:400px;" border="0" cellpadding="0" cellspacing="0" width="100%"> \
<tbody id="yui_3_16_0_ym19_1_1510251141194_6776"><tr \
id="yui_3_16_0_ym19_1_1510251141194_6775"> <td \
id="yui_3_16_0_ym19_1_1510251141194_6774" class="card-primary-image-cell" \
style="background:#000 \
url('https://s.yimg.com/vv//api/res/1.2/pRVwB03of2_aX_BXhn7gdQ--~A/YXBwaWQ9bWFpbDtmaT1 \
maWxsO2g9MjAwO3c9NDAw/https://cfl.dropboxstatic.com/static/images/logo_catalog/glyph_m1@2x-vflA6lTFZ.png.cf.jpg') \
no-repeat center center;background-size:cover;height:200px;position:relative;" \
background="https://s.yimg.com/vv//api/res/1.2/pRVwB03of2_aX_BXhn7gdQ--~A/YXBwaWQ9bWFp \
bDtmaT1maWxsO2g9MjAwO3c9NDAw/https://cfl.dropboxstatic.com/static/images/logo_catalog/glyph_m1@2x-vflA6lTFZ.png.cf.jpg" \
bgcolor="#000000" valign="top"> <!--[if gte mso 9]><v:rect fill="true" stroke="false" \
style="width:400px;height:218px;position:absolute;top:0;left:0;"><v:fill type="frame" \
color="#000000" src="https://s.yimg.com/vv//api/res/1.2/pRVwB03of2_aX_BXhn7gdQ--~A/YXB \
waWQ9bWFpbDtmaT1maWxsO2g9MjAwO3c9NDAw/https://cfl.dropboxstatic.com/static/images/logo_catalog/glyph_m1@2x-vflA6lTFZ.png.cf.jpg"/></v:rect><![endif]--> \
<table id="yui_3_16_0_ym19_1_1510251141194_6773" class="yahoo-ignore-table" \
valign="top" style="width:100%;" border="0" cellpadding="0" cellspacing="0"> <tbody \
id="yui_3_16_0_ym19_1_1510251141194_6772"><tr \
id="yui_3_16_0_ym19_1_1510251141194_6771"> <td \
id="yui_3_16_0_ym19_1_1510251141194_6770" style="background:transparent \
url('https://s.yimg.com/nq/storm/assets/enhancrV2/12/overlay-tile.png') repeat left \
top;height:200px;" background="https://s.yimg.com/nq/storm/assets/enhancrV2/12/overlay-tile.png" \
bgcolor="transparent" valign="top"> <!--[if gte mso 9]><v:rect fill="true" \
stroke="false" style="width:400px;height:218px;position:absolute;top:-18px;left:0;"><v:fill \
type="pattern" color="#000000" \
src="https://s.yimg.com/nq/storm/assets/enhancrV2/12/overlay-tile.png"/><v:textbox \
inset="0,0,20px,0"><![endif]--> <table id="yui_3_16_0_ym19_1_1510251141194_6769" \
class="yahoo-ignore-table" style="width:100%;height:185px;min-height:185px;" \
height="185"> <tbody id="yui_3_16_0_ym19_1_1510251141194_6768"><tr \
id="yui_3_16_0_ym19_1_1510251141194_6767"> <td \
id="yui_3_16_0_ym19_1_1510251141194_6821" class="card-richInfo2" \
style="text-align:left;text-align:left;padding:15px 0 0 15px;vertical-align:top;">  \
</td> <td id="yui_3_16_0_ym19_1_1510251141194_6766" class="card-actions" \
style="text-align:right;padding:15px 15px 0 0;vertical-align:top;"> <div \
class="card-share-container"></div> </td> </tr> </tbody></table> <!--[if gte mso \
9]></v:textbox></v:rect><![endif]--> </td> </tr> </tbody></table> </td> </tr> <tr \
id="yui_3_16_0_ym19_1_1510251141194_6831"> <td \
id="yui_3_16_0_ym19_1_1510251141194_6830"> <table \
id="yui_3_16_0_ym19_1_1510251141194_6829" class="card-info yahoo-ignore-table" \
style="background:#fff;position:relative;z-index:2;width:95%;max-width:380px;border:1px \
solid #e0e4e9;border-bottom:3px solid \
#007ee4;margin-top:-40px;margin-left:auto;margin-right:auto;" align="center" \
border="0" cellpadding="0" cellspacing="0"> <tbody \
id="yui_3_16_0_ym19_1_1510251141194_6828"><tr \
id="yui_3_16_0_ym19_1_1510251141194_6827"> <td \
style="background-color:#ffffff;padding:16px 0 16px 12px;vertical-align:top;"> <img \
data-id="0a754ae0-9d93-414c-e0fd-617ebe678283" class="card-object-1 \
yahoo-ignore-inline-image ymail-preserve-class" \
src="https://s.yimg.com/nq/storm/assets/enhancrV2/23/logos/dropbox.png" \
style="min-width:32px;border:1px solid #e0e4e9;margin-top:3px;" height="32">  </td> \
<td id="yui_3_16_0_ym19_1_1510251141194_6826" \
style="vertical-align:middle;padding:16px 12px;width:99%;"> <h2 class="card-title" \
style="font-size: 16px; line-height:19px; margin:0 0 4px 0;font-family:'Helvetica \
Neue', Helvetica, Arial, \
sans-serif;word-break:break-word;">vtkUGridToFemMesh.zip</h2>  <div \
class="card-description" \
style="font-size:11px;line-height:15px;color:#999;word-break:break-word;">Shared with \
Dropbox</div> </td> <td style="text-align:right;padding:16px 12px 16px 0;">  </td> \
</tr> </tbody></table> </td> </tr> </tbody></table> </td> </tr> </tbody></table> \
</a></div><div id="yui_3_16_0_ym19_1_1510251141194_6824"><br></div><div \
id="yui_3_16_0_ym19_1_1510251141194_6810" dir="ltr"><font \
id="yui_3_16_0_ym19_1_1510251141194_6825" size="2">Thank you for any \
help!<br></font></div></div></body></html>



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://www.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-users


[prev in list] [next in list] [prev in thread] [next in thread] 

Configure | About | News | Add a list | Sponsored by KoreLogic