Transcript Real Time Deformation of Volumetric Models for Surgery Simulation Andrew B. Raij ([email protected]) Comp 259 Class Presentation Department of Computer Science, UNC-Chapel Hill.

Real Time Deformation of
Volumetric Models for
Surgery Simulation
Andrew B. Raij ([email protected])
Comp 259 Class Presentation
Department of Computer Science, UNC-Chapel Hill
Surgery Simulation


Realistic surgery simulators could improve
surgeon training, surgery planning
Requirements for Realism







Volumetric models of the human body that capture
the internals of organs, etc
Bio-mechanical descriptions of the properties of
different body parts
Interactive, real-time, realistic deformation of
volumetric models
Cutting, modification of model geometry
Haptic Interaction
Collision detection (haptic device to organs, organ to
organ)
Real-time rendering of volumetric objects
Volumetric Deformation

I will discuss two approaches



An FEM-based approach from Morten BroNielsen and Stephane Cotin at INRIA
3D Chainmail approach from Sarah Frisken at
MERL
See my references for information about
other approaches
Bro-Nielsen, Cotin

Requirements




Speed Speed Speed
Try to save time with offline computation no
matter how small the gain
Elastic behavior only has to seem correct,
doesn’t have to be physically correct
Although not implemented in the paper, they
want enough flexibility to add cutting later
An Object




 - a linear elastic solid
x - 3D point in the object
u - displacement of x due to
deformation
G - Surface of the object



G = G0  G1, where G0  G1 = 
A force f(x) is applied at G1
G0 has a fixed displacement u(x)
The Linear Elastic Model

Strain energy of  is




is the strain vector
 is the stress vector
Stress and strain are related by Hooke’s Law
where D is a 6x6 symmetric material stiffness matrix
For isotropic, homogeneous materials D is
where ,  are lamé
material constants
The Linear Elastic Model (cont’d)

If the elements of
follows
are defined as

Then we can define a matrix
such that
The Linear Elastic Model (cont’d)





The potential energy of the system is the strain energy
minus the work done by the external
forces applied to G1
Therefore, using the variables defined on previous
slides, the potential energy is
A body’s desired deformation is found when its potential
energy is minimum.
If we had an equation describing the displacement over
, we could take the derivative of E, set it to 0, and solve
for the displacement
We don’t have this so we use FEM
Apply FEM

Steps




Discretize model into tetrahedral finite
elements
Find the displacement for each element
Find the potential energy for each element
Build a system describing the displacement
of all the elements and solve the system
Apply FEM (cont’d)

Displacement of a tetrahedron



Li are the natural coordinates of the tetrahedron
Ve is the volume of the element
ai, bi, ci, di, are found by inverting the equation on the left below


xi, yi, zi, are the global coordinates of each node
Apply FEM (cont’d)

If we apply our equation for the potential energy
of a body to this element and take its derivative
we get


where
is the discretized force applied to the
element
To find Be, we remember that B for the whole body is
Apply FEM (cont’d)

Be is B for each node in the element applied
to the element’s displacement already.
Apply FEM (cont’d)


Everything inside the integral is constant
so we can rearrange and integrate to get
where
Ke is called the stiffness matrix. Build the
global stiffness matrix from the stiffness
matrices for each element to get the global
system:
Condensation





Don’t need to display displacement of internal nodes
during simulation, only surface node displacements.
If the nodes in our system
are ordered such that
the surface nodes are first and internal nodes are second
then our system can be rewritten in block matrix form.
describes the stiffness between the surface nodes
and
describes the stiffness between the internal
nodes.
and
describe the stiffness where the internal and
surface nodes meet.
and
are the displacements for the surface and
internal nodes, respectively.
and
are the forces
applied to the surface and internal nodes, respectively.
Condensation (cont’d)
So the surface displacement is


The system is condensed to the complexity of a
surface model but still conveys the behavior of
the volumetric model
The condensed system is typically dense but
this doesn’t matter (you’ll see why next).
Explicit Matrix Inversion

Typically we use something like the Conjugate
Gradient algorithm to implicitly solve a system


It’s hard to guess how long this will take, which is bad
when you need to display at a set frame rate
Bro-Nielsen, Cotin explicitly invert the stiffness
matrix




If the stiffness matrix doesn’t change (i.e. no cutting),
we can do the inversion offline and save it
Solving for displacement becomes a cheap matrixvector multiply
But we lose precision due to numerical error
We need more memory to store the large, dense,
inverted stiffness matrix
Dynamic System




If we want a more physically correct model, we need to
add mass and damping
Use Lagrangian equation of motion
M = mass
C = damping
K = stiffness
Masses are only at nodes so
the mass and damping matrices
are block diagonal
For the condensed
system:
V = element volume
 = mass density
 = scaling factor
Dynamic System (cont’d)

If we plug finite difference estimates into the Lagrangian
equation we get
Or

where
If forces are only applied at a few surface nodes, will be
very sparse

Can’t take advantage of that since the other components of
not be sparse.
will
Static System, SMMV

Use simpler, less physically realistic system
instead of the Lagrangian equation to take
advantage of the sparseness of the force vector.



When solving for displacement, skip elements in force
vector that are 0 since they will not contribute to the
result.
The authors call this Selective Matrix Vector
Multiplication (SMVM)
Static system and dynamic system give similar
results for slow changes in forces
3D Chainmail

Basic Idea – Chain Analogy



Elements of a volumetric object are chained together
When a link in the chain is pushed or pulled,
neighboring links move along with it
If a link is compressed or expanded to its limit, motion
is transferred to its neighbors
Deformation Constraints


Elements have deformation constraints controlling
where they can be relative to their neighbors
Stretching, Contraction


Elements must be in the
horizontal range [minDx,
maxDx] and the vertical
range [minDy, maxDy] wrt
to its four neigbors
Shearing

Element must be within +/maxHorizDy of left and right
neighbor’s y and +/maxVertDx of top and
bottom neighbor’s x
Algorithm Data Structures

A list of pointers to elements that moved
on the last time step and the elements’ old
positions


For turning back time if there’s a collision
between objects
Deformation Candidate Lists


Contains elements that might have to move
on the next time step
There is one global candidate list for each
direction, i.e. top, bottom, etc
The Algorithm

When an element is moved



It is added to the list of moved elements along with its
original position
Its closest neighbors are added to the proper
deformation candidate list, i.e. the element to the left of
the moved element is added to the left candidate list
Then the candidate lists are processed in order
until the lists are empty or the system is not
permissible (constraints not met, collision)

When the system is not permissible, time is rolled back
and a smaller time step is tried
The Algorithm (cont’d)

EXAMPLE: For each element in the right
candidate list



Check shear and stretch constraints against the
neighbor that sponsored it – the sponsor must be the
candidate’s left element
If the element’s constraints are not met, move it a
minimum distance until they are met
If the element is moved, sponsor its neighbors for
motion by putting them in the appropriate candidate
list.


Note the left candidate is not sponsored since it sponsored
the motion that lead to this element’s motion
The other candidate lists are processed similarly
Elastic Relaxation




The Chainmail algorithm alone does not minimize
the overall potential energy of the system
After every application of the Chainmail algorithm
and when extra computation time is available, the
elements are adjusted to satisfy energy
constraints.
Define system energy in terms of the relative
positions of elements to their neighbors.
Result: the object deforms approximately, and then
eases into a more correct shape over time
Analysis


Very simple
Fast






Elements are considered at most once per
deformation
Elements only check themselves against their
sponsoring neighbor to decide if/how they should
move
Deformation propagates outward and ends as soon
as it can
Efficient since it uses extra time to refine the
deformation
Not physically-based like FEM
Not obvious how to convert biomechanics
information about body parts to element
deformation constraints
References

M. Bro-Nielsen and S. Cotin, Real-time
volumetric deformable models for surgery
simulation using finite elements and
condensation, Computer Graphics Forum,
15(3):57-66, 1996. Eurographics '96.

H. Delingette, Towards realistic soft tissue
modeling in medical simulation, In
Proceedings of the IEEE : Special Issue
on Surgery Simulation, pages 512--523,
Apr. 1998.
References (cont’d)

S. Frisken, et al, Simulating Surgery using
Volumetric Object Representations, RealTime Volume Rendering and Haptic
Feedback, Feb 1997, Not Published.

S. Gibson, 3D Chainmail: a Fast Algorithm
for Deforming Volumetric Objects, In Proc.
Symp. on Interactive 3D Graphics, 149154, 1997.