Advanced Solutions in Object‐Oriented Mechatronic Simulation

(1)

Budapest University of Technology and Economics Dept. of Control Engineering and Information Technology

Advanced Solutions in Object‐Oriented Mechatronic Simulation

Summary of PhD Theses

Tamás Juhász

Supervisor: Dr. László Vajta, associate professor,

Dept. of Control Engineering and Information Technology

Budapest, 2008

(2)

List of theses

Thesis 1... 6

Thesis 2.1... 11

Thesis 2.2... 14

Thesis 3... 17

(3)

1. Introduction

The virtual engineering methods applied in the computer‐aided product development can yield to a synergic outcome: the planning, analysing and optimizing of a product in virtual reality can be carried out faster in multiple engineering domains at the same time. This is emphasized as approximately 70‐80 percent of the final production costs are determined during the product development. During the first phase there is a possibility to analyse the behaviour and optimize various parameters of a mechatronic system using its functional model, before any physical prototype would be manufactured.

Nevertheless, today’s widely spread CAD‐construction systems do not allow integrated testing of the product being developed by any means of multi‐domain simulation. Based on large amount of data the functional simulation model of the product can be created manually: the more detailed this model is, the better it can approximate the real world behaviour. So far such modelling tasks demanded tremendous amount of manual calculations from a human expert.

In order to describe and simulate a functional model of a system, a big number of descriptive languages and simulation‐environments exist already. Such a simulation model has to be both syntactically and semantically correct, which requires extreme attention if it is created manually. As the development process is usually iterative (e.g.: changes in the requirements, varying some parameters during optimization) the continuous tuning of the model multiplies the risk that the human modeller can fail. This explains that there is a large demand for a computerized solution for the translation from a construction model to a multibody simulation model. With this help all the relevant physical parameters of a construction can be embedded largely automated in the simulation.

In majority of the simulation cases the contacts between moving mechanical parts and / or the environment has to be detected and their physical effect must be taken into account. The simulation of industrial manipulators or mobile robots lays stress on the performance and accuracy of this contact processing task.

One practical application of mechatronic simulation is the planning of a teleoperation [7] mission of a mobile platform by means of virtual reality techniques (e.g.: for training purposes). In such cases the proper visual feedback is very important for the human operator. During the real‐time visualization of a mechatronic simulation it can be taken into account that the relative motion and the changing occlusions in the sight have already a strong depth clue, which can assist the operator to navigate in the virtual reality quite well. The binocular techniques using auxiliary devices for stereoscopic visualization can also help, but these have usually more disadvantages: wearing extra glasses can be inconvenient and the simulator‐sickness (because of the eyes’ lack of ability to focus on any virtual object) can occur already after 10‐15 minutes of usage.

(4)

4/24 Summary of Theses

2. Aims of research

The dissertation emphasizes three scopes of virtual engineering problems introduced in the previous section:

• A solution is needed for automating the task of creating a functional simulation model directly out of a mechanical CAD construction plan and offering the possibility to extend the model with other mechatronic elements.

• The imperfection of the chosen multi‐domain simulation environment – having lack of support for contact processing – has to be solved and the virtual contact model needs to be compared with real world contact problems.

• A mechatronic simulation in an interactive virtual environment can be used for virtual teleoperation trainings. It is known that in general the communication bandwidth between the remote system and the operator can be rather low, limiting the (visual) information being transmitted. The usage of monocular visualization (demanding less bandwidth than the binocular stereoscopic methods) can be an alternative way, in case we assure that it presents enough depth information for a human operator navigating in the 3D environment.

3. Chosen methods and tools

During my research I followed the sequence of literature overviews, problem definition, program‐, system‐ and experiment design, early implementation tests, approval tests, evaluation and case studies comparing real system behaviour.

The Pro/ENGINEER CAD system and the Dymola simulator are the two software tools that define the start and end of my presented automated model translation task. These are chosen for multiple good reasons.

3.1. The Pro/ENGINEER CAD system

The PTC (Parametric Technology Corporation [13]) is a partner of such major international automotive manufacturers as Audi, Volkswagen, Skoda, Seat, Maserati, Porsche or Toyota.

The Pro/ENGINEER system is very popular in Hungary with both the local affiliates of multinational companies and the national supply chain enterprises.

Audi Hungaria Motor Kft., the company with the largest export production in Hungary, and the third largest engine manufacturer in the world is one of the biggest Hungarian Pro/ENGINEER users. They are planning further investment and will continue to increase their weight in Hungary, as a result of their recently announced joint research and development project with the BME University. As Pro/E is a trend‐setter in the Hungarian industry, I chose it as the source CAD system for the model translation process.

(5)

3.2. Modelica, the object‐oriented physical modelling language

The synergic outcome that was mentioned in the introduction is promoted by the object‐oriented modelling paradigm. It was infiltrated from the software engineering also to the modelling of physical system behaviour. This concept allows reusing and integration of acquired knowledge of previous projects and modelling domains. By choosing the modelling language (Modelica) and its supporting simulation environment (Dymola) I focused on the complexity and diversity of mechatronic systems I would like to simulate.

The Modelica standard is developed continuously for almost 10 years by an international non‐profit association [14]. This modern, object‐oriented modelling language allows describing a complex mechatronic system over multiple engineering domains. It is based on bi‐directional mathematic equations, whereby the causal relations are only determined by the direction of information flow.

Besides process and control system‐oriented components it supports modelling of mechanical‐, electronic‐, hydraulic, pneumatic and thermal subsystems of a complex mechatronic system.

The Modelica description is based on ordinary differential equations (ODE) and differential‐algebraic equations and it can also contain state machines or Petri‐

nets [15]. The building blocks in Modelica can connect to each other with uni‐ or bi‐directional connectors, too. The latter types of connectors implicitly define the Kirchhoff‐laws on connecting nodes (the value of potential‐variables are the same and the sum of flow‐variables equals to zero in each connected node).

The bi‐directional information‐flow between the components and the declarative mathematical description of Modelica is emphasized: in comparison to that, Matlab / Simulink system supports causal description with assignments only.

3.3. The Dymola simulator

The Dymola (Dynamic Modelling Laboratory) simulator is a highly developed commercial environment for Modelica‐based multi‐domain simulation [16]. It utilizes the complete Modelica standard, thus allowing mechatronic simulation based on a general mathematical description of mechanical, electrical, thermodynamic, hydraulic, pneumatic, thermal, power and control components.

As Dymola automatically manipulates the equation system of a complex model, the number of equations and the rank of the system can be significantly reduced, thus the non‐linear solver can simulate the model in real‐time, too.

Moreover there is a possibility for embedding a complete Modelica model into a Simulink block, integrating the advantages of both simulation environments.

(6)

4. Scientific outcome

According to the presented motivation, my research can be divided into three main parts: I summarize these in the following main sections.

4.1. Automated model generation based on a CAD construction

In order to analyse the dynamic behaviour of macroscopic multi‐body systems in advance, there is a need for functional simulation models. These can be used already in the early product development phase in order to detect construction problems, thus can spare costs before any physical prototype is manufactured.

At this time the creation of such models based on a CAD construction is a challenging manual problem for modelling engineers. My first statement relates to the automation of the model translation process:

Thesis 1 [1][2][3] – Chapter 2 of the dissertation

a) I have created a new, automated method that translates a CAD construction plan into a consistent mechanical model in Modelica language. During this method I build a special, hierarchical graph structure: on the one hand it eases handling kinematical loops, and on the other hand it can also be used later in the effective visualization of the hierarchical model.

b) In my extended method active components can be chosen from a multi‐domain model library interactively, and these parametric models can be assigned to the nodes of my mechanical model’s hierarchical graph. This extended method can serve the complete functional mechatronic model of the CAD construction. □

4.1.1. Source of information

In Pro/ENGINEER the geometry information of the original K construction (hierarchy, transformations and shapes) stay also in VRML format at our disposal.

The kinematical (initial poses, joint axes) and dynamic (centres of mass, masses, inertias) parameters are temporally stored in a hierarchical XML file. These files are the source of information during the model translation process.

4.1.2. Building a tree structure

In my new method I create a new F tree graph, the nodes of which refer to the original components of K (they are either joints or rigid bodies). The edges of F represent the relations between these components (including the hierarchical structure). All of these relations are extracted directly out of the XML data.

Preserving the original kinematical structure of K asked for attention. As the F graph is a tree, it cannot contain a cycle. However in K there could be more kinematical loops: these have to be cut along special joints (refer to cutJoint in the thesis) into multiple branches of F. The cutJoints are located always in leaf‐nodes

(7)

of F and are detected during the construction of F. These joints have a reference to their second original neighbour node in the graph with a special virtual edge. In F these virtual edges close the original kinematical loops of K.

The following example can help to understand the cutJoint concept:

Figure 1: Example construction (8 bodies + 9 joints) and its F graph

Two kinematical loops in this example articulated structure are cut by the „D‐‐

C” and „F‐‐G” revolute joints forming a tree with two main branches and introducing two virtual edges.

4.1.3. Assigning geometry to bodies in F

For rigid body models there is an obvious need for geometries in order to visualize the construction during a simulation. They are also essential when contact processing (including collision detection) is also considered (refer to section 4.2).

The F scene graph – that I interposed in the translation process – yields to a spanning tree of the components in K. The T subset of nodes in F refers to the rigid bodies in the K construction. So far these are point‐mass objects only. Besides the XML description, there are VRML files originated from the CAD construction (section 4.1.1), the structure of which follow exactly the hierarchy of CAD assemblies and parts.

Storing multi‐body information in such a tree structure has a performance advantage during visualization (section 4.3), as the most graphical APIs support hierarchical transformations on the objects being visualized. For example it can be enough to rotate a revolute joint’s axis to transform all child‐segments connecting to that certain joint.

As it is described in my work, I can assign VRML nodes to the respective point‐

mass bodies in T. As the VRML description of geometries is also hierarchical, such an assigned geometrical node might also have subordinate child geometries. If at that node the VRML hierarchy is deeper than the XML hierarchy, there can be multiple geometrical shapes associated to the same body, as well.

This multiple geometry‐association possibility carries also a great advantage: as there is a possibility to model a concave shape with multiple, adjoined convex ones (refer to concaveÆconvex partitioning in [24]), the geometrically more robust, highly developed convex collision detection algorithms can be used in the contact processing (see section 4.2), allowing a more stable simulation.

(8)

4.1.4. Creating a Modelica model based on the scene graph structure

In case of the target Modelica language description each node of F is transformed to a respective model instance in the M mechanical model. All the necessary parameters for these instances are originated from the XML description.

The example graph on Figure 1 has 17 nodes, which translate to 8 rigid body instances and 9 abstract revolute joint instances.

The original topology (kinematical structure) of K must be represented by flat point‐to‐point connect equations in M. In other words: while processing the edge list of F (also including the virtual edges), for both nodes at the endpoint of each edge the associated Modelica model instances’ frame connectors must be connected with an ideal rigid junction element’s two frame connectors. A connect(X, Y) equation operator in Modelica is used to define Kirchhoff laws on both related X and Y connectors implicitly: the potential variables are equal, while the flow variables sum to zero in connectors.

By extending the F tree with its virtual edges, you can get the F’ graph, which is isomorph to the topology of K. Using a breadth‐first search of F’, the C set of kinematical loops of K can be found: each loop in C has to be analysed. Each cycle in F’ must contain a cutJoint, as there is no cycle in the F tree, and – by definition – a cutJoint cuts a cycle in F’ into two branches in F.

Although the equation system in the background has an algebraic loop in case of structures having kinematic loops, Dymola is able to resolve this problem in the background (without the user’s intervention), as long as these loops are not planar ones. In a planar kinematic loop (e.g.: consisting 4 revolute joints, where all joint axes are parallel to each other) the contained joints have only a single degree of freedom. Within mechanical systems containing planar kinematic loops there is no longer a unique mathematical solution and the symbolic algorithms of Dymola will fail. At such a case one has to use a special planarCutJoint instance explicitly at exactly one joint in the planar loop in order to support the solver of Dymola.

My cutJoint detection rule negotiates between the normal and the planarCut joint instances in a kinematic loop (based on its planar status), thus there is no need to correct the final model manually. My new, completely automated model translation procedure allows building a consistent M dynamic model in Modelica language, which can represent the functional behaviour of the original K construction.

4.1.5. Extending the model with subsequent active components

It is obvious that the passive M mechanical model by itself cannot be analysed remarkably: in order to inspect its dynamic behaviour, there must be an external source of energy that sets the system in motion. This can be achieved when the components of the pure mechanical model (more precisely the contained joint instances) are extended with an A set of active elements of other engineering domains (e.g.: with electromagnetic actuators). A superior control system can

(9)

influence the behaviour of these new active elements. By adding the A set to M, the derived result of my method is the final M’ mechatronic model. M’ will be the multi‐domain model that is simulated in Dymola environment.

As the M’ Modelica model inherits the base M mechanical model, a containing relation exists between their set of equations: H_M ⊂H_M_' The inheritance allows the M’ and M models to keep some level of independence: in my thesis I distinguished between two major conditions to represent this phenomenon:

• As long as K remains unchanged (both the structure and the names of all components are constant), the elements in the A set can be freely replaceable and M’ can be recreated without affecting M (thus keeping it intact).

• If a parameter of K is varied only, it influences M exclusively. As long as the topology of F’ (the structure of K) stays the same, there is no need to recreate M’.

Concrete examples:

• By designing servo‐drives the fine‐tuning of electrical parameters (changes in A) won’t require to recreate M, as only M’ will be updated (via A).

• Varying some material parameters in K relates to changes of dynamic parameters (masses, inertias). If A is kept constant, the changes in K trigger only the recreation of M, so the old M’ can be reused with the new M, as they remain

“compatible”.

The conditional independency of M and M’ in my extended method is emphasized:

this way the individual domains of a mechatronic system can be optimized separately.

(10)

4.2. Extending Modelica models with contact processing

The Modelica.Multibody library [17] contains the models that are used in the standard Modelica implementation of mechanical multi‐body simulation. At the present time only the models of ideal rigid bodies and joints with permanent degrees of freedom stay at our disposal. The new elements I had to create for the model translation process (which is related to my Thesis 1) are also based on these standard implementations.

Surprisingly there is no standard support in the library for collision processing between the mechanical models. According to the lack of any retaining forces the standard bodies can interpenetrate freely into each other. There have been multiple suggestions to solve contact problems in special cases in a Modelica environment (e.g.: for tyre + road contact), but these neither give a general solution, nor build on standard models. The standard support for collision response needs a solution for many years.

Three years ago a proposal appeared in the literature for solving the general problem [18], but it was not mature enough. The solution I present in my thesis is different from all previously discussed ones.

4.2.1. Mechanical interactions

There are two major groups of possible interactions within a mechanical system:

• The mechanical joints define continuous constraints between the relative movements of their coupled mechanical elements.

• A mechanical contact occurs when surfaces of multiple bodies are touching. Two phenomena can be distinguished according to the duration of interaction:

o During a resting / sliding contact static and dynamic friction forces arise, the duration of which can be the whole inspected time

o A collision contact is usually only a temporal reaction between two colliding bodies

In majority of the dynamic multi‐body simulation cases the handling of mechanical contacts (both contact detection and contact response calculation) between the virtual bodies is a key importance task. Between the performance and accuracy of a simulation a golden middle way has to be found, especially when real‐time simulation and interactivity is demanded (e.g.: in virtual teleoperation, refer to section 4.3, Thesis 3).

(11)

4.2.2. Tasks of contact processing

The general handling of mechanical contacts in a virtual environment bases on multiple subtasks:

• The solution has to be backwards compatible with the standard Modelica description of rigid body mechanics, in order to give models of older projects a fair chance to introduce it by easily replacing some basic elements with their extended counterparts.

• A central task is to trace the pose of all possible colliding mechanical parts (including their motion) and in case of each overlapping contact, extract all necessary information for the contact model (contact points, penetration depth and ‐direction, relative velocities, cross‐section area or ‐volume, etc.). For determination of most of these contact features there are mathematical libraries (both free and commercial ones) available on the market.

• Based on the available contact information, a response reaction has to be introduced in the system. According to the actual material properties of all pairs of bodies in contact, the magnitude and direction of corrective forces must be calculated in order to prevent further interpenetration.

Thesis 2.1 [2][3] – Chapter 3 of the dissertation

In Modelica environment I separated the description of the standard rigid body model from multiple optional collision shapes. Joining upon Thesis 1, this new method allows the decomposition of complex shapes (being extracted from a CAD system) into multiple convex geometries, creating a better basis for a more stable – geometrically more robust – collision processing. □

4.2.3. Determining the contact parameters

There are a few articles about detecting the status of overlapping between spatial shapes with various accuracy, however the majority of which is not related to Modelica at all. There are some analytical results [20][23], but mostly only approximations are considered [21][22][25], as they can improve simulation performance a lot. As a result of these researches there are multiple free and commercial software libraries on the market that support collision detection. They offer well‐defined C or C++ interfaces, thus I could integrate some of those into Modelica environment.

In order to unify the calculation of contact features in Modelica environment I created a common wrapper interface that allows reusing software packages of various researchers.

(12)

4.2.4. Contact response calculation

I discussed two major well‐known methods for handling mechanical collisions:

the impulse‐ and the force‐based theories. After introducing these methods and pointing out the disadvantages of the impulse‐based technique, I chose the force‐

based calculation.

Following the force‐based technique keeps the physical background of mechanical contacts: a body conveys energy to its colliding pair, during a finite time of compression and decompression phases. Meanwhile the work of the retraining F collision force – that prevents the bodies from further penetration – consumes energy and produces heat. There are some (boundary) conditions on F:

• Before and after the contact the magnitude of the vector F is zero.

• During the whole contact time F is directed to separate the bodies (no adhesion is modelled)

• After the collision the relative speed of the bodies is decreased: the overall energy of the system decreases after each contact event, because some amount of kinetic energy is transformed to heat.

• When two point‐mass bodies are colliding ideally, the impulse‐preservation law must be valid

• If a T body is resting on a horizontal plane, the system balances its mg weight with an opposite direction, but same magnitude F=‐mg reaction force.

The method I discuss in my dissertation applies virtual, interim spring and damper elements between the contacting bodies:

Figure 2: Contact model based on spring and damper forces

In the following paragraphs I describe the collision case of two elastic, isotropic, homogenous bodies (A and B), which are separated by an F responsive force having a normal and a tangential component.

Let’s take a contact point pair (p^A és p^B) between bodies A and B. Note that when there are more contact point pairs between those, the simple vector‐sum of the calculated F_i forces yield to the result F force. There is a p penetration vector between these contact points (pA–pB), the magnitude of which relates to the p penetration depth: the shortest distance A and B can be separated along that direction. The collision detection step I mentioned before can serve this vector for various primitive and triangular mesh geometries.

(13)

Figure 3: Determining the normal direction

The direction of p matches the normal direction of F: it will act in the pB point on body B and -F will act in the pA point on body A, respectively (see Figure 3).

The magnitude of F needs to be calculated.

As I investigated the function F(p), I proved that p (the penetration normal) alone is not enough contact feature information to determine F. By introducing a common ε restitution factor (which depends on the materials of both A and B) I could confirm that besides the penetration depth and normal, the relative speed of the contact points is also essential to calculate the normal component of F.

The pA and pB contact points have a relative velocity vrel, its projection in the tangential plane is:

(

rel

)

rel t

rel v e e v

v = − _n⋅ _n• (1)

where eⁿ is the (unit length) normal direction vector of the penetrating contact:

it is defined at t=0 to be parallel to the vector p^A–p^B. The tangential contact speed

t

vrelis the magnitude of the (bold written) velocity vector:

t rel t

vrel = v (2)

If this speed is greater than zero, a frictional force arises in the contact’s tangential direction: along an e^t unit length vector in Cartesian space:

t rel t rel

t v

e = v (3)

In my model the static and dynamic friction forces are switched at a relative tangential speed limit:v_rel^t^,^lim, which is dependent on the contacting materials. A force acts along the e^t tangential direction as either static or dynamic friction. Until the tangential relative speed limitvrel^t^,^lim, a virtual spring from a fixed anchor point causes an increasing static friction arising between the colliding bodies.

Figure 4: Static and dynamic friction effects with interim springs

If the colliding bodies in the tangential plane have larger relative velocity thanv_rel^t^,^lim, a dynamic friction force will replace the static one, and the “anchor”

point is sliding together with the actual contact point in the tangential plane.

(14)

Thesis 2.2 [5] – Chapter 3 of the dissertation

For calculating the normal‐ and tangential reaction forces in simulated mechanical contacts in Modelica environment, I introduced a coupled spring + damper material model based on new material properties:

Normal contact force:

⎪⎩

⎪⎨

⎧

≤

⎟⎟ >

⎠

⎜⎜ ⎞

⎝

⎛ + −

=

0 ) ( , 0

0 ) ( , ) ( )

1 ( ) 1

( 0

t x

t x t x S

t t x

F_NORMAL ν ε & ^CONTACT

ε

(4)

where ν0 is the normal component of the relative collision impact speed at the start of a contact, S^CONTACT is a common spring stiffness factor [N/m], ε is the restitution factor (both depend on the materials of A and B), x is the penetration depth and x& is the normal component of the relative velocity of the contact points.

Tangential contact forces:

A static friction force is caused by a virtual spring having stiffness SFRIC that is anchored in the contact point and extends to the following length:

τ

≤

=

∫

^v ^t ^dt ^t

t l

t t

rel( ) , 0 )

(

0

( )

otherwise v v F if

t l t S

F

t rel t NORMAL rel

static static FRIC

FRICTION

lim

0 ,

, 0

),

( ) min

( < <

⎩⎨

⎧ ⋅ ⋅

= μ

(5)

(6) At larger relative tangential velocities the dynamic friction force replaces the static one.

Using a material‐dependent μdynamic dynamic friction coefficient:

otherwise v F v

t F

t rel t NORMAL rel dynamic

dynamic FRICTION

lim ,

, 0

) ,

( >

⎩⎨

⎧ ⋅

= μ

(7) The final contact response is the vector sum of the normal component along en and the two (alternative) tangential forces along et, which also acts as torques on the (cA, cB) centres of mass:

( )

⎩⎨

⎧

⋅

>

+

⋅ +

= ⋅

otherwise t

F

v t

F t F

t t F

NORMAL

t rel dynamic

FRICTION static

FRICTION t

NORMAL CONTACT

, )

(

0 ,

) ( )

( )

) ( (

n n

e

e F e

( )

(

_B _B

)

CONTACT B

A A CONTACT

A

t t

CONTACT CONTACT

p c F

τ

p c F

τ

−

×

=

−

×

−

=

) ( )

(

) ( )

( □

(8)

(9)

(15)

4.2.5. A case study

A complex case study of my thesis is the six‐legged mobile robot – Anton, which is developed in the RobotsLab [19] of the Otto‐von‐Guericke University of Magdeburg together with the Fraunhofer Institute IFF.

I created the complete mechatronic model of the robot using my automated process stated by Thesis 1 (on page 6). On the left side of next figure there is a photorealistic rendered image of Anton’s construction in the CAD system, while on the right side a screenshot of Dymola can be seen, presenting the virtual contact forces during a simulation:

Figure 5: The construction and the simulation model of Anton

In order to analyse the motions of the robot I had to implement the onboard control system in Modelica environment, including the motion‐planning algorithms with inverse kinematics – the details of these tasks can be read in the dissertation.

One of the inspected motions was a forward walk with tripod‐gait mode on a horizontal plane, started from the stable base stance. I compared the results of the simulation with the data of the real contact force sensors in the feet. The next figure shows the normal component of both the virtual and the real contact forces, respectively:

(16)

Figure 6: Simulation and real sensor data of a tripod gait of Anton

On this figure it can be observed that the weight of the “ant‐robot” Anton is balanced in a way that the front and middle legs carry a little bit more load than the rear ones, because on the head there are two (stereoscopic) cameras mounted.

As the Start command was executed after 0.5 seconds, it can be seen that the load on each feet started to change periodically. The robot was walking forwards in a tripod‐gait mode: triple‐support way of insect‐like walking, using the following way: while three legs are pushing the body forward they are supporting the complete weight of the robot (approx. 7,8 kg) the other legs (lifted in the air) can swing forwards rapidly. Because of the relocation of the robot’s centre of mass, the load on the rear feet is always decreasing while it is increasing on the front ones. As the other three legs accelerate during swinging forward, small peaks can be observed in the force values. During the swing phases of the observed three feet (the front / rear left and the middle right ones) there is no relevant, measurable contact force, of course.

Figure 6 compares my simulation results with the real robot’s behaviour, validating my contact model that I introduced in Modelica environment with Thesis 2.2, showing a very good match to the reality.

0 ⁵ ¹⁰ ¹⁵ ²⁰ ²⁵

-5 0

5

10

15

20

25

30

35

40

45 ANTON - Simulation

time [s]

Force [N]

0 ⁵ ¹⁰ ¹⁵ ²⁰ ²⁵

-5 0

5

10

15

20

25

30

35

40

45 ANTON - Sensors

time [s]

Force [N]

Front Left Middle Right Rear Left

(17)

4.3. Spatial visual feedback in teleoperation

Generally the phrase “teleoperation” can indicate circumstances in which a device or machine is operated by a person from a distance. Before a teleoperation mission involving expensive robot‐systems would start, it is often that the operating personal does some training in a simulated environment. This is a high‐

level application of mechatronic simulation (refer to Thesis 1 and 2), which must be also real‐time and interactive.

Nevertheless the 3D visual feedback is also very important in teleoperation.

The human operator must be able to navigate in the virtual reality by interpreting a large amount of visual information. For navigation the correct localisation of possible obstacles demands significant visual depth information, as well.

The human visual depth perception has multiple components that can be grouped in two sensory categories of either monocular or binocular effects:

• The monocular effects can be perceived also by a single eye: linear perspective, size‐consistency, occlusion / overlapping, aerial perspective, shading / shadows and motion parallax.

• All binocular clues are based however on the retinal disparities of both concurrently stimulated eyes.

The superposition of multiple monocular and binocular effects leads to a complex depth perception process in our brain. The relative strength of some binocular and some static monocular clues have already been analysed [28][29].

4.3.1. Posing the problem

The amount of visual information that can be sent during a teleoperation mission can be limited by the capacity of the transmission channel. In comparison to the binocular methods, the monocular ones need definitely fewer resources.

Therefore the use of pure monocular visualization can be considered, if it can present adequate depth information that is necessary for the task.

The most people are able to navigate a vehicle in a dynamic virtual environment by using monocular information exclusively. The motion‐parallax effect – meaning a steady change of angular position of two observations of the same object projected to the retina while the observer and the environment have a relative velocity – is established to present a relatively strong depth clue:

Thesis 3 [4][6] – Chapter 4 of the dissertation

Proceeding from the important role of motions during the interactive simulation of mobile mechatronic systems I recognized that besides the monocular visualization requires less resource compared to the binocular stereoscopic techniques, it has also enough strength in creating adequate depth impressions in the operator. □

(18)

During my literature research I did not find any relevant articles about the quantitative analysis of the intrinsic depth effect caused by the motion. For the applicability in mobile teleoperation a measurement was demanded: thus it can be assured whether the monocular visualization has enough strength for that purpose.

In Chapter 4 of my thesis I present a series of psychophysical experiments that were done at the department. These experiments involved 44 voluntary student of the University, and were organized by our 3DMR laboratory with the help of my colleagues, Tamás Urbancsek and Dr. Ferenc Vajda. The purpose of the experiments was the measurement of the intrinsic (objective) quality of depth information presented by the conventional rendering technique including the optional usage of auxiliary stereoscopic devices:

• Anaglyph spectacles

• LCD shutter glasses

• Polarization goggles

I combined all of these enumerated techniques with two known camera models and an optional relative motion (introducing the motion‐parallax effect). For each possible combination all voluntary students are prompted to give subjective scores during a 10 minute experiment (pro person) judging its strength to reproduce depth sensation.

Using a full‐factorial analysis I could extract a quantitative measure for the strength of motion‐parallax in 3D visualization:

4.3.2. Design of the experiment

As depth perception is a subjective quantity, one can have the question whether it is measurable at all or the measurements delivered by two different scenes are really comparable. Hereafter it is shown that it is possible to build an appropriate measurement model.

The first problem is that the attendants cannot quantize exactly their depth experiences even on an exactly defined scale, thus it can happen that one estimates the same experience differently at a later time. A bigger issue is that the same scene generates different depth sensation level to the observer under the same circumstances but at a later time. The first impression with a novel technique may cause over‐estimation of the experience.

Therefore, the measured value depends on the preceding scenes viewed by the experimental person. This effect can be modelled with an internal state of the observer person.

Def.: Constellation: the measurement arrangement that includes the visualization technique, the presented 3D scenario and the environment (related later to the x input vector). The observer person is separated from that.

(19)

Figure 7: Model of a series of experiments

• x denotes the parameters of the scene arrangement (which technique is used, is the virtual camera moving, depth clues).

• c denotes all the other parameters of the measurement arrangement, which affect the intrinsic depth perception level and had to be set to constant during the series of experiments.

• I is the intrinsic depth perception level. I want to estimate this value at the different scene arrangement parameter values (x).

• u influences the initial q0 state of the observers (e.g.: they were always told the same a‐priori information), it can be taken also constant during an experiment.

• z values are the direct answers of the observer (absolute and relative depth sensation levels).

• Every person has his internal habits: a personal bias (b, the liability to use higher numbers on a scale) and a personal interpretation of scale (a, a tendency to exaggerate).

The expected values of the a and b factors are 1 and zero, respectively. It was assumed that these factors are constant for each individual observer. These individual properties are treated as a linear transformation of the intrinsic depth sensation. The x input vector of the constellations has three factors:

• Applied 3D technique (x¹)

o Spectral separation (anaglyph): A

o Time‐shared technique (with LC shutters): S o Polarization‐based: D

o Conventional technique: M

• Camera model (x2)

o Monocular, central: 0 o Binocular, convergent: T o Binocular, parallel: P

x input factors

Constellation

I = F(x, c) I intrinsic

effect

Observer (O) z = a(I + ξ) + b

a(O), b(O) ξ (O, q(t),u(t))

Internal state: q(t)

z measured effect uexternal factors

c constant parameters

ξ noise

(20)

x³

N R

x¹

A S D M

0 T

P x²

0 1

2 3

4 5

6 8 7

9

10 11

12 13

14 15

16 17

18 19

Figure 8: The GX constellation‐graph

• Dynamism (x3)

o Moving scene: R o Static scene: N

Four combinations of x have no meaning, thus there are 20 possible constellations that are represented by the nodes of the GX graph on Figure 8. The adjacencies between the nodes are due to the following requirement: each person must evaluate all constellations once, and between two adjacent cases only a single factor of the x input should change. An edge in GX

shows that the connected constellations differ only in either x1 or x2 or x3 factor.

The requirement of adjacent constellation sequences during the experiments demanded Hamiltonian paths (HPs) in GX: such a path contains 19 edges and passes each node in the graph only once, by definition. The start node of each path had to be selected randomly for a given observer, thus a heuristic algorithm had to be created that searches HPs in GX. It can be proven that there are less than enough independent HPs in the graph, thus I allowed my algorithm to skip one edge, resulting a quasi Hamiltonian path that uses only 18 real edges of GX.

Using an electronic questionnaire each person had to score two depth effect levels for each of the 20 constellations on a scale in the [1, 9] interval: the absolute zabs subjective and the zrel relative quality of depth (comparing always to the previous scenario).

4.3.3. Results

The answers of each observer have been logged by the application I created to conduct these experiments. I did the pre‐processing and evaluation of the results using the statistical toolbox of Matlab for the two main questions:

Absolute Z values:

First I eliminated the distortion of the individual factors (a and b) from the absolute answers (zabs) using a linear transformation (normalization).

In the next step I determined the I=F(x,c) function for each node of the GX graph (x⁰…x¹⁹) by fitting normal distribution functions for the normalized answers. By determining the mean values (presented in Table 1), I chose confidence intervals at 90% confidence level:

(21)

Table 1: Intrinsic depth perception levels F(X1,X2,N) 0 T P F(X1,X2,R) 0 T P

A 3.08 4.24 4.40 A 5.23 5.82 5.93 S 3.32 5.29 5.83 S 5.99 6.78 7.17 D 3.87 6.44 6.64 D 5.38 7.13 8.10

M 2.97 ‐ ‐ M 4.79 ‐ ‐

On a scale from 1 to 9 the motion (M0NÆM0R) improves +1.8 units in the depth perception levels compared to a still scene. Using auxiliary goggles the F(X1,X2,R) – F(X1,X2,N) differences lay at +1,8±0,9 units.

Due to the bigger distance between the observers and the canvas, the parallel camera model was more suitable: the difference caused by binocular information using the P projection model (F(X1,0,N) – F(X1,P,N)) is at +2,5±1,2 units. Thus it can be seen that the visual depth information presented by the motion‐parallax technique is in many cases as big as that of the binocular visualization methods.

It is also obtained that by simply wearing a “placebo” auxiliary spectacle (in A, S or D modes) the observer might have an illusion he / she perceives binocular information, although only a central 0 camera model was used in six constellations covering this psychic effect (e.g.: S0N or A0R). These constellations gained in average +0.4 units more with “placebo” goggles than without one.

Relative answers:

The zrel answers gave me a balanced set of information about the edges of GX.

Although my heuristic (quasi)Hamiltonian path‐search algorithm used random start nodes, each observer could score only 18 real edges in GX, while it has 52 edges altogether. Nevertheless the comparison of adjacent constellations was directly possible using these answers, as each edge of GX relates to a single factor change in the x input vector of constellations. The mean value of the scores of the 10 edges related to the change of the x³ factor (NÆR) gave the result that the motion‐parallax effect scored +2,05 (± 0,39 @ 90%) more units in the quality improvement level of depth perception. Similarly, by analysing the x1 factor, the usage of binocular visualization improves +1,91 (± 0,56 @ 90%) units in the quality.

4.3.4. Conclusion

Using statistical analysis of measurements gained from a series of psychophysical experiments I determined the quantitative strength of the motion among the depth clues of various visualization techniques.

Consequently, in a dynamic environment the motion parallax effect can compare to the effect of real stereoscopic visualization methods in creating depth impressions, thus in the application field of teleoperation the support of this monocular technique can be satisfactory.

(22)

5. Summary

I put the guide‐line of my thesis along a straight path of virtual engineering. It starts from the subtask of automated translation of a mechanical Pro/Engineer construction to functional mechatronic simulation models described in Modelica language (Thesis 1), which is extended by handling mechanical contacts (Thesis 2), and leads to a sensory assessment of visualization techniques during virtual teleoperation mission trainings (Thesis 3). I created an integrated mechatronic model authoring and visualization tool to encompass all of my research results under a common software environment. My system allows real‐time visualization of an interactive mechatronic simulation.

In our days the digital engineering becomes more important in the lifecycle of a product, including computer‐aided design, ‐production planning and

‐manufacturing, allowing a broad spectrum of analysis and optimization

techniques in these areas. The competition induced by the digital revolution challenges the small‐ and medium sized enterprises (SME), as well. It is rather often that these companies don’t have enough expertise on modelling and simulation of complex mechatronic systems, though they are usually contractors with big enterprises that can demand such knowledge from the suppliers.

The solutions discussed in my work can make the workflow between the design and the optimization more effective for the SME, thus they can be more competitive on the digital market of today. For example the smaller companies in the supply chain can have better interoperability with the larger enterprises, because they can send electronically a verified / optimized virtual model much earlier. All the supplied components can be virtually integrated into the final digital product, which can also be optimized along given criteria before any physical production would start.

6. Acknowledgement

Hereby I would like to express my sincere gratitude and appreciation to my director of studies and adviser, Dr. László Vajta, for his expert guidance and mentorship, his corrective suggestions for this dissertation and for his encouragement and continuous support at all levels since my undergraduate years and during my PhD studies.

I would like to thank my colleagues at the Department of Control Engineering and Informatics: Dr. István Loványi, for generously providing the advanced 3D visualization hardwares and Dr. Ferenc Vajda and Tamás Urbancsek for their kind assistance in organizing and conducting the psychophysical experiment and nonetheless for Miklós Vogel, who was my adviser for the master’s thesis in informatics for his suggestions on the direction of research on mobile robotic simulator systems.

(23)

Moreover I am very thankful for my Marie Curie programme‐supervisor, Prof.

Ulrich Schmucker that he offered me the opportunity to gather great experience on the research field of multi‐domain modelling and simulation as a fellow in the Early‐Stage Research Training mobility programme of the European Commission at the Virtual Engineering Expert Group of the Fraunhofer Institute for Factory Operation and ‐Automation in Magdeburg.

Related publications

[1] Juhász, T.; Schmucker, U.: “From Construction CAD to a Modelica model: Structural Manipulations throughout the Translation Process”; Proc. of 11^th IFF Science Days, 25^th‐26^th June 2008, Magdeburg, Germany

[2] Juhász, T.; Schmucker, U.: “CAD to SIM: CAD Model Conversion for Dymola‐based Mechatronic Simulation”; Proc. of 10^th International Conference on Computer Modelling and Simulation, 1^st‐3^rdApril 2008, Cambridge, England, ISBN 0‐7695‐3114‐8

[3] Juhász, T.; Schmucker, U.: “Automatic Model Conversion to Modelica for Dymola‐based Mechatronic Simulation”; Proc. of 6^th International Modelica Conference, 3^rd‐4^th March, 2008, Bielefeld, Germany, Vol. 2, pp. 719‐726.

[4] Vajta, L.; Urbancsek, T.; Vajda, F.; Juhász, T.: “Comparison of Different 3D (Stereo) Visualization Methods ‐ Experimental Study”; Proc. of 6^th Eurosim Congress on Modelling and Simulation, 9^th‐13^th September 2007, Ljubljana, Slovenia, Vol‐1 (abstracts), p. 351.

[5] Juhász, T.; Konyev, M.; Rusin, V.; Schmucker, U.: “Contact Processing in the Simulation of CLAWAR”; Proc. of 10^th CLAWAR International Conference, 16^th‐18^th July 2007, Singapore, pp. 583‐590.

[6] Vajta, L.; Juhász, T.: “The Role of 3D Simulation in the Advanced Robotic Design, Test and Control”; International Journal of Advanced Robotic Systems – Cutting Edge Robotics 2005, ISBN 3‐86611‐038‐3, pp. 47‐61.

Other publications

[7] Juhász, T.: “Áttekintés a telerobotikáról és a dinamikus robot modell paramétereinek identifikációjáról”; MicroCAD 2004 Nemzetközi Tudományos Konferencia, 2004. március 17‐

18, ISBN 9‐636‐61619‐1, Miskolc, Hungary, Vol. K, pp. 211‐216.

[8] Juhász, T.; Urbancsek, T.: “Beyond the limits of kinematics in planning keyframed biped locomotion”; Proc. of 6^th Eurosim Congress on Modelling and Simulation, 9^th‐13^th September 2007, Ljubljana, Slovenia, Vol‐1 (abstracts), p. 78.

[9] Juhász, T.: “Behavioral Simulation for the Fast and Flexible Manufacturing Enterprise”;

Production Systems and Information Engineering 2004, Volume 2, University of Miskolc, Hungary, pp. 121‐129.

[10] Juhász, T.: “Graphics Acceleration Techniques for a Mobile Robot Simulator”; Proc. of 7^th Central European Seminar on Computer Graphics, CESCG ’03, 22^nd‐24^th April, 2003, Budmerice, Slovakia, pp. 147‐156.

Advanced Solutions in Object‐Oriented Mechatronic Simulation