Dynamic control of an embedded cavity resonator

1 July 2002

Optics Communications 208 (2002) 145–153 www.elsevier.com/locate/optcom

Dynamic control of an embedded cavity resonator
n a,*, J.V. Jos
e b A. Antillo a b

Centro de Ciencias F
ısicas, Universidad Nacional, Aut
onoma de M
exico, Apartado Postal 48-3, 62251 Cuernavaca, Morelos, Mexico Physics Department, Center for the Interdisciplinary Research on Complex Systems, Northeastern University, Boston, MA 02115, USA Received 29 November 2001; received in revised form 2 April 2002; accepted 17 April 2002

Abstract We propose two embedded elliptic billiards as a possible cavity resonator. We analyze the ray and wave optic properties of the resonators, mostly in 2D but also in 3D. We show that this system can produce special whispering gallery modes that can be controlled by changing the parameters of the billiards. Our results indicate that this configuration may work well as a short wavelength practical table-top cavity resonator. Ó 2002 Elsevier Science B.V. All rights reserved. PACS: 05.45; 95.10

There is significant interest in producing short wavelength cavity resonators for basic research and applications, for example, in photoelectron spectroscopy, microscopy, holography and lithography. There are, for example, X-ray sources which use synchrotron radiation [1], high harmonic generations [2], intense laser–atom interactions [3] and plasma lasers [4]. Table-top low power lasers are also promising radiation sources [5], where laser cavities play an important role [6,7]. No significant practical use of optical cavity resonators has been, however, implemented (see [8]). Here we suggest a new optical embedded cavity configuration scheme that uses special whispering

*

Corresponding author. Tel.: +52-777-3291727; fax: +52777-3291775. E-mail address: armando@fis.unam.mx (A. Antill
on).

gallery modes (WhgM) that can potentially be used to produce short wavelength radiation resonators. A free electron laser (FEL) resonator has been considered in [9]. There the active medium was a ring of relativistic electrons rotating near a cylindrical surface opposite to a beam of intense pumped light that travels along the cylinder’s surface in the WhgM regime. Braud and Hagelstein (BH) [10] have proposed a WhgM laser cavity resonators with a plasma as the active medium. There the coupling between the amplifier and the mirrors is problematic since each time the rays reflect from a mirror the beam broadens by the time it reaches the amplifier. The net result is a net power loss and beam degradation with each pass through the amplifier. BH tried to reduce this problem by using WhgM in cavity resonators with non-circular curvature. They concluded that an elongated geometry could greatly reduce the beam divergence problems. This configuration did not

0030-4018/02/$ - see front matter Ó 2002 Elsevier Science B.V. All rights reserved. PII: S 0 0 3 0 - 4 0 1 8 ( 0 2 ) 0 1 4 9 5 - 5

146

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

completely eliminate the broadening because mode mixing is still present in the geometries they considered. Here we propose reflecting embedded cavity geometries that allow to carefully control the grazing angles of reflection with the boundary. These billiard geometries can be easily changed so that the ray orbits show chaotic and resonant trajectories that allow to significantly improve and control the resonator properties. The embedded billiards do not have to be elongated, but they do need to have components with small radii of curvature. The design allows to place straight mirror sections in the active media. These sections can be chosen according to the experimental requirements, such as the distance from the mirrors to the active media, or the circumference of the resonator to improve gain. Adding inner mirrors to the cavity allows the selection of the unique wanted resonator mode, which is stable even including the errors produced by the approximate location of the mirrors. Our embedded billiard proposal can also be implemented in 3D, as we discuss below. We start our analysis considering the two embedded elliptic billiards shown in Fig. 1, in the short wavelength ray dynamic limit. In this case the light rays bounce off the boundaries elastically, they are confined to the plane and bounded by the region between the two elliptical boundaries. It is the complexity of the resulting dynamics that plays an important role in this paper. Billiards have been widely studied in classical and quantum mechanics [11,12], and there are well-developed techniques for their analysis. The two concentric billiards considered here (see Fig. 1) are defined by their semi-axes: for the internal ellipse by (a1 , b1 ), and by (a2 , b2 ) for the external one. The potential is infinite on both elliptical boundaries, leading to full specular ray reflections. The WhgM are produced by the reflections of light circulating between the two elliptic boundaries. The ray dynamics can be analyzed in terms of the Poincar
e–Birkhoff (PB) phase space canonical variables h and sinðvÞ [11]. Here h is the polar angle, and v is the ray angle of incidence with respect to the normal at the boundary, as shown in Fig. 1. The ray dynamics can have periodic, quasi-periodic and chaotic trajectories.

Fig. 1. Here we show the two embedded elliptic billiard resonator considered in this paper. h is the polar angle and v the reflection angle with respect to the normal to the surface. The Poincar
e–Birkhoff (PB) surface of section is defined by the canonical conjugate variables (h; sin v). The light rays are restricted to move on the plane defined by the shaded region between the two concentrically rigid ellipses. Light propagates indefinitely satisfying the law of total specular reflection at collisions with the boundary of the resonator.

Here we present some results that provide the basis for the resonator model proposed in this paper. To graphically illustrate our results, we begin by choosing typical parameters for the billiards; ða1 ; b1 Þ ¼ ð0:89; 0:85Þ and ða2 ; b2 Þ ¼ ð1:0; 1:0Þ. In Fig. 2(a) we show the PB phase space for 60 ray trajectories with initial condition x0 ¼ 0, and for several y0 ½0:85; 0:91. The initial h orientation is along the x-axis. We have sin v ¼constant and h varying when the initial y0 values are within the interval [0.89,0.91], the resulting trajectories never hit the internal elliptic billiard. A typical trajectory for this case is shown in Fig. 2(b), for 200 reflections, and for y0 ¼ 0:905. There are some trajectories that just start to hit the internal billiard when y0 is below 0.89. The orbits are trapped into elliptic resonances in the PB phase space, usually of high order. One of these orbits is shown in Fig. 2(d) for y0 ¼ 0:88733. Smaller values of y0 give chaotic trajectories embedded in a chaotic sea in the PB phase space. Fig. 2(c) shows a chaotic trajectory, with no regularity in both sin v or h (here y0 ¼ 0:875). Inside the chaotic region, how-

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

(a)

0.895 sin χ

1.1

y

147

(b)

0.0

θ

0.835 0.0 1.1 y

-1.1 3.1 -1.1 1.1 y (c)

0.0

-1.1 -1.1

x 0.0

1.1 (d)

0.0

x 0.0

x -1.1 1.1 -1.1

0.0

1.1

Fig. 2. In (a) we show the PB phase space for 60 ray trajectories with initial condition x0 ¼ 0 and for several y0 . Rays with y0 [0.89,0.91] have sin v ¼ constant. A typical trajectory for this case is shown in (b), for 200 reflections, and for y0 ¼ 0:905. Some trajectories with initial conditions y0 below 0:89 begin to hit the internal billiard and then are trapped into elliptic resonances, as the one shown in (d) for y0 ¼ 0:88733. Other values give chaotic trajectories and a chaotic sea in phase space, as shown in (c) for y0 ¼ 0:875.

ever, there are resonances of smaller order, as the eight order resonance that appears in the chaotic region for y0 [0.86,0.76]. We show this resonance in configuration space in Fig. 3(a) for y0 ¼ 0:8719361, which is close to a fixed point. The grazing angle in this case is close to 30°. We now compare our scheme to specific resonator values and design goals of [10] (if we disregard tight grazing angles). Taking the parameter values used by BH [10], we can scale Fig. 3(a) by a factor of 16 to get Fig. 3(b). The trajectory in the later figure has y0 ¼ 13:96853368, which is quite close to an elliptic fixed point. This result shows that beam stability can be achieved not only with a continuous boundary billiard, but also with some concave and a few convex properly separated mirror segments (at least for the case of large grazing angles.) The inset in Fig. 3(b) corresponds to one of the eight resonant dynamic islands that appear in the PB phase space of Fig. 2(a). In Fig. 3(c) we amplify the top horizontal segment of Fig. 3(b) for the elliptic fixed point at y0 ¼ 13:95097768. Two more rays are 25 lm above and below that value, y0 ¼ 13:95100268 and y0 ¼ 13:95095268 to

cover the plasma height of [10]. Fig. 3(c) shows a light beam divergence in configuration space that arises when the height dimension of 50 lm is covered by the plasma active media. The segments in the figure give approximately the dimension of the plasma of size 2 cm  50 lm. In Fig. 3(d) we again show one of the eight islands of the PB phase space for the rays that appear in the previous figure. Note the sin v scale in the vertical coordinate. This figure will look just as a point in the inset of Fig. 3(b), which shows that the area of stability of this kind of resonances is quite broad. The stability shown here is due to the combination of a focusing effect by the concave mirrors and a defocussing effect by the convex mirrors. This result resembles the particle accelerator approach, where a magnetic focusing–defocussing system is used to control the charged particle trajectories [13]. A material’s photoabsorption cross-section gets larger at short wavelengths and, consequently, its reflectivity diminishes. Good reflectivities can be found, however, in some materials for specific wavelengths and small grazing angles. This means that the boundary between the chaotic and regular

148

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

1.1

(a)

y

16

y 0.873

0

0.0

sin χ

(b)

0.872 θ 0.870 0.90 1.05 1.20

-16 -16 1.1 0.871938 sin χ

x

-1.1 -1.1 13.9525 y

0.0 (c)

0

x 16 (d)

13.9510

13.9495

x -8

-4

0

4

θ 0.871934 8 1.0588

1.0595

Fig. 3. In (a) we show how resonances of smaller order are excited inside the chaotic region. Here x0 ¼ 0, y0 ¼ 0:8719361, and the grazing angle is close to 30°. The parameters are scaled in (a) by a factor of 16 to get (b). The trajectory in (b) has y0 ¼ 13:96853368, and it shows a beam that is stable, for large grazing angles, with only some concave and few convex mirror segments. The inset in (b) is one of the eight elliptic islands appearing in the phase space of Fig. 2(a). In (c) we show the light beam divergence in configuration space that arises when the height dimension of 50 lm of the plasma active media is covered [10]. The segments in the figure give approximately the dimensions of the plasma of size 2 cm  50 lm. The elliptic fixed point is at y0 ¼ 13:95097768. Two extra rays are at 25 lm with y0 ¼ 13:95100268 and y0 ¼ 13:95095268. In (d) we show again one of the eight islands that appear in (c). This figure will look just as a point in the inset of (b).

regions in the PB phase space of Fig. 2(a) needs to be moved up. In this case the number of reflections increases and the total return transmittance to the active medium will be lower. For grazing angles of 1°, for example, the number of reflections is of the order of 100 (the number of reflections for a concave surface is p/grazing angle). If for a given material and wavelength the reflectivity is 99%, then a 360° ray deviation gives a total reflectivity of 37%. For wavelengths as low as hard X-ray, typical values of grazing angles are of the order of 0.26°, which implies that sin v  0:99999 (see resonances of the type shown in Figs. 2(a) and (d). The resonances shown in Fig. 3(a) can also be used to design a cavity with less constrained requirements and with longer straight sections. They have smaller values of sin v but with fewer reflections and where the active medium can amplify the residual radiation with a gain larger than the loss. The resonator geometries proposed here can be optimized to specific experimental situations. Consider, for example, the resonator suggested by

Vinogradov et al. [14]. This resonator can be obtained by cutting and rigidly separating the two semi-circles and the two semi-ellipses from the resonator of Fig. 3(a) (as in the embedded stadium billiard). Under such conditions the elliptic fixed points remain in the same place in the PB phase space, but the beam divergence will be smaller (see Fig. 3(c)). The rearrangement of the internal semiellipses can thus be modified as needed to optimize the resonator properties. A fast buildup of gain saturation can be achieved provided that the single-pass gain of the ondulator (e.g., for a FEL) should be greater than the inverse of the total single-pass transmittance of the resonator. Assume that the transmittance is given by Rp=a , where R is the mirrors reflectivity, ‘‘a’’ is the grazing angle and p=a is the number of reflections in the resonator. If damage to the mirrors is not a problem, as may happen with low power lasers, smaller reflectivities can allow larger incident angles, as those shown in Fig. 3. In the high power laser case, as in the Linac Coherent Light Source (LCLS)

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

[15], things change and we may have to go to very small grazing angles to avoid mirror damage. In Fig. 4 we show other possible resonances that can help improve the light beam stability. In Fig. 4(a) the PB phase space is shown in the vicinity of 0.23° angles of incidence. We chose two of these resonances and plotted the configuration space for the places where the ray trajectory hits the internal mirrors. In Fig. 4(b) we show the two mirror locations needed to achieve beam stability. In Fig. 4(c) we show the ray trajectory that hits four mirrors. In Fig. 4(d) we show, for comparison, the whole internal and external elliptical boundaries (not shown) that are needed to confine a chaotic trajectory. Next we consider defects in the resonator mirrors that can affect the life time of the resonances. The duration of the gain in the resonator will also be important to determine if the dwell time is long enough to observe the resonances. To mimic the presence of defects or rough surfaces, we chose a Gaussian distribution of random errors, with standard deviation r about the normal angle. For a trajectory amplitude shown in Fig. 3(c), we took r ¼ 5:8  106 , which is 1°/3000. The results are shown in Figs. 5(a) and (b), for 1:5  105 mirror

0.99999205

collisions. The phase space grows inside the boundary of the resonance. As shown in Fig. 5(b), however, the motion oscillates staying bounded. This does not happen, for example, when r gives 1°/1000, where the motion becomes eventually chaotic. In some cases the trajectory in phase space returns to the resonances and then goes back again to the chaotic region. Other resonances can also be excited after a chaotic transient. To achieve stability there is a compromise between the rough mirrors, beam width and dwell time for the rays to stay inside the resonator. In the example shown in Fig. 5(b), the ray turns 105 times about the resonator. Fig. 5(c) shows the number N of collisions with the mirrors as a function of r. N is related to the lifetime of the resonance, and it is defined as the time after which the ray dynamics become chaotic. The number of turns about the resonator is N =14, with each point in the figure obtained from an average over ensembles of 10 random angles. Note that the decay appears to be exponential with N. Fig. 5(d) shows that a trajectory inside the resonance, with r ¼ 1:74  105 , has N independent of the amplitude y0 ðx0 ¼ 0Þ, except for small statistical errors. Amplitudes close to the last

0.015

(a)

sin χ

0

θ

0.015 0.1 0.9998 0.015 (c) y

0

0.015 y

0

0.015 0.9998

(b)

y

0.99999185

0.99999165

149

mirror

x 1 (d)

0

mirror

x

0.015 1 0.9998

x 1

Fig. 4. Here we show possible resonances that can help improve the light beam stability for a high power laser. In (a) the PB phase space is shown in the vicinity of angles of incidence of 0.23°. In (b) only two (one in each side) internal mirrors are needed, and in (c) we show the resonance for a ray that hits four internal mirrors. In (d) a chaotic trajectory is shown for comparison.

150

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

Fig. 5. (a) shows the PB phase space for 1:5  105 hits with mirrors for the trajectory used in Fig. 3(c). The normal to the boundary has been perturbed by random errors given by a Gaussian distribution with r ¼ 5:8  106 , which is 1°/3000. (b) shows the oscillatory but bounded behavior of h as a function of time for the same trajectory as in (a). The time units used here allows the ray to turn around the resonator for 105 times. (c) shows the variation of the number N of collisions with the mirrors vs the r of the distribution for the same trajectory. N gives the time after which the ray dynamics goes chaotic. N =14 is the number of turns around the resonator. Each point in the figure is obtained as an average over 10 ensembles of random angles. The N dependence appears to be exponential. In (d) we plot N as function of the amplitude y0 (x0 ¼ 0) for the trajectory inside the resonance, for r ¼ 1:74  105 . Note that, except for statistical errors, N is independent of the amplitude if it is not too close to the last invariant torus. Each point is an average over 15 ensembles of random angles. A, B and C are, respectively, the fixed point, the amplitudes of trajectories used in Figs. 3(c) and (b) for the resonator scaled as in Fig. 3(a). (e) shows the statistical results for systematic errors representing the mirrors misplacement around the curved geometry for the resonance of Fig. 3, at the collision points of the fixed point orbit shown in Fig. 3(a). We considered a random rotation of the curved mirrors by an angle given by a Gaussian distribution. The calculation was done again for the amplitudes corresponding to Fig. 3(c). (f) displays the behavior of the number of collisions as a function of the amplitude of the trajectory for r ¼ 9  104 . Each point in (e)–(f) is an average over 10 ensembles of randomly rotated mirrors. Both random and systematic perturbations have a similar qualitative effects, but the latter is dynamically more stable.

invariant torus fall quickly, however, into the chaotic region. Each point in this figure represents an average over 15 ensembles of random angles.

Here A, B and C are respectively, the fixed point, the amplitudes of trajectories used in Figs. 3(c) and (b) for the resonator scaled as in Fig. 3(a). We

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

expect that resonances for which sin v 1 have tight tolerances since their widths are smaller. For very small grazing angles, (sin v 0:99999), the r’s that allow for an N  1  104 , are about 108 . The misplacement of the mirrors around the curved geometry are simulated by systematic errors. Note that the dynamics is more robust against systematic than the random errors considered before. This is shown in Fig. 5(e), where we considered the resonance of Fig. 3. We made random rotations of the curved mirrors by an angle chosen from a Gaussian distribution around the collision points corresponding to the fixed point orbit shown in Fig. 3(a). At each reflection point the ray sees a fixed but tilted mirror. The calculation was done for the amplitudes of Fig. 3(c). In Fig. 5(f) we plot the number of collisions as a function of the amplitude of the trajectory for r ¼ 9  104 , which is larger than the one used in Fig. 5(d). Each point in Figs. 5(e) and (f) is obtained from an ensemble averages over 10 randomly rotated mirrors. Both random and systematic perturbations have similar qualitative effects, but the latter is more robust dynamically.

151

To show the existence of corresponding wave modes, one could use physical optics as in BH, or the WKB approximation. One could also use results from Quantum Chaos theory to show a direct correlation between classical configuration space orbits and quantum eigenfunctions. For the 2D billiards we can use the similarity between the Schr€ odinger and Maxwell equations. In Fig. 6 we show wave results for a quarter of the resonator. Here we used the finite element to integrate the Helmholtz equation. The results exhibit the correspondence between the optical mode and the classical resonance in Fig. 3. The calculation involved 67,584 triangles for the eigenvalue 70,421. The corresponding wavelength is large, but scaling to higher eigenvalues is presently computationally very CPU demanding. However, the general trend will improve the sharpness of the mode as we go to smaller wavelengths, approaching the semi-classical limit. In the same figure we superpose the fixed point classical orbit corresponding to the wave like solution. In this paper we are mainly interested in solutions very close to the plane, however, we will also

Fig. 6. There is a direct correlation between classical orbits and quantum eigenfunctions in planar billiards [12]. Here we show a classical periodic orbit superposed to its corresponding wave solution, obtained by using the finite element method (see text for a more detailed discussion). We display only a quarter of the resonator. The calculation involved 67,584 triangles for an eigenvalue of 70,421.

152

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

briefly discuss the 3D case. In 3D the diffraction losses transverse to the plane will be important, unless appropriate transverse curvature mirrors are included to minimize this effect. We have simulated a 3D resonator with two component walls. One is the external face of a circular boundary and the other an elliptic boundary for the internal one. The vertical mirrors in the external boundary have a concave curvature, both in the plane and transverse to the plane direction. The inner boundary will have convex curvatures. The boundaries in phase space will be controlled by using focussing and defocussing elements mirrors. This is shown in Fig. 7(a), where drawn are only the collision points with the mirrors. In Fig. 7(b) we plot a space where h is the angle made by the projection of the position vector of the collision point in the plane, and the x-axis. v is the ray angle of incidence with respect to the normal to the surface in 3D. The radii of curvature in the vertical direction were taken as 0.35. The stability is sen-

sitive to the parameter values. When including the convex curvature in the inner mirrors there is, however, a wide range of radii parameters to achieve stability. In conclusion, in this paper we have mainly considered the complex ray dynamics trapped inside 2D embedded elliptic billiards. The analysis was briefly extended to the 3D case, where stability can still be produced and controlled. Error tolerances on mirror surface design were also considered. Although the errors limit the resonances lifetimes, we expect that statistically one can still have resonant modes for appropriate mirror perturbation strengths. We found that both the geometry and the dynamics can be controlled, and there can be regular and chaotic behavior that we have analyzed in detail. This model billiard considered here is an extension of the model proposed by Braud and Hagelstein for whispering gallery modes, where they used only one boundary resonator with variable curvature [10]. Here, we have

Fig. 7. In this figure we show how to achieve stability in 3D embedded billiards. We used a resonator with two wall components. A circular torus defines the external boundary (with concave vertical curvature of mirrors, both in the plane and transverse to the plane direction). The external face of an elliptical torus for the internal one (convex curvature). In (a) we only show the collision points with the mirrors drawn for radii of curvature in the vertical direction of 0.35. The stability seems to be slightly sensitive to this value. A wide area of stability is expected in the space of radii parameters when introducing the convex curvatures in the inner mirrors. (b) we plot an extension of phase space. Here h is the angle between the x-axis and the projection of the position vector of the point of collision into the plane. v is the angle of the ray incidence with respect to the normal to the surface in 3D.

A. Antill
on, J.V. Jos
e / Optics Communications 208 (2002) 145–153

added an internal boundary, that plays a crucial role in our analysis. Using this model system, we have shown how to achieve ray stability in a grazing angle resonator, by controlling the beam size and its divergence. Using the finite element method to solve the Helhmoltz equation, we obtained the appropriate optical mode that corresponds to the periodic classical solution. Further research is needed to achieve a successful short wavelength resonator design. Specifically, to find the laser power, tunability, and available materials for high reflectivity mirrors. For a more detailed analysis of the geometric ray dynamics of this and related billiards see [16].

Acknowledgements We wish to thank J.C. Gallardo from Brookhaven National Laboratory for his comments and for giving to us some relevant references to this work, and to M. Cornacchia from Stanford Linear Accelerator Center for useful suggestions. The work by J.V.J. was supported in part by the NSF.

References [1] R. Coisson, Appl. Opt. 34 (1995) 904. [2] P. Salieres et al., Adv. At. Mol. Opt. Phys. 41 (1998) 83.

153

[3] J.C. Kieffer et al., Phys. Fluids B 5 (1993) 2676; A. Rousse et al., Phys. Rev. E 50 (1994) 2200. [4] D.C. Eder et al., in: S. Svanberg, C.G. Wahlstrom (Eds.), X-ray Lasers 1996, Institute of Physics Conference Series, vol. 151, Institute of Physics, Bristol, UK, 1996, p. 136. [5] J.J. Rocca, Reviews of Scientific Instruments 70 (1999) 3799. [6] N.M. Ceglio, D.G. Stearns, D.P. Gaines, A.M. Hawryluk, J.E. Trebes, Opt. Lett. 13 (1988) 108. [7] N. Yamaguchi, T. Hara, T. Ohchi, Ch. Fujikawa, T. Sata, Jpn. J. Appl. Phys. 38 (1999) 5114. [8] R. Tatchyn, ‘‘Optical cavity vs. electron beam requirements LCLS in a regenerative for the operation of a 1.5 A amplifier mode’’, SLAC-PUB-95-7036 (1995). Available from: http://www.slac.stanford.edu/pubs/slacpubs/8000/ slac-pub-7036.html. [9] I.N. Bukreeva, I.V. Kozhevnikov, A.V. Vinogradov, J. Xray Sci. Technol. 5 (1995) 396. [10] J.P. Braud, P.L. Hagelstein, IEEE J. Quantum Electron. 27 (1991) 1069; IEEE J. Quantum Electron. 28 (1992) 254. [11] M.V. Berry, Eur. J. Phys. 2 (1981) 91; M.C. Gutzwiller, Chaos in Classical and Quantum Mechanics, Springer, Berlin, 1990. [12] H.J. St€ ockman, Quantum Chaos an Introduction, Cambridge University Press, Cambridge, UK, 1999, and references therein. [13] E.D. Courant, H.S. Snyder, Ann. Phys. 3 (1958) 1. [14] A.V. Vinogradov et al., Zh. Tekh. Fiz. 55 (1985) 567, Sov. Phys. Tech. Phys. 30 (1985) 335. [15] M. Cornacchia, The LCLS X-ray FEL at SLAC SLACPUB-8053 (1999). Available from . [16] A. Antill
on, J.V. Jos
e (work in progress).