Extractors

"Extractors" have been added to the model. They allow the user to read particular results such as the equivalent Von-Mises stress at given positions in the mesh after the FE analysis. The values are displayed in the console at the end of the calculation (if only 1 node is concerned). They are also saved to disc

• in a .tsv (tab-separated values) file when fossils is used,
• or in a .ascii file when Metafor is used.

The input syntax is documented [here in the wiki].(https://gitlab.uliege.be/rboman/fossils/-/wikis/input-file-for-boneload-simulations#extractors).

Some extractors have been added to models/dolicorhynchops/dolicorhynchops_10k.py:

    p['extractors'] = [
        {
            'name': 'tooth1',
            'nodes': [-10.20362, -17.46838, -229.9061],
            'variables': [ 'sigvm' ]
        },
        {
            'name': 'tooth2',
            'nodes': [-11.92466, 26.3042, -229.5354],
            'variables': [ 'sigxx', 'sigyy', 'sigzz', 'sigvm' ]
        },
        {
            'name': 'Lmuscle',
            'nodes': f'{path}/Lmuscle.stl',
            'variables': [ 'x0', 'y0', 'z0', 'sigvm' ]
        }
    ]

The values are calculated at the classical extractor list:

extracting results...
        internal energy = 95.219 N.mm
        contact_pts_Fx = 17.915 N
        [...]
        tooth1_sigvm = 4.507 MPa
        tooth2_sigxx = -1.000 MPa
        tooth2_sigyy = -1.025 MPa
        tooth2_sigzz = -0.687 MPa
        tooth2_sigvm = 3.803 MPa
        Lmuscle_x0 = (95 values in "Lmuscle_x0.tsv")
        Lmuscle_y0 = (95 values in "Lmuscle_y0.tsv")
        Lmuscle_z0 = (95 values in "Lmuscle_z0.tsv")
        Lmuscle_sigvm = (95 values in "Lmuscle_sigvm.tsv")
Total Execution time: 2.2s

cxxfem/src/femPost.cpp

@@ -24,11 +24,45 @@ Post::~Post()
         std::cout << "~Post()\n";
 }
 
-/// rebuild solution vectors for gmsh and add views
+/// rebuild solution vectors for gmsh and add views:
+///
+///  X
+///  Y
+///  Z
+///  displacement_vector
+///  force_vector
+///  stress_tensor
+///      stress_xx
+///      stress_yy
+///      stress_zz
+///      stress_xy
+///      stress_xz
+///      stress_yz
+///  strain_tensor
+///      strain_xx
+///      strain_yy
+///      strain_zz
+///      strain_xy
+///      strain_xz
+///      strain_yz
+///  smooth_stress_tensor
+///      smooth_stress_xx
+///      smooth_stress_yy
+///      smooth_stress_zz
+///      smooth_stress_xy
+///      smooth_stress_xz
+///      smooth_stress_yz
+///  smooth_strain_tensor
+///      smooth_strain_xx
+///      smooth_strain_yy
+///      smooth_strain_zz
+///      smooth_strain_xy
+///      smooth_strain_xz
+///      smooth_strain_yz
 
 void Post::build_views()
 {
-    std::cout << "creating views...\n";
+    std::cout << "creating views..." << std::endl;
     // nodal dofs
     add_dofs_view(solver.pbl, solver.dofs, solver.X, solver.dofPart, solver.dofIdx);
 
@@ -330,6 +364,19 @@ Post::probe(std::string const &name, std::string const &grpname)
     // get physical group from grpname
     Group *grp = solver.pbl.groups_by_name.at(grpname); // TODO: better handling of errors
 
+    std::vector<std::size_t> nodeTags;
+    std::vector<double> coord;
+    gmsh::model::mesh::getNodesForPhysicalGroup(grp->dim, grp->tag, nodeTags, coord);
+    // j'ai essayé de boucler sur les entités du groupe et faire un getNode de l'entité 
+    // mais ça ne me retourne rien (?)
+    // std::cout << "nodeTags.size()=" << nodeTags.size() << '\n';
+
+    // retourne les coordonnées des noeuds du groupe si demandé
+    if(name=="coords")
+    {   
+        return coord;
+    }
+
     // get the view
     auto it = views.find(name);
     if (it == views.end())
@@ -339,14 +386,6 @@ Post::probe(std::string const &name, std::string const &grpname)
     // probe values
     std::vector<double> all_values;
 
-    std::vector<std::size_t> nodeTags;
-    std::vector<double> coord;
-    gmsh::model::mesh::getNodesForPhysicalGroup(grp->dim, grp->tag, nodeTags, coord);
-    // j'ai essayé de boucler sur les entités du groupe et faire un getNode de l'entité 
-    // mais ça ne me retourne rien (?)
-
-    // std::cout << "nodeTags.size()=" << nodeTags.size() << '\n';
-
     for (size_t i = 0; i < nodeTags.size(); ++i)
     {
         std::vector<double> values;
@@ -365,6 +404,21 @@ Post::probe(std::string const &name, std::string const &grpname)
         // std::cout << "\tgot from gmsh: " << values << '\n';
         // std::cout << "\tdistance = " << distance << std::endl;
 
+        if(values.size() == 9) // tensor -> compute j2
+        {
+            // xx, xy, xz, yx, yy, yz, zx, zy, zz
+            //  0,  1,  2,  3,  4,  5,  6,  7,  8
+            double xx = values[0];
+            double yy = values[4];
+            double zz = values[8];
+            double xy = values[1];
+            double zx = values[6];
+            double yz = values[5];
+            double j2 = sqrt( ((xx-yy)*(xx-yy) + (yy-zz)*(yy-zz) + (zz-xx)*(zz-xx) )/2 + 3*(xy*xy+yz*yz+zx*zx) );
+            values.clear();
+            values.push_back(j2);
+        }
+
         all_values.insert(all_values.end(), values.begin(), values.end());
     }