by
Rebecca McWilliam
B.Eng., McMaster University, 2015
A Thesis Submitted in Partial Fulfillment of the Requirements for the Degree of
MASTER OF APPLIED SCIENCE
in the Department Mechanical Engineering
c
Rebecca McWilliam, 2017 University of Victoria
All rights reserved. This thesis may not be reproduced in whole or in part, by photocopying or other means, without the permission of the author.
Local Model Interpolation for Stable Haptic Interaction
by
Rebecca McWilliam
B.Eng., McMaster University, 2015
Supervisory Committee
Dr. D. Constantinescu, Co-Supervisor (Mechanical Engineering)
Dr. D. Damian, Co-Supervisor (Computer Science)
Supervisory Committee Dr. D. Constantinescu, Co-Supervisor (Mechanical Engineering) Dr. D. Damian, Co-Supervisor (Computer Science) ABSTRACT
Haptic, or force, interaction in virtual environments can support the development of physical manipulation skills if users feel realistic forces in response to their motion in the virtual environment. The realism of the forces felt by users depends on: (i) how accurately the virtual environment simulates a real life situation such as surgery; and (ii) how faithfully the haptic controller renders the simulated interactions to users. Accurate simulations of real life situations such as surgery run at variable frequencies of the order of 20 − 100 Hz. However, the haptic controller needs updated stiffness and direction of contact at 1000 Hz to faithfully convey the shape and hardness of the virtual objects to the user. This thesis proposes to bridge the gap between the required fast haptic control rate and the slower virtual environment updates through a passive local model of interaction. This model comprises an approximation of the shape and stiffness of the virtual world in the area near the point of interaction. It also monitors its exchange of energy with the user to ensure its own passivity and thus, the stability of the haptic system. Lastly, the local model eliminates the spurious discontinuities that arise in contact direction at model updates by interpolating the contact normal before rendering it to the user.
Contents
Supervisory Committee ii
Abstract iii
Table of Contents iv
List of Tables vi
List of Figures vii
1 Introduction 1
1.1 Surgical Training . . . 1
1.2 Challenges in Haptic interaction . . . 2
1.2.1 Slow Virtual Environment . . . 4
1.2.2 Mismatch in Update Rates . . . 6
1.2.3 Instability . . . 10
1.3 Objectives of the Thesis . . . 13
1.4 Overview of the Thesis . . . 13
2 Methodology 15 2.1 System Overview . . . 15
2.1.1 Virtual Environment . . . 16
2.1.2 Local Model . . . 17
2.2 Stability . . . 17
2.2.1 Overview of Force Bounding Approach . . . 17
2.2.2 Force Bounding Approach applied to Deformable VE . . . 21
3 Experiments 30 3.1 Simulation Models . . . 30
3.2.1 Conclusion . . . 33 3.3 3 DOF experiments . . . 33 3.3.1 Conclusions . . . 38
4 Conclusion 39
A 41
A.1 Collision Detection extension . . . 41 A.2 Finite Element Model . . . 42
List of Tables
Table 1.1 Summary of how challenges with haptic interactions are addressed 4 Table 2.1 Notation used in this thesis . . . 15
List of Figures
Figure 1.1 Ambiguous as to proper collision response without prior knowl-edge of the system. When the collision a) occurs, it is unclear whether b) or c) should be the proper response . . . 3 Figure 1.2 During collision detection, the object in the virtual environment
is replaced by a plane tangential to the closest point of the object 6 Figure 1.3 Coarse top mesh layer is divided to have a finer child mesh in
each subregion . . . 7 Figure 1.4 Norton equivalence applied at the shared nodes to represent the
full object dynamics to the subregion mesh . . . 8 Figure 1.5 Multiple points of contact requiring separate local models for
each interaction point . . . 9 Figure 1.6 Difference in distance between proxy and haptic device position
for stiff or soft objects . . . 9 Figure 1.7 Non linear virtual coupling . . . 11 Figure 1.8 Hysteresis from different pressing and releasing path for
position-force behaviour . . . 12 Figure 2.1 Interconnection and data exchange between the user, haptic
de-vice, local model, and virtual environment . . . 16 Figure 2.2 Overall configuration of haptic system as described in [1] . . . . 18 Figure 2.3 Projection of the vector between the haptic device position and
a point on the plane onto the normal of the plane . . . 23 Figure 2.4 Projection of haptic device position onto local plane used to
de-termine if it lies within the other three nodes . . . 24 Figure 2.5 Interpolation between the current local nodes and the virtual
environment local nodes . . . 27 Figure 2.6 Interpolation between the current local nodes and the virtual
Figure 2.7 Flow chart of process in the local model . . . 29 Figure 3.1 Haptic interaction with a deformable virtual environment
with-out force bounding: positions of the haptic device (dash-dotted red line) and of the virtual constraint in the local model (solid black line), as updated by the virtual environment. . . 31 Figure 3.2 Force bounding on constraint updates and when moving out of
the constraint: positions of the haptic device (dash-dotted red line), and of the virtual constraint in the local model (solid black line) and in the virtual environment (dashed black line, mostly overlapped with the solid black line). . . 32 Figure 3.3 Detail of Fig. 3.2: positions of the haptic device (dash-dotted
red line) and of the virtual constraint in the local model (solid black line) and in the virtual environment (dashed black line). Note the interpolation in the local model between environment updates when force bounding occurs. . . 32 Figure 3.4 Virtual object composed of 11 tetrahedral elements. . . 33 Figure 3.5 Jumps are seen in the position of the haptic device at times
corresponding to updates coming from the virtual environment. 34 Figure 3.6 Large changes to the normal sent from the virtual environment
after a deformation has occurred. . . 35 Figure 3.7 The position of the haptic device changes smoothly as the user
tries to move along the z axis. . . 36 Figure 3.8 Updates to the normal vector after a large deformation from the
virtual environment and how it is treated by the local model . . 37 (a) Virtual Environment . . . 37 (b) Local Model . . . 37 Figure A.1 Simple FEM model with two spring elements. Element 1 between
Introduction
In order to create a virtual reality simulation, visual feedback is extremely important. Equaly important in bringing a real immersive quality to the simulation is to integrate haptics into the simulation. Haptics is the use of force to recreate a touch sensation. It can be used to give a user the perception that they are interacting with real objects that are in fact virtual. A haptic device is one that can receive inputs from a user as to how they would like to interact with a virtual environment and can provide a force to the user as an output . Haptics can be used in a variety of different fields including, military training, creative painting, surgical training and teaching, remote surgery, and in video games or entertainment, among others [2].
1.1
Surgical Training
Various models have be used in order to train for live surgery. The closer the model is to the real surgery the greater the benefit. Many of the current methods have some serious drawbacks. The use of human cadavers gives good anatomical accuracy, the tissue response is not that close to live tissue. There is also limited availability and a high cost associated with storage. When surgery is practised on these subjects, there is little indication as to the success of the attempt and the lack of readily available new cadavers can cause the surgeons in training to have a greater fear of making mistakes. An alternative that can be used in order to remedy some of the issue of the lack of realistic tissue response is to use live animals. Unfortunately this comes with additional ethical issues and the anatomy may be quite different than for humans. Using surgical dummies can alleviate ethical concerns and allow for greater
repeatability, however there is quite a large upfront cost and storage cost necessary. These dummies might also be limited in the types of surgery that can be practised on them.
The solution presented to remedy the concerns of previous surgical training tools, is to create a virtual surgery training simulation with realistic haptic feedback. This would require a low equipment cost and allow for extremely high degree of repeata-bility and diverse surgery options. There is even the possirepeata-bility to create a custom surgery based on existing patient information. The simulation may also be able to provide feedback as to the student surgeon’s skill and areas requiring improvement. With improved haptic feedback, virtual simulations can become very close to having the feeling of live tissue.
After training with these methods the next step is to move on to supervised live surgery. In this case it is important that the student surgeon comes into the environment having been well practised and with confidence in order to reduce patient harm. Haptic simulations could be used not only for surgical practice but could also be used to develop new surgical techniques. In addition they may be used by a trained surgeon to plan and practice an upcoming surgery that may be unique or especially difficult.
1.2
Challenges in Haptic interaction
A virtual environment can be used to implement new designs and ideas cheaply without the need to use a machine shop. However much of the important information to be gained comes from interacting with the design/implementation making haptic interactions that allow one to physically interact with a virtual environment crucial. A virtual environment without interaction can be useful in developing procedure but not skill [3]. Interactions require that any action performed by the operator lead to an instantaneous response from the simulated organ, with both force and texture feedback being important [4].
Haptics need to operate at a very fast sampling rate in order to give a realistic sensation to the user of approximately 1000 Hz. On the other hand, the virtual environment is much slower in updating as it is computationally expensive. The updates from the virtual environment are also usually variable rather than fixed as is needed for haptics. When the virtual environment is very stiff this can also lead to stability problems in haptics.
An important use of haptics is to provide a user with the experience of having a tool interact with an environment [5]. In [3], the major challenge in this stems from the unilateral constraints occurring during the collision of rigid bodies. These must be computed in real time, and the simulation must be stable. Collision detection and response must occur in real time in order to generate reaction forces. In order to accurately perform collision response it is often necessary to have knowledge of the system from previous time steps, as is the case in Fig. 1.1.
Figure 1.1: Ambiguous as to proper collision response without prior knowledge of the system. When the collision a) occurs, it is unclear whether b) or c) should be the proper response
Colgate also defines the Z-width in [3] as the range of impedances that a particular haptic device can display without loss of stability. To achieve higher impedance physical dampers are required, while the introduction of negative virtual damping to cancel the physical dampers can be used to display low impedance in the same device. In complex virtual environments, the apparent impedance is greatly affected by the geometry.
Challenge Approach Used to Address Papers Slow virtual environment Model reduction [6] [19] [7] [9] Model interpolation [8] [11] [4] Pre computations [12] Mismatch in update rates Nested/hierarchical model [22] [10]
Local model & Multi-rate [13] [14] [15] [16] [21] [17] [18]
Multi-rate only [20]
Instability
Human operator modelling [23]
System analysis framework [24]
Controller/bounding [25] [26] [27] [28] [29] [1]
Table 1.1: Summary of how challenges with haptic interactions are addressed Multiple approaches have been implemented to try to address various challenges associated with haptic interactions. These approaches are summarized in Table 1.1.
1.2.1
Slow Virtual Environment
In order to deal with the instability issues in the haptic interaction with virtual en-vironments, some papers have relied either in whole or in part on the speeding up of the virtual environment with pre-computations and model reductions. Virtual envi-ronments, especially those that must undergo physical deformations large enough to be visible are known to required computationally expensive simulations [6]. Physical based models are important to give the most realistic characteristics to a virtual envi-ronment. Finite element models (FEM), for example can give a good representation of soft tissue, even with multiple non-homogeneous layers. Though these models are intensive both on CPU and memory [7]. FEM becomes increasingly complex with more elements, which are required for larger and more complex environments.
A solution is proposed in [6] to use dimensional model reduction for deformable models that have geometric non-linearities, yet are comprised of linear materials. This model reduction can lead to fast and pre-computable approximations.
In order to display visual real time environment, [7] uses Spectral Lanczos De-composition method to approximate tissue characteristics.
The Multiplicativce Jacobian Energy Decomposition method is implented by [19] to simulate soft tissue dynamics. Using this method, non-linear hyperelastic
mate-rials are discretized on linear tetrahedral meshes. This is done in order to speed up assembly of the stiffness matrix. An approximation of the compliance matrix is used as an asynchronous pre conditioner that is updated at a low rate to give a close ap-proximation of the mechanical coupling interaction making the assumption that only small perturbations can occur within the pre conditioner updates. Additionally the contact problem is linearised.
In [9] it is noted that the use of surface models to simplify computations are not sufficient as in surgical simulations there will undoubtedly be a need for layering, par-ticularly in the case of cutting surgeries. Instead the suggestion is made to condense the finite element model to remove interior nodes in intermediary steps. In addition performing necessary matrix inversions in the pre calculation phase.
In the case of complicated non-linear heterogeneous soft tissue, [8] measures phys-ical deformation in objects in the real world and then creates a strain space with these finding. This avoids the material parameter selection process. Once the strain space has been made, the behaviour can be interpolated to model the full non-linear behaviour.
Similarly, [4] performs pre computations of elementary deformations to all free (i.e. unconstrained) nodes on a surface of mesh for all DOF. The displacements of all other free nodes in the mesh, including interior nodes, are computed for each of these nodes and stored in a tensor. These elementary deformations can be superimposed during real-time simulation to quickly provide a reasonable force to the user. Though with more complicated meshes the pre computations required to accomplish this can take multiple hours.
An interpolation-extrapolation data structure known as a force grid is designed by [11]. This grid is defined over a virtual environment workspace, and each vertex of the grid contains a force. The force grid is filled with the new object shaped at a rate of 30 Hz. The point of contact with the haptic device is assumed to act as a light source, such that it has a larger influence closer to the interaction point. The influence becomes negligible as the distance increases, so it is bounded in order to reduce computation cost. It is found that this bounding distance should be approximately 2-3 times the distance that may be travelled between time steps. This data structure is to be used for haptic interactions with continuous contact.
Green’s functions are used in [12] as a basis for describing linear elastostatic model deformations. Green’s functions allow for fast computation of boundary value prob-lems. The simulation can be sped up by only boundary data needing to be computed,
which is sufficient for surface deformations and solving contact forces. A capacitance matrix is also used to determine the linear combination of the pre-computed Green’s functions needed to solve the boundary value problem, and lets the local force re-sponse be decoupled from the global deformation.
1.2.2
Mismatch in Update Rates
A haptic interaction in general consists of three stages. First the operator’s movement is tracked, then collision detection is performed, finally a mechanical impedance for the device is generated. The stage of collision detection is one that can significantly slow down a simulations, especially in the case of complex objects. In order to remedy this, [13] represents the virtual object as a virtual plane during collision detection see Fig. 1.2. The plane is tangent to the nearest point on the virtual object. Provided that the update rate is not too slow, this can feel smooth, otherwise the object may feel lumpy, with many edges.
Figure 1.2: During collision detection, the object in the virtual environment is re-placed by a plane tangential to the closest point of the object
To speed up collision detection in the virtual environment, [22] uses a point-based representation of one virtual object and a signed distance field for the other. A nested point shell is used in order to prevent point by point checking. A each haptic cycle, the hierarchy of the point shell is traversed, looking up the respective signed value in the signed distance field. If the value is greater than the radius of the point the
corresponding sub tree is eliminated from the search. If the value is within the radius the sub tree is looked into for the contact point, which occurs when the distance field value is negative. When a negative value is found a penalty force is assigned to that node. In order to further speed up the search, at every tree node processing an estimate of the earliest possible time that node can collide is determined and stored so that it may be skipped until that time. This method allows for multiple points of contact. The quality of the haptic force feedback is greatly dependant on the quality of the point shell and with this method it is difficult to render the feeling of sharp corners.
In [14] the user is represented by a probe in the virtual environment, and by a virtual proxy in a local model. The local model is a very simple non-deformable physical model. This is where collision detection and force computation occurs. The local model is updated based on the position of the probe in the virtual environment and its partial derivatives.
A proposal is made in [10] to use a multiple layered methods. This can reduce the computation time in areas away from the point of interaction. The model is comprised of a top mesh that is coarse and a child mesh for each subregion of the coarse mesh which is much finer as seen in Fig. 1.3.
Figure 1.3: Coarse top mesh layer is divided to have a finer child mesh in each subregion
The meshes are decoupled by utilizing Norton equivalence at the shared nodes between the coarse and fine mesh in Fig. 1.4. Instead of voltage and current as
typical seen in Norton and Thevenin equivalence in electrical circuits, displacement and force are used respectively to find equivalent impedance.
Figure 1.4: Norton equivalence applied at the shared nodes to represent the full object dynamics to the subregion mesh
Due to the decoupling, the two meshes can be updated independently, with the coarse mesh being solved 10 times slower than the fine mesh (1 Hz vs 10 Hz). At each slower update of the coarse mesh, the appropriate subregion for the finer mesh is selected based on the haptic device position. In addition any displacements seen to the shared nodes during the faster fine mesh are used to update the coarse mesh. The case of multiple user interaction is explored by [15] and [16], with a set of deformable objects. A local model, computed at a higher rate than the virtual environment, is needed in this case for each interaction point shown in Fig. 1.5.
Figure 1.5: Multiple points of contact requiring separate local models for each inter-action point
A PI controller is also implemented in order to reduce the error between the proxy position and the haptic device position, by controlling the deformation and force back to the user. For a softer object deformation can be more heavily relied on then for a stiff object as seen in Fig. 1.6
Figure 1.6: Difference in distance between proxy and haptic device position for stiff or soft objects
To model intersample behaviour, [17] and [18] use a low order linear approxima-tion. Balance model reduction is also performed on the linear approximation as it is still too computationally expensive to model the entire virtual environment. This lead to the conclusion that the balanced model reduction produced almost the same reduced model as simply choosing a local area around the point of interaction. It is noted however that this would not necessarily be the case in a object that had
inhomogeneous material, especially with a object that has a stiff outside layer and a softer inside.
A discussion is presented in [20] of the issues related to interaction between a soft object and a stiff object. Since the soft object is sufficiently modelled at a lower frequency than the stiff object. Collision response is performed at both frequencies using data structures to share constraints between the loops.
The division of the fast and slow rate loops is approached in [21] with a client-server based architecture. The client-server performs the computationally expensive object deformation and send and intermediate representation of the virtual environment to the client. The client uses this to render the haptic force feedback. In order to compensate for delay and its effects on the generated force, adaptive stiffness adjusting is done. The stiffness is also modulated and a moving average filter applied to avoid sudden stiffness jumps.
1.2.3
Instability
A haptic system has difficult stability characteristics due to the human operator being introduced to the system. It is not sufficient for the haptic interaction to be stable in isolation as the dynamics change with a human operator. The notion of passivity is therefore useful in ensuring the stability of these systems [3]. A virtual environment may be caused to be non-passive or unstable due to sampled and hold operators, sampling delay, large virtual impedance value, quantization error [23].
Hybrid switching systems are analyzed by [24] that have a controller to switch between operating modes in order to account for time delays and react to user actions. Passivity of a system can be guaranteed if a global Lyapunov function is found. It is found that it is sufficient for hybrid systems to have each operating mode as well as the switching behaviour of the system be passive. This is the case because the combination of passive systems is also passive.
The problem of the uncertainties associated with having a human operator in the system are addressed in [23] by developing a human response model. The human arm impedance is modelled as a linear time-invariant 2nd order model. The effective force from the operator is broken down into two components: the residual force and the reaction force. The reaction force takes into account how the operator responds, adapting to the simulated force. There is also damping included in the model, however this is quite variant between individuals. It is found that by including the human
operator in the model less conservative stability constraints can be applied than those from passivity alone.
Linear virtual coupling is also expanded on by using non-linear virtual coupling. In this case the stiffness in the coupling is increased over a number of time steps shown in Fig. 1.7 to increase the stiffness that can be handled by the device and human.
Figure 1.7: Non linear virtual coupling
A passive observer is introduced by [25] which is used to monitor the energy flowing in and out of subsystems. In addition, a passivity controller is put in place as dissipative element to adapt to changes in the energy flow and absorb net output that has been measured by the observer. The controller has two methods of responding, either by conserving velocity or force. This can run into the problem that the force that would be required to dissipate the energy is greater than actuator limits. Due to this the energy may need to be dissipated over the course of a number of time steps. When this must be done, instability can initially occur though this instability is bounded, it is still too much to be pleasant for the operator. Another suggestion noted would be to add a rule to reset the observer to zero, potentially in the case where free motion is detection.
It is noted in [26] that because passivity only represents a portion of the whole sta-ble region it suffers from intrinsic conservatism. Inspired by hysteresis non-linearities, a less conservative input-to-state stable approach is proposed. Input-to-state stability required that for a large sampling time, the size of the state must be bounded by a function of the inputs. A system is input-to-state stable if and only if it is dissipative.
A sufficient condition for dissipativity is that the available energy in the system be finite for every initial condition.
The virtual environment can be regarded as a one port network. The use of both discrete and continuous time domains as the user interacts with the virtual environment is a major source of energy leak. This energy leak results in the pressing path, as the user contacts a virtual object, being different from the releasing path, as the user loses contact, for position-force behaviour. The behaviour in Fig. 1.8 displays counter-clockwise hysteresis.
Figure 1.8: Hysteresis from different pressing and releasing path for position-force behaviour
The property of hysteresis that he time derivative of input-output pairs being sector bounded is used. A structure that allows a finite extractable energy to pass through the system is proposed, by guaranteeing a lower bound to make the system dissipative.
Passive-set-position-modulation is used in [27] to allow a continuous time robot to track the position to track a slow updating set position signal. The position of a virtual mass moves as close to the desired position as possible subject to passivity, via a series of spring coupling with damping injection. At each interval where the spring coupling is fixed, the system behaves passively, but upon switching, the spring energy jumps. Either enforcing or breaking passivity in doing so. The modulation strategy is applied during the jumps and allows the jump energy to be that which could be dissipated by the damping of the haptic device in the previous time step.
In [28] and [29] an energy bounding algorithm is introduced to restrict the energy generated by the zero order hold. The bound is found based on the amount of energy that can be dissipated by the physical damping in the haptic device. The bound is implemented by placing a restriction on the change in the actuator force between time steps. An important difference between this method and virtual coupling is that the virtual environment is not assumed to be passive, which may happen in certain surgical scenarios (i.e. beating heart).
The idea of bounding the force is extended in [1] by applying it to the 3 DOF case, still with non-deformable virtual objects. Two approaches are presented, one more conservative than the other. A more in-depth look at into the force bounding methodology will be done further in this thesis in Section 2.2.1.
1.3
Objectives of the Thesis
The objective of this thesis is to address some of the problems that still remain in haptic interactions for surgical simulation by using a local model implemented a multi-rate simulation with the use of a bounding algorithm to maintain stability. Specifically, the need for a local model that can ensure passivity of a haptic interaction and that can be implemented universally to any virtual environment regardless of the physical model used or its complexity. The force bounding technique from [1] is used as the basis to calculate the bound on the force feedback to the user. In addition the thesis also addresses problems felt by the user when having a multiple DOF virtual environment. When updates to a virtual environment are slow, in addition to creating instability due to the magnitude of the force feedback, the direction of the force can also cause the feedback to the user to jump and be less smooth than desired.
1.4
Overview of the Thesis
The remainder of the thesis is comprised of the following sections:
• Chapter 2: Methodology - An overview of the general system and how the local model and force bounding are implemented. Extension of previous force bounding method being applied to a deformable virtual environment. Interpo-lation is used in the local model for both bounding the force magnitude and
change in direction.
• Chapter 3: Experiments - Experiments done using a haptic device demon-strating the implementation in a 1-DOF and 3-DOF virtual environment. In the 1-DOF experiments it is found that applying the force bounding approach in the local model greatly decreases the oscillations felt by the user. While in the 3-DOF case, the importance of limiting the change in the direction of the force between time steps is noted by large jumps in haptic device position seen by the user when the limitation is not in place.
• Chapter 4: Conclusion - Summary of the work done in this thesis as well as problems encountered and opportunities for further development
Chapter 2
Methodology
2.1
System Overview
Looking at the overall system to be examined, there is a user that is using a haptic device to interact with a deformable virtual environment. In order to provide realistic force feedback to the user, the haptic device must be updated by a force signal at a constant rate of approximately 1 KHz. The virtual environment that will be modelled must also be deformable, and due to this will have a much slower update rate, dependant on complexity, from 1-100 Hz. It is also desirable that this deformable virtual environment is able to be visualized and thus a rate of at least 20 Hz is desired. In order to bridge the gap between the different update rates, an intermediate buffer will be used, referred to hereafter as the local model. This local model will represent a small area of the virtual environment that is approximated with simpler dynamics as to speed up calculations. In Fig. 2.1 the interconnection of the system is displayed, where Fh is the force applied by the human operator, Fd is the force feedback to the
Category Notation Example
Vector bold lower case font b
Matrix bold upper case font B
Scalar regular font b or B
Directional vector ending point/starting point (as subscript) rA/P
Further description subscript bd
Time step parenthesis containing the index of the time step b(n)
user by means of an actuator force applied by the haptic device, xd is the haptic
device position, xV E is the geometry information from the virtual environment, and
K is the object stiffness received from the virtual environment.
Figure 2.1: Interconnection and data exchange between the user, haptic device, local model, and virtual environment
Since this buffer is asynchronously updated by the virtual environment, it acts as a switching system. It is important that the system is stable to prevent harm to the user and have a realistic feeling interaction. Since passivity is a sufficient condition for stability, we can ensure that the system is passive. For a switching system to be passive, each mode as well as the switching must be passive. The energy in the system must be tracked in order to ensure that an excess of energy is not being transferred to the user.
2.1.1
Virtual Environment
The virtual environment receives and an input from the rest of the system the current position of the haptic device and the force being applied by the haptic device and the force being applied onto/by the environment, as estimated in the local model. The virtual environment must take these inputs and perform a deformation calculation so that the local model can be updated accurately. In order to perform this deformation the virtual environment must be able to detect where collision is occurring between the haptic device and the virtual environment. The environment force must then be applied to the object at this location. Then finally the deformation calculation can occur to find the new shape of the virtual object. In order to accomplish these things, an appropriate physics-based model is required for the virtual environment. In this case, a finite element model is used, where the virtual object is represented by elements and nodes.
2.1.2
Local Model
The local model is a separate thread from the virtual environment that runs at a much faster rate. It receives asynchronous updates from the virtual environment related to the local geometry. Additionally, the local model is updated with position information from the haptic device. By separating the local model from the virtual environment, the environment force is able to be computed at a haptic rate of 1 kHz due to the dynamics approximations used to simplify calculations. The local model receives the nodes from the surface element closest to the haptic device as well as the normal to that surface. This information is able to be used for simplified collision detection. Based on the haptic device position, penetration depth is able to be determined, which is used in calculating the environment force. The velocity of the haptic device is also useful in determining if the motion of the haptic device is entering or leaving the virtual object.
The overall goal of the local model is to synchronize the local model surface position to the virtual environment. When an update is received from the virtual environment after a deformation, the synchronization must be done passively without injecting energy into the loop. Therefore, the local model is also used to monitor the energy flow during each mode (i.e. when there is no deformation change) as well as during a switching instance (i.e. an update is received from the virtual environment).
2.2
Stability
2.2.1
Overview of Force Bounding Approach
In order to maintain stability, in the local model, an extension of the force bounding approach in [1] is implemented. In their paper, the force bounding approach is used to maintain passivity of a non-deformable virtual environment that loses information through a sample and hold operator. For the system to be passive, the energy flow into the haptic system (HSCV in Fig. 2.2) must not be negative; there must not be energy flow out of the system to the human operator.
Figure 2.2: Overall configuration of haptic system as described in [1]
For the case of a single-DOF haptic device with mass and damping elements, the energy flow for 0 ≤ t < nT , where T is the size of the time step, can be written:
E(n) = Z nT 0 m ˙vd(t)vd(t)dt + Z nT 0 bvd2(t)dt + Z nT 0 Fdh(t)vd(t)dt ≥ 0, (2.1)
where E(n) is the energy flow into the system at the n-th time step, m is the mass of the haptic device, vd is the velocity of the haptic device, b is the physical damping
of the haptic device, Fh
d is the hold actuator force.
The first term in Equation (2.1) represents the energy stored via inertia. If it is assumed that the initial velocity is zero, this value can be neglected for passivity considerations as it is always finite and non-negative as demonstrated below:
Z nT 0 m ˙vd(t)vd(t)dt = 1 2mv 2 d(n) ≥ 0 (2.2)
The second term from Equation (2.1) represents the energy that is dissipated due to physical damping in the haptic device. Using the Cauchy-Schwarz inequality 2 vectors u and v with an inner product space Equation (2.3), the second term from Equation (2.1) can be written as shown in Equation (2.4).
(2.4) Z nT 0 bvd2(t)dt = n−1 X k=0 b Z (k+1)T kT v2d(t)dt ≥ n−1 X k=0 b T " Z (k+1)T kT vd(t)dt #2 = n−1 X k=0 B∆x2d(k + 1) ∆xd(k + 1) = xd(k + 1) − xd(k), (2.5)
where B = b/T representing the damping present over a single time-step.
The last term from Equation (2.1) represents the energy flow into the discrete portion of the system, which can be written as Equation (2.6) when a zero-order hold operator is used. (2.6) Z nT 0 Fdh(t)vd(t)dt = n−1 X k=0 ( Z (k+1)T kT Fdh(t)vd(t)dt ) = n−1 X k=0 ( Fd(k) Z (k+1)T kT vd(t)dt ) = n−1 X k=0 Fd(k)∆xd(k + 1)
Now, the sufficient condition for passivity can be stated as:
E(n) ≥ E1(n) ≡ n−1 X k=0 B∆x2d(k + 1) + n−1 X k=0 Fd(k)∆xd(k + 1) ≥ 0, (2.7)
where E(n) is the energy flow into the system at the n-th time-step, E1(n) is the
reduced energy flow term at the the n-th time-step, B is the physical damping over one time step, ∆xd(k + 1) is the change in haptic device position from the current
time step to the next time step, Fd(k) is the actuator force of the haptic device at
the current time step.
Rewritten for the multi-DOF case, the condition becomes:
E1(n) ≡ n−1 X k=0 ∆xd(k + 1)TB∆xd(k + 1) + n−1 X k=0 Fd(k)T∆xd(k + 1) ≥ 0, (2.8)
where B is a diagonal matrix with entries B along the diagonal. Fd(k) is the
three-dimensional interaction force of the k-th control time step, and xd(k + 1 is the
three-dimensional displacement of the haptic device at the k + 1-th control time step. Since the energy flow into the haptic system should be assessed at each time step, Equation (2.8) can be calculated recursively:
E1(n + 1) = ∆xd(n + 1)TB∆xd(n + 1) + Fd(n)T∆xd(n + 1) + E1(n) ≥ 0 (2.9)
Now we have an equation that can be used to determine the feedback force needed at a given time step to ensure passivity until the next time step. However, ∆xd(n + 1)
is unknown and thus, the value chosen for Fd(n) must satisfy the passivity condition
for all possible ∆xd(n + 1). Using Cholesky decomposition such that B = CCT. The
passivity condition can be again rewritten:
(2.10) E1(n + 1) = ∆xd(n + 1)TCCT∆xd(n + 1) + 1 2Fd(k) T C−1C∆xd(n + 1) +1 2Fd(k) TC−1 C∆xd(n + 1) + E1(n) ≥ 0 (2.11) E1(n + 1) = ∆xd(n + 1)TCCT(k)∆xd(n + 1) + 1 2Fd(n) TC−1 C∆xd(n + 1) +1 2Fd(n) T C−1C∆xd(n + 1) + E1(n) +1 4Fd(n) T(CTC)−1 Fd(n) − 1 4Fd(n) T(CTC)−1 Fd(n) ≥ 0 1 4y(n) Ty(n) + E 1(n) − 1 4Fd(n) TB−1 Fd(n) ≥ 0 (2.12) where y(n) = 2CT∆xd(n + 1) + C−1Fd(n) (2.13)
Since y(n)Ty(n) is non-negative for all possible ∆xd(n+1) then if Equation (2.14)
is true, Equation (2.12) is also true. E1(n) −
1 4Fd(n)
T
To make the haptic system passive, the condition in Equation (2.14) must be satisfied for all n. Thus the magnitude of the feedback force should be upper-bounded by the total energy accumulated up until the previous step. If the desired feedback force does not satisfy this condition, then the bound in Equation (2.15) needs to be imposed on the magnitude of the actuator force.
ud(n)TB−1ud(n)
4 fd(n)
2 ≤ E
1(n), (2.15)
where fd(n) and ud(n) are the magnitude and unit vector of the feedback force
re-spectively, i.e., Fd(n) = fd(n)ud(n).
The new bounded energy flow into the haptic system equation becomes:
E1(n) = ∆xd(n)TB∆xd(n) + Fd(n − 1)T∆xd(n) + E1(n − 1) (2.16) fd,max(n) = s 4E1(n) ud(n)TB−1ud(n) , (2.17)
where fd,max(n) is the maximum allowable feedback force magnitude at the n-th time
step based on the passivity requirements.
2.2.2
Force Bounding Approach applied to Deformable VE
In the force bounding approach presented in Section 2.2.1 there is no local model needed as the virtual environment is non-deformable. For a deformable virtual envi-ronment, the force bounding will be implemented within the local model, and there are additional consideration regarding passivity and user feeling when an update is received from the virtual environment. Since Equation (2.18) is used to calculate the feedback force to the user in the local model, this force can be reduced either by reducing the stiffness or by reducing the penetration depth. In order to be able to include visualization between virtual environment time steps and for a smooth transition, the penetration depth is chosen to be reduced.
Fd(i) = K(I)dLM(i)nLM(i) (2.18)
dLM(i) =
(xd(i) − xLM(i)) · nLM(i)
||nLM(i)||2 , (2.19)
where i is the time index for the time step in the local model, I is the time index for the time step of quantities received as updates from the virtual environment, Fd
is the force feedback to the user, K is the object stiffness received from the virtual environment, dLM is the penetration depth of the haptic device into the local model
constraint, nLM is the normal to the constraint in the local model, xd is the haptic
device position, and xLM is some point on the local model constraint plane
To maintain the convention that the energy flowing into the haptic system is positive; ∆x(i), the change in haptic device position between the current and previous time step in the local model, is defined to be negative when the haptic device is moving in the direction opposite to the normal of the plane of interaction, i.e., into the virtual object, and positive when moving in the same direction as the normal, i.e., out of the virtual object. Since the actuator force is in the same direction as the normal, the second term in Equation (2.16) will be negative when the haptic device is moving out of a virtual object. If this is the case, the actuator force will should be subject to the force bound when moving out of the virtual object. However when moving into the virtual object all terms are positive so it is guaranteed to be passive.
At each local model time step, the energy flow into the haptic system is updated by Equation (2.16).
Then the position of the local nodes, A(i), B(i), C(i) and the associated normal are used to determine if the haptic device is in contact with the local wall. A vector, rA/P(i) is formed between between one of the nodes A(i) and the haptic device position
P (i). The side of the plane the haptic device is on is determined in Equation (2.20) by checking the sign of the dot product between the normal of the plane (positive facing out of virtual object) and rA/P(i):
rA/P(i) · nLM(i) =
≤ 0 on same side of plane as the object
> 0 not on the same side of plane as the object
(2.20)
The projection of rA/P onto the normal as in Fig. 2.3 is used to find the penetration
Figure 2.3: Projection of the vector between the haptic device position and a point on the plane onto the normal of the plane
Once it has been determined that the haptic device position is on the interior side of the plane, the position is projected onto the plane as P0. The projected point is checked using Equations (2.21)-(2.23) to ensure it is on the right side of each boundary as shown in Fig. 2.4.
(rA/B × rA/P) · (rA/B × rA/C) ≥ 0 (2.21)
(rB/C× rB/P0) · (rB/C × rB/A) ≥ 0 (2.22)
Figure 2.4: Projection of haptic device position onto local plane used to determine if it lies within the other three nodes
If in addition to being on the right side of the plane, Equations (2.21)-(2.23) are all satisfied, the haptic device is deemed to be in contact with the local model. If the nodes received from the virtual environment are different than the local nodes, then these steps are repeated to find the penetration depth with the virtual environment, and if there is contact. If the nodes are the same then the same values are used without additional calculations.
Once this information has been gathered it becomes possible to check if the force bound will be satisfied and what adjustments need to be made in the local model to interpolate between its current position and the required updated position. First after the desired actuator force is calculated Equation (2.24), a bound on the magnitude of the force is found in Equation (2.26) and the corresponding allowable penetration depth calculated in Equation (2.27) is checked against the penetration depth in the virtual environment.
Fd(i) = K(I)dV E(i)nV E(I) (2.24)
dV E(i) =
(xd(i) − xV E(I)) · nV E(I)
||nV E(I)||2 , (2.25)
where Fd is force feedback to the user, K is the object stiffness received from the
virtual environment, nV E is the normal to the plane of interaction from the virtual
environment, xd is the haptic device position, and xV E is some point on the plane of
interaction received from the virtual environment.
fmax(i) =
s
4E1(i)
nV E(I)TB−1nV E(I)
, (2.26)
where fmax is the maximum allowable magnitude of the feedback force to the user to
ensure passivity, E1(n) is the reduced energy flow into the system at the n-th time
step, nV E is the normal to the plane of interaction from the virtual environment,
B is a diagonal matrix with components along the diagonal indicating the damping between local model time steps.
dmax(i) =
fmax(i)
K(I) , (2.27)
where dmax is the maximum allowable penetration depth of the haptic device based
on the restrictions imposed by the passivity condition.
If the penetration depth from the virtual environment is too large, and there is also currently contact with the local model, then interpolation should be performed between the local model and virtual environment nodes to satisfy the force bound.
In addition to interpolating the penetration depth, for a multi-DOF virtual envi-ronment the force can also be smoothed by interpolating between the local normal and the updated normal from the virtual environment . A small allowable change in angle φ is selected. The angle between the local model normal nLM(i) and the normal
of the updated plane from the virtual environment nV E(I) is calculated as follows:
θ(i) = cos−1(nV E(I)TnLM(i)), (2.28)
where θ is the angle between the local model normal vector and the virtual environ-ment normal vector, nV E is the normal to the plane of interaction in from the virtual
environment, nLM is the current normal being used to represent the constraint in the
local model.
If the angle θ is greater than the allowable angle φ, then interpolation of the normal is needed. This is done by multiplying the local model normal by the rotation matrix in Equations (2.30)-(2.31), where the axis of rotation, µ calculated in Equation (2.29), is the cross product between the local normal and virtual environment normal, and
the angle of rotation is φ.
µ(i) = nLM(i) × nV E(I) ||nLM(i) × nV E(I)|| (2.29) Rµ(i),φ(i) = cos(φ) + µ2 xα µxµyα − µzsin(φ) µxµzα + µysin(φ)
µyµxα + µzsin(φ) cos(φ) + µy2α µyµzα − µxsin(φ)
µzµxα − µysin(φ) µzµyα + µxsin(φ) cos(φ) + µ2zα
(2.30)
α = (1 − cos(φ)), (2.31)
where φ is the maximum allowable angle for the normal to change between time steps, µ is the axis of rotation and subscripts the projection of that vector onto the x, y, or z axis.
The new local normal becomes:
nLM(i + 1) = Rµ(i),φ(i)nLM(i) (2.32)
Once the new local normal has been established, the new local model node posi-tions can be found. The interpolation is shown in Fig. 2.5. The new local model plane is located the maximum allowable penetration depth away from the haptic device po-sition and has the new normal found in the previous step. The location of the new nodes A(i + 1),B(i + 1), and C(i + 1) are found by finding where the line between the corresponding nodes (i.e. A(i), the current local node and A(I) the updated virtual environment node) intersects with the new plane.
Figure 2.5: Interpolation between the current local nodes and the virtual environment local nodes
In the case where there is no current contact with the local model but is contact with the update from the virtual environment; interpolation of the normal is not required, and the normal corresponding to the nodes from the virtual environment can be used directly. The new local model node positions are found by moving each of the nodes from the virtual environment moved opposite to the normal in Equation (2.33) such that the maximum allowable penetration depth is satisfied.
A(i + 1) B(i + 1) C(i + 1) = A(I) B(I) C(I)
− (dV E(i) − dmax(i)) 1 1 1 nV E(I) T, (2.33)
where A(I), B(I) and C(I) are the node positions provided by the update from the virtual environment, positions A(i + 1), B(i + 1), and C(i + 1) are the new positions of the nodes in the local model, dV E is the penetration depth of the haptic device
with the virtual environment plane of interaction, dmax is the maximum allowable
Figure 2.6: Interpolation between the current local nodes and the virtual environment local nodes
If there is neither contact with the virtual environment nor with the current local model, then the local model nodes and normal are directly updated to be the same as those received from the virtual environment.
In summary, the local model behaves according to Fig. 2.7, where dLM is the local
model penetration depth, dV E is the virtual environment penetration depth, dmax is
the maximum penetration depth, i is the current time step in the local model, I is the time step of the most recent update from the virtual environment, E1 is the reduce
energy flow into the system, fmax is the maximum magnitude of the actuator force, θ
is the angle between the local model normal and the virtual environment normal, φ is the maximum allowable change in angle to the local model normal between time steps, nLM is the normal to the constraint in the local model, nV Eis the normal to the plane
of interaction in received from the virtual environment, µ the normal vector to nLM
and nV E, Rµ,φ is a rotation matrix, Fdis the actuator force, [A(i + 1)B(i + 1)C(i + 1)]
are the new node coordinates for the local model to update, [A(I)B(I)C(I)] are the node coordinates received from the virtual environment, K is the virtual environment stiffness.
Chapter 3
Experiments
3.1
Simulation Models
A Simulink model is used in order to conduct experiments to validate experimentally the passivity of haptic interaction with a deformable virtual environment through force bounding in a local model. This model consists of multiple S-functions used to separate the local model calculations from the virtual environment. A finite element model is used in the virtual environment in order to compute the deformations to the virtual objects.
A Novint Falcon 3-DOF haptic device is used to implement the experiments. Experiments where conducted first with only 1 DOF then expanded to use 3 DOF.
3.2
1 DOF experiments
In the experiments, two threads run in parallel: a local model thread that computes the force bounding feedback to the users at a fixed 1 kHz rate; and a virtual en-vironment thread that, for the purpose of these experiments, consists of a 1-DOF finite element model artificially slowed down at a variable rate of 20 − 50 Hz. In the finite element model, the virtual constraint is modelled as a linear spring with stiffness K = 10 kN/m. This stiffness is chosen sufficiently high to destabilize the haptic interaction in the absence of force bounding even when constraint deformation is not considered. A physical damping of b = 0.03 Ns/m is assumed for the Novint Falcon. The virtual wall is located at x = 0 m.
the results of successive experiments comparable. In all experiments, the device is first placed at the fully retracted position x = −0.061 m, with its end effector unobstructed. It is then released and a constant force of gravity moves its end effector down into the horizontal virtual environment. In addition to making them repeatable, this manner of controlling the experiments ensures that the physical damping comes only from the Novint Falcon. Adding a user into the haptic system will add physical damping and will only increase the stability of the interaction.
Fig. 3.1 depicts the results obtained when force bounding is not used in the local model and the local constraint is directly updated by the virtual environment. In other words, the constraint in the local model follows exactly the constraint in the virtual environment, with only a fixed delay equal to one time step of the local model 0.001 s at each update of the virtual environment. Note the significant oscillations that occur, with a peek to peek amplitude of roughly 5 mm. Further, the deformation of the virtual environment never settles. However, when a user pushes the device into the environment with a constant force as in this experiment, the device should settle into static contact after an initial transient.
Figure 3.1: Haptic interaction with a deformable virtual environment without force bounding: positions of the haptic device (dash-dotted red line) and of the virtual constraint in the local model (solid black line), as updated by the virtual environment.
Fig. 3.2 illustrates that the oscillations of the haptic device upon contact with the virtual environment die out quickly when force bounding is applied both on the con-straint updates and when the device moves away from the local concon-straint. The lack of oscillations increases the realism of the interaction with the virtual environment.
Figure 3.2: Force bounding on constraint updates and when moving out of the straint: positions of the haptic device (dash-dotted red line), and of the virtual con-straint in the local model (solid black line) and in the virtual environment (dashed black line, mostly overlapped with the solid black line).
Fig. 3.3 plots the behaviour of the virtual constraint in the local model in greater detail. As expected, the local constraint follows the update from the virtual environ-ment only in as much as permitted by the passivity requireenviron-ment when the update leads to contact with, or increased penetration into, the virtual environment of the haptic device.
Figure 3.3: Detail of Fig. 3.2: positions of the haptic device (dash-dotted red line) and of the virtual constraint in the local model (solid black line) and in the virtual environment (dashed black line). Note the interpolation in the local model between environment updates when force bounding occurs.
3.2.1
Conclusion
The experiment done in 1 DOF has shown the importance of force bounding when interacting with a stiff deformable virtual wall. Without the use of force bounding, the system is unstable and oscillations felt by the user are quite intense. It has also shown that it is possible to implement the force bounding in the local model while allowing the virtual environment to remain at a slower update rate.
3.3
3 DOF experiments
Instead of a simple linear spring model for the virtual environment. The complexity of the experiment is increased by instead using a sphere formed by 16 tetrahedral elements as the virtual object, Fig. 3.4. This results in 11 nodes and 16 outer surface triangles. In order to prevent the virtual object from moving out of the area of interaction available to the haptic device, 2 of the interior nodes are fixed to have zero displacement.
Figure 3.4: Virtual object composed of 11 tetrahedral elements.
This increased complexity add a need for more rigorous collision detection. A symmetrical object is used in the experiment however in order to keep the collision
straightforward. The surface triangle on which collision takes place is determined by checking the normal vectors of each surface triangle and determining which most closely matches the haptic device position relative to the center of the object. In these experiments, in addition to the force bounding done in the 1 DOF experiments, experiments are also done to validate the use of interpolation of the local normal in order to have a smooth feeling change in the force direction. In order to show larger deformation in the 3 DOF case the experiment is done by having a user start outside of the virtual object and attempt to move along one axis (z-axis) towards the center of the virtual object. Similarly to in the 1 DOF experiments, when looking at the penetration depth, the difference between the local model constraint position and the constraint from the virtual environment can be assessed.
When the updates from the virtual environment come slowly or it is a very de-formable virtual object, the new constraint normal updated from the virtual environ-ment can be large enough to cause a severe disturbance to the haptic device position. In the case of Fig. 3.5, the updates come from the virtual environment every 0.5 s. This results in large jumps to the position of the haptic device when there are large changes in the updated normal in Fig. 3.6.
Figure 3.5: Jumps are seen in the position of the haptic device at times corresponding to updates coming from the virtual environment.
Figure 3.6: Large changes to the normal sent from the virtual environment after a deformation has occurred.
In order to remedy this issue, interpolation of the normal vector in the local model is done to smooth out the change in force direction felt by the user. When normal interpolation is utilized in the local model, a maximum allowable change to the angle of the normal direction can be put in place to prevent these jumps.
Figure 3.7: The position of the haptic device changes smoothly as the user tries to move along the z axis.
Looking at the change in the updated normal from the virtual environment in Fig. 3.8a there are still large changes to the normal vector that occur when receiving an update from the virtual environment. However due to how this jump is treated in the local model Fig. 3.8b the user only experiences small changes to the normal between time steps. (Note: the plots are viewed in x-y coordinates for simplicity to the viewer, but they also have a z component which causes the angle between local model updates to be inconsistent).
(a) Virtual Environment
(b) Local Model
Figure 3.8: Updates to the normal vector after a large deformation from the virtual environment and how it is treated by the local model
3.3.1
Conclusions
From the 3 DOF experiments, it is found that not only is it important for user experience to bound the magnitude of the actuator force, but it is also useful to interpolate between normal vectors to improve the quality of the haptic interaction. When the normal is allowed to be used directly from the results of the deformation in the virtual environment, significant jumps are experience in the haptic device position, due to sudden changes in direction. These large jumps are likely caused by the user applying an input force to counteract the current direction of the actuator force then suddenly receiving the actuator force from a different direction. This issue can be solved by limiting the change in the normal direction between time steps so that the user does not perceive a large shift.
Chapter 4
Conclusion
This thesis is devoted to implementing a force bounding approach and local model interpolation to guarantee passivity and improve smoothness of haptic interactions with arbitrary deformable virtual environments. The different time steps needed for deformation calculations and for haptic force feedback are addressed through the use of a local model. This local model is specifically designed to be applied with any arbitrary virtual environment. Regardless of the physically based model used for deformation, the local model simply requires information about the geometry near the haptic interaction point (normal to constraint and outer boundary contours) as well as local stiffness information. As such the virtual environment creator has no need to worry about passivity issues in their model as these can be addressed in the local model.
The first part of this work focussed on applying the force bounding techniques developed in [1] to an arbitrary virtual environment that may be deformable. This was done through a local model that limited the penetration depth of the haptic device in the virtual object based on the energy flow in the system. In the next portion the problems seen in the direction of the force when the virtual object becomes 3 DOF are addressed. This is done by limiting the change in angle of the force between local model time steps. This prevents large changes from the virtual environment from causing the haptic device to jump in position.
The experimental results suggest that by implementing these techniques, the hap-tic interaction remains passive and is more smooth than without. It was also seen however that there are some limitations in the ability to create a reasonably fast, deformable, virtual environment using Matlab/Simulink.
environment such as SOFA. Current SOFA documentation is limited however making it difficult to integrate with. There may be other deformable model that could be integrated with as well. This work can also be extended to deal with multiple points of interaction. If this multiple interaction points are on the same tool, the effects of including torque in the energy calculations must be thoroughly addressed.
Appendix A
A.1
Collision Detection extension
The first step in determining where a collision occurs is to find on which surface plane of the virtual environment the interaction takes place. This can be a very slow process and difficult, especially for complex object. In order to simplify this, only uniform objects centered around some point will be considered.
A vector, rP /C indicates the displacement from the center of mass, C, to the
current haptic device position, P . This vector is compared to the normal vector of each surface plane normal. The plane which has a normal closest to the vector rP /C
will be considered the plane of interaction. Once the appropriate plane is determined this will be saved to send back to the local model.
In the case of a more simple 1 DOF virtual environment with a single wall, this first step can be skipped as the plane of interaction is known.
The next step is to find if the haptic device position is inside or outside of the virtual object. The haptic device position must be on the side of the interaction plane opposite to the normal (sufficient condition for a 1 DOF environment) as well as being between the nodes that form the surface plane.
Once it is determined that the haptic device position is inside the object. An avatar position can be determined (i.e. the closest position on the surface of the object to the haptic device position) to use as the point where a force will be applied. This is done by finding the interception point between the plane and a line from the haptic device position in the direction of the surface normal.
A.2
Finite Element Model
As a quick overview of the Finite element method, the simple case of a 1 DOF object is presented. This may be done by having multiple elements, such as springs, in series as shown in Fig. A.1.
Figure A.1: Simple FEM model with two spring elements. Element 1 between nodes 1 and 2, element 2 between nodes 2 and 3
Each spring element has a stiffness matrix Equation (A.1), where K is the stiffness coefficient of the spring for that element.
"
K −K
−K K
#
(A.1)
These element stiffness matrices must then be appropriately assembled into a single matrix to represent the entire object in Equation (A.2). In this case the matrix will overlap at the 2nd row, 2nd column as node 2 is shared between elements.
K1 −K1 0 −K1 K1+ K2 −K2 0 −K2 K2 (A.2)
In order to solve for the node we need to use Equations (A.3)-(A.4). The displace-ment at node 1 and 2 are x1 and x2 respectively. At node 3 the wall is fixed therefore
the displacement is known to be 0. Similarly there is no force at the second node, as the force from element 1 on element 2 is equal and opposite to the force from element 2 on element 1. There is only a force applied at node 1 and the reaction force at node 3.
K1 −K1 0 −K1 K1+ K2 −K2 0 −K2 K2 x1 x2 0 = F1 0 F3 (A.4)
In order to solve for the displacements in this equation, the finite element method step of elimination is performed. First the known displacement is used to eliminate the 3rd column, then the row of the force which has no information (row 3) is eliminated in order to create a reduced stiffness matrix than can be inverted Equation (A.5).
" K1 −K1 −K1 K1+ K2 #−1" F1 0 # = " x1 x2 # (A.5)
Returning to our virtual environment model, it is possible to assemble our stiffness matrix and compute the inverse of the reduced stiffness matrix off-line before the real-time interaction begins. This is important as 3D objects can have much more complex stiffness matrices. Though it is still needed to determine on which nodes force is applied in order to determine the displacement at each node.
Since the environment force and the collision location have already been deter-mined, this is simply a matter of determining what portion of the force is applied to each node. If there is only 1 DOF, where there is only one ”surface” node, or if the collision point is directly on a node this is simple as the entire force is applied at that node.
The when there is a plane with multiple nodes associated with it the percentage of the total environment force to be applied at each node is determined by comparing the distances from the collision point to each node. Such that the node that is closest to the collision point will have the highest portion of the force applied. The direction of the force will be applied in the direction opposite to the normal vector of the surface plane of interaction.
Then the displacement at each node is calculated. In order to find the new node positions these displacements are added to the current node positions.
Once this has completed, the virtual environment sends back to the local model, the new geometry (nodes and normal associated with the plane of interaction) and any relevant local stiffness information if it is variable in the virtual environment.
Bibliography
[1] J. Kim, S. Y. Baek, and J. Ryu, “A force bounding approach for multi-degree-of-freedom haptic interaction,” IEEE/ASME Transactions on Mechatronics, vol. 20, no. 3, pp. 1193–1203, June 2015.
[2] C. Basdogan and A. A. Srinivasan, “Haptic rendering in virtual environments,” in in Stanney, K. (Ed.), Handbook of Virtual Environments, Lawrence Erlbaum, Inc, 2002, pp. 117–134.
[3] J. Colgate, M. Stanley, and J. Brown, “Issues in the haptic display of tool use,” in Intelligent Robots and Systems 95. ’Human Robot Interaction and Cooperative Robots’, Proceedings. 1995 IEEE/RSJ International Conference on, vol. 3, Aug 1995, pp. 140–145 vol.3.
[4] S. Cotin, H. Delingette, and N. Ayache, “Real-time elastic deformations of soft tissues for surgery simulation,” Visualization and Computer Graphics, IEEE Transactions on, vol. 5, no. 1, pp. 62–73, Jan 1999.
[5] M. Mahvash and V. Hayward, “High-fidelity haptic synthesis of contact with deformable bodies,” Computer Graphics and Applications, IEEE, vol. 24, no. 2, pp. 48–55, March 2004.
[6] J. Barbiˇc and D. L. James, “Real-time subspace integration for st. venant-kirchhoff deformable models,” in ACM SIGGRAPH 2005 Papers, ser. SIGGRAPH ’05. New York, NY, USA: ACM, 2005, pp. 982–990. [Online]. Available: http://doi.acm.org/10.1145/1186822.1073300
[7] C. B. P. D, “Real-time simulation of dynamically deformable finite element mod-els using modal analysis and spectral lanczos decomposition methods,” 2001.
[8] B. Bickel, M. B¨acher, M. A. Otaduy, W. Matusik, H. Pfister, and M. Gross, “Capture and modeling of non-linear heterogeneous soft tissue,” in ACM SIG-GRAPH 2009 Papers, ser. SIGSIG-GRAPH ’09. New York, NY, USA: ACM, 2009, pp. 89:1–89:9. [Online]. Available: http://doi.acm.org/10.1145/1576246.1531395 [9] M. Bro-nielsen and S. Cotin, “Real-time volumetric deformable models for surgery simulation using finite elements and condensation,” in Computer Graph-ics Forum, 1996, pp. 57–66.
[10] O. R. Astley and V. Hayward, “Real-time finite-elements simulation of general visco-elastic materials for haptic presentation,” 1997.
[11] F. Mazzella, K. Montgomery, and J.-C. Latombe, “The forcegrid: a buffer struc-ture for haptic interaction with virtual elastic objects,” in Robotics and Automa-tion, 2002. Proceedings. ICRA ’02. IEEE International Conference on, vol. 1, 2002, pp. 939–946 vol.1.
[12] D. L. James and D. K. Pai, “A unified treatment of elastostatic contact simulation for real time haptics,” in ACM SIGGRAPH 2005 Courses, ser. SIGGRAPH ’05. New York, NY, USA: ACM, 2005. [Online]. Available: http://doi.acm.org/10.1145/1198555.1198614
[13] Y. Adachi, T. Kumano, and K. Ogino, “Intermediate representation for stiff virtual objects,” in Virtual Reality Annual International Symposium, 1995. Pro-ceedings., Mar 1995, pp. 203–210.
[14] R. Balaniuk, “Using fast local modeling to buffer haptic data,” in The Fourth PHANTOM Users Group Workshop -PUG99, 1999.
[15] F. Barbagli, D. Prattichizzo, and K. Salisbury, “A multirate approach to hap-tic interaction with deformable objects single and multipoint,” Contacts, The International Journal of Robotics Research, pp. 703–715, 2005.
[16] F. Barbagli, K. Salisbury, and D. Prattichizzo, “Dynamic local models for stable multi-contact haptic interaction with deformable objects,” in Haptic Interfaces for Virtual Environment and Teleoperator Systems, 2003. HAPTICS 2003. Pro-ceedings. 11th Symposium on, March 2003, pp. 109–116.
