svVascularize

svv.tree.Tree

Core Class v1.0.0

The Tree class defines a branching tree structure that is used to abstract the physical representation of the generated vascular network. It implements stochastic constructive optimization algorithms for growing vascular trees within defined domains.

Trees grow iteratively by adding vessels using optimized bifurcation placement, managing collision detection, and maintaining physiological constraints. The class provides export capabilities for 3D models, centerlines, and simulation-ready geometries.

Quick Example

from svv.tree.tree import Tree
from svv.domain.domain import Domain
import pyvista as pv

# Create domain from mesh
mesh = pv.read('tissue_region.stl')
domain = Domain(mesh)
domain.create()
domain.solve()
domain.build()

# Initialize tree
tree = Tree()
tree.set_domain(domain)

# Set root position
tree.set_root([0, 0, 0])  # Starting point

# Grow tree with n terminals
tree.n_add(100)  # Add 100 vessels

# Visualize
tree.show()

# Export results
solid_model = tree.export_solid()
centerlines, polys = tree.export_centerlines()

Constructor

Tree(*, parameters=None, unit_system=None, preallocation_step=4000000)

Initialize an empty tree structure. The tree must be configured with a domain and root before growth.

Parameters

Parameter Type Default Description
parameters TreeParameters None Pre-configured parameters object. If None, creates default parameters.
unit_system UnitSystem None Unit system for physical quantities (CGS by default). Ignored if parameters is provided.
preallocation_step int 4,000,000 Number of vessel segments to preallocate. Larger values reduce reallocations during growth but use more memory.

Constructor Examples

from svv.tree.tree import Tree
from svv.tree.data.data import TreeParameters
from svv.tree.data.units import UnitSystem

# Default tree (CGS units)
tree = Tree()

# Tree with SI units
si_units = UnitSystem(length='m', mass='kg', time='s')
tree_si = Tree(unit_system=si_units)

# Tree with custom parameters
params = TreeParameters()
params.murray_exponent = 2.7
params.terminal_pressure = 60.0  # mmHg equivalent in CGS
tree_custom = Tree(parameters=params)

# Memory-optimized tree for smaller networks
tree_small = Tree(preallocation_step=10000)

Attributes

Core Data Structures

Attribute Type Description
data TreeData Main data structure containing vessel geometry and properties
parameters TreeParameters Configuration parameters for tree growth
vessel_map TreeMap Mapping of vessel connectivity and hierarchy
domain Domain Spatial domain constraining tree growth
connectivity np.ndarray Connectivity matrix for vessel segments

Tree Properties

Attribute Type Description
n_terminals int Current number of terminal vessels
segment_count int Total number of vessel segments
tree_scale float Total tree volume scale factor
characteristic_length float Characteristic length from domain
physical_clearance float Minimum physical spacing between vessels
convex bool Whether domain is treated as convex

Internal Structures

Attribute Type Description
kdtm KDTreeManager KD-tree for spatial queries
hnsw_tree USearchTree HNSW index for fast nearest neighbor search
times dict Performance timing statistics

TreeParameters

The TreeParameters class holds physical and algorithmic settings that steer tree generation. Parameters are stored in the unit system referenced by unit_system. By default, the centimetre-gram-second (CGS) convention is used.

TreeParameters(*, unit_system=None)

Attributes

Attribute Type Default (CGS) Description
unit_system UnitSystem CGS Active unit system governing numerical values
kinematic_viscosity float 0.036 cm²/s Blood kinematic viscosity
fluid_density float 1.06 g/cm³ Blood mass density
murray_exponent float 3.0 Murray's law exponent for radius scaling at bifurcations
radius_exponent float 2.0 Exponent for radius contribution to cost function
length_exponent float 1.0 Exponent for length contribution to cost function
terminal_pressure float 60 mmHg Pressure at terminal vessels (capillary bed)
root_pressure float 100 mmHg Pressure at the root vessel (inlet)
terminal_flow float 6.25e-5 cm³/s Volumetric flow per terminal vessel
root_flow float or None None Total flow at root (computed during growth if None)
max_nonconvex_count int 100 Max attempts before switching to convex mode in non-convex domains

Methods

set(parameter, value, *, unit=None)

Update a named parameter with optional unit conversion.

Parameters

  • parameter (str): Parameter name (e.g., 'kinematic_viscosity', 'murray_exponent')
  • value: New value
  • unit (str or UnitSystem, optional): Source unit for conversion
set_unit_system(unit_system, *, convert_existing=True)

Switch to a new unit system, optionally converting existing values.

Parameters

  • unit_system (UnitSystem): Target unit system
  • convert_existing (bool): If True, convert stored values to new units

TreeParameters Examples

from svv.tree.data.data import TreeParameters
from svv.tree.data.units import UnitSystem

# Create parameters with defaults
params = TreeParameters()
print(params)
# Tree Parameters:
# ----------------
# Kinematic Viscosity: 0.036 cm^2/s
# Fluid Density: 1.06 g/cm^3
# Murray Exponent: 3.0
# ...

# Modify parameters
params.murray_exponent = 2.7  # Direct assignment
params.set('terminal_pressure', 80.0, unit='mmHg')  # With unit conversion

# Create SI parameters
si_system = UnitSystem(length='m', mass='kg', time='s', pressure='Pa')
params_si = TreeParameters(unit_system=si_system)

# Convert existing parameters to SI
params.set_unit_system(si_system, convert_existing=True)
Murray's Law: The murray_exponent (typically 2.7-3.0) controls how parent vessel radius relates to daughter radii at bifurcations. A value of 3.0 minimizes power dissipation; lower values (2.0-2.5) are observed in some vascular beds.

Setup Methods

set_domain(domain, convexity_tolerance=1e-2)

Set the spatial domain for tree growth.

Parameters

  • domain (Domain): Domain object defining growth region
  • convexity_tolerance (float): Tolerance for convexity check

Notes

The domain must have its implicit function built before setting. The method automatically determines if the domain is convex, which affects collision detection strategies.

set_root(*args, inplace=True, **kwargs)

Set the root vessel of the tree.

Parameters

  • *args: Variable arguments
    • Single argument: starting point (3D coordinates)
    • Two arguments: starting point and direction vector
  • inplace (bool): Modify tree in place
  • **kwargs: Additional parameters for root configuration

Notes

When direction is provided, the root becomes "clamped" with fixed orientation.

Growth Methods

add(inplace=True, **kwargs)

Add a single vessel to the tree using optimized bifurcation placement.

Parameters

  • inplace (bool): Modify tree in place
  • decay_probability (float, default=0.9): Probability decay for revisiting regions
  • **kwargs: Additional growth parameters

Returns

  • If inplace=True: None
  • If inplace=False: Tuple of modification data

Algorithm

The method uses constructive constrained optimization (CCO) to find optimal bifurcation points that minimize total tree volume while maintaining flow constraints.

n_add(n, **kwargs)

Add multiple vessels to the tree.

Parameters

  • n (int): Number of vessels to add
  • **kwargs: Parameters passed to each add() call

Notes

Shows a progress bar during growth. Each vessel is added sequentially with the tree structure updated after each addition.

Export Methods

export_solid(outdir=None, shell_thickness=0.0, watertight=False, **kwargs)

Export tree as a 3D solid model.

Parameters

  • outdir (str, optional): Output directory path
  • shell_thickness (float): Wall thickness for hollow vessels
  • watertight (bool): Create watertight mesh with smooth junctions

Returns

  • model (pv.PolyData): 3D mesh of the vascular tree
export_centerlines(outdir=None, **kwargs)

Export tree centerlines for simulation or analysis.

Returns

  • centerlines (pv.PolyData): Centerline geometry
  • polys: Connectivity information

Persistence (Save/Load)

Trees can be saved to and loaded from .tree files, enabling checkpointing during long growth processes and sharing of generated vascular structures.

save(path, include_timing=False)

Save the tree to a .tree file.

Parameters

  • path (str): Output filename. If no extension is provided, .tree is appended.
  • include_timing (bool): Include generation timing data for profiling/debugging. Default is False.

Notes

The saved file contains all tree data, parameters, and connectivity information needed to reconstruct the tree. The domain is NOT saved and must be set separately after loading via set_domain().

Tree.load(path) classmethod

Load a tree from a .tree file.

Parameters

  • path (str): Path to a .tree file.

Returns

  • Tree: Loaded tree instance.

Notes

The loaded tree will NOT have a domain set. You must call set_domain() after loading to enable domain-dependent operations like collision detection.

Save/Load Example

from svv.tree.tree import Tree
from svv.domain.domain import Domain
import pyvista as pv

# Create and grow a tree
domain = Domain(pv.Cube())
domain.create()
domain.solve()
domain.build()

tree = Tree()
tree.set_domain(domain)
tree.set_root()
tree.n_add(100)

# Save the tree
tree.save("my_vascular_tree.tree")

# Later, load it back
loaded_tree = Tree.load("my_vascular_tree.tree")

# Re-attach the domain for further operations
loaded_tree.set_domain(domain)

# Continue growing or export
loaded_tree.n_add(50)
loaded_tree.show()
Domain Not Saved: The domain implicit function can be expensive to compute and is not serialized with the tree. Save your domain separately using domain.save() if you need to preserve the complete workspace.

Visualization

show(**kwargs)

Display interactive 3D visualization of the tree.

Parameters

  • color (str): Vessel color (default: 'red')
  • plot_domain (bool): Display domain boundary
  • return_plotter (bool): Return the PyVista plotter instead of showing
  • **kwargs: Additional PyVista plotting parameters

Examples

Basic Tree Generation

from svv.tree.tree import Tree
from svv.domain.domain import Domain
import numpy as np
import pyvista as pv

# Create spherical domain
sphere = pv.Sphere(radius=10.0, center=(0, 0, 0))
domain = Domain(sphere)
domain.create()
domain.solve()
domain.build(resolution=30)

# Initialize and configure tree
tree = Tree()
tree.set_domain(domain)
tree.parameters.radius_exponent = 2.0  # Murray's law
tree.parameters.length_exponent = 1.0

# Set root at center
tree.set_root([0, 0, 0])

# Grow tree
tree.n_add(50)

# Show statistics
print(f"Terminals: {tree.n_terminals}")
print(f"Segments: {tree.segment_count}")
print(f"Tree scale: {tree.tree_scale:.3f}")

# Visualize
tree.show(color='royalblue', plot_domain=True)

Constrained Growth with Physical Clearance

# Set minimum vessel spacing
tree = Tree()
tree.physical_clearance = 0.5  # mm minimum spacing
tree.set_domain(domain)

# Root with fixed direction
start_point = [0, 0, -5]
direction = [0, 0, 1]
tree.set_root(start_point, direction)

# Grow with custom parameters
for i in range(100):
    tree.add(decay_probability=0.95)

    # Check growth progress
    if i % 20 == 0:
        print(f"Step {i}: {tree.n_terminals} terminals")

# Export results
solid = tree.export_solid(watertight=True)
solid.save('vascular_tree.stl')

Export for Simulation

import numpy as np
# Generate tree
tree = Tree()
tree.set_domain(domain)
tree.set_root([0, 0, 0])
tree.n_add(100)

# Export centerlines with radius data
centerlines, connectivity = tree.export_centerlines()

# Add flow properties
centerlines['flow_rate'] = np.zeros(centerlines.n_points)
centerlines['pressure'] = np.ones(centerlines.n_points) * 80  # mmHg

# Save for simulation
centerlines.save('tree_centerlines.vtp')

# Export solid model
solid = tree.export_solid(
    watertight=True,
    shell_thickness=0.1  # Create hollow vessels
)
solid.save('tree_solid.stl')

# Access raw data for custom processing
vessel_data = tree.data  # TreeData object
vessel_map = tree.vessel_map  # Connectivity information

Performance Considerations

Spatial Indexing: The Tree class uses multiple spatial data structures for efficient collision detection: KD-trees for exact nearest neighbor queries and HNSW (Hierarchical Navigable Small World) graphs for approximate nearest neighbor search in high-dimensional spaces.

Memory Management: Trees preallocate memory in chunks (default 4 million segments) to avoid repeated allocations during growth. This significantly improves performance for large trees.

Convexity: Non-convex domains require additional collision checks which can slow growth. The tree automatically switches to convex mode after repeated failures in non-convex regions (controlled by parameters.max_nonconvex_count).

Timing Analysis

Access detailed timing information through the times attribute:

# After growing tree
import numpy as np
import pandas as pd

# Convert timing data to DataFrame
timing_df = pd.DataFrame(tree.times)

# Analyze performance bottlenecks
print("Average times (ms):")
for key in ['vessels', 'collision', 'local_optimization']:
    if key in tree.times and tree.times[key]:
        avg_time = np.mean(tree.times[key]) * 1000
        print(f"{key}: {avg_time:.2f}")

See Also