Towards a Unified Physical Model for Virtual Environments

Peter M Chapman and Derek P M Wills

Virtual Environment Research Centre

Department of Computer Science

University of Hull

Hull

HU6 7RXTel:

(01482) 466348Fax:

(01482) 466666e-mail:

Peter.Chapman@dcs.hull.ac.uk

To simulate the behaviour apparent in the real world, virtual environments must include a detailed physical model based on the concepts of Newtonian mechanics. Although considerable research in computer graphics has led to the application of physically-based models in animation, the need for real-time performance in a fully interactive environment has hampered their progress in virtual environments. This paper will outline the development of a unified physically-based model suitable for virtual environments. Focusing on the technique ofModal Analysis(Essaet al., 1992), a series of experiments are discussed which demonstrate the applicability of this method for virtual environments.

Keywords:Physically-based modelling, Rigid and Non-rigid Dynamics, Finite Element Method, Modal Analysis

**1. Physically-Based Modelling
in Virtual Environments**

In a physically-based virtual environment
objects interact and behave in the same manner as their
counterparts in the real world. This ‘realism’ is
brought about through the application of the physical laws that
govern motion. However, the computational expense of implementing
these laws prohibits their inclusion in most current virtual
environment systems. The high degree of realism required by
applications such as surgical training systems (Logan *et al.*,
1996) demands the incorporation of detailed physical models.

This paper describes an investigation into
a unified physical model suitable for virtual environments. The
finite element method provides a detailed model of the rigid and
deformable behaviour of objects. Using a modal solution (Cook *et
al.*, 1989) it is possible to dramatically reduce the
complexity of this model without affecting its accuracy. This
paper will demonstrate how this technique, known as *Modal
Analysis* (Essa *et al*., 1992), could make the finite
element model suitable for use in virtual environment
applications. Experiments are undertaken to assess how the
complexity of the model may be reduced, along with the
limitations and effects on accuracy of these efficiency savings.
Through this analysis the feasibility of using modal analysis for
future virtual environments will become clearer.

*1.1 Other Physically-Based
Animation Models*

A unified physical model for virtual
environments must simulate the behaviour of rigid and non-rigid
motion. Due to the complexity of deformable behaviour it is
commonly disregarded (for example Keller *et al.*, 1995), in
order to provide an efficient model for virtual environments.
However, the resultant rigid body models have limited versatility
and cannot be extended to model general environments, containing
objects of different materials.

On the other hand non-rigid dynamic models,
based on the properties of elasticity, are typically very complex
and computationally intensive (for example Terzopoulos and
Metaxas, 1991), thus prohibiting their direct application in
virtual environments. Various abstractions of complete deformable
behaviour have been made in the past, in particular the models of
cloth are frequently based on approximations such as 2D
spring-mass meshes (Ng and Grimsdale, 1996). Alternatively,
deformable behaviour can be modelled through the manipulation of
object geometry alone (for example Wyvill *et al.*, 1986),
or through the incorporation of dynamic behaviour to geometric
models (for example Qin and Terzopoulos, 1996). Whilst such
models are more computationally efficient they are not based on
the physical laws of elasticity and they cannot be expanded to
formulate equations for rigid models.

Finite element equations of motion describe
the complete (rigid and non-rigid) motion of an object in a
single system of equations. Typically these equations are very
complex and their solution computationally expensive. However, by
performing a modal analysis (Essa *et al.*, 1992), the
complexity of the equations is reduced without a comparative loss
in accuracy. The following sections describe how savings are made
and investigate the use of this technique for virtual environment
applications.

**2. Modal Analysis**

Modal analysis is an implementation of the
finite element method for transient dynamic analysis of object
motion. Although an overview of the method follows, it is beyond
the scope of this paper to give a detailed description. The
interested reader will find rigorous derivations of the method in
any finite element text such as Cook *et al.* (1989).

Finite element analysis is based on a discrete representation of continuous behaviour. Each object is divided into a mesh of volume elements connected at nodes (Figure 1).

**Figure 1 finite
element representation of object**

A distribution of stress-strain forces act
across the individual elements, linking the displacement of
element nodes to the elasticity of the material. From the
stress-strain distribution a system of second order equations of
motion is derived and expressed as equation *(1)*.

*(1)*

where: **R** is the load vector of
forces acting on the nodes,** U** is the vector of nodal
displacements, **M** is the mass matrix, **D** is the
damping matrix and **K** the stiffness matrix which describes
the material properties of the object. This system of equations
describes the complete motion, both rigid and non-rigid, of an
object in response to the application of external forces.

Numerical techniques such as Gaussian
Elimination (Press *et al.*, 1992) may be used to solve
equation *(1)* for the displacements of each node in
response to applied forces. However this is very computationally
expensive, *O(n*^{3}*)* for *n*
degrees of freedom. An alternative solution is to decouple the
motion equations by solving the eigen problem equation *(2)*,
that results from the free undamped behaviour of equation *(1)*.

*(2)*

The eigenvalues are the frequencies of vibration of the
object and the corresponding eigenvectors are the modes
of vibration (the displacements of each node of the mesh
corresponding to the frequency of vibration). Transforming
equation *(1)* by the eigenvectors results in a system of *n*
distinct equations *(3)*, these relate the vibration modes
of the object to the frequencies of vibration, and represent a
change of basis from the original equation of motion.

*(3)*

with:

where: *q*_{i}(*t*)
is the modal displacement, _{i} is the modal damping factor, _{ij} is
the Kronecker delta and *i,j* = 1-> *n*. The final
stage of the analysis determines the actual nodal displacements
of the mesh by completing the superposition of equation *(4)*.

*(4)*

The solution of the eigen problem is
computationally expensive, however this is a pre-processing
stage. Only the simple damped motion equations *(3)* and the
summation *(4)* are computed at each time-step, which is a
vast improvement upon solving the coupled system of equations *(1)*.

*2.1 Applicability of Modal
Analysis*

On first inspection there appears no real benefit to modal analysis over the direct solution of the finite element model. However, in vibration analysis it can be observed that the high frequency components, which model only small local deformations, are heavily damped and contribute little to the overall displacement. It is therefore possible to ignore the high frequency modes without a substantial loss of accuracy. The effect of this is to reduce the number of calculations required at each time-step and hence improve efficiency.

The simulation time-step is dependent upon
the highest frequency mode used in the analysis. If the higher
frequencies are dropped then the time-step can be increased
without risk of loosing accuracy. Furthermore, this method is
simpler and more robust than other techniques (such as the *hybrid
model* of Terzopoulos and Metaxas, 1991), that model the rigid
and deformable behaviour separately before combining them to
create a general system of equations for motion.

The range of shapes that can be suitably meshed and analysed by this method is limited. The initial object shape is based on a implicit surface description and is limited to simple volumes such as cuboids and ellipsoids. The analyses and comments reported in this paper are only valid for these simple volume shapes.

**3 Analysis of Method**

In order to assess the effects of simplifying the full dynamic model a series of experiments have been undertaken using the finite element analysis software ANSYSŪ. The dynamic simulation of object motion under stress from external forces for a full model is compared against the simulation for simplified models. From this investigation it is possible to determine how the complexity of the finite element model may be reduced and still produce (visually) meaningful results.

There are certain well established areas
where savings can be made, for example reducing the mesh density
(that is, reducing the number of elements and degrees of freedom
in the mesh), or reducing the number of modes used in the
superposition *(4)*. In addition to these experiments the
following sections also describe how the choice of material and
element type can effect the simulation.

*3.1 Mesh Density*

It is known that high density meshes provide more accurate models of behaviour than less dense meshes. To ascertain how the mesh density effects the accuracy of the model, the nodal displacements for each mode of vibration were investigated. Consider the beam in Figure 2.

**Figure 2 beam to
undergo modal analysis**

**Figure 3 tenth
mode of vibration for three meshes of different density**

The element nodes are displaced by slightly different amounts for meshes of different density. The node displacements of the finest mesh are used as a metric to assess the accuracy of the other (coarser) meshes. If a node displacement for a coarse mesh is compared against the node displacement of the fine mesh (Figure 4), the relative difference of the node displacements can be calculated. Taking the average of this relative difference across a series of nodes will give a measure of the accuracy of the coarse mesh. Figure 5 shows the results for a set of these comparisons.

**Figure 4 comparison of node
displacements for meshes of different density**

**Figure 5 graph
showing the relative difference in node displacements for
hexahedral meshes of density (different numbers of degrees of
freedom)**

As can be seen from Figure 5, the coarser the mesh the less accurate the finite element analysis. Furthermore, it can be noted from the graph that the loss of accuracy becomes more appreciable when the mesh has less than 1000-1500 degrees of freedom. There is little difference in accuracy between meshes with 5000 degrees of freedom and those with 1000 degrees of freedom (although this would represent a huge reduction in computation costs). From these results it is suggested that a finite element mesh with around 1000 degrees of freedom produces sufficient accuracy for a physical simulation. This will not apply in all situations with all meshes and object shapes, however, it is a good guide to expected behaviour. As stated previously, this analysis was undertaken on simple volumetric shapes and is not valid for complex shapes.

*3.2 Element Types and Shapes*

There are two main types of element used in
volumetric finite element meshes, hexahedral and tetrahedral
shaped elements. Generally, the greater the number of nodes in
one element the better that element will model continuous
behaviour, although the equations describing the element
behaviour also dictate the accuracy (see Cook *et al.*, 1989
for more details).

The analysis described above for hexahedral meshes of different density was also undertaken with meshes of tetrahedral elements. Figure 6 shows the same mode of vibration for two different density tetrahedral meshes. Comparing these meshes with those of the hexahedral meshes in Figure 3, it can be noted that in the case of the tetrahedral meshes the mode for the less dense mesh is appreciably different in shape to the finer mesh. In general the (8-node) hexahedral meshes are more stable than the (10-node) tetrahedral meshes. This result bears out accepted wisdom (NAFEMS, 1986), which suggests that whilst it is easier to mesh a general object using tetrahedral elements, hexahedral elements produce more accurate results during analysis.

**Figure 6
deformed shape of the twentieth mode for two different density
tetrahedral meshes**

*3.2 Material Properties*

The fall-off in accuracy with a reduction
in mesh density is also dependent upon the type of material used.
The Modulus of Elasticity determines the overall flexibility of
the object, so this will govern the actual amount of deformation
for each mode of vibration. However, this does not effect the
relative displacements for meshes with different densities, this
is dominated by the effect of Poisson’s ratio (*v*).

When a bar is loaded in tension, the width
of the bar becomes smaller as the length increases.
Poisson’s ratio expresses this relationship, between the
strain in the lateral direction and that in the longitudinal
direction. This ratio ranges between 0.2 and 0.5 for most
materials (see Gere and Timoshenko, 1991 for examples). Lower
values of Poisson’s ratio signify the body will not stretch
as much under tension, accordingly the value for rubber (*v*
= 0.45ƒ 0.5) is higher than for aluminium (*v* =
0.33).

Figure 5 demonstrates that coarse meshes of
rubber are poorer approximations than the equivalent meshes of
aluminium or nylon (*v* = 0.4), this is the effect of
Poisson’s ratio. Consequently, in terms of physical
simulation, materials with low Poisson’s ratio are more
accurately modelled with coarse meshes. Although this situation
does not constitute a limitation on materials that can modelled,
it is worth noting during the building of a virtual environment
application as it may effect the choice of mesh density used.

*3.3 Number of Modes in
Superposition*

As mentioned previously reducing the number of modes used in the superposition will vastly reduce the computational costs. How many modes can be cut out of the superposition without a great loss in accuracy? This section focuses on results from an experiment designed to discover the answer to this question.

A dynamic transient analysis was performed using the finite element analysis software ANSYSŪ. A small rubber beam was fixed at one end and forces were applied at free nodes. A transient analysis of this situation demonstrates how the beam deforms under the action of the forces (Figure 7 shows the forces and the resulting deformation). The simulation can be performed with any number of modes used in the superposition.

**Case 1**

**Case 2**

**Figure 7 beam
with forces and displacement constraints and the resultant
deformation due to the forcing functions. Case1 is a simple
bending deformation, Case 2 slightly twists the beam**

To measure the loss in accuracy for different numbers of modes in the superposition, the nodal displacements of a series of mode superpositions were compared against an analysis using 1000 modes. Figure 8 shows the average results for this kind of comparison.

**Figure 8 graph
showing, for a rubber beam, the average relative difference in
node
displacements of a series of mode superpositions with different
numbers of modes.
This is for the two cases of loading shown in Figure 7.**

It can be seen that there is not a great fall off in accuracy for analyses using between 1000 and 100-200 modes. However, if less than 100 modes are used the results become less reliable. A reduction in the number of modes used from 1000 to 100 would represent a 90% reduction in the number of equations solved at each time-step. This efficiency saving produces only a comparatively small loss in accuracy.

The two cases shown in Figure 7 emphasise the effects of high frequency modes upon object deformation. Figure 8 shows that, when the number of modes is reduced, the simulation for Case 1 is more accurate than the simulation for Case 2. The twisting effects of deformation are the result of higher frequency modes of vibration, therefore when the higher frequency modes are eliminated the twisting effects are removed. Case 1 depends less on the high frequency modes than Case 2 and thus remains more accurate when the high frequency modes are removed.

**4 Conclusion**

The results of the experiments described in this paper, suggest the computational cost of a finite element analysis, can be significantly reduced without a great loss in accuracy. This reduction in complexity increases the potential of using Modal Analysis in a unified physical model for virtual environment applications.

The use of relatively coarse meshes reduces the number of calculations during the generation of the finite element equations and the modes of vibration. It has been demonstrated that relatively coarse meshes can be used without a noticeable effect upon the accuracy. Also, hexahedral elements provide a more stable finite element solution than tetrahedral elements. An interesting point raised from the experiments is the effect of Poisson’s ratio upon the accuracy of simulation. Materials with low Poisson’s ratio maintain accuracy longer when the number of degrees of freedom is reduced. This could affect the choice of mesh density in some simulations.

A considerable saving on processing time can be made by reducing the number of modes in the superposition. It is demonstrated that accuracy is maintained when only a fraction of the available modes are used in the simulation.

Overall, the results here suggest more effort should be directed at evaluating this method in a purpose built virtual environment. Timings for dynamic simulation must be compared against rendering times and other operations in the virtual environment (such as collision detection), in order to demonstrate the use of this model in future systems.

**Acknowledgements**

This work is supported by the Engineering and Physical Sciences Research Council (EPSRC) of the United Kingdom. The ANSYSŪ software product is furnished and licensed by ANSYS, Inc.

**References**

Cook, R. D., Malkus, D. S. and Plesha, M.
E., (1989), *Concepts and Applications of Finite Element
Analysis*., 3rd Edition, John-Wiley and Sons, New York, USA.

Essa, I. A., Scarloff, S. and Pentland, A.,
(1992), ‘A Unified Approach for Physical and Geometric
Modeling for Graphics and Animation’, *Computer Graphics
Forum*, **11**, (3), C130-C138.

Gere, J. M. and Timoshenko, S. P., (1991), *Mechanics
of Materials*, 3rd SI Edition, Chapman & Hall, London, UK,
777-782.

Keller, H., Stolz, H., Ziegler, A., Braunl,
T., (1995), ‘Virtual Mechanics: Simulation and Animation of
Rigid Body Systems with AERO’, *Simulation*, **65**,
(1), 74-79.

Logan, I. P., Wills, D. P. M., Avis, N. J.,
Mohsen, A. M. M. A. and Sherman, K. P., (1996), ‘Virtual
Environment Knee Athroscopy Training System’, *Society for
Computer Simulation, Simulation Series*, **28**, (4),
17-22.

NAFEMS (National Agency for Finite Element
Methods and Standards), (1986), *A Finite Element Primer*,
Department of Trade and Industry, UK.

Ng, H. N. and Grimsdale, R. L., (1996),
‘Computer Graphics Techniques for Modelling Cloth’, *IEEE
Computer Graphics and Applications*, **16**, (5), 28-41.

Press, W. H., Flannery, S. A., Teukolsky,
S. A. and Vetterling, W. T., (1992), *Numerical Recipes in C:
The Art of Scientific Computing*, 2nd Edition, Cambridge
University Press, Cambridge, England.

Qin, H. and Terzopoulos, D., (1996),
‘D-NURBS: A Physics Based Framework for Geometric
Design’, *IEEE Transactions on Visualization and Computer
Graphics*, **2**, (1), 85-96.

Terzopoulos, D. and Metaxas, D., (1991),
‘Dynamic 3D Models with Local and Global Deformations:
Deformable Superquadrics’, *IEEE Transactions on Pattern
Analysis and Machine Intelligence*, **13**, (7), 703-714.

Wyvill, G., McPheeters, C. and Wyvill, B.,
(1986) ‘Animating Soft Objects’, *Visual Computer*,
**2**, 235-242.