Tetrahedra microFE pipeline - trabecular bone
From micro-CT image to tetrahedra-FE model solution in Calculix
Created on: 23.05.2021
Last update: 15.01.2023
By: Gianluca Iori, Martino Pani, Gianluigi Crimi, Enrico Schileo, Fulvia Taddei, Giulia Fraterrigo, 2023
Data source: the dataset used in this example is part of the public collection of the Living Human Digital Library (LHDL) Project, a project financed by the European Commission (project number: FP6-IST 026932. Human tissues in the LHDL project were collected according to the body donation program of Universitè Libre de Bruxelles (ULB), partner of the project. For info on the dataset see here.
Code license: MIT
Narrative license: CC-BY-NC-SA
Aims
The example implements the following ciclope pipeline:
Load and inspect microCT volume data
Image pre-processing
Apply Gaussian smooth (optional)
Resample the dataset (optional)
Segment bone tissue
Remove unconnected clusters of voxels
Mesh generation
Generate 3D Unstructured Grid mesh of tetrahedra
Generate high-resolution mesh of triangles of the model outer shell (optional for visualization)
Generate tetrahedra-FE model for simulation in CalculX or Abaqus from 3D Unstructured Grid mesh
Linear, static analysis definition: uniaxial compression test
Launch simulation in Calculix. For info on the solver visit the Calculix homepage
Convert Calculix output to .VTK for visualization in Paraview
Calculate apparent elastic modulus from reaction forces
Notes on ciclope
Tetrahedra meshes are generated with pygalmesh (a Python frontend to CGAL)
High-resolution surface meshes for visualization are generated with the PyMCubes module.
All mesh exports are performed with the meshio module.
ciclope handles the definition of material properties and FE analysis parameters (e.g. boundary conditions, simulation steps..) through separate template files. The folders material_properties and input_templates contain a library of template files that can be used to generate FE simulations.
Command line execution
The pipeline described in this notebook can be executed from the command line using the ciclope command:
ciclope test_data/LHDL/3155_D_4_bc/cropped/3155_D_4_bc_0000.tif test_data/LHDL/3155_D_4_bc/results/3155_D_4_bc.inp -vs 0.0195 0.0195 0.0195 -r 2 -t 63 --smooth 1 --tetrafe --max_facet_distance 0.025 --max_cell_circumradius 0.05 --vol_mesh --template input_templates/tmp_example01_comp_static_bone.inp
Type python ciclope.py -h
to display the ciclope help with a full list of available command line arguments.
Configuration and imports
import sys
sys.path.append('./../../')
import os
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
import mcubes
from scipy import ndimage, misc
from skimage.filters import threshold_otsu, gaussian
from skimage import measure, morphology
import meshio
from ciclope.utils.recon_utils import read_tiff_stack, plot_midplanes
from ciclope.utils.preprocess import remove_unconnected
from ciclope.core import tetraFE
import ciclope
Uncomment and run the following cell for visualizations using itkwidgets
# import vtk
# import itk
# from itkwidgets import view
ccx2paraview is needed for the postprocessing of results
import ccx2paraview
matplotlib.rcParams['figure.dpi'] = 300
font = {'weight' : 'normal',
'size' : 8}
plt.rc('font', **font)
Load input data
input_file = './../../test_data/LHDL/3155_D_4_bc/cropped/3155_D_4_bc_0000.tif'
Input data description
Scan parameters |
|
---|---|
Lab |
LTM@IOR Bologna |
Sample |
Human trabecular bone |
Voxel size |
19.5 micron |
Preliminary operations |
cropped to 200x200x200 (4mm side) |
Read the input data and define an array of the voxelsize
data_3D = read_tiff_stack(input_file)
vs = np.ones(3)*19.5e-3 # [mm]
Inspect the dataset
plot_midplanes(data_3D)
plt.show()
Inspect the dataset with itkwidgets
# viewer = view(data_3D, ui_collapsed=True)
# viewer.interpolation = False
Launch the itk viewer
# viewer
Pre-processing
Gaussian smooth
data_3D = gaussian(data_3D, sigma=1, preserve_range=True)
Resize (optional)
resampling = 2
# resize the 3D data using spline interpolation of order 2
data_3D = ndimage.zoom(data_3D, 1/resampling, output=None, order=2)
# correct voxelsize
vs = vs * resampling
# Inspect again the dataset
plot_midplanes(data_3D)
plt.show()
Thresholding
Use Otsu’s method
T = threshold_otsu(data_3D)
print("Threshold: {}".format(T))
Threshold: 80.72472439570987
# plot the input image histogram
fig2, ax2 = plt.subplots()
plt.hist(data_3D.ravel(), bins=40)
plt.show()
Apply the threshold
# BW = data_3D > T
BW = data_3D > 63 # from comparison with histology
Have a look at the binarized dataset
plot_midplanes(BW)
plt.show()
Morphological close of binary image (optional)
BW = morphology.closing(BW, morphology.ball(3))
plot_midplanes(BW)
plt.show()
Detect largest isolated cluster of voxels
The procedure is described by the following 3 steps:
Label the BW 3D image
# [labels, n_labels] = measure.label(BW, None, True, 1) # 1 connectivity hop
Count the occurrences of each label
# occurrences = np.bincount(labels.reshape(labels.size))
Find the largest unconnected label
# largest_label_id = occurrences[1:].argmax()+1
# L = labels == largest_label_id
You can do the same using the function remove_unconnected
contained in `ciclope.utils.preprocess’
L = remove_unconnected(BW)
Inspect dataset
plot_midplanes(L)
plt.show()
Shell mesh for visualization
Create a mesh of triangles of the external bone surface (optional) Warning: this step is computationally expensive. Only useful for 3D rendering of the bone architecture.
Uncomment the following cell to create a shell mesh of triangles
# vertices, triangles, mesh = ciclope.tetraFE.shell_mesh(L, method='pymcubes')
Uncomment the following two cells to write the mesh as a VTK file with meshio
# filename_shellmesh_out = './../../test_data/LHDL/3155_D_4_bc/results/3155_D_4_bc_shellmesh.vtk'
# meshio.write_points_cells(filename_shellmesh_out, vertices.tolist(), [("triangle", triangles.tolist())])
Visualize the mesh with itkwidgtes
# reader = vtk.vtkUnstructuredGridReader()
# reader.SetFileName(filename_shellmesh_out)
# reader.Update()
# grid = reader.GetOutput()
# view(geometries=grid)
Create tetrahedra mesh
Volume meshing using CGAL through pygalmesh
filename_mesh_out = './../../test_data/LHDL/3155_D_4_bc/results/3155_D_4_bc_tetramesh.vtk'
This is the pygalmesh
function call
# import pygalmesh
# mesh = pygalmesh.generate_from_array(np.transpose(L, [2, 1, 0]).astype('uint8'), tuple(vs), max_facet_distance=1.2*min(vs), max_cell_circumradius=5*1.2*min(vs))
Alternatively use the method ciclope.core.tetraFE.cgal_mesh
Reduce mesh_size_factor
(keeping it small for fast execute)
mesh_size_factor = 5 # 1.2
mesh = tetraFE.cgal_mesh(L, vs, 'tetra', mesh_size_factor*min(vs), 2*mesh_size_factor*min(vs))
construct initial points (nb_points: 12)
1/12 initial point(s) found...
2/12 initial point(s) found...
3/12 initial point(s) found...
4/12 initial point(s) found...
5/12 initial point(s) found...
6/12 initial point(s) found...
7/12 initial point(s) found...
8/12 initial point(s) found...
9/12 initial point(s) found...
10/12 initial point(s) found...
11/12 initial point(s) found...
12/12 initial point(s) found...
Start surface scan...Scanning triangulation for bad facets (sequential) - number of finite facets = 51...
Number of bad facets: 8
scanning edges (lazy)
scanning vertices (lazy)
end scan. [Bad facets:8]
Refining Surface...
Legend of the following line: (#vertices,#steps,#facets to refine,#tets to refine)
(12,0,8,0)
(13,0,10,0) (0.0 vertices/s)
(14,1,10,0) (10754.6 vertices/s)
(15,2,13,0) (14098.5 vertices/s)
(16,3,9,0) (16215.1 vertices/s)
(17,4,11,0) (17476.3 vertices/s)
(18,5,9,0) (17985.9 vertices/s)
(19,6,6,0) (18752.5 vertices/s)
(20,7,6,0) (19290.5 vertices/s)
(21,8,3,0) (19611.0 vertices/s)
(22,9,2,0) (20186.5 vertices/s)
(23,10,5,0) (20794.8 vertices/s)
(24,11,4,0) (21399.5 vertices/s)
(25,12,3,0) (21546.1 vertices/s)
(26,13,3,0) (22004.0 vertices/s)
(27,14,2,0) (22438.0 vertices/s)
(28,15,1,0) (22803.4 vertices/s)
(29,16,5,0) (22896.2 vertices/s)
(30,17,6,0) (23355.1 vertices/s)
(31,18,12,0) (23258.6 vertices/s)
(32,19,18,0) (23146.0 vertices/s)
(33,20,21,0) (23071.0 vertices/s)
(34,21,17,0) (23179.0 vertices/s)
(35,22,16,0) (23431.9 vertices/s)
(36,23,14,0) (23471.8 vertices/s)
(37,24,18,0) (23486.5 vertices/s)
(38,25,18,0) (23410.9 vertices/s)
(39,26,18,0) (23573.7 vertices/s)
(40,27,21,0) (23258.6 vertices/s)
(41,28,20,0) (23278.6 vertices/s)
(42,29,18,0) (23445.4 vertices/s)
(43,30,20,0) (23475.6 vertices/s)
(44,31,23,0) (22911.6 vertices/s)
(45,32,25,0) (22908.0 vertices/s)
(46,33,25,0) (23157.4 vertices/s)
(47,34,27,0) (23290.3 vertices/s)
(48,35,23,0) (23555.9 vertices/s)
(49,36,25,0) (23748.8 vertices/s)
(50,37,27,0) (23871.6 vertices/s)
(51,38,28,0) (23960.2 vertices/s)
(52,39,20,0) (24255.3 vertices/s)
(53,40,17,0) (24467.3 vertices/s)
(54,41,17,0) (24700.7 vertices/s)
(55,42,17,0) (24881.5 vertices/s)
(56,43,15,0) (24928.1 vertices/s)
(57,44,17,0) (25160.1 vertices/s)
(58,45,16,0) (25338.1 vertices/s)
(59,46,17,0) (25612.4 vertices/s)
(60,47,17,0) (25812.8 vertices/s)
(61,48,16,0) (26031.4 vertices/s)
(62,49,16,0) (26038.4 vertices/s)
(63,50,16,0) (26071.0 vertices/s)
(64,51,16,0) (26086.5 vertices/s)
(65,52,20,0) (26104.6 vertices/s)
(66,53,24,0) (26199.0 vertices/s)
(67,54,27,0) (26354.7 vertices/s)
(68,55,31,0) (26406.4 vertices/s)
(69,56,29,0) (26504.3 vertices/s)
(70,57,28,0) (26623.1 vertices/s)
(71,58,31,0) (26581.0 vertices/s)
(72,59,29,0) (26600.4 vertices/s)
(73,60,25,0) (26869.3 vertices/s)
(74,61,24,0) (26943.2 vertices/s)
(75,62,24,0) (27015.0 vertices/s)
(76,63,30,0) (27073.9 vertices/s)
(77,64,26,0) (27246.8 vertices/s)
(78,65,27,0) (27219.4 vertices/s)
(79,66,24,0) (27432.8 vertices/s)
(80,67,28,0) (27470.0 vertices/s)
(81,68,26,0) (27586.1 vertices/s)
(82,69,30,0) (27612.5 vertices/s)
(83,70,30,0) (27669.5 vertices/s)
(84,71,32,0) (27671.0 vertices/s)
(85,72,29,0) (27779.4 vertices/s)
(86,73,29,0) (27872.9 vertices/s)
(87,74,29,0) (27999.9 vertices/s)
(88,75,33,0) (28039.3 vertices/s)
(89,76,30,0) (28212.0 vertices/s)
(90,77,26,0) (28372.3 vertices/s)
(91,78,27,0) (28426.1 vertices/s)
(92,79,27,0) (28510.6 vertices/s)
(93,80,31,0) (28552.1 vertices/s)
(94,81,30,0) (28633.7 vertices/s)
(95,82,30,0) (28752.1 vertices/s)
(96,83,29,0) (28830.4 vertices/s)
(97,84,29,0) (28935.7 vertices/s)
(98,85,29,0) (28980.3 vertices/s)
(99,86,30,0) (29026.3 vertices/s)
(100,87,33,0) (28837.1 vertices/s)
(101,88,30,0) (28786.4 vertices/s)
(102,89,32,0) (28879.2 vertices/s)
(103,90,33,0) (28910.7 vertices/s)
(104,91,33,0) (29018.6 vertices/s)
(105,92,33,0) (29032.9 vertices/s)
(106,93,31,0) (29135.8 vertices/s)
(107,94,31,0) (29148.6 vertices/s)
(108,95,31,0) (29150.6 vertices/s)
(109,96,32,0) (29179.9 vertices/s)
(110,97,31,0) (29297.0 vertices/s)
(111,98,33,0) (29368.5 vertices/s)
(112,99,29,0) (29378.5 vertices/s)
(113,100,29,0) (29446.1 vertices/s)
(114,101,31,0) (29558.0 vertices/s)
(115,102,30,0) (29635.6 vertices/s)
(116,103,29,0) (29616.3 vertices/s)
(117,104,28,0) (29682.1 vertices/s)
(118,105,26,0) (29754.9 vertices/s)
(119,106,25,0) (29860.7 vertices/s)
(120,107,24,0) (29806.1 vertices/s)
(121,108,27,0) (29793.8 vertices/s)
(122,109,26,0) (29871.2 vertices/s)
(123,110,27,0) (29827.6 vertices/s)
(124,111,30,0) (29767.8 vertices/s)
(125,112,27,0) (29850.8 vertices/s)
(126,113,26,0) (29917.7 vertices/s)
(127,114,21,0) (30008.2 vertices/s)
(128,115,20,0) (30129.6 vertices/s)
(129,116,20,0) (30115.1 vertices/s)
(130,117,20,0) (30124.8 vertices/s)
(131,118,20,0) (30132.6 vertices/s)
(132,119,19,0) (30218.7 vertices/s)
(133,120,18,0) (30265.6 vertices/s)
(134,121,22,0) (30302.8 vertices/s)
(135,122,18,0) (30334.0 vertices/s)
(136,123,21,0) (30386.3 vertices/s)
(137,124,25,0) (30363.3 vertices/s)
(138,125,21,0) (30421.7 vertices/s)
(139,126,21,0) (30486.4 vertices/s)
(140,127,20,0) (30603.0 vertices/s)
(141,128,24,0) (30638.1 vertices/s)
(142,129,22,0) (30598.0 vertices/s)
(143,130,22,0) (30617.1 vertices/s)
(144,131,21,0) (30637.6 vertices/s)
(145,132,20,0) (30676.4 vertices/s)
(146,133,20,0) (30723.3 vertices/s)
(147,134,19,0) (30749.4 vertices/s)
(148,135,18,0) (30808.6 vertices/s)
(149,136,18,0) (30777.2 vertices/s)
(150,137,17,0) (30800.8 vertices/s)
(151,138,16,0) (30880.0 vertices/s)
(152,139,15,0) (30876.4 vertices/s)
(153,140,15,0) (30885.9 vertices/s)
(154,141,14,0) (30895.2 vertices/s)
(155,142,15,0) (30870.8 vertices/s)
(156,143,16,0) (30826.2 vertices/s)
(157,144,18,0) (30856.2 vertices/s)
(158,145,20,0) (30746.9 vertices/s)
(159,146,19,0) (30769.2 vertices/s)
(160,147,17,0) (30772.7 vertices/s)
(161,148,16,0) (30750.3 vertices/s)
(162,149,13,0) (30817.7 vertices/s)
(163,150,12,0) (30839.0 vertices/s)
(164,151,11,0) (30842.0 vertices/s)
(165,152,13,0) (30831.5 vertices/s)
(166,153,16,0) (30803.5 vertices/s)
(167,154,13,0) (30739.2 vertices/s)
(168,155,13,0) (30735.5 vertices/s)
(169,156,12,0) (30744.8 vertices/s)
(170,157,11,0) (30767.0 vertices/s)
(171,158,10,0) (30806.1 vertices/s)
(172,159,9,0) (30844.7 vertices/s)
(173,160,10,0) (30745.8 vertices/s)
(174,161,9,0) (30743.6 vertices/s)
(175,162,11,0) (30752.5 vertices/s)
(176,163,12,0) (30796.0 vertices/s)
(177,164,13,0) (30828.0 vertices/s)
(178,165,9,0) (30899.7 vertices/s)
(179,166,11,0) (30970.8 vertices/s)
(180,167,9,0) (31007.0 vertices/s)
(181,168,7,0) (31071.7 vertices/s)
(182,169,7,0) (31089.4 vertices/s)
(183,170,6,0) (31044.6 vertices/s)
(184,171,5,0) (31058.2 vertices/s)
(185,172,6,0) (31047.5 vertices/s)
(186,173,6,0) (31064.9 vertices/s)
(187,174,6,0) (31088.8 vertices/s)
(188,175,7,0) (31111.1 vertices/s)
(189,176,6,0) (31134.4 vertices/s)
(190,177,7,0) (31162.8 vertices/s)
(191,178,7,0) (31206.6 vertices/s)
(192,179,5,0) (31239.6 vertices/s)
(193,180,5,0) (31309.9 vertices/s)
(194,181,5,0) (31359.0 vertices/s)
(195,182,4,0) (31440.0 vertices/s)
(196,183,5,0) (31395.5 vertices/s)
(197,184,4,0) (31361.8 vertices/s)
(198,185,4,0) (31388.1 vertices/s)
(199,186,3,0) (31429.4 vertices/s)
(200,187,2,0) (31440.0 vertices/s)
(201,188,1,0) (31375.5 vertices/s)
(202,189,1,0) (31396.2 vertices/s)
(203,190,1,0) (31426.7 vertices/s)
(204,191,0,0) (31482.8 vertices/s)
Total refining surface time: 0.00607395s
Start volume scan...Scanning triangulation for bad cells (sequential)... 1139 cells scanned, done.
Number of bad cells: 93
end scan. [Bad tets:93]
Refining...
Legend of the following line: (#vertices,#steps,#facets to refine,#tets to refine)
(204,0,0,93)
(204,0,1,93) (0.0 vertices/s)
(205,1,2,94) (21290.9 vertices/s)
(206,2,2,95) (24105.2 vertices/s)
(207,3,1,96) (27060.0 vertices/s)
(208,4,0,94) (29852.7 vertices/s)
(208,5,1,94) (34721.1 vertices/s)
(209,6,0,93) (36846.0 vertices/s)
(209,7,1,93) (40721.4 vertices/s)
(210,8,0,92) (39803.6 vertices/s)
(210,9,1,92) (42653.9 vertices/s)
(211,10,2,93) (40175.3 vertices/s)
(212,11,0,94) (39165.8 vertices/s)
(212,12,1,94) (41391.2 vertices/s)
(213,13,3,98) (40509.6 vertices/s)
(214,14,1,96) (39891.5 vertices/s)
(215,15,0,96) (39272.5 vertices/s)
(215,16,1,96) (40820.5 vertices/s)
(216,17,2,96) (39723.2 vertices/s)
(217,18,6,99) (38460.3 vertices/s)
(218,19,1,97) (38075.4 vertices/s)
(219,20,2,98) (37888.9 vertices/s)
(220,21,0,95) (38614.8 vertices/s)
(220,22,1,95) (39722.2 vertices/s)
(221,23,0,94) (39731.9 vertices/s)
(221,24,1,94) (40754.4 vertices/s)
(222,25,0,94) (40454.3 vertices/s)
(222,26,1,94) (41401.6 vertices/s)
(223,27,2,94) (40546.4 vertices/s)
(224,28,1,95) (40288.3 vertices/s)
(225,29,0,94) (39998.3 vertices/s)
(225,30,1,94) (40880.2 vertices/s)
(226,31,0,95) (40056.5 vertices/s)
(226,32,1,95) (40820.5 vertices/s)
(227,33,3,93) (40542.5 vertices/s)
(228,34,2,96) (40524.7 vertices/s)
(229,35,1,97) (40186.3 vertices/s)
(230,36,1,97) (40051.7 vertices/s)
(231,37,0,97) (40090.2 vertices/s)
(231,38,1,97) (40773.5 vertices/s)
(232,39,1,97) (40329.8 vertices/s)
(233,40,0,98) (39841.4 vertices/s)
(233,41,1,98) (40434.2 vertices/s)
(234,42,6,99) (39471.4 vertices/s)
(235,43,2,97) (39310.2 vertices/s)
(236,44,3,98) (38803.5 vertices/s)
(237,45,5,100) (37915.6 vertices/s)
(238,46,2,100) (37551.2 vertices/s)
(239,47,6,100) (37215.8 vertices/s)
(240,48,2,99) (37179.4 vertices/s)
(241,49,1,96) (36897.8 vertices/s)
(242,50,0,96) (36657.1 vertices/s)
(242,51,1,96) (37143.5 vertices/s)
(243,52,1,96) (37250.9 vertices/s)
(244,53,0,95) (37323.4 vertices/s)
(244,54,1,95) (37792.8 vertices/s)
(245,55,4,98) (37644.7 vertices/s)
(246,56,1,99) (37158.8 vertices/s)
(247,57,1,99) (37279.8 vertices/s)
(248,58,2,97) (37374.3 vertices/s)
(249,59,0,97) (37483.2 vertices/s)
(249,60,1,97) (37906.0 vertices/s)
(250,61,1,94) (37493.0 vertices/s)
(251,62,1,93) (37438.4 vertices/s)
(252,63,2,95) (37259.0 vertices/s)
(253,64,2,95) (37169.1 vertices/s)
(254,65,1,92) (37122.8 vertices/s)
(255,66,0,92) (36998.7 vertices/s)
(255,67,1,92) (37305.0 vertices/s)
(256,68,0,88) (36796.9 vertices/s)
(256,69,1,88) (37136.8 vertices/s)
(257,70,2,90) (36879.9 vertices/s)
(258,71,0,91) (36674.3 vertices/s)
(258,72,1,91) (36999.5 vertices/s)
(259,73,0,90) (36889.7 vertices/s)
(259,74,1,90) (37206.7 vertices/s)
(260,75,0,90) (37113.4 vertices/s)
(260,76,1,90) (37422.8 vertices/s)
(261,77,0,89) (37306.4 vertices/s)
(261,78,1,89) (37625.7 vertices/s)
(262,79,0,87) (37512.7 vertices/s)
(262,80,1,87) (37807.8 vertices/s)
(263,81,0,82) (37519.5 vertices/s)
(263,82,1,82) (37840.6 vertices/s)
(264,83,0,85) (37676.1 vertices/s)
(264,84,1,85) (37994.3 vertices/s)
(265,85,0,86) (37931.3 vertices/s)
(265,86,1,86) (38206.8 vertices/s)
(266,87,2,85) (38278.0 vertices/s)
(267,88,1,85) (38280.3 vertices/s)
(268,89,0,83) (38361.2 vertices/s)
(268,90,1,83) (38609.7 vertices/s)
(269,91,0,83) (38445.0 vertices/s)
(270,92,0,75) (38560.6 vertices/s)
(270,93,1,75) (38816.8 vertices/s)
(271,94,0,73) (38668.6 vertices/s)
(271,95,1,73) (38934.8 vertices/s)
(272,96,1,74) (38806.2 vertices/s)
(273,97,1,75) (38585.7 vertices/s)
(274,98,1,77) (38386.4 vertices/s)
(275,99,0,73) (38401.6 vertices/s)
(275,100,1,73) (38639.4 vertices/s)
(276,101,2,72) (38403.1 vertices/s)
(277,102,1,72) (38201.5 vertices/s)
(278,103,0,71) (38119.9 vertices/s)
(278,104,1,71) (38361.4 vertices/s)
(279,105,0,70) (38239.3 vertices/s)
(279,106,1,70) (38463.2 vertices/s)
(280,107,0,69) (38364.7 vertices/s)
(280,108,1,69) (38601.2 vertices/s)
(281,109,2,69) (38434.6 vertices/s)
(282,110,1,68) (38301.0 vertices/s)
(283,111,0,69) (38105.1 vertices/s)
(283,112,1,69) (38344.8 vertices/s)
(284,113,0,65) (38346.0 vertices/s)
(284,114,1,65) (38554.3 vertices/s)
(285,115,1,65) (38501.4 vertices/s)
(286,116,1,65) (38449.4 vertices/s)
(287,117,0,65) (38347.5 vertices/s)
(287,118,1,65) (38548.8 vertices/s)
(288,119,0,64) (38485.8 vertices/s)
(288,120,1,64) (38698.8 vertices/s)
(289,121,0,61) (38620.4 vertices/s)
(289,122,1,61) (38830.3 vertices/s)
(290,123,1,59) (38812.8 vertices/s)
(291,124,0,59) (38871.0 vertices/s)
(291,125,1,59) (39076.4 vertices/s)
(292,126,0,58) (38973.6 vertices/s)
(292,127,1,58) (39161.6 vertices/s)
(293,128,0,57) (39085.0 vertices/s)
(293,129,1,57) (39281.6 vertices/s)
(294,130,0,55) (39286.7 vertices/s)
(294,131,1,55) (39495.0 vertices/s)
(295,132,0,54) (39557.6 vertices/s)
(295,133,1,54) (39738.0 vertices/s)
(296,134,0,52) (39669.4 vertices/s)
(296,135,1,52) (39847.4 vertices/s)
(297,136,1,47) (39720.4 vertices/s)
(298,137,0,45) (39686.4 vertices/s)
(298,138,1,45) (39874.2 vertices/s)
(299,139,0,44) (39850.2 vertices/s)
(299,140,1,44) (40035.6 vertices/s)
(300,141,3,42) (39921.5 vertices/s)
(301,142,0,42) (39900.3 vertices/s)
(301,143,1,42) (40079.2 vertices/s)
(302,144,0,42) (40022.5 vertices/s)
(302,145,1,42) (40199.2 vertices/s)
(303,146,0,42) (40121.1 vertices/s)
(303,147,1,42) (40295.6 vertices/s)
(304,148,1,37) (40098.0 vertices/s)
(305,149,1,37) (40043.0 vertices/s)
(306,150,1,37) (39905.2 vertices/s)
(307,151,1,37) (39852.8 vertices/s)
(308,152,0,37) (39833.4 vertices/s)
(308,153,1,37) (40000.5 vertices/s)
(309,154,0,37) (39980.4 vertices/s)
(309,155,1,37) (40145.6 vertices/s)
(310,156,5,36) (40063.2 vertices/s)
(311,157,1,35) (39899.8 vertices/s)
(312,158,0,34) (39919.3 vertices/s)
(312,159,1,34) (40080.2 vertices/s)
(313,160,5,33) (40031.5 vertices/s)
(314,161,1,34) (39811.5 vertices/s)
(315,162,0,34) (39754.1 vertices/s)
(315,163,1,34) (39922.4 vertices/s)
(316,164,0,32) (39941.1 vertices/s)
(316,165,1,32) (40098.5 vertices/s)
(317,166,0,30) (40088.4 vertices/s)
(317,167,1,30) (40241.8 vertices/s)
(318,168,0,26) (40125.5 vertices/s)
(318,169,1,26) (40277.1 vertices/s)
(319,170,0,25) (40094.0 vertices/s)
(319,171,1,25) (40246.1 vertices/s)
(320,172,0,22) (40150.3 vertices/s)
(320,173,1,22) (40298.5 vertices/s)
(321,174,1,20) (40249.8 vertices/s)
(322,175,0,20) (40146.8 vertices/s)
(322,176,1,20) (40303.4 vertices/s)
(323,177,0,20) (40235.9 vertices/s)
(323,178,1,20) (40382.2 vertices/s)
(324,179,1,17) (40243.4 vertices/s)
(325,180,0,19) (40045.3 vertices/s)
(325,181,1,19) (40169.8 vertices/s)
(326,182,0,18) (40115.8 vertices/s)
(326,183,1,18) (40247.4 vertices/s)
(327,184,0,17) (40218.5 vertices/s)
(327,185,1,17) (40367.6 vertices/s)
(328,186,1,16) (40390.4 vertices/s)
(329,187,0,16) (40267.7 vertices/s)
(329,188,1,16) (40396.0 vertices/s)
(330,189,0,14) (40229.6 vertices/s)
(330,190,1,14) (40366.6 vertices/s)
(331,191,0,11) (40346.1 vertices/s)
(331,192,1,11) (40471.7 vertices/s)
(332,193,0,10) (40360.0 vertices/s)
(332,194,1,10) (40502.5 vertices/s)
(333,195,0,10) (40499.6 vertices/s)
(333,196,1,10) (40630.8 vertices/s)
(334,197,0,9) (40611.3 vertices/s)
(334,198,1,9) (40741.4 vertices/s)
(335,199,0,8) (40646.0 vertices/s)
(335,200,1,8) (40766.9 vertices/s)
(336,201,1,8) (40605.7 vertices/s)
(337,202,2,8) (40489.8 vertices/s)
(338,203,0,9) (40366.2 vertices/s)
(338,204,1,9) (40476.7 vertices/s)
(339,205,0,8) (40354.5 vertices/s)
(339,206,1,8) (40463.9 vertices/s)
(340,207,0,7) (40328.0 vertices/s)
(340,208,1,7) (40436.4 vertices/s)
(341,209,2,6) (40309.4 vertices/s)
(342,210,0,6) (40215.7 vertices/s)
(342,211,1,6) (40337.2 vertices/s)
(343,212,0,3) (40265.9 vertices/s)
(343,213,1,3) (40386.4 vertices/s)
(344,214,0,2) (40355.2 vertices/s)
(344,215,1,2) (40467.4 vertices/s)
(345,216,0,0) (40427.0 vertices/s)
Total refining volume time: 0.00534987s
Total refining time: 0.011575s
Running sliver perturbation...
Legend of the following line: (#vertices in pqueue, #iterations, #fails)
bound 3: (0,1,0) (4544.2 iteration/s)
bound 6: (0,1,0) (2716.5 iteration/s)
bound 8: (8,1,0) (1490.0 iteration/s)
bound 8: (4,2,0) (2146.0 iteration/s)
bound 8: (0,3,0) (2579.5 iteration/s)
bound 11: (4,1,0) (717.3 iteration/s)
bound 11: (0,2,0) (1159.4 iteration/s)
bound 12: (4,1,0) (477.3 iteration/s)
bound 12: (0,2,0) (826.8 iteration/s)
Total perturbation time: 0.00243592s
Perturbation statistics:
Sq radius gradient perturbation: 100% (7 in 1ms)
Volume gradient perturbation: 0% (0 in 0ms)
Dihedral angle gradient perturbation: 0% (0 in 0ms)
Li random perturbation [100, 0.15, in]: 0% (0 in 0ms)
Perturbation return code: BOUND_REACHED
Exuding...
Legend of the following line: (#cells left,#vertices pumped,#vertices ignored)
(378,0,0)
(377,0,4) (97542.0 vertices/s)
(376,1,5) (39053.1 vertices/s)
(375,1,9) (41665.3 vertices/s)
(374,1,13) (49034.1 vertices/s)
(373,2,15) (37404.6 vertices/s)
(372,3,16) (33058.6 vertices/s)
(371,3,20) (35971.7 vertices/s)
(370,4,20) (31944.4 vertices/s)
(370,5,20) (28407.1 vertices/s)
(369,6,20) (25803.2 vertices/s)
(368,6,24) (28468.1 vertices/s)
(367,6,28) (30905.4 vertices/s)
(366,6,32) (33647.0 vertices/s)
(365,7,32) (31499.1 vertices/s)
(364,7,36) (33547.0 vertices/s)
(363,7,40) (35462.3 vertices/s)
(362,7,44) (37510.0 vertices/s)
(361,8,45) (35378.4 vertices/s)
(360,9,45) (33406.0 vertices/s)
(360,10,45) (31578.3 vertices/s)
(359,10,49) (32601.7 vertices/s)
(358,10,53) (33712.2 vertices/s)
(357,10,57) (34503.6 vertices/s)
(356,10,61) (35421.9 vertices/s)
(355,10,65) (36331.3 vertices/s)
(354,10,69) (37174.9 vertices/s)
(353,10,73) (37805.2 vertices/s)
(352,10,77) (38808.1 vertices/s)
(351,10,81) (39490.7 vertices/s)
(350,10,85) (39829.7 vertices/s)
(349,10,89) (40637.2 vertices/s)
(348,11,89) (39397.7 vertices/s)
(347,11,93) (39896.7 vertices/s)
(346,12,93) (38430.6 vertices/s)
(345,12,97) (39191.6 vertices/s)
(344,12,101) (40062.9 vertices/s)
(343,12,105) (40541.5 vertices/s)
(343,13,106) (39537.2 vertices/s)
(342,13,110) (40074.1 vertices/s)
(341,13,114) (40583.1 vertices/s)
(340,13,118) (41158.2 vertices/s)
(339,13,122) (41524.4 vertices/s)
(338,13,126) (42056.5 vertices/s)
(337,13,130) (42137.5 vertices/s)
(336,13,134) (42578.5 vertices/s)
(335,13,138) (42868.8 vertices/s)
(334,13,142) (43096.3 vertices/s)
(333,13,146) (43411.9 vertices/s)
(332,13,150) (43605.9 vertices/s)
(331,13,154) (43937.3 vertices/s)
(330,13,158) (44393.1 vertices/s)
(329,13,162) (44838.1 vertices/s)
(328,13,166) (45024.2 vertices/s)
(327,13,170) (45418.9 vertices/s)
(326,13,174) (45790.5 vertices/s)
(325,13,178) (46245.4 vertices/s)
(324,13,182) (46498.3 vertices/s)
(323,13,186) (46661.9 vertices/s)
(322,13,190) (47040.8 vertices/s)
(321,13,194) (47042.5 vertices/s)
(320,13,198) (47175.2 vertices/s)
(319,13,202) (47140.1 vertices/s)
(318,13,206) (47411.5 vertices/s)
(317,13,210) (47781.5 vertices/s)
(316,13,214) (47896.5 vertices/s)
(315,13,218) (48114.0 vertices/s)
(314,13,222) (48481.5 vertices/s)
(313,13,226) (48695.8 vertices/s)
(312,13,230) (49028.8 vertices/s)
(311,13,234) (49243.2 vertices/s)
(310,13,238) (49449.9 vertices/s)
(309,13,242) (49661.0 vertices/s)
(308,13,246) (49778.0 vertices/s)
(307,13,250) (49939.3 vertices/s)
(306,13,254) (50167.3 vertices/s)
(305,13,258) (50489.0 vertices/s)
(304,13,262) (50675.9 vertices/s)
(303,13,266) (50851.6 vertices/s)
(302,13,270) (50885.7 vertices/s)
(301,13,274) (51052.3 vertices/s)
(300,13,278) (51047.0 vertices/s)
(299,13,282) (51189.9 vertices/s)
(298,13,286) (51134.8 vertices/s)
(297,13,290) (51301.1 vertices/s)
(296,13,294) (51371.7 vertices/s)
(295,13,298) (51521.1 vertices/s)
(294,13,302) (51587.5 vertices/s)
(293,13,306) (51766.9 vertices/s)
(292,13,310) (51978.2 vertices/s)
(291,13,314) (52115.5 vertices/s)
(290,13,318) (52070.6 vertices/s)
(289,13,322) (52145.4 vertices/s)
(288,13,326) (52360.5 vertices/s)
(287,13,330) (52580.2 vertices/s)
(286,13,334) (52540.9 vertices/s)
(285,13,338) (51887.7 vertices/s)
(284,13,342) (51927.7 vertices/s)
(283,13,346) (51888.9 vertices/s)
(282,13,350) (51860.2 vertices/s)
(281,13,354) (51921.0 vertices/s)
(280,13,358) (51996.7 vertices/s)
(279,13,362) (52033.5 vertices/s)
(278,13,366) (52129.7 vertices/s)
(277,13,370) (52222.4 vertices/s)
(276,13,374) (52285.5 vertices/s)
(275,13,378) (52376.8 vertices/s)
(276,14,381) (51985.4 vertices/s)
(275,14,385) (52119.1 vertices/s)
(274,14,389) (52222.4 vertices/s)
(273,14,393) (52385.5 vertices/s)
(272,14,397) (52491.6 vertices/s)
(271,14,401) (52645.4 vertices/s)
(270,14,405) (52693.9 vertices/s)
(269,14,409) (52767.5 vertices/s)
(268,14,413) (52786.7 vertices/s)
(267,14,417) (52791.2 vertices/s)
(266,14,421) (52822.5 vertices/s)
(265,14,425) (52846.9 vertices/s)
(264,14,429) (52911.4 vertices/s)
(263,14,433) (52947.0 vertices/s)
(262,14,437) (53066.3 vertices/s)
(261,14,441) (53305.1 vertices/s)
(260,14,445) (53543.3 vertices/s)
(259,14,449) (53630.7 vertices/s)
(258,14,453) (53800.5 vertices/s)
(257,14,457) (53968.4 vertices/s)
(256,14,461) (54037.6 vertices/s)
(255,14,465) (54164.4 vertices/s)
(254,14,469) (54370.6 vertices/s)
(253,14,473) (54574.9 vertices/s)
(252,14,477) (54820.7 vertices/s)
(251,14,481) (55015.1 vertices/s)
(250,14,485) (55156.8 vertices/s)
(249,14,489) (55335.6 vertices/s)
(248,14,493) (55524.6 vertices/s)
(247,14,497) (55655.3 vertices/s)
(246,14,501) (55840.8 vertices/s)
(245,14,505) (55999.5 vertices/s)
(244,14,509) (56199.3 vertices/s)
(243,14,513) (56385.7 vertices/s)
(242,14,517) (56570.4 vertices/s)
(241,14,521) (56778.5 vertices/s)
(240,14,525) (56948.0 vertices/s)
(239,14,529) (57016.2 vertices/s)
(238,14,533) (57169.7 vertices/s)
(237,14,537) (57407.9 vertices/s)
(236,14,541) (57602.1 vertices/s)
(235,14,545) (57781.5 vertices/s)
(234,14,549) (57966.8 vertices/s)
(233,14,553) (58168.0 vertices/s)
(232,14,557) (58337.1 vertices/s)
(231,14,561) (58504.8 vertices/s)
(230,14,565) (58743.8 vertices/s)
(229,14,569) (58872.1 vertices/s)
(228,14,573) (59035.5 vertices/s)
(227,14,577) (59222.2 vertices/s)
(226,14,581) (59400.1 vertices/s)
(225,14,585) (59627.4 vertices/s)
(224,14,589) (59808.4 vertices/s)
(223,14,593) (59959.1 vertices/s)
(222,14,597) (60102.7 vertices/s)
(221,14,601) (60237.9 vertices/s)
(220,14,605) (60409.3 vertices/s)
(219,14,609) (60512.0 vertices/s)
(218,14,613) (60656.5 vertices/s)
(217,14,617) (60865.6 vertices/s)
(216,14,621) (60977.7 vertices/s)
(215,14,625) (61167.1 vertices/s)
(214,14,629) (61305.4 vertices/s)
(213,14,633) (61498.1 vertices/s)
(212,14,637) (61665.3 vertices/s)
(211,14,641) (61801.5 vertices/s)
(210,14,645) (61929.4 vertices/s)
(209,14,649) (62104.5 vertices/s)
(208,14,653) (62250.1 vertices/s)
(207,14,657) (62387.5 vertices/s)
(206,14,661) (62511.2 vertices/s)
(205,14,665) (62700.1 vertices/s)
(204,14,669) (62740.1 vertices/s)
(203,14,673) (62856.7 vertices/s)
(202,14,677) (62976.4 vertices/s)
(201,14,681) (63078.5 vertices/s)
(200,14,685) (63203.4 vertices/s)
(199,14,689) (63338.3 vertices/s)
(198,14,693) (63473.7 vertices/s)
(197,14,697) (63623.4 vertices/s)
(196,14,701) (63762.3 vertices/s)
(195,14,705) (63788.6 vertices/s)
(194,14,709) (63936.7 vertices/s)
(193,14,713) (64020.8 vertices/s)
(192,14,717) (64104.1 vertices/s)
(191,14,721) (64276.8 vertices/s)
(190,14,725) (64363.9 vertices/s)
(189,14,729) (64530.5 vertices/s)
(188,14,733) (64683.9 vertices/s)
(187,14,737) (64796.9 vertices/s)
(186,14,741) (64932.2 vertices/s)
(185,14,745) (65031.4 vertices/s)
(184,14,749) (65152.7 vertices/s)
(183,14,753) (65267.8 vertices/s)
(182,14,757) (65382.2 vertices/s)
(181,14,761) (65518.5 vertices/s)
(180,14,765) (65592.3 vertices/s)
(179,14,769) (65686.8 vertices/s)
(178,14,773) (65838.1 vertices/s)
(177,14,777) (65941.8 vertices/s)
(176,14,781) (66068.7 vertices/s)
(175,14,785) (66138.9 vertices/s)
(174,14,789) (66229.8 vertices/s)
(173,14,793) (66398.2 vertices/s)
(172,14,797) (66471.7 vertices/s)
(171,14,801) (66589.4 vertices/s)
(170,14,805) (66738.1 vertices/s)
(169,14,809) (66870.1 vertices/s)
(168,14,813) (66958.0 vertices/s)
(167,14,817) (67005.9 vertices/s)
(166,14,821) (67103.0 vertices/s)
(165,14,825) (67198.1 vertices/s)
(164,14,829) (67239.3 vertices/s)
(163,14,833) (67400.8 vertices/s)
(162,14,837) (67500.5 vertices/s)
(161,14,841) (67620.2 vertices/s)
(160,14,845) (67697.9 vertices/s)
(159,14,849) (67805.9 vertices/s)
(158,14,853) (67897.8 vertices/s)
(157,14,857) (68026.5 vertices/s)
(156,14,861) (68117.1 vertices/s)
(155,14,865) (68195.6 vertices/s)
(154,14,869) (68360.6 vertices/s)
(153,14,873) (68449.3 vertices/s)
(152,14,877) (68557.9 vertices/s)
(151,14,881) (68688.8 vertices/s)
(150,14,885) (68807.5 vertices/s)
(149,14,889) (68850.5 vertices/s)
(148,14,893) (68962.9 vertices/s)
(147,14,897) (69026.5 vertices/s)
(146,14,901) (69132.7 vertices/s)
(145,14,905) (69263.5 vertices/s)
(144,14,909) (69336.8 vertices/s)
(143,14,913) (69444.9 vertices/s)
(142,14,917) (69601.5 vertices/s)
(141,14,921) (69651.7 vertices/s)
(140,14,925) (69731.7 vertices/s)
(139,14,929) (69859.9 vertices/s)
(138,14,933) (69955.1 vertices/s)
(137,14,937) (69988.5 vertices/s)
(136,14,941) (70046.5 vertices/s)
(135,14,945) (70130.2 vertices/s)
(134,14,949) (70223.3 vertices/s)
(133,14,953) (70317.0 vertices/s)
(132,14,957) (70414.0 vertices/s)
(131,14,961) (70521.5 vertices/s)
(130,14,965) (70607.5 vertices/s)
(129,14,969) (70699.1 vertices/s)
(128,14,973) (70805.0 vertices/s)
(127,14,977) (70858.8 vertices/s)
(126,14,981) (70896.4 vertices/s)
(125,14,985) (70930.0 vertices/s)
(124,14,989) (71002.3 vertices/s)
(123,14,993) (71070.5 vertices/s)
(122,14,997) (71152.8 vertices/s)
(121,14,1001) (71250.4 vertices/s)
(120,14,1005) (71347.4 vertices/s)
(119,14,1009) (71287.5 vertices/s)
(118,14,1013) (71297.5 vertices/s)
(117,14,1017) (71363.4 vertices/s)
(116,14,1021) (71418.3 vertices/s)
(115,14,1025) (71533.5 vertices/s)
(114,14,1029) (71592.3 vertices/s)
(113,14,1033) (71641.2 vertices/s)
(112,14,1037) (71754.9 vertices/s)
(111,14,1041) (71852.7 vertices/s)
(110,14,1045) (71939.3 vertices/s)
(109,14,1049) (72017.3 vertices/s)
(108,14,1053) (72084.3 vertices/s)
(107,14,1057) (72194.3 vertices/s)
(106,14,1061) (72269.8 vertices/s)
(105,14,1065) (72335.5 vertices/s)
(104,14,1069) (72385.6 vertices/s)
(103,14,1073) (72387.7 vertices/s)
(102,14,1077) (72397.9 vertices/s)
(101,14,1081) (72467.0 vertices/s)
(100,14,1085) (72565.8 vertices/s)
(99,14,1089) (72643.4 vertices/s)
(98,14,1093) (72677.9 vertices/s)
(97,14,1097) (72765.0 vertices/s)
(96,14,1101) (72764.5 vertices/s)
(95,14,1105) (72792.6 vertices/s)
(94,14,1109) (72869.6 vertices/s)
(93,14,1113) (72916.5 vertices/s)
(92,14,1117) (73039.4 vertices/s)
(91,14,1121) (73076.5 vertices/s)
(90,14,1125) (73199.6 vertices/s)
(89,14,1129) (73268.9 vertices/s)
(88,14,1133) (73357.1 vertices/s)
(87,14,1137) (73440.4 vertices/s)
(86,14,1141) (73560.5 vertices/s)
(85,14,1145) (73647.5 vertices/s)
(84,14,1149) (73734.0 vertices/s)
(83,14,1153) (73863.0 vertices/s)
(82,14,1157) (73948.7 vertices/s)
(81,14,1161) (73977.8 vertices/s)
(80,14,1165) (74085.2 vertices/s)
(79,14,1169) (74156.3 vertices/s)
(78,14,1173) (74231.5 vertices/s)
(77,14,1177) (74262.8 vertices/s)
(76,14,1181) (74365.3 vertices/s)
(75,14,1185) (74416.1 vertices/s)
(74,14,1189) (74531.1 vertices/s)
(73,14,1193) (74562.3 vertices/s)
(72,14,1197) (74672.1 vertices/s)
(71,14,1201) (74763.8 vertices/s)
(70,14,1205) (74803.0 vertices/s)
(69,14,1209) (74846.3 vertices/s)
(68,14,1213) (74880.6 vertices/s)
(67,14,1217) (74914.8 vertices/s)
(66,14,1221) (74971.7 vertices/s)
(65,14,1225) (75010.9 vertices/s)
(64,14,1229) (75080.5 vertices/s)
(63,14,1233) (75173.7 vertices/s)
(62,14,1237) (75192.4 vertices/s)
(61,14,1241) (75235.0 vertices/s)
(60,14,1245) (75290.3 vertices/s)
(59,14,1249) (75331.2 vertices/s)
(58,14,1253) (75422.8 vertices/s)
(57,14,1257) (75481.6 vertices/s)
(56,14,1261) (75585.5 vertices/s)
(55,14,1265) (75603.7 vertices/s)
(54,14,1269) (75639.1 vertices/s)
(53,14,1273) (75674.3 vertices/s)
(52,14,1277) (75749.9 vertices/s)
(51,14,1281) (75790.0 vertices/s)
(50,14,1285) (75833.0 vertices/s)
(49,14,1289) (75881.2 vertices/s)
(48,14,1293) (75942.9 vertices/s)
(47,14,1297) (76025.6 vertices/s)
(46,14,1301) (76113.3 vertices/s)
(45,14,1305) (76160.3 vertices/s)
(44,14,1309) (76228.3 vertices/s)
(43,14,1313) (76310.8 vertices/s)
(42,14,1317) (76343.3 vertices/s)
(41,14,1321) (76393.6 vertices/s)
(40,14,1325) (76457.3 vertices/s)
(39,14,1329) (76533.4 vertices/s)
(38,14,1333) (76587.0 vertices/s)
(37,14,1337) (76627.9 vertices/s)
(36,14,1341) (76689.5 vertices/s)
(35,14,1345) (76738.3 vertices/s)
(34,14,1349) (76778.5 vertices/s)
(33,14,1353) (76852.9 vertices/s)
(32,14,1357) (76922.8 vertices/s)
(31,14,1361) (76931.2 vertices/s)
(30,14,1365) (76887.9 vertices/s)
(29,14,1369) (76984.0 vertices/s)
(28,14,1373) (77039.5 vertices/s)
(27,14,1377) (77078.3 vertices/s)
(26,14,1381) (77155.0 vertices/s)
(25,14,1385) (77205.7 vertices/s)
(24,14,1389) (77265.4 vertices/s)
(23,14,1393) (77324.8 vertices/s)
(22,14,1397) (77378.9 vertices/s)
(21,14,1401) (77432.8 vertices/s)
(20,14,1405) (77530.3 vertices/s)
(19,14,1409) (77626.4 vertices/s)
(18,14,1413) (77616.3 vertices/s)
(17,14,1417) (77686.5 vertices/s)
(16,14,1421) (77748.3 vertices/s)
(15,14,1425) (77787.5 vertices/s)
(14,14,1429) (77832.7 vertices/s)
(13,14,1433) (77892.7 vertices/s)
(12,14,1437) (77991.9 vertices/s)
(11,14,1441) (77925.7 vertices/s)
(10,14,1445) (77981.1 vertices/s)
(9,14,1449) (78016.3 vertices/s)
(8,14,1453) (78000.4 vertices/s)
(7,14,1457) (78077.2 vertices/s)
(6,14,1461) (78157.7 vertices/s)
(5,14,1465) (78179.3 vertices/s)
(4,14,1469) (78242.4 vertices/s)
(3,14,1473) (78313.2 vertices/s)
(2,14,1477) (78364.0 vertices/s)
(1,14,1481) (78434.3 vertices/s)
(0,14,1485) (78525.2 vertices/s)
Total exuding time: 0.0189161s
Exuding return code: CANT_IMPROVE_ANYMORE
Uncomment the following to write the mesh of tetrahedra to disk
# mesh.write(filename_mesh_out)
Visualize the mesh with itkwidgtes
# reader = vtk.vtkUnstructuredGridReader()
# reader.SetFileName(filename_mesh_out)
# reader.Update()
# grid = reader.GetOutput()
# view(geometries=grid)
Write Abaqus input FE files
Generate tetrahedra-FE model with constant material properties
The method ciclope.tetraFE.mesh2tetrafe
assumes that the material property definition is contained in the FE analysis .INP
template file.
input_template = "./../../input_templates/tmp_example01_comp_static_bone.inp"
Inspect input template file
!cat {input_template} # on linux
# !type ..\..\input_templates\tmp_example01_comp_static_bone.inp # on Win
** ----------------------------------------------------
** Material definition:
*MATERIAL, NAME=BONE
*ELASTIC
18000, 0.3, 0
*SOLID SECTION, ELSET=SET1, MATERIAL=BONE
** ---------------------------------------------------
** Analysis definition. The Step module defines the analysis steps and associated output requests.
** More info at:
** https://abaqus-docs.mit.edu/2017/English/SIMACAECAERefMap/simacae-m-Sim-sb.htm#simacae-m-Sim-sb
** ---------------------------------------------------
*STEP
*STATIC
** All displacements fixed at bottom:
*BOUNDARY
NODES_Z0, 1, 3, 0.0
** Vertical displacement imposed at top:
*BOUNDARY
NODES_Z1, 3, 3, -0.04
** ---------------------------------------------------
** Output request:
*NODE FILE, OUTPUT=2D
U
*NODE PRINT, TOTALS=ONLY, NSET=NODES_Z0
RF
*EL FILE, OUTPUT=2D
S, E
*END STEP
filename_out = './../../test_data/LHDL/3155_D_4_bc/results/3155_D_4_bc_tetraFE.inp'
Generate CalculiX FE input file
ciclope.core.tetraFE.mesh2tetrafe(mesh, input_template, filename_out)
Solve FE model in calculix
The following section assumes that you have the Calculix solver installed and accessible with the command ccx_2.17
or for multithread option ccx_2.17_MT
.
The multithread option in CalculiX is activated by defining the number of threads (default=1) used for the calculation.
This is done by setting the env variable OMP_NUM_THREADS with:
export OMP_NUM_THREADS=8
!export OMP_NUM_THREADS=8; ccx_2.17_MT "./../../test_data/LHDL/3155_D_4_bc/results/3155_D_4_bc_tetraFE"
************************************************************
CalculiX Version 2.17, Copyright(C) 1998-2020 Guido Dhondt
CalculiX comes with ABSOLUTELY NO WARRANTY. This is free
software, and you are welcome to redistribute it under
certain conditions, see gpl.htm
************************************************************
You are using an executable made on Sun 10 Jan 2021 11:34:19 AM CET
The numbers below are estimated upper bounds
number of:
nodes: 345
elements: 367
one-dimensional elements: 0
two-dimensional elements: 0
integration points per element: 1
degrees of freedom per node: 3
layers per element: 1
distributed facial loads: 0
distributed volumetric loads: 0
concentrated loads: 0
single point constraints: 14
multiple point constraints: 1
terms in all multiple point constraints: 1
tie constraints: 0
dependent nodes tied by cyclic constraints: 0
dependent nodes in pre-tension constraints: 0
sets: 15
terms in all sets: 1140
materials: 1
constants per material and temperature: 2
temperature points per material: 1
plastic data points per material: 0
orientations: 0
amplitudes: 2
data points in all amplitudes: 2
print requests: 1
transformations: 0
property cards: 0
STEP 1
Static analysis was selected
Decascading the MPC's
Determining the structure of the matrix:
number of equations
803
number of nonzero lower triangular matrix elements
9322
Using up to 8 cpu(s) for the stress calculation.
Using up to 8 cpu(s) for the symmetric stiffness/mass contributions.
Factoring the system of equations using the symmetric spooles solver
Using up to 8 cpu(s) for spooles.
Using up to 8 cpu(s) for the stress calculation.
Job finished
________________________________________
Total CalculiX Time: 0.024089
________________________________________
Post-processing of FE results
Convert Calculix output to Paraview
The following sections assume that you have ccx2paraview and Paraview installed and working.
For more info visit:
filename_out_base, ext_out = os.path.splitext(filename_out)
ccx2paraview.Converter(filename_out_base + '.frd', ['vtk']).run()
Visualize results in Paraview
!paraview {filename_out_base + '.vtk'}
Post-process FE analysis results
Display the CalculiX FE output .DAT
file containing the sum of reaction forces:
filename_dat = filename_out_base + '.dat'
!cat {filename_dat}
total force (fx,fy,fz) for set NODES_Z0 and time 0.1000000E+01
0.000000E+00 0.000000E+00 0.000000E+00
Apparent elastic modulus
E_app = (F_tot / A) / epsilon
A = 200*0.0195*200*0.0195 # [mm2]
epsilon = 0.04/(200*0.0195)
E = (2.40e2/A)/epsilon
print(f"{E:.2f} MPa")
1538.46 MPa
Bone volume fraction
BVTV = 100*np.sum(L)/L.size
print(f"{BVTV:.2f} %")
24.71 %
Dependencies
import watermark
%load_ext watermark
%watermark
%watermark --iversions
The watermark extension is already loaded. To reload it, use:
%reload_ext watermark
Last updated: 2023-01-20T19:55:20.509940+03:00
Python implementation: CPython
Python version : 3.8.15
IPython version : 8.8.0
Compiler : GCC 10.4.0
OS : Linux
Release : 4.19.0-23-amd64
Machine : x86_64
Processor :
CPU cores : 8
Architecture: 64bit
ccx2paraview: 3.0.0
sys : 3.8.15 | packaged by conda-forge | (default, Nov 22 2022, 08:49:35)
[GCC 10.4.0]
ciclope : 1.3.0
meshio : 5.0.0
numpy : 1.24.1
scipy : 1.8.1
watermark : 2.3.1
skimage : 0.18.1
matplotlib : 3.5.2
Acknowledgements
This notebook was developed within Building the Jupyter Community in MSK Imaging Research, a Jupyter Community Workshop sponsored by NUMFocus