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

models/bonemodel2.py

@@ -64,6 +64,7 @@ def parms(d={}):
     ]
 
     p['loads'] = [] # additional loads
+    p['extractors'] = [] # additional extractors
 
     # material properties
     p['Young'] = 17000.      # [MPa] elastic modulus - bone: 17-20 GPa
@@ -135,6 +136,10 @@ def solve(p={}):
     for fix in p['fixations']:
         create_group(fix['nodes'], nods_no, nods_pos, groups, fix['name'])
 
+    # create groups of additional extractors
+    for extr in p['extractors']:
+        create_group(extr['nodes'], nods_no, nods_pos, groups, extr['name'])
+
     # --------------------------------------------------------------------------
     # print('== ADDING MEDIUM...')
     print('')
@@ -189,6 +194,7 @@ def solve(p={}):
 
     print('extracting results...')
     print(f'\tinternal energy = {solver.int_energy:.3f} N.mm')
+    # fixations (reaction forces)
     for fix in p['fixations']:
         grpname = fix['name']
         v = np.array(post.probe('force_vector', grpname))
@@ -198,6 +204,7 @@ def solve(p={}):
         print(f'\t{grpname}_Fy = {v[1]:.3f} N')
         print(f'\t{grpname}_Fz = {v[2]:.3f} N')
 
+    # muscle groups
     for name, group in mgroups.items():
         v=np.array( [0.0, 0.0, 0.0] )
         for tag in group.ntags:
@@ -206,6 +213,46 @@ def solve(p={}):
         print(f'\t{name}_Fy = {v[1]:.3f} N')
         print(f'\t{name}_Fz = {v[2]:.3f} N')
 
+    # additional extractors
+
+    var2field = { 
+        'dx': ('X', 'mm'), 
+        'dy': ('Y', 'mm'), 
+        'dz': ('Z', 'mm'),  
+        'sigxx': ('smooth_stress_xx', 'MPa'),  
+        'sigyy': ('smooth_stress_yy', 'MPa'), 
+        'sigzz': ('smooth_stress_zz', 'MPa'), 
+        'sigxy': ('smooth_stress_xy', 'MPa'), 
+        'sigxz': ('smooth_stress_xz', 'MPa'), 
+        'sigyz': ('smooth_stress_yz', 'MPa'), 
+        'epsxx': ('strain_xx', ''), 
+        'epsyy': ('strain_yy', ''), 
+        'epszz': ('strain_zz', ''), 
+        'epsxy': ('strain_xy', ''), 
+        'epsxz': ('strain_xz', ''), 
+        'epsyz': ('strain_yz', '')                          
+    }
+
+    for extr in p['extractors']:
+        # if only one variable is given, convert it to a list
+        if isinstance(extr['variables'], str):
+            extr['variables'] = [extr['variables']]
+
+        for var in extr['variables']:
+            fname = f'{extr["name"]}_{var}.tsv'
+            grpname = extr['name']
+            if var.lower() in var2field:
+                field = var2field[var.lower()][0]
+                units = var2field[var.lower()][1]
+                v = np.array(post.probe(field, grpname))
+                np.savetxt(fname, v, delimiter='\t', header='value')
+                if len(v) == 1:
+                    print(f'\t{grpname}_{var} = {v[0]:.3f} {units}')
+                else:
+                    print(f'\t{grpname}_{var} = ({len(v)} values in "{fname}")')
+            else:
+                print(f'unknown variable: "{var}"')
+
     end_time = time.perf_counter()
     elapsed_time = end_time - start_time
     print(f'Total Execution time: {elapsed_time:.1f}s')