**Exploring Quantum, Classical and Semiclassical
Chaos in the Stadium Billiard
**Chris King

Mathematics Department – University of Auckland – 9 Jul 2009

**You can now run the wave function CAs in my open source application CA2D for Mac!**

*Abstract: This paper explores quantum and classical chaos in the stadium billiard using Matlab simulations to investigate the behavior of wave functions in the stadium and the corresponding classical orbits believed to underlie wave function scarring. The simulations use three complementary methods. The quantum wave functions are modeled using a cellular automaton (CA) simulating a Hamiltonian wave function with discrete (square pixel) boundary conditions approaching the stadium in the classical limit. The classical orbits are computed by solving the reflection equations at the classical boundary thus giving direct insights into the wave functions and eigenstates of the quantum stadium. **Finally a simplified semi-classical algorithm is developed to show the comparison between this and the quantum wave function method.*

* *

Fig 1: Evidence of scarring: Superposition of the first 3300 iterations of a wave function initiated by two different generating functions (a and c in fig 4) with wavelength chosen to be a fraction of the associated classical repelling orbit. Upper two figures show superposition of
phase amplitudes and lower two the corresponding superimposed probability
distributions of the squared amplitudes, showing evident scarring. Generally
the probability distributions give better resolution of the scarring.

* *

**Introduction:**

Estimation of the quantum wave functions in
chaotic systems, from the nucleus of atoms from helium to uranium to the
eigenfunctions of a wave-particle in the Bunimovich [i]stadium
(the classical Coliseum stadium shape consisting of a rectangle or square
capped off with semi-circular discs) proved to be initially intractable and
even in the post-modern era of revived semi-classical methods computationally
complex to the point of being counter-intuitive. Part of the purpose of this paper is to simplify
understanding these quantum solutions so their dynamics can be more directly
appreciated and understood.

A full quantum solution to a quantum
chaotic wave function, given its boundary constraints proved illusive for the
founders of quantum theory in the case of the helium nucleus. In 1970 [ii]
Gutzwiller [iii]
developed a semi-classical method using a trace formula, which integrates the
GreenÕs function of wave propagation over all coordinates so that it can be
expressed in terms of the repelling orbits hidden within the classically
chaotic version of the quantum system.

The trace formula for a particle of mass *m* and momentum inside a box with arbitrarily shaped walls is

[1]

where
is a smooth
function giving the mean density of states, summed over all orbits
of length
with
the phase shift
counting the focal points and twice the number of reflections off the walls and
over all retracings *k* of these orbits. The stability matrix
records the sensitivity of the orbit to changes in initial
conditions. From the left hand side of the equation it can be seen that if the
sum converges (this is generally possible only with some difficulty and depends
on reordering the sum terms), the eigenvalues will appear as singularities in
the sum. The general similarity of the trace formula to the Riemann zeta
function
[iv]
and this functionÕs GUE type statistics shared by quantum chaotic systems (fig
2e) has led to extensive attempts to solve quantum chaotic systems through the
paradigm of the zeta zeros.

Both quantum and semi-classical approaches
have been used to model systems such as the stadium billiard. In the quantum
approach a wave packet is propagated throughout the potential well using a time
dependent Schršdinger equation via Fourier transform techniques. This wave
propagation approach is discretely modeled in the CA below. Alternatively the
semi-classical approach uses the least repelling periodic orbits in the
classical stadium to generate amplitudes based on the path lengths of all classical
trajectories connecting a given location to another in the region.

Heller and
Tomsovic [v]
[vi]
among others have used the semi-classical orbit-based approach to develop the
theory of scarring of the quantum wave [vii].
The authors comment: ÒAccurate excited-state eigenvalues have been computed
from knowledge of relatively few periodic orbits. However, deep questions about
the convergence of the trace formula and its modifications remain unanswered. A
particularly complete expose of the methods involved is provided in their
Physical Review paper [viii].

Fig 2: Quantum chaos: (a) [ix],
the stadium is densely filled with repelling periodic orbits, three of which
are shown in black in (d). The quantum solution of the stadium potential well
(b) [x]
and (d) [xi]
shows ÔscarringÕ of the wave function along these repelling orbits, thus
repressing the classical chaos, through probabilities clumping on the repelling
orbits. A semi-classical simulation (c) shows why this is so. The quantum
solution is scarred on precisely these orbits (d). This causes to coincide with the eigenfunctions of the
repelling periodic orbits, just as the orbital waves of an atom constructively
interfere with themselves, in completing an orbit to form a standing wave, like
that of a plucked string. Over time, although the behavior may be transiently
chaotic, the quantum system eventually settles into a periodic solution or a
linear combination of these. Experimental realizations, such as a magnetized
electron in a quantum dot (e) [xii]
((*a*) and (*b*) give GOE and GUE statistics respectively), or the scanning tunneling
view of an electron on a copper sheet bounded by a stadium of carefully-placed
iron atoms (f) [xiii], confirm
the general picture, although quantum tunneling in the latter case leaked the
wave function outside too much to demonstrate proper scarring. The
semi-classical approach matches closely the quantum calculation (g).

To get a more direct picture of the
process of wave spreading and wave function scarring, this simulation depends
on a simple wave propagation formula for a cellular automaton which preserves
the essential nature of Hamiltonian dynamics and the GreenÕs function governing
transmission of wave amplitudes from point to point. Somewhat paradoxically, although the local rules are defined on a rectangular grid with only N, S E and W neighbors, the process well models wave transmission at a variety of angles and under situations where the curved surfaces of the potential well boundaries refocus the wave fronts continuously in various ways. Again this is possible although the local rules at the boundary provide for exclusively horizontal and vertical reflection.

In effect, the wave period and implicit frequency of the generating wave gives the energy and momentum of the transmitted wave-particle and the individual pixels on the boundary provide a simulation of the smaller quantum features at the atomic level occurring in any real physical implementation of the quantum stadium, as exemplified in even more rudimentary discrete form, in the atomic sized stadium of fig 2(f).

Fig 3: The wave function CA iterates in two stages. First in a > b white positions move using black neighbour values. Then in b>c black positions move. Finally in c>d white moves again. This enables a discrete simulation of a Hamiltonian process [xiv].

The cellular automation iterates the wave
function according to the rule:

[2]

As in figure 3, black positions whose four neighbors are white are all iterated in one step, followed by the white positions iterated in terms of their newly updated black neighbors. This provides a discrete iteration modeling a Hamiltonian process resembling action-angle variables with an angle shift of
.

In keeping with the quantized nature of the
system, the boundary of the region consists of a discrete series of vertical
and horizontal edges corresponding to a discrete region of square cells
approaching the curved stadium shape in the infinite limit. Vertical and
horizontal edges have a reflecting rule illustrated by the following two edge
and corner situations:

[3]

The process is also capable of handling
applied forces by integrating their effect:

[4]

The simulation uses a variety of
generating function to set up excitations whose period corresponds to a
fraction of the orbit length of a classical repelling orbit in the stadium.
Broadly they are of the type:

[5]

combining one or more Ôplane cosine waves
of appropriate period lying along an orbit path clipped being multiplied by a
circular Gaussian, although circular and uni-directional traveling wave packets
can also be established by clipping a sector of an emerging circular wave.

Six of the scars resulting from
superimposing successive probabilities measured as the square of the amplitude
are shown in fig 4. The predominant repelling orbits clearly show as
probability peaks coinciding with the corresponding repelling orbits in just
the matter of the semi-classical simulations already noted in the references
and discussion. In fig 5 are shown the corresponding generating functions for
each of the images in fig 4.

Fig 4: Scars of the wave stadium function illustrated
corresponding to 3,300 iterations of 6 distinct repelling periodic orbits in
the classical system. The wave function is iterated from the corresponding generating function in fig 6, whose wavelength and direction is chosen to correspond to a fraction of the orbit length. The scars are formed by superimposing the squared wave functions at each complete revolution of the phase angle (corresponding to 11 steps of iteration because of the choice of stadium size and initial generator). The rectangular orbit (g) scar also shows resonance with the bow tie (b). Its orbit-length-dividing period of 1.2 relative to (a) is close to the bow tie's 1.29 and it closely shares the vertical section of its orbit with (b).(Bottom) (h) A vertical wave generator corresponding to the neutral vertical orbit, produces an ordered wave function.

The images in fig 4 were generated using the first 3300 iterations of the generating function, however when the scars (b) and (f) in fig 4 were iterated on from 3300 to 6600 iterations the scar remained stable , as shown in fig 8, verifying it is not just a transient feature of the expanding generating functions.

Fig 5: The generating functions for the scars in fig 4
consisted of radial Gaussian wave packets of ÔplaneÕ cosine waves aligned along
the classical repelling orbits shown with wave periods a fixed fraction of the
orbit length.

Fig 6: (Left):The corresponding non-chaotic wave functions
for generators a, b in fig 4 in a rectangular potential well do not show
scarring, but form regular ordered solutions. (Right): Non-chaotic wave function corresponding to a
quasi-periodic orbit in a lemon-shaped potential well with the quasi-periodic
orbit superimposed below.

These patterns of scarring can also be
compared with other potentials wells such as a rectangular one which leads to
ordered wave functions and the lemon shaped potential well, which provides an
example of so-called Ôsoft chaosÕ in which ordered quasi-periodic solutions
coexist with chaotic ones. In fig 6 the generating functions of scars (a) and
(b) in fig 4 are applied again for 3300 iterations in the rectangular well
resulting in ordered wave functions, and a vertical generating function in the
lemon well produces an ordered wave function corresponding to the
quasi-periodic stable orbit superimposed below.

Fig 7: The cellular automaton also well models the
semi-classical view of the spreading of a traveling wave packet. A small
traveling wavelet, generated by clipping sector of a circular ripple with a
Gaussian envelope, bounces back and forth, ultimately forming a periodic wave
pattern, due to internal resonances generated by its characteristic frequency
and wave spreading. Wave spreading
compensates for the exponential spreading of the classical repelling orbit and
the wave function can form standing wave constructive interference when its
energy and consequently its frequency and period corresponds to a fraction of
the length of the periodic orbit.

Significantly the CA is able to well
handle oblique generating functions with periods not reducible directly to an
integer number of iterations. In
fig 1 are shown both the probability superposition scars and the amplitude
superposition scars, confirming that even when an oblique wave with an
irrational period is modeled, the CA is able to successfully superimpose
amplitudes in phase over thousands of iterations.

Fig 8: The long-term stability of the scars in the CA simulation is confirmed when the cases (b) and (f) in fig 4 are run for a further 3,300 iterations and the initial transient states of the generator are thus omitted entirely. Top: Iterations 1-3300 superimposed. Bottom: Iterations 3300-6600 superimposed. In these figures the absolute wave amplitudes have been added rather than the squared amplitudes to illustrate various methods of scar highlighting.

The CA can also be used to model the simulations of traveling wave packets illustrated from other research. In fig 7 is shown a series of images of iterations of a right-moving wave packet giving a similar evolution profile to the one generated by Heller and Tomsovic in fig 2(c).

Fig 8b: A more detailed long-term evolution of scarring of a wave packet with twice the frequency and half the Gaussian variance

for the same sequence of orbits as in fig 4, superimposing stages: 0-3300, 3300-6600, 6600-9900 9900-13200 for (f) and (g)).

The long-term evolution suggests that the systems are converging to a linear combination of the principal orbits due to the initial waveform

exciting more than one orbital solution. Hence (b) and (e), being two of the 'strongest' orbits, retain their excitation almost unchanged

through 9900 iterations, while for example (g) already shows features of (b) in its ititial excitations.

The CA can also well-model other wave experiments in physics such as the double-slit experiment illustrated in fig 9.

Fig 9: A modification of the rectangular potential well consists of a discrete model of a double slit interference experiment, showing the CA well models circular wave fronts emerging from a double slit and forming constructive and destructive interference fringes.

It is also possible to compare the
effects of different boundaries directly by choosing a circular wave generator
of period a fraction of the overall dimensions (both dimensions for stadium and
rectangle and the vertical dimension for the lemon). In fig 10 the three
boundaries are compared and the scarring characteristic of the stadium shows up
clearly and quite elegantly suggestive of several superimposed orbit solutions.

Fig 10: Probability distributions for a circular wave
packet of period a fraction of the overall dimensions for 3300 iterations for
the stadium, rectangle and lemon display the differences between chaotic and ordered
wave forms.

The CA can also be easily interrogated to
produce an autocorrelation function [6] by
multiplying the states of the generating wave function by the
propagated wave function and integrating
over all cells in the array. In fig 11 is shown the autocorrelation function
for the bow tie generator up to 3300 iterations consisting of about 7 sweeps
across the length of the stadium and the fast Fourier transform, giving back
the energy/frequency/period of the generated wave function.

Fig 11: (Left) Autocorrelation function of the
generated Ôbow tieÕ scar in fig 4(b). (Right) FFT giving frequency domain.

To complement the wave CA and draw
superimposed classical orbits, a ray drawing algorithm was developed for the
same matrix of elements used in the wave CA algorithm. This also enables
portrayal of several features of the associated classical chaotic system.

Fig 12: Classical orbits in the stadium display the classic
features of chaotic billiards. Top left: Sensitive dependence on initial
conditions illustrated for two closely related initial conditions (red and
green) show exponentiating divergence with each reflection. Chaos also causes
the space to be densely permeated with repelling periodic orbits, two of which
are illustrated below. Because
they are repelling, neighboring orbits are thrown further away, rather than
being attracted into a stable periodic orbit as in an ordered system. Consequently almost all orbits are
recurrent and densely fill phase space (top right).

In fig 12 are shown how sensitive
dependence on initial conditions (in this case a mall angle difference in the
initial orbit leads on multiple reflection in the curved sections of the stadium
to exponentially incr3asing angular divergences between the orbits, eventually
ending in essentially a random relationship with one another. Letting an orbit
run for a large number of reflections also shows that most orbits are not
periodic, but rather densely permeate the phase space of the stadium. Although
the repelling periodic points are rare, they are dense which means there are
periodic orbits arbitrarily close to any trajectory. Most orbits are recurrent
so that they eventually come arbitrarily close to a give trajectory, as shown
in fig 15

Fig 13: A collection of classical 16-step close
recurrences approximating periodic repelling orbits, some with reverse orbit
symmetrization to 32 steps. Top row symmetric by horizontal reflection, second
row by a rotation by 180^{o}, third by vertical reflection and the
fourth row asymmetric.

* *

To more closely model the classical
chaotic stadium, a second version of the ray-drawing algorithm was developed
which rather than operating on a matrix of pixels, solves the reflection
equations at the boundary on a classical continuous stadium. In turn this
algorithm was refined to search for and display the closed orbits through a
given point for reflection periods up to 20. Because these are testing
numerically for a closed orbit by close approximation, they are really
detecting close recurrence in a finite number of reflections, rather than an
exact periodic solution.

Thus is a dramatic increase is made in
the number of reflections, such an orbit will eventually break stability as the
periodic solution is repelling, and may over time recur closely to other
repelling orbits. In fig 15, the computed solution to the Ôbow tieÕ 4-period in
fig 4(b) is iterated for 200 steps in each direction, showing noticeable
recurrences to another period 4 skewed bow tie and a rotationally symmetric
period 8 orbit.

In fig 13 sample orbits displaying
reflective symmetry in the axes, rotational symmetry by half a revolution and
asymmetric orbits are shown illustrating the variety of closed orbits and their
symmetries.

The algorithm was set to compute all
close recurrences back to the starting trajectory for 2 up to 20
reflections. The algorithm
initially uses 80000 steps of angle through a sample right angle of directions
from a given point and tests for local minima of a metric combining the
distance from the initial vertex and divergence of angular direction. The local
minima are then scanned down a further four orders of magnitude to find as
close a recurrence to a closed orbit as possible given the numerical resolution
of the floating point process.

These solutions are then plotted as
closed orbits and their lengths and sensitivity to initial conditions are
determined as an initial step towards calculations based on the trace formula [1].

The distribution of the lengths was also investigated as a comparison with the prime numbers whose separations increase and Riemann zeta zero distributions whose separations decrease with increasing imaginary value and whose step differences have the statistics of a Gaussian unitary ensemble (GUE).

Fig 14: (Above) Distribution of lengths of 6444 close
orbit recurrences of up to 20 reflections generated from positions (0,0), (2,0)
and (0,1), having rotational and reflective symmetry shaded to give relative
densities. (Below) Corresponding distribution of orbit length differences of
the ranked length array. The right hand peak is an artifact caused by missing orbits
with higher reflection numbers.

The orbit length profile shows some
regions that are sparsely filled, and others which are more densely filled.
Certain small finite values of length, such as the perimeter distance of the
stadium and its multiples have an infinite number of orbits of length tending
to these values. The statistic of gap distances shows a descending power law
trend, which is repeated internally, for example when the orbits from 3000 to
5000 are chosen, which, from fig 14, is a relatively consistent central section
of the orbit lengths not suffering from anomalies of small value divergences or
missing large orbits due to cutting off the search at 20 reflections.

Nevertheless when the inverses of orbit lengths are taken to give a measure of energy as inversely proportional to wavelength, a GUE-like statistical distribution does arise, suggesting the GUE statistics seen in experiments on the quantum stadium do result from a collection of eigenvalues each corresponding to a closed orbit.

Fig 15: (left) Distribution of length differences by frequency shows a descending law in which larger gaps are rarer. This pattern is repeated in the inset when the outliers are removed by choosing only the values from 2000 to 5000. This distribution differs markedly from the distribution f primes and is unlikely to be an artifact of the orbit selection because it is clear there are sequences of orbits with orbit lengths converging in the limit, such as the oneÕs exemplified top left in fig 11, which tend asymptotically to the circumference of the stadium. According to the trace formula [1] the phase shifts and sensitivity of the orbits also have to be included before a GUE distribution is arrived at, however, as shown (right) a roughly GUE distribution, with prominent outliers corresponding to the shortest groups of orbits, does arise when the inverses of the orbit lengths are taken as a measure of energy.

The perimeter-angle coordinate system described in fig 16 makes it
possible to investigate the classical billiard as a discrete dynamical system.
In this system orbits over a finite number of reflections become homologous in
a symbolic dynamical representation in terms only of which curved and straight
sections the orbits reflect off (e.g. rtblltrrb É). This makes it possible to reduce summations over all orbits, in forming classical or semiclassical correlation functions and wave plots, to summations over a single member from each family.

Fig 16: (Top left) Extended orbit of 200 reflections in
both directions from a close homoclinic recurrence generating the repelling bow
tie orbit representation (bottom left black) later shows heteroclinic
recurrence to other orbits, a related period 4 skewed bow tie (blue) and a
period 8 orbit (red). Because the dynamics is conservative and area in phase
space is conserved, the chaotic process consists of horseshoe saddles and
repulsion in one eigenspace is compensated by attraction on the other. (Top
right) Single directed orbit of 500 reflections from (0,1) with direction
(1,-3) can be transformed by the coordinate transformation shown into a
Ômomentum-positionÕ representation in which *P* is the cosine of the angle of the ray to the tangent
and *Q* is the perimeter distance e.g. from the mid-point of the right circular
arc. This 2D coordinate system results in the discrete dynamic below for 50000
iterations, in which phase space is ergodically (pseudo-randomly) filled except
for a diminishing pair of central regions corresponding to the neighborhood of
vertical stable orbits.

The dynamic can thus be described in terms of the Markov partition of phase space determined by the orbits which strike or emanate from the four junctions between curved and straight sections at (±1, ±1), labelled in red, magenta, green and
blue dots in fig 16.

Fig 17: (Top left) The
one-step partition of phase space induced by colouring all the orbits fanning
off (±1, ±1), colour-coded as in the discussion and fig 16. (Top right) the
corresponding partition for four steps. (Bottom left) The region mapped from
the red ellipse in two reflection steps becomes fractally recurrent after six
steps (bottom right).

The semi-classical approach retrieves very comparable results to a full
quantum representation, fig 2(g), by tracing a key member of each family of
orbits, in the partition of the classical billiard, whose momentum divides the
length by the time
and applying an oscillating wave propagator to the member of the form ^{[viii]}:

[7]

where *S* is the action, with *j* weighted according to the product of the orthogonal stability matrices:

[8]

along the orbit, giving a correlation function for a single orbital family member as follows:

To explore a simplified semi-classical process, illustrating some of the
principles of the above approach, we developed a Matlab program superimposing
successive reflections of a 360^{o} fan of 2^{12} – 2^{18} rays from the centre of the stadium with an oscillating propagator having period a fraction of the linear dimensions of the stadium. Cells are summed by superimposing the amplitude of each propagator at the point its trajectory crosses a given cell, weighted by the stability factor as above. In theory this should have a similar effect to the selection of the orbit with the fractional momentum arriving at a given cell at a given time, since, as we are adding the cumulative effects over successive times, orbits with an integral number of wavelengths over a given path will interfere constructively and those not maintaining phase will effectively cancel one another out.

A more refined method would be to consider the heteroclinic orbits passing between two Gaussian distributions around the initial and final states, by selecting from the fractal intersections of the forward orbits from the initial state and the backward orbits from the final state those with a fractional momentum/period. One can then integrate this process over all the final states to give a time evolving distribution.

In fig 18 the semi-classical method for a propagator having wavelength 1/60^{th} of the width of a 600 pixel stadium is compared with the CA wave function method with a circular wave packet of wavelength a 1/54th fraction of the 594 cell width of the stadium. Although the profiles of the two simplified methods and their scarring patterns are not identical, both display similar features of scarring by wave coherence, illustrating how the two methods can highlight the quantum suppression of chaos. The semi-classical version also shows prominent effects of caustics caused by the curved ends as verified in fig 19.

Fig 18: (Top) Semiclassical
superimposed periodic ray trace of the stadium from a fan of rays at the centre
with propagator having period a fraction of the stadium dimensions for 4
reflections. (Bottom) The wave function CA, with a circular wave packet
generator, having a similar fractional period, relative to the stadium width,
to the semiclassical example above.

As shown in fig 19, because it is developed from a single fan of rays, the corresponding plots for successive reflection paths using superimposed absolute (or squared) amplitudes (a-e) do not manifest closed orbit scarring, but rather display periodic versions of the caustics found when classical rays are superimposed (f). We have thus made a complete transition from classical chaos to caustic order.

Fig 19: Semiclassical
superimposed periodic ray trace of the stadium from a circular wave generator
at the centre with period a fraction of the stadium dimensions for (a) 2, (b)
3, (c) 4, (d) 5, and (e) 8 reflections. Rather than highlighting the scars
around existing periodic orbits, the method highlights periodic versions of the
successive caustics caused by repeated reflection also evident in (f), the
classical ray trace for 4 reflections, shown for comparison.

Matlab Files

Source matlab
files for the project can be downloaded from:

http://www.dhushara.com/DarkHeart/QStad/QStad.zip

Further matlab
files for the zeta function and discrete chaos are available at:

http://www.dhushara.com/DarkHeart/RH2/RH.htm

http://www.dhushara.com/DarkHeart/DarkHeart.htm

References

[i] Bunimovich L (1974)
Funct. Anal. Appl. 8 254, (1979) Commun. Math. Phys. 65, 295.

[ii] M. C. Gutzwiller,
(1971) J. Math. Phys. 12, 343.

[iii] Gutzwiller Martin (1992) Quantum Chaos
Scientific American Jan. 266, 78 – 84.

http://www.scientificamerican.com/article.cfm?id=quantum-chaos-subatomic-worlds

For an
illustrated version with additions, see: http://www.dhushara.com/book/quantcos/qchao/quantc.htm

[v] Heller Eric, Tomsovic Steven (1993) *Postmodern
Quantum Mechanics* Physics
Today July 38-46.

[vi] S. Tomsovic, E. J.
Heller, (1993) Semiclassical Construction of Quantum Eigenstates Phys. Rev.
Lett. 70 1405-8.

[vii] Peterson Ivars (1991) *Back to the
Quantum Future* Science News
140 282-285.

[viii] Tomsovic S, Heller E
(1993) Long-term dynamics of Chaos: The Stadium Billiard Physical Review E 47
282-299.

[ix] Alisa Bokulich (2008)
*Can Classical Structures Explain Quantum Phenomena?* Brit. J. Phil. Sci. **59** 217–235

[x] Gutzwiller, M.C.
(1992). *Quantum Chaos*. Scientific American **266** 78 - 84.

[xi] Heller E, Tomsovic S
(1993) *Postmodern Quantum Mechanics* Physics Today **July** 38-46.

[xii] Liu B, Zhang G, Dai J Zhang H (1998)
Eigenvalues and Eigenfunctions of a Stadium-Shaped Quantum Dot Subjected to a
Perpendicular Magnetic Field Chin. Phys. Lett. 15/9 628.

[xiii] http://www.aip.org/png/html/corral.htm

[xiv]
Kwon Y & Hosoglu S (2008) *Application of lattice Boltzmann method, finite element method,
and cellular automata and their coupling to wave propagation problems* Computers and Structures 86 663-670.