text:
# References

I. Introduction

With the advent of modern computer technology, robust numerical calculations of time-dependent phenomena in quantum mechanics are now common, as is the visualization of the resulting effects . The time-evolution of wave packet solutions for many scattering geometries as well as for bound state systems  are discussed with increasing frequency in the pedagogical literature, illustrating not only such familiar aspects as wave packet spreading, but also extending student experience to more novel phenomena such as wave packet revivals. Reviews of wave packet revivals (initially highly localized quantum wave packet solutions which exhibit quasi-classical behavior, then disperse or spread in time to a so-called collapsed phase, only to reform later to something very much like it’s initial state) in many familiar model one-dimensional quantum mechanical systems such as the harmonic oscillator and infinite well have appeared  providing students with accessible examples of an important quantum effect which can be probed experimentally in atomic systems.

Besides studies of such effects in the two ’classic’ quantum mechanical model systems mentioned above, wave packet propagation, and especially the structure of revivals, has also recently been discussed in the quantum version of another familiar classical system, the so-called quantum bouncer , . This system is the bound state version of the classical ’falling object’, and is defined by the potential energy

$$V(z)=\{\begin{array}{cc}\mathrm{}\hfill & \text{for }z<0\hfill \\ Fz\hfill & \text{for }0<z\hfill \end{array}$$
(1)
corresponding to a constant downwards force acting on a particle above an impenetrable flat surface at $`z=0`$. Generalizing on many early papers  which discuss the time-independent solutions for this problem, the author of Ref.  focuses attention on deriving expressions for the collapse and revival times for initially Gaussian wave packet solutions, using a mixture of numerical and analytic techniques. In this note, we will revisit this problem, along much the same lines as in Ref. , by examining, in detail, the short-term (quasi-classical) and long-term (revival) time-development of such solutions in terms of their position and momentum expectation values and uncertainties. This type of expectation value analysis, coupled with existing analytic, numerical, and visualization studies can then help form a more complete picture of the highly non-trivial time-development possible in one of the ’classic’ one-dimensional model systems in introductory quantum mechanics.

As often happens, the more familiar and tractable cases of the harmonic oscillator and the infinite well examples exhibit special features, in this case, in terms of their time-dependence, compared to a more general form such as in Eqn. 1. For that reason, in Sec. II, we first very briefly review Gaussian wave packet solutions in these two more familiar systems. The bouncer system shares with the infinite well the feature of the ’bounce’ at the impenetrable wall, but some observables, such as the time-dependent wave packet spread, exhibit a much more typical cyclic structure than for the infinite well, where only the free-particle spreading time  is of relevance. Since the numerical studies presented here can also make contact with analytic solutions to the problem of a particle undergoing uniform acceleration, we also mention those solutions in Sec. III and then proceed to detailed results for the quantum bouncer.

II. Gaussian wave packet solutions for the harmonic oscillator and infinite well

Because of the special nature of the energy level structure of the harmonic oscillator, any time-dependent solution of this problem, $`\psi (x,t)`$, will be exactly periodic and return precisely to its initial state after the classical period, $`T_{cl}=2\pi /\omega `$. Using standard propagator techniques or other methods, it is, in fact, easy to construct closed form Gaussian wave packet solutions  such as

$`|\psi (x,t)|^2`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{\pi }L(t)}}e^{[xx_0\mathrm{cos}(\omega t)]^2/L^2(t)}`$ (2)
$`|\varphi (p,t)|^2`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{\pi }p_L(t)}}e^{[p+m\omega x_0\mathrm{sin}(\omega t)]^2/p_L^2(t)}`$ (3)
which have expectation values which satisfy the classical equations of motion, namely

$$x_t=x_0\mathrm{cos}(\omega t)\text{and}p_t=m\omega x_0\mathrm{sin}(\omega t).$$
(4)
The fact that the expectation values for any wave packet solutions for the harmonic oscillator satisfy the classical equations of motion does not depend on these specific Gaussian forms, but can be shown to be true analytically in a quite general way . The solutions in Eqns. (2) and (3) have time-dependent position and momentum uncertainties given by $`\mathrm{\Delta }x_t=L(t)/\sqrt{2}`$ and $`\mathrm{\Delta }p_t=p_L(t)/\sqrt{2}`$ where

$`[L(t)]^2`$ $`=`$ $`L^2\mathrm{cos}^2(\omega t)+\left({\displaystyle \frac{\mathrm{}}{m\omega L}}\right)^2\mathrm{sin}^2(\omega t)`$ (5)
$`[p_L(t)]^2`$ $`=`$ $`\left({\displaystyle \frac{\mathrm{}}{L}}\right)^2\mathrm{cos}^2(\omega t)+(Lm\omega )^2\mathrm{sin}^2(\omega t)`$ (6)
and the parameter $`L`$ sets the scale for both the position- and momentum-space spreads. Unless $`L=\sqrt{\mathrm{}/m\omega }`$, these position- and momentum-uncertainties oscillate with a period twice that of the classical motion , regaining the initial values of $`\mathrm{\Delta }x_0`$ and $`\mathrm{\Delta }p_0`$ at two opposing points in phase space. Once again, these results are not specific to the Gaussian solution, but have been derived quite generally in Ref. : they are also sometimes rediscovered ,  in these specific cases.

Another tool which is standardly used to probe the wave packet’s approximate return to the initial state is the auto-correlation function , defined by

$$A(t)=_{\mathrm{}}^+\mathrm{}\psi ^{}(z,t)\psi (z,0)𝑑z=_{\mathrm{}}^+\mathrm{}\varphi ^{}(p,t)\varphi (p,0)𝑑p=\underset{n=1}{\overset{\mathrm{}}{}}|c_n|^2e^{iE_nt/\mathrm{}}$$
(7)
which measures the overlap of the position- or momentum-space wavefunction at later times with the initial state. For the oscillator, one can argue that because the wave packets are exactly periodic and never collapse, there are no revivals.

On the other hand, for the infinite well, because the energy levels are integral multiples of a common value (but not equally spaced), there are exact revivals , but initial Gaussian-type wave packets do undergo dispersion into a collapsed phase. For example, wave packets constructed from energy eigenstates via

$$\psi (x,t)=\underset{n=0}{\overset{\mathrm{}}{}}c_nu_n(x)e^{iE_nt/\mathrm{}}$$
(8)
where $`u_n(x)=\mathrm{sin}(n\pi x/L)/\sqrt{L}`$ for the infinite well, with quasi-Gaussian expansion coefficients of the form

$$c_n=\sqrt{\frac{\alpha \mathrm{}\sqrt{\pi }}{L}}e^{\alpha ^2(p_np_0)^2}e^{ip_nx_0/\mathrm{}}$$
(9)
(where $`p_n=n\mathrm{}\pi /L`$) give an initial momentum distribution very close to the standard one for Gaussian free-particle packets, namely

$$\varphi _0(p)=\frac{\alpha }{\sqrt{\pi }}e^{\alpha ^2(pp_0)^2/2}e^{ipx_0/\mathrm{}}$$
(10)
and yield localized Gaussian-like position wave packets which are initially very close to the free-particle form

$$|\psi (x,t)|^2=\frac{1}{\beta _t\sqrt{\pi }}e^{(x[x_0+p_0t/m])^2/\beta _t^2}$$
(11)
where $`\beta _t=\mathrm{}\alpha \sqrt{1+t^2/t_0^2}`$ with $`t_0m\mathrm{}\alpha ^2`$. Such packets undergo quasi-classical motion, bouncing back and forth between the two infinite walls, with an increasing width given roughly by an envelope defined by the time-dependent free-particle width, $`\mathrm{\Delta }x_t=\beta _t/\sqrt{2}`$ . During the impulsive collisions with the walls, the position-space width temporarily decreases during the time that the momentum-space distribution is ‘flipping’ from positive to negative values.

As mentioned above, since the quantum bouncer has features in common with both systems, namely one infinite wall boundary at which ’bounces’ will occur and a second, smoother potential barrier, we will keep both these cases in mind as we examine the short- and long-term behavior of expectation values for the quantum bouncer.

III. Expectation values for Gaussian wave packets for the quantum bouncer

Before turning to wave packet solutions of the quantum bouncer problem, we first recall some results from the related problem of a quantum particle undergoing uniform acceleration, namely subject to a constant force, $`F`$, with a potential given by $`V(x)=Fx`$ everywhere in space. Transforming the resulting Schrödinger equation into momentum-space, one can construct arbitrary solutions  of the form

$$\varphi (p,t)=\varphi _0(pFt)e^{i[(pFt)^3p^3]/6mF\mathrm{}}$$
(12)
where $`\varphi _0(p)`$ is any initial momentum distribution. This result already implies that the momentum distribution for this case simply translates uniformly in time, with no change in shape, since

$$|\varphi (p,t)|^2=|\varphi _0(pFt)|^2$$
(13)
so that $`\mathrm{\Delta }p_t=\mathrm{\Delta }p_0`$. (We expect to see this behavior initially in the quantum bouncer case, until the wave packet nears the wall, as well as after the collision, at least for the first few classical periods.) Using this momentum-space solution, we can also construct closed-form position-space solutions; for example using the distribution

$$\varphi _0(p)=\sqrt{\frac{\alpha }{\sqrt{\pi }}}e^{\alpha ^2(pp_0)^2/2}e^{ipx_o/\mathrm{}}$$
(14)
we find a general Gaussian wave packet for the accelerating particle

$$|\psi (x,t)|^2=\frac{1}{\beta _t\sqrt{\pi }}e^{(x[x_0+p_0t/m+Ft^2/2m])^2/\beta _t^2}$$
(15)
with $`\beta _t=\alpha \mathrm{}\sqrt{1+t^2/t_0^2}`$. (The same results can be obtained in a variety of ways, using operator  or propagator  methods or other techniques .) This solution gives

$$x_t=x_0+\frac{p_0t}{m}+\frac{F}{2m}t^2\text{and}\mathrm{\Delta }x_t=\frac{1}{\sqrt{2}}\beta _t=\frac{\alpha \mathrm{}}{\sqrt{2}}\sqrt{1+t^2/t_0^2}.$$
(16)
Thus, the expectation value tracks the classical trajectory, while the wave packet spreads in exactly the same way as the standard free-particle Gaussian and we expect to see remnants of this behavior as well.

The author of Ref.  has constructed quasi-Gaussian wave packet solutions for the quantum bouncer using a combination of numerical and semi-analytic techniques, focusing on the derivation of quantities such as approximate revival time, $`T_{rev}`$. (See Ref.  for a very general derivation of both the classical period and quantum revival time for arbitrary bound state energy spectra). For an initial wave packet which is constructed so as to represent a particle released from rest at a height $`z_0`$, the classical period of motion, $`T_{cl}`$, and the revival time, $`T_{rev}`$, are given by

$$T_{cl}=2\sqrt{\frac{2mz_0}{F}}\text{and}T_{rev}=\frac{4}{\pi }\left(\frac{2mz_0^2}{\mathrm{}}\right)$$
(17)
respectively. For simplicity, following Ref. , we take units such that

$$\mathrm{}=1,2m=1,\text{and}F=1$$
(18)
and we choose particular initial values, namely

$$z_0=25\text{and}\mathrm{\Delta }z_0=1.$$
(19)
With these values, the classical period and revival time are given by

$$T_{cl}=10\text{and}T_{rev}=\frac{4}{\pi }(5)^2=795.78.$$
(20)
The collapse time, the time it takes for the initially localized wave packet to de-phase and become spread over the entire well, is given, in these dimensionless units as well as for our specific wavepacket parameters, as

$$T_{coll}=\frac{T_{cl}^3}{8\mathrm{\Delta }z_0}\text{or}T_{coll}=125$$
(21)
respectively.

Using these parameters, we can evaluate both the position-space and momentum-space wavefunctions in a manner which is similar to that in Ref. . We generate the expansion coefficients, $`c_n`$, in Eqn. (8) by explicit numerical calculation of the overlap integrals of the initial Gaussian form with the normalized energy eigenstates obtained by numerical solution of the Schrödinger equation; the eigenstates, normalization constants, and expansion coefficients are then compared to the Airy function solutions used in Ref. , supplemented by the more exact analytic results in Ref. , as a cross-check. (We truncate the expansion in eigenstates at a level consistent with the double-precision accuracy of our numerical evaluation, typically at values of the $`c_n10^6`$; this is to be be compared to the maximum values of $`c_n0.5`$. The resulting wavefunctions are then found to be appropriately normalized to unity to within $`10^6`$ or better over all time intervals considered here.)

We use this more numerical procedure partly because this technique runs more efficiently in our computing environment, but also because we intend to eventually extend our investigations of wave packet behavior to other, more general power law potentials of the form $`V_{(k)}(x)=V_0|x/a|^k`$ for comparison to recent work on the revival times in such potentials . This class of potentials contains the harmonic oscillator ($`k=2`$) and the infinite well ($`k=\mathrm{}`$) as special cases and can also be analyzed with an infinite wall added at $`x=0`$ to give the quantum bouncer case (where $`k=1`$): this fact also helps motivate our brief review of those two special cases.

We then plot the corresponding probability distributions over the first classical period in Fig. 1. On the left, the position-space probability is seen to spread in a manner which is numerically consistent with Eqn. (15), while the calculated position value $`z_t`$ (solid curve) agrees well with the classical expectation for the trajectory (dashed) except, of course, for the cusp at the ’bounce’. The packet exhibits the standard ‘interference’ pattern during the collision with the wall , , at the ’bounce’, and then reforms into something like the initial packet (compare to the dotted initial packet superimposed on the $`t=10`$ case), only wider.

For the momentum-space distributions (shown on the right of Fig. 1), we also see features of both the classical motion and the uniformly accelerated wave packet. The expectation value of momentum $`p_t`$, calculated from $`|\varphi (p,t)|^2`$ and plotted as the solid curve, is once again consistent with the classical trajectory (dashed curve), except near the discontinuous, impulsive change in momentum values at the ‘bounce’. The shape of the momentum-distribution follows the form expected from Eqn. (13), namely uniform translation with no change in shape, from $`t=1t=3`$ and then again from $`t=8t=10`$, that is, during the time when it is not in collision with the wall, but with a definite final change in shape, compared to the initial $`|\varphi _0(p)|^2`$ superimposed on the $`t=10`$ result. The dotted vertical lines indicate the values of $`p=0`$, and also the classically expected minimum and maximum values of momentum given by $`\pm p_M=\pm \sqrt{2mE}=\pm \sqrt{2mFz_0}`$.

The asymmetric shape of $`|\varphi (p,t)|^2`$ which is obvious after one period can be understood using purely classical arguments. For example, we can examine the classical trajectories of particles released from the same height, $`H`$, but with varying initial velocities (or momenta); the resulting time-dependent momenta, $`p(t)`$ versus $`t`$, are shown in Fig. 2. We can see the same asymmetry arising at $`t=T_{cl}`$ as for the quantum mechanical $`|\varphi (p,t)|^2`$ in Fig. 1; trajectories with momentum values $`p\stackrel{<}{_{}}0`$ at $`t=0`$ are shifted upwards towards $`p0`$ at $`t=T_{cl}`$, while those with $`p\stackrel{>}{_{}}0`$ are also moved to higher values to form the observed high-momentum ’tail’. (A very similar classical explanation of the time-dependent behavior of the momentum-space probability density for an otherwise free particle undergoing a ’bounce’ from an infinite wall is discussed in Ref. ; in that case, the observed asymmetry in $`|\varphi (p,t)|^2`$ during the ’bounce’ is due to the fact that the fastest (highest momentum) components of the wavepacket strike the wall first, and hence are reversed in direction, before the slow components.)

We next focus on the time-development of quantum mechanical expectation values over the first eight classical periods and plot $`z_t`$ and $`p_t`$ over this interval in Fig. 3. Once again, the solid curves are the quantum mechanical expectation values, while the dashed curves are the classical trajectories. Clearly, the quantum results generally track the classical values, but become increasingly smoothed out as the packet spreads in time. We show in the same figure the time-development of the autocorrelation function, $`|A(t)|^2`$ versus $`t`$, over the same time interval which shows the (increasingly) approximate nature of the return to the initial state at multiples of the classical period: these results look very similar to data for wave packets in the infinite well ,  as it approaches the collapsed state, while the same plot for the oscillator gives a completely periodic structure with $`A(t)`$ returning exactly to $`1`$ at integral multiples of $`\tau =2\pi /\omega `$.

In Fig. 4, we plot the time-dependent uncertainties in the position ($`\mathrm{\Delta }z_t`$, top) and momentum ($`\mathrm{\Delta }p_t`$, bottom) over the first eight classical periods. We note that the behavior of $`\mathrm{\Delta }z_t`$ is somewhat similar to that for the oscillator case: it initially increases in a manner consistent with Eqn. (16) for much of the first half period, but then shrinks dramatically near $`T_{cl}/2`$ during the collision with the wall (as in Ref. ), and finally shrinks more slowly during the second half period, returning to something close to its initial spread. This behavior is different than the case of the infinite well where the wavepacket spreading is mostly determined by the free-particle expansion and more similar to the harmonic oscillator case in Eqn. (6) and general bound state systems.

In Fig. 4, the expectation value during each period is superimposed with two expressions for the spread (shown as pairs of dashed curves); specifically, for the first period we use

$$\mathrm{\Delta }z_t=\mathrm{\Delta }z_0\sqrt{1+\left(\frac{t}{t_0}\right)^2}\text{and}\mathrm{\Delta }z_t=\mathrm{\Delta }z_0\sqrt{1+\left(\frac{T_{cl}t}{t_0}\right)^2}$$
(22)
and we extend these cyclically to later periods as shown in the figure. These agree fairly well with the observed time-dependence during the first cycle except during the collision phase, but become increasing bad representations as the wave packet spreads into its collapsed phase. Thus, over at least the first few periods, the return to the initial spread is similar to that in Eqn. (6).

In the same figure, the spread in momentum is also illustrative of both classical motion and results for accelerated wave packets. Over the first cycle, when the packet is not colliding with the wall, the spread in momentum remains constant, as expected from Eqn. (13), while it increases dramatically during the collision time, just as with the ‘bouncing wavepacket’ in Ref. . After each bounce, however, the spread in momentum has grown slightly so that $`\mathrm{\Delta }p_t`$ increases in a quasi-stepwise manner between collision, with flat plateaus between each impulsive ‘spike’.

We now turn our attention to the long-term behavior of these solutions and plot in Fig. 5 the same quantities, $`z_t`$, $`p_t`$, and $`|A(t)|^2`$, as in Fig. 3, but now over an interval containing two expected revivals; these are indicated by the bold vertical dashed lines. We also indicate the collapse time, given by Eqn. (21), after $`t=0`$ as well as one unit of collapse time around each expected revival. The reformation of the wavepacket at the revivals in clearly apparent in all three variables, as is the approach to the more constant values during the collapsed phases.

The horizontal dashed lines around which the quantum mechanical expectation values cluster during the collapsed phases are given by the average values derived from purely classical probability densities. Using standard arguments about how much time a particle spends in a small interval of position, we can find the normalized position-space classical probability as

$$P_{CL}(z)=\frac{1}{2\sqrt{A(Az)}}\text{for }0<z<A$$
(23)
and which vanishes for all other values of $`z`$. The classical turning point, $`A`$, is given by $`E=FA`$, and we naturally associate $`A`$ with the initial position $`z_0`$ in our quantum analysis. With this identification, the classical average values are

$`z`$ $`=`$ $`{\displaystyle _0^{z_0}}zP_{CL}(z)𝑑z={\displaystyle \frac{2}{3}}z_0`$
$`z^2`$ $`=`$ $`{\displaystyle _0^{z_0}}z^2P_{CL}(z)𝑑z={\displaystyle \frac{8}{15}}z_0^2`$ (24)
$`\mathrm{\Delta }z`$ $`=`$ $`\sqrt{z^2z^2}={\displaystyle \frac{2}{\sqrt{45}}}z_0.`$
For the momentum-space probability distribution, we note that for a constant force, the probability density is uniform over the allowed space of values, namely

$$P_{CL}(p)=\frac{1}{2p_M}\text{for}p_M<p<+p_M$$
(25)
and vanishing for $`|p|>p_M`$: in this expression we have $`p_M=\sqrt{2mE}=\sqrt{2mFA}=\sqrt{2mFz_0}`$. The classical averages are then given by

$`p`$ $`=`$ $`0`$
$`p^2`$ $`=`$ $`{\displaystyle \frac{1}{3}}p_M^2`$ (26)
$`\mathrm{\Delta }p`$ $`=`$ $`\sqrt{p^2p^2}={\displaystyle \frac{1}{\sqrt{3}}}p_M.`$
These results are known to agree well (after local averaging) with the position- and momentum-space probability densities corresponding to time-independent energy eigenstates for this type of potential . We would expect that the quantum wave packet expectation values to closely match these classical results (again, after local averaging) during the collapsed phase as well, since the wave function is closer to being an incoherent sum of many such stationary states rather than the highly coherent, well-localized, quasi-periodic state close to the initial time and later revivals. The time-dependent values of $`z_t`$ and $`p_t`$ in Fig. 5 do collapse to just these values.

For the autocorrelation function we have no classical analog, but we can formally evaluate $`|A(t)|^2`$ in terms of the expansion coefficients, $`c_n`$, of the wave packet. In the limit that the various time-dependent components of the packet get out of phase, we have the expectation that during the collapsed interval

$$|A(t)|^2=\left|\underset{n=1}{\overset{\mathrm{}}{}}|c_n|^2e^{iE_nt/\mathrm{}}\right|^2=\underset{n,m=1}{}|c_n|^2|c_m|^2e^{i(E_nE_m)t/\mathrm{}}\underset{n=1}{}|c_n|^4$$
(27)
as the ‘off-diagonal’ terms tend to cancel as they have lose the phase coherence built into the initial state. This limiting value, calculated using the numerically determined $`c_n`$, is plotted in Fig. 5 (as the bold horizontal dashed line) and compares well to the locally averaged value of $`|A(t)|^2`$ during much of the collapsed period.

In Fig. 6, we plot the position- and momentum-spreads, as in Fig. 4, but once again over the longer time interval, and superimpose the classical (constant) expectation values for $`\mathrm{\Delta }z`$ and $`\mathrm{\Delta }p`$ as bold horizontal dashed lines. In both plots, the quantum wave packet expectation values and uncertainties cluster appropriately around the classical results of Eqns. (24) and (26) during the collapsed phase as expected.

Finally, one of the most interesting results of Ref.  is that the wave packet revivals happen in such a way that the probability densities are almost completely out of phase with the classical motion. In order to visualize this effect in the expectation value approach followed here, we plot again the same quantities as in Fig. 3, namely $`z_t`$, $`p_t`$, and $`|A(t)|^2`$, also over a time interval corresponding to 8 classical periods, but now centered around the expected revival time. In each case we indeed see that the quantum expectation values (solid) are approximately half a cycle out of phase with the corresponding classical trajectory values (dashed) and that the autocorrelation function is peaked at regular intervals, but offset from integral multiples of the classical period by roughly $`T_{cl}/2`$.

In conclusion, we have reexamined the interesting problem of quasi-Gaussian wave packets for the quantum bouncer, focusing on the numerical calculation and visualization of the results through the expectation values and uncertainties of position and momentum. We have focused on a single classical period to illustrate the important similarities and significant differences between this case and well-studied periodic wave packet solutions to other familiar problems such as the harmonic oscillator and infinite well as well as to the closed-form solutions for the constant acceleration case in Eqns. (15) and (16). We have also examined these quantities over the first few periods, well before the collapsed phase and compared them to classical trajectory expectations. Finally, we have confirmed, numerically and through visualization, that the quantum expectation values approach (in a locally averaged sense) the purely classical results from the time-independent probability distributions in Eqns. (24) and (26) during the collapsed phase as well as illustrating the expected phase relationships between the classical trajectories and quantum expectation values near the revivals.

Acknowledgments

The work of MAD was supported, in part, by the Commonwealth College of The Pennsylvania State University under a Research Development Grant (RDG); the work of RR was supported, in part, by NSF grant DUE-9950702.

Figure Captions

* Gaussian wave packet solutions for the quantum bouncer in position-space ($`|\psi (z,t)|^2`$ versus $`z`$, left) and momentum-space ($`|\varphi (p,t)|^2`$ versus $`p`$, right) for various times during the first classical period. The solid curves represent the time-dependent expectation values of position ($`z_t`$, left) and momentum ($`p_t`$, right) for these solutions. The similar dashed curves are the classical trajectories, $`z(t)`$ (left) and $`p(t)`$ (right), superimposed. The wave packet parameters in Eqn. (19) are used. For the momentum-space figure, the vertical dotted lines represent the values $`p=0`$ and the classical extremal values of momentum, $`\pm p_M=\pm \sqrt{2mFz_0}`$.
* Plots of the classical momentum versus time, $`p(t)`$ versus $`t`$, for objects released from a height $`H`$ with initial momenta given by $`p(0)=0`$ (solid), $`p(0)=\pm \mathrm{\Delta }p`$ (dotted), and $`p(0)=\pm 2\mathrm{\Delta }p`$ (dashed). After one classical period, the distribution of momentum values is asymmetric about $`p=0`$, in just the same way as the $`|\varphi (p,t)|^2`$ plot in Fig. 1; compare the momentum values defined by the $`t=0`$, $`(+2\mathrm{\Delta }p,2\mathrm{\Delta }p)`$ ’band’ to the corresponding values at $`t=T_{cl}`$.
* The solid curves represent the expectation values of position ($`z_t`$, top) and momentum ($`p_t`$, middle) and the autocorrelation function ($`|A(t)|^2`$, bottom) versus time for the first eight classical periods of motion for the bouncing wave packet. The wave packet parameters in Eqn. (19) are used. The dashed lines are the classical trajectories for comparison. The horizontal line in the top figure corresponds to the average value of position over one classical period, namely $`z=2z_0/3`$, from Eqn. (24).
* Same as Fig. 3 except that the time-dependent values of the uncertainty in position ($`\mathrm{\Delta }z_t`$, top) and momentum ($`\mathrm{\Delta }p_t`$, bottom) are shown over the first eight periods. The dashed curves in the top picture correspond to the spreading (anti-spreading) of a Gaussian wave packet as it ’slides’ downhill (uphill) as in Eqn. (22); these expressions agree well with the observed values of $`\mathrm{\Delta }z_t`$ during the first period, except during those times when the packet is ‘bouncing’ against the infinite wall at $`z=0`$, near $`T_{cl}/2=5`$. The dashed line in the $`\mathrm{\Delta }p_t`$ figure (bottom) indicates the classical spread in momentum defined in Eqn. (26), $`\mathrm{\Delta }p=p_M/\sqrt{3}`$.
* Same as Fig. 3, but now over a longer time scale, containing two revivals. The expected locations of the revivals, according to Eqn. (17), are shown as vertical dashed lines, and times within one collapse time, $`T_{coll}`$, of the initial time and these two revival times are also shown as vertical dotted lines. The classical average values of position and momentum ($`z=2z_0/3`$ and $`p=0`$) from Eqns. (24) and (26) are shown as bold horizontal dashed lines, as is the expected value of the autocorrelation function during the collapsed phase, $`|A(t)|^2=_n|c_n|^4`$, from Eqn. (27), using the numerically determined values of the $`c_n`$.
* Same as Fig. 3, expect for the long term values of $`\mathrm{\Delta }z_t`$ and $`\mathrm{\Delta }p_t`$. The classical average values of these quantities ($`\mathrm{\Delta }z=2z_0/\sqrt{45}`$ and $`\mathrm{\Delta }p=p_M/\sqrt{3}`$) from Eqns. (24) and (26) are shown as bold horizontal dashed lines.
* Same as Fig. 3, except for an eight classical period time window around the first revival which is indicated by the bold vertical dashed line. Note how the quantum mechanical averages, $`z_t`$ and $`p_t`$ (solid curves), are out of phase with the classical trajectory values (dashed curves). In the bottom figure, the dashed curve corresponds to the time-development of the auto-correlation function, $`|A(t)|^2`$, over the first eight periods of the motion, for comparison.


format: text

text:
# Inclusion relations among separability criteria

## I Introduction

The development of criteria for entanglement and separability is one aspect of the current research efforts in quantum information theory that is receiving, and certainly deserves, considerable attention . Indeed, much progress has recently been made in consolidating such a cornerstone of the theory of quantum entanglement . The relevant state-space here is of a high dimensionality, already 15 dimensions in the simplest instance of two-qubit systems. The systematic exploration of these spaces can provide us with valuable insight into some of the theoretical questions extant.

As a matter of fact, important steps have been recently made towards a systematic exploration of the space of arbitrary (pure or mixed) states of composite quantum systems  in order to determine the typical features exhibited by these states with regards to the phenomenon of quantum entanglement . Entanglement is one of the most fundamental and non-classical features exhibited by quantum systems , that lies at the basis of some of the most important processes studied by quantum information theory  such as quantum cryptographic key distribution , quantum teleportation , superdense coding , and quantum computation .

It is well known  that, for a composite quantum system, a state described by the density matrix $`\rho `$ is called “entangled” if it can not be represented as a mixture of factorizable pure states. Otherwise, the state is called separable. The above definition is physically meaningful because entangled states (unlike separable states) cannot be prepared locally by acting on each subsystem individually .

The separability question has quite interesting echoes in information theory and its associate information measures or entropies. When one deals with a classical composite system described by a suitable probability distribution defined over the concomitant phase space, the entropy of any of its subsystems is always equal or smaller than the entropy characterizing the whole system. This is also the case for separable states of a composite quantum system . In contrast, a subsystem of a quantum system described by an entangled state may have an entropy greater than the entropy of the whole system. Indeed, the von Neumann entropy of either of the subsystems of a bipartite quantum system described (as a whole) by a pure state provides a natural measure of the amount of entanglement of such state. Thus, a pure state (which has vanishing entropy) is entangled if and only if its subsystems have an entropy larger than the one associated with the system as a whole.

Regrettably enough, the situation is more complex when the composite system is described by a mixed state: there are entangled mixed states such that the entropy of the complete system is smaller than the entropy of one of its subsystems. Alas, entangled mixed states such that the entropy of the system as a whole is larger than the entropy of either of its subsystems do exist as well. Consequently, the classical inequalities relating the entropy of the whole system with the entropy of its subsystems provide only necessary, but not sufficient, conditions for quantum separability. There are several entropic (or information) measures that can be used in order to implement these criteria for separability. Considerable attention has been paid, in this regard, to the $`q`$-entropies , which incorporate both Rényi’s  and Tsallis’  families of information measures as special instances (both admitting, in turn, Shannon’s measure as the particular case associated with the limit $`q1`$). The reader is referred to Appendix A for a brief review on $`q`$-entropies.

The early motivation for the studies reported in  was the development of practical separability criteria for density matrices. The discovery by Peres of the partial transpose criteria, which for two-qubits and qubit-qutrit systems turned out to be both necessary and sufficient, rendered that original motivation somewhat outmoded. In point of fact, it is not possible to find a necessary and sufficient criterion for separability based solely upon the eigenvalue spectra of the three density matrices $`\rho _{AB},\rho _A=Tr_B[\rho _{AB}]`$, and $`\rho _B=Tr_A[\rho _{AB}]`$ associated with a composite system $`AB`$ .

Interesting concepts that revolve around the separability issue have been developed over the years. A beautiful account is given in Terhal in . Among them we find criteria like the so-called Majorization, Reduction and Positive Partial Transpose (PPT) together with the concept of distillability . Quantum entanglement is a fundamental aspect of quantum physics that deserves to be investigated in full detail from all possible points of view. The chain of implications, and the related inclusion relation, among the different separability criteria is certainly a vantage point worth of detailed scrutiny. It is our purpose here to revisit, with such a goal in mind, the separability question by means of an exhaustive Monte Carlo exploration involving the whole space of pure and mixed states. Such an effort should shed some light on the inclusion issues that interest us here. Concrete numerical evidence will thus be provided on the relations among the separability criteria. We will then be able to quantify, for a bipartite system of arbitrary dimension, the proportion of states $`\rho `$ that can be distilled according to a definite criterion. This numerical exploration could be viewed as a complement on the review paper by Terhal , because some questions posed by her will receive an answer in this work.

The paper is organized as follows. We sketch in Section II the different separability criteria to be investigated and discuss some mathematical and numerical techniques used in our survey in Section III. Our results are reported in Section IV, and some conclusions are drawn in Section V. For the sake of completeness, we include an Appendix on $`q`$-entropies.

## II Brief sketch on separability criteria

From a historic viewpoint, the first separability criterion is that of Bell (see  and references therein). For every pure entangled state there is a Bell inequality that is violated. It is not known, however, whether in the case of many entangled mixed states, violations exist. There does exist a witness for every entangled state though . It was shown by Horodecki et al. that a density matrix $`\rho \rho _{AB}`$ is entangled if and only if there exists an entanglement witness (a hermitian super-operator $`\widehat{W}=\widehat{W}^{}`$) such that

$`Tr\widehat{W}\rho `$ $``$ $`0,\mathrm{while}`$ (1)
$`Tr\widehat{W}\rho `$ $``$ $`0,\mathrm{for}\mathrm{all}\mathrm{separable}\mathrm{states}.`$ (2)
A special, but quite important LOCC operational separability criterion, necessary but not sufficient, is provided by the positive partial transpose (PPT) one. Let $`T`$ stand for matrix transposition. The PPT requires that

$$[\widehat{1}\widehat{T}](\rho )0.$$
(3)
Another operational criterion is called the reduction criterion, that is satisfied, for a given state $`\rho \rho _{AB}`$, when both

$`\widehat{1}\rho _B\rho `$ $``$ $`0`$ (4)
$`\rho _A\widehat{1}\rho `$ $``$ $`0.`$ (5)
Intuitively, the distillable entanglement is the maximum asymptotic yield of singleton states that can be obtained, via LOCC, from a given mixed state. Horodecki et al.  demonstrated that any entangled mixed state of two qubits can be distilled to obtain the singleton. This is not true in general. There are entangled mixed states of two qutrits, for instance, that cannot be distilled, so that they are useless for quantum communication. In our scenario an important fact is that all states that violate the reduction criterion are distillable .

Entanglement witnesses completely characterize the set of separable states. Alas, they are not usually associated to a simple computational treatment, except in the PPT instance. Thus, in order to decide whether a given state $`\rho `$ is entangled one needs additional criteria, functional separability ones . One of them associates PPT to the rank of a matrix. Consider two subsystems $`A`$, $`B`$ whose description is made, respectively, in the Hilbert spaces $`_n`$ and $`_m`$. Focus attention now in the density matrix $`\rho \rho _{AB}`$ for the associated composite system. If

1. $`\rho `$ has PPT, and
2. its rank $``$ is such that $`\mathrm{max}(n,m)`$,

then, as was proved in , $`\rho `$ is separable. The above referred to entropic criteria are also functional separability ones. Still another one is majorization.

Let $`\{\lambda _i\}`$ be the set of eigenvalues of the matrix $`\xi _1`$ and $`\{\gamma _i\}`$ be the set of eigenvalues of the matrix $`\xi _2`$. We assert that the ordered set of eigenvalues $`\stackrel{}{\lambda }`$ of $`\xi _1`$ majorizes the ordered set of eigenvalues $`\stackrel{}{\gamma }`$ of $`\xi _2`$ (and writes $`\stackrel{}{\lambda }\stackrel{}{\gamma }`$) when $`_{i=1}^k\lambda _i_{i=1}^k\gamma _i`$ for all $`k`$. It has been shown  that, for all separable states $`\rho _{AB}\rho `$,

$`\stackrel{}{\lambda }_{\rho _A}`$ $``$ $`\stackrel{}{\lambda }_\rho ,\mathrm{and}`$ (6)
$`\stackrel{}{\lambda }_{\rho _B}`$ $``$ $`\stackrel{}{\lambda }_\rho .`$ (7)
There is an intimate relation between this majorization criterion and entropic inequalities, as discussed in .

## III Separability probabilities: exploring the whole state space

We promised in the Introduction to perform a systematic numerical survey of the properties of arbitrary (pure and mixed) states of a given quantum system by recourse to an exhaustive exploration of the concomitant state-space $`𝒮`$. To such an end it is necessary to introduce an appropriate measure $`\mu `$ on this space. Such a measure is needed to compute volumes within $`𝒮`$, as well as to determine what is to be understood by a uniform distribution of states on $`𝒮`$. The natural measure that we are going to adopt here is taken from the work of Zyczkowski et al. . An arbitrary (pure or mixed) state $`\rho `$ of a quantum system described by an $`N`$-dimensional Hilbert space can always be expressed as the product of three matrices,

$$\rho =UD[\{\lambda _i\}]U^{}.$$
(8)
Here $`U`$ is an $`N\times N`$ unitary matrix and $`D[\{\lambda _i\}]`$ is an $`N\times N`$ diagonal matrix whose diagonal elements are $`\{\lambda _1,\mathrm{},\lambda _N\}`$, with $`0\lambda _i1`$, and $`_i\lambda _i=1`$. The group of unitary matrices $`U(N)`$ is endowed with a unique, uniform measure: the Haar measure $`\nu `$ . On the other hand, the $`N`$-simplex $`\mathrm{\Delta }`$, consisting of all the real $`N`$-uples $`\{\lambda _1,\mathrm{},\lambda _N\}`$ appearing in (8), is a subset of a $`(N1)`$-dimensional hyperplane of $`^N`$. Consequently, the standard normalized Lebesgue measure $`_{N1}`$ on $`^{N1}`$ provides a natural measure for $`\mathrm{\Delta }`$. The aforementioned measures on $`U(N)`$ and $`\mathrm{\Delta }`$ lead then to a natural measure $`\mu `$ on the set $`𝒮`$ of all the states of our quantum system , namely,

$$\mu =\nu \times _{N1}.$$
(9)
All our present considerations are based on the assumption that the uniform distribution of states of a quantum system is the one determined by the measure (9). Thus, in our numerical computations we are going to randomly generate states according to the measure (9).

## IV Survey’s results

### A The overall scenario

An overall picture of the situation we encounter is sketched in Fig. 1, that is to be compared to Fig. 3 of . Notice that our numerical exploration allows us to dispense with Terhal’s interrogation signs. This constitutes part of the original content of the present Communication.

The set of all mixed states presents an onion-like shape, as conjectured by Terhal . Which among these states are separable? As reviewed above, several criteria are available. We start with the $`q`$-entropic one (see the Appendix and ). By using a definite value of $`q`$, namely $`q=\mathrm{}`$, and the sign of the associated, conditional $`q`$-entropy, we are able to define a closed sub-region, whose states are supposedly separable. This region has a definite border, that separates it from the sub-region of states entangled according to this criterion. What we see now is that, if we use now other separability-criteria, the associated sub-regions shrink in a manner prescribed by the particular criterion one employs. The shrinking process ends when one reaches the sub-region defined by the Positive Partial Transpose (PPT) criterion, which is a necessary and sufficient separability condition for $`2\times 2`$ and $`2\times 3`$ systems, being only necessary for higher dimensions.

Summing up, the volume of states which are separable according to different criteria diminish as we use stronger and stronger criteria. There is a first shrinking stage associated to entropic criteria, from its Von Neumann ($`q=1`$) size, as $`q`$ grows, to the limit case $`q\mathrm{}`$ . A second stage involves majorization, reduction, and, finally positive partial transpose (PPT) .

### B PPT and Reduction

We report now on our state-space exploration with regards to the probability of finding a state with positive partial transpose. The results are depicted in Fig. 2. The solid line corresponds to states with dimension $`N=2\times N_2`$, while the dashed line corresponds to $`N=3\times N_2`$ states. Note how similar are the pertinent values in both cases. The tiny difference between them can be inspected in the inset (a semi-logarithmic plot). To a good approximation, our PPT probabilities decrease exponentially.

Fig. 3 deals instead with the probability of finding a state which obeys the strictures of the reduction criterion, for $`N=2\times N_2`$ (solid line) and $`N=3\times N_2`$ (dashed line). As a matter of fact, PPT and reduction coincide for $`N=2\times N_2`$. It is known that if a state satisfies PPT, it automatically verifies the reduction criterion . Here we have demonstrated that, at least in the $`N=2\times N_2`$-instance, the converse is also true. However, in the $`N=3\times N_2`$-case, it is much more likely to encounter a state that verifies reduction than one that verifies PPT.

### C Entropic criteria and Majorization

We begin with a brief recapitulation of former $`q`$-entropic results. The situation encountered in  was that the “best” result within the framework of the “classical $`q`$-entropic inequalities” as a separability criterion was reached using the limit case $`q\mathrm{}`$, but considerably less attention was paid to other values of $`q`$. This was remedied in , where the question of $`q`$-entropic inequalities for finite $`q`$-values was extensively discussed. It was there re-confirmed that the above mentioned limit case does indeed the better job as far as separability questions are concerned . For such a reason, this limit $`q`$-value is the only one to be employed below. See the Appendix for more details on $`q`$-entropies.

In Fig. 4 we depict the probability of finding a state which, for $`q\mathrm{}`$, has its two relative $`q`$-entropies positive (dashed curves). In view of the intimate relation of entropic inequalities with majorization , we also analyze in Fig. 4 the probability that a state is completely majorized by both of their subsystems (solid line). It is shown in  that, if $`\rho _{AB}`$ satisfies the reduction criterion, its two associated relative $`q`$\- entropies are non-negative as well.

In the same work the authors assert that majorization is not implied by the relative entropy criteria. Our results confirm this assessment. In Fig. 4, the lower curves correspond to states $`\rho `$ with $`N=2\times N_2`$, while the upper curves have $`N=3\times N_2`$. Majorization results and $`q`$-entropic do coincide for two-qubits systems ($`N_1=N_2=2`$). More generally, majorization probabilities are a lower bound to probabilities for relative $`q`$-entropic positivity, an interesting new result, as far as we know. Notice also that the two approaches yield quite similar results in the $`N=3\times N_2`$ case.

### D Comparing more than two criteria together

We compare now the reduction criterion to the PPT one. The former is implied by the latter but is nonetheless a significant condition since its violation implies the possibility of recovering entanglement by distillation, which is as yet unclear for states that violate PPT . Fig. 5 a) depicts the probability that state $`\rho `$ with $`N=3\times N_2`$ either:

1. has a positive partial transpose and does not violate the reduction criterion, or
2. has a non positive partial transpose and violates reduction.

Remember that in the case $`N=2\times N_2`$, the two criteria always coincide . For $`3\times N_2`$ the agreement between the two criteria becomes better and better as $`N_2`$ augments.

Of more interest is to compare the relations among PPT, majorization, and the entropic criteria (Fig. 5b), since it is not yet known how the majorization criterion is related to other separability criteria like PPT, undistillability, and reduction . In this vein, Fig. 5 b) plots the “coincidence-probability” between, respectively,

1. PPT and majorization (solid line), and
2. PPT and the $`q`$-entropic criterion (dashed line).

The curves on the top correspond to $`N=2\times N_2`$, while those at the bottom to $`N=3\times N_2`$. In this last case the two curves agree with each other quite well.

The conclusion here is that, as $`N_2`$ augments, the probability of coincidence among the three criteria, and in particular between majorization and PPT (our main concern), rapidly diminishes at first, and stabilizes itself afterwards. For two qubits the three criteria do agree with each other to a large extent.

Fig. 6 a) depicts the probability that, for a given state $`\rho `$,

1. reduction and majorization (solid line) and
2. reduction and the $`q`$-entropic criterion (dashed line)

yield the same conclusion as regards separability. Without PPT in the game, and opposite to what we encountered in Fig. 5, we find better coincidence for $`N=3\times N_2`$ systems (top) than for $`N=2\times N_2`$ (bottom). The deterioration of the degree of agreement as $`N_2`$ grows is similar to that of Fig. 5, though.

Fig. 6 b) represents the probability that a state, for $`q\mathrm{}`$, either:

1. has both positive relative $`q`$-entropies and satisfies the majorization criterion, or
2. has a negative relative $`q`$-entropy and is majorized by both of their subsystems.

The solid line corresponds to the case $`N=2\times N_2`$, while the dashed lines corresponds to the $`N=3\times N_2`$ instance. These results together with those of Figs. 4-5 could be read as implying that majorization and the $`q`$-entropic criteria provide almost the same answer for dimensions greater or equal than $`N=3\times N_2`$.

Finally, in Fig. 7 we look for the probability $`P_{agree}`$ that all criteria considered in the present work do lead to the same conclusion on the separability issue. $`P_{agree}`$ is plotted as a function of the total dimension $`N=N_1\times N_2`$, with $`N_1=2`$ (solid line) and $`N_1=3`$ (dashed line). The agreement is quite good for two qubits, deteriorates first as $`N_2`$ grows, and rapidly stabilizes itself around a value of 0.26 for $`N_1=2`$ and of 0.1 for $`N_1=3`$.

### E Distilling

Let us at now consider the results plotted in Fig. 8. We ask first for the relative number of states that violate the reduction criterion and are thus distillable  (solid line), and appreciate the fact that, as $`N`$ grows, so does the probability of finding distillable states. On the other hand, the probability of encountering states that violate the majorization criterion, represented by dashed lines, is much lower than that associated to distillation.

For both criteria, the upper solid line corresponds to the case $`N=2\times N_2`$, and the lower one to $`N=3\times N_2`$. The dashed curve with crosses represents the case $`N=2\times N_2`$, while the one with squares indicates the $`N=3\times N_2`$ instance. The dependence with $`N_2`$ of the dashed curves (majorization violation) is not so strong as that of the solid ones (distillability). Our results are lower bounds to the total volume of states that can be destilled.

## V Conclusions

We have performed a systematic numerical survey of the space of pure and mixed states of bipartite systems of dimension $`2\times N_2`$ and $`3\times N_2`$ in order to investigate the relationships ensuing among different separability criteria. Our main results are

* Regarding the line of separability implication, see our graph in Fig. 1 and compare with the similar one of Terhal’s (her Fig. 3). Her interrogation signs receive an answer in our Fig.1
* It is known that if a state satisfies PPT, it automatically verifies the reduction criterion . In the present work we show that in the $`N=2\times N_2`$-instance, the converse is also true. In the $`N=3\times N_2`$-case, it is much more likely to encounter a state that verifies reduction than one that verifies PPT.
* We have numerically verified the assertion made in  that majorization is not implied by the relative entropic criteria. Majorization results and $`q`$-entropic criteria coincide for two-qubits systems. In general, majorization probabilities constitutes lower bound for relative $`q`$-entropic positivity.
* Regarding the relation between majorization and PPT, the agreement between the criteria deteriorates as $`N_2`$ grows.
* For dimensions $`3\times N_2`$, as illustrated by Figs. 4-5, majorization and the $`q`$-entropic criteria provide almost the same answers.

The present authors believe that the results of this numerical exploration shed some light on the intricacies of the separability issue.

###### Acknowledgements.

This work was partially supported by the MCyT grants BFM2002-0341 and SAB2001-006(Spain), by the Gouvernment of Balearic Islandsand by CONICET (Argentine Agency).

FIGURE CAPTIONS

Fig. 1 - Schematics of the inclusion relations among separbility criteria as given by the volume occupied by states $`\rho `$ for a given dimension $`N`$ which fulfill them.

Fig. 2 - Probability of finding a state with positive partial transpose. The solid line corresponds to states with dimension $`N=2\times N_2`$, while the dashed line corresponds to $`N=3\times N_2`$ states. The difference between these curves can be appreciated in the inset (semi-logarithmic plot). Our probabilities decrease, to a good approximation, in exponential fashion.

Fig. 3 - Probability of finding a state fulfilling the reduction criterion for $`N=2\times N_2`$ (solid line) and $`N=3\times N_2`$ (dashed line). The two probabilities coincide for $`N=2\times N_2`$.

Fig. 4 - Probability of finding a state whose two relative $`q`$-entropies are positive for $`q\mathrm{}`$ (dashed curves). The probability that a state be completely majorized by both of their subsystems is represented by the solid line. Bottom: curves correspond to states $`\rho `$ with $`N=2\times N_2`$. Top: $`N=3\times N_2`$.

Fig. 5 a) Probability that the state $`\rho `$ with $`N=3\times N_2`$ either has i) a positive partial transpose and does not violate the reduction criterion, or ii) has a non positive partial transpose and violates reduction. In the case $`N=2\times N_2`$ the outcome is always unity. Fig. 5 b) Probability that

1. PPT and majorization (solid line) and,
2. PPT and the $`q`$-entropic criterion (dashed line)

lead to the same conclusion regarding separability. Top: $`N=2\times N_2`$. Bottom: $`N=3\times N_2`$.

Fig. 6 a) Probability that reduction and majorization (solid line) and reduction and the $`q`$-entropic criterion (dashed line) yield the same conclusion reagrding separability. Top: $`N=3\times N_2`$. Bottom: $`N=2\times N_2`$ (lower curves). Fig. 6 b) Probability that a state, for $`q\mathrm{}`$, either i) has both positive relative $`q`$-entropies and fulfills majorization, or ii) has a negative relative $`q`$-entropy and is majorized by both of their subsystems. The solid line corresponds to the case $`N=2\times N_2`$, while the dashed line corresponds to $`N=3\times N_2`$.

Fig. 7 - Total probability that all criteria considered in the present work lead to the same conclusion regarding separability. Probabilities are plotted as a function of the total dimension $`N=N_1\times N_2`$, with $`N_1=2`$ (solid line) and $`N_1=3`$ (dashed line).

Fig. 8 - Solid line: probability that a state violates the reduction criterion. Dashed line: the same for violation of the majorization criterion. Top: $`N=2\times N_2`$. Bottom: $`N=3\times N_2`$. The dashed curve with crosses represents the case $`N=2\times N_2`$, while the one with squares indicates the $`N=3\times N_2`$ instance.

## A q-Information measures and the issue of quantum separability

There are several useful entropic (or information) measures for the investigation of a quite important subject: the violation of classical entropic inequalities by quantum entangled states. The von Neumann measure

$$S_1=Tr\left(\rho \mathrm{ln}\rho \right),$$
(A1)
is important because of its relationship with the thermodynamic entropy. The $`q`$-entropy, which is a function of the quantity

$$\omega _q=Tr\left(\rho ^q\right),$$
(A2)
provides one with a whole family of entropic measures. In the limit $`q1`$ these measures incorporate (A1) as a particular instance. Most of the applications of $`q`$-entropies to physics involve either the Rényi entropy ,

$$S_q^{(R)}=\frac{1}{1q}\mathrm{ln}\left(\omega _q\right),$$
(A3)
or the Tsallis entropy

$$S_q^{(T)}=\frac{1}{q1}\left(1\omega _q\right).$$
(A4)
We reiterate that the von Neumann measure (A1) constitutes a particular instance of both Rényi’s and Tsallis’ entropies, obtained in the limit $`q1`$. The most distinctive single property of Tsallis’ entropy is its nonextensivity. The Tsallis entropy of a composite system $`AB`$ whose state is described by a factorizable density matrix, $`\rho _{AB}=\rho _A\rho _B`$, is given by Tsallis’ $`q`$-additivity law,

$$S_q^{(T)}(\rho _{AB})=S_q^{(T)}(\rho _A)+S_q^{(T)}(\rho _B)+(1q)S_q^{(T)}(\rho _A)S_q^{(T)}(\rho _B).$$
(A5)
In contrast, Rényi’s entropy is extensive. That is, if $`\rho _{AB}=\rho _A\rho _B`$,

$$S_q^{(R)}(\rho _{AB})=S_q^{(R)}(\rho _A)+S_q^{(R)}(\rho _B).$$
(A6)
Tsallis’ and Rényi’s measures are related through

$$S_q^{(T)}=F(S_q^{(R)}),$$
(A7)
where the function $`F`$ is given by

$$F(x)=\frac{1}{1q}\left\{e^{(1q)x}1\right\}.$$
(A8)
An immediate consequence of equations (A7-A8) is that, for all non vanishing values of $`q`$, Tsallis’ measure $`S_q^{(T)}`$ is a monotonic increasing function of Rényi’s measure $`S_q^{(R)}`$.

Considerably attention has been recently paid to a relative entropic measure based upon Tsallis’ functional defined as

$$S_q^{(T)}(A|B)=\frac{S_q^{(T)}(\rho _{AB})S_q^{(T)}(\rho _B)}{1+(1q)S_q^{(T)}(\rho _B)}.$$
(A9)
Here $`\rho _{AB}`$ designs an arbitrary quantum state of the composite system $`AB`$, not necessarily factorizable nor separable, and $`\rho _B=Tr_A(\rho _{AB})`$. The relative $`q`$-entropy $`S_q^{(T)}(B|A)`$ is defined in a similar way as (A9), replacing $`\rho _B`$ by $`\rho _A=Tr_B(\rho _{AB})`$. The relative $`q`$-entropy (A9) has been recently studied in connection with the separability of density matrices describing composite quantum systems . For separable states, we have

$`S_q^{(T)}(A|B)`$ $``$ $`0,`$ (A10)
$`S_q^{(T)}(B|A)`$ $``$ $`0.`$ (A11)
On the contrary, there are entangled states that have negative relative $`q`$-entropies. That is, for some entangled states one (or both) of the inequalities (A10) are not verified.

Notice that the denominator in (A9),

$$1+(1q)S_q^{(T)}=w_q>\mathrm{\hspace{0.17em}0}.$$
(A12)
is always positive. Consequently, as far as the sign of the relative entropy is concerned, the denominator in (A9) can be ignored. Besides, since Tsallis’ entropy is a monotonous increasing function of Rényi’s (see Equations (A7-A8)), it is plain that (A9) has always the same sign as

$$S_q^{(R)}(A|B)=S_q^{(R)}(\rho _{AB})S_q^{(R)}(\rho _B).$$
(A13)


text:
# Octupole response and stability of spherical shape in heavy nuclei.

## I INTRODUCTION

It is well known that there is an intimate connection between the shell structure in quantum systems like nuclei and metallic clusters and the properties of classical trajectories within these systems (see e.g. , p. 579). In particular, for nuclei it has been argued that closed orbits of triangular shape might lead to an instability of the spherical shape against octupole-type deformations in the region beyond 208Pb , p.560. Here we would like to investigate in detail this possibility by using a semiclassical theory of nuclear response based on the linearized Vlasov equation . Instabilities are expected to show up as some kind of pathological behaviour in response functions (vanishing eigenfrequencies and diverging response) and the semiclassical theory of  and , that has already been shown to give good qualitative results for lower multipolarities , is an ideal tool to study the role played by classical trajectories in determining the response of large quantum systems. Of course quantum corrections are expected to modify the results of this theory , especially at very low energy, however having a clear picture of what should be expected already at the classical level might help in making progress. We study the isoscalar octupole response of a sample “nucleus” of $`A=208`$ nucleons contained in a square-well potential of radius $`R=r_0A^{\frac{1}{3}}`$, with $`r_0=1.2\mathrm{fm}`$. This is done with different degrees of approximation: in a first approximation, discussed in Sect. II, the interaction between nucleons is neglected, while in Sect. III this interaction is taken into account in a separable approximation. Finally in Sect. IV the effect of surface vibrations on the octupole response function is included.

Because of an interesting scaling property of the zero-order semiclassical response, the $`A`$-dependence can be factorized and becomes trivial. Thus, although we perform explicit calculations only for $`A=208`$, our results on stability apply also to the region beyond 208Pb.

## II Single-particle response (fixed surface)

Following Ref. , the connection between shell structure and classical trajectories is most easily illustrated by considering the quantum energy levels of a particle in a spherical potential (no spin-orbit interaction). In this case the energy levels depend on the radial and angular momentum quantum numbers $`n_r`$ and $`l`$ only. For sufficiently large values of these quantum numbers, the energy difference between two neighbouring levels can be approximated as

$$ϵ(n_r^{},l^{})ϵ(n_r,l)(n_r^{}n_r)\frac{ϵ}{n_r}+(l^{}l)\frac{ϵ}{l}.$$
(1)
Now the important point is that the derivatives $`ϵ/n_r`$ and $`ϵ/l`$ are essentially the derivatives of the classical Hamiltonian with respect to the corresponding action variable, hence they can be recognized as the radial $`(\omega _r)`$ and angular $`(\omega _\phi )`$ frequencies of the classical orbits respectively. Thus

$$ϵ(n_r^{},l^{})ϵ(n_r,l)(n_r^{}n_r)\omega _r+(l^{}l)\omega _\phi ,$$
(2)
(we use units such that $`\mathrm{}=c=1`$). This combination of frequencies is exactly that appearing in the denominator of the zero-order response function obtained in  from the solution of the linearized Vlasov equation with fixed-surface boundary conditions. Moreover, from a comparison of the zero-order Vlasov response function with the analogous quantum propagator , it can be seen that in the Vlasov propagator the quantum matrix elements of the excitation operator are replaced by appropriate Fourier coefficients that can be evaluated as integrals along the classical trajectories. Since here we are interested in the octupole response, we report the propagator describing the response of a nucleus to an external field of the type

$$Q(𝐫)=r^3Y_{3M}(\widehat{𝐫})$$
(3)
in the single-particle approximation. This propagator is given by

$$_{L=3}^0(s)=\frac{9A}{8\pi }\frac{R^6}{ϵ_F}\underset{n=\mathrm{}}{\overset{+\mathrm{}}{}}\underset{N=\pm 1,\pm 3}{}C_{3N}^2_0^1𝑑xx^2s_{nN}(x)\frac{\left(Q_{nN}^{(3)}(x)/R^3\right)^2}{s+i\epsilon s_{nN}(x)}.$$
(4)
Instead of the frequency $`\omega `$, as independent variable we have used the dimensionless quantity $`s=\omega /(v_F/R)`$, as a consequence the $`A`$-dependence of this propagator is factorized as $`AR^6A^3`$. The eigenfrequencies (2) are accordingly replaced by the functions (for a square-well potential)

$$s_{nN}(x)=\frac{n\pi +N\alpha (x)}{x}.$$
(5)
The variable $`x`$ is related to the classical angular momentum $`\lambda `$ of a nucleon. The relation is $`x=\mathrm{sin}\alpha `$, where $`\alpha `$ is the angle spanned by the radial vector when the particle moves from the inner to the outer turning point. Clearly, for a square-well potential one has $`\mathrm{cos}\alpha =\lambda /\overline{\lambda }`$, where $`\overline{\lambda }`$ is the maximum particle angular momentum $`\overline{\lambda }=p_FR`$, $`p_F`$ is the Fermi momentum, while $`v_F`$ and $`ϵ_F`$ are the corresponding velocity and energy. The sum over discrete angular momentum values of the quantum propagator is replaced here by the integration over $`x`$, thus the particle angular momentum is treated as a continuous variable.

The quantities $`C_{3N}`$ in Eq. (4) are classical limits of the Clebsh-Gordan coefficients coming from the angular integration of the quantum matrix elements . Their explicit value is $`C_{3\pm 1}=\frac{\sqrt{3}}{4}`$ and $`C_{3\pm 3}=\frac{\sqrt{5}}{4}`$. In principle the integer $`N`$ takes values between $`L`$ and $`L`$, however only the coefficients $`C_{LN}`$ where $`N`$ has the same parity as $`L`$ are nonvanishing. The coefficients $`Q_{nN}^{(3)}(x)`$ appearing in the numerator of Eq. (4) have been defined in Ref. , they are essentially the classical limit of the radial matrix elements of the octupole operator $`r^3`$ and can be evaluated explicitly as

$$Q_{nN}^{(3)}(x)=()^nR^3\frac{3}{s_{nN}^2(x)}\left(1+\frac{4}{3}N\frac{\sqrt{1x^2}}{s_{nN}(x)}\frac{2}{s_{nN}^2(x)}+4(|N|1)\frac{1x^2}{s_{nN}^2(x)}\right).$$
(6)
For $`N=\pm 1`$ this expression coincides with that appearing in the (compression) dipole response (cf. Eq.(A.5) of Ref. ), however for the octupole response we also need terms with $`N=\pm 3`$. These new modes have an interesting property since the associated eigenfrequencies $`s_{nN}(x)`$ can vanish in the interval $`0<\alpha (x)<\frac{\pi }{2}`$ \[the equation $`(n\pi +N\alpha )=0`$ has a solution for $`\alpha =\frac{\pi }{3}`$, corresponding to closed triangular orbits\]. In Ref. it has been pointed out that the vanishing of this eigenfrequency might give rise to a possible instability against octupole-type deformations in nuclei heavier than 208Pb. Although at first sight it might seem that the coefficients (6) would diverge when $`s_{nN}(x)0`$, it is actually possible to check that

$$\underset{x\frac{\sqrt{3}}{2}}{lim}Q_{1\pm 3}^{(3)}(x)=\frac{1}{4}R^3.$$
(7)
The very fact that this limit is finite is important for our discussion about the role of triangular nucleon orbits. The linear response function in the single-particle approximation is well-behaved because

$$\underset{s0}{lim}\left(\mathrm{Im}_3^0(s)\right)=0,$$
(8)
and

$$\underset{s0}{lim}\left(\mathrm{Re}_3^0(s)\right)=\frac{r_0^6}{8\pi ϵ_F}A^3,$$
(9)
a result that we shall need later. Since the zero-order octupole strength function is proportional to the imaginary part of the response function (4), we can see that the contribution of triangular trajectories to the octupole strength function does not give any pathology, at least at the single-particle level. Of course instabilities could still arise in the collective response, so we study also this response by using the same model, however, before doing this we look more in detail to the single-particle octupole strength given by the Vlasov theory.

The zero-order octupole strength function $`S_{L=3}(E)`$, given by ($`E=\mathrm{}\omega `$)

$$S_{L=3}(E)=\frac{1}{\pi }\mathrm{Im}_3^0(\mathrm{E}),$$
(10)
is shown in Fig.1. It can be seen that the single-particle octupole strength is concentrated in two regions around 8 and 24 MeV. As pointed out already in , in this respect our semiclassical response is strikingly similar to the quantum response, which is concentrated in the $`1\mathrm{}\omega `$ and $`3\mathrm{}\omega `$ regions. The modes that contribute most are those with $`(n,N)=(0,1),(1,1),(1,3),(2,3),`$ for the $`1\mathrm{}\omega `$ region and with $`(n,N)=(0,3),(1,1),(2,1),(3,3),`$ for the $`3\mathrm{}\omega `$ region. This concentration of strength is quite remarkable because our static distribution, which is taken to be of the Thomas-Fermi type $`f_0(𝐫,𝐩)\theta (ϵ_Fh_0(𝐫,𝐩))`$, does not include any shell effect, however, because of the close connection between shell structure and classical trajectories expressed by Eqs. (2) and (6), we still obtain a stregth distribution that is very similar to the one usually attributed to shell effects. Clearly the integers $`n`$ and $`N`$ correspond to the difference of radial and angular momentum quantum numbers in Eq. (2).

## III Collective response (fixed surface)

The zero-order response function (4) gives only a first approximation to the nuclear response. When the residual interaction between nucleons is taken into account, a collective response function can be obtained by solving the Vlasov equation with appropriate boundary conditions . If the interaction is assumed to be of the octupole-octupole type,

$$V(𝐫_1,𝐫_2)=\kappa _3r_1^3r_2^3\underset{M}{}Y_{3M}^{}(\widehat{𝐫}_1)Y_{3M}(\widehat{𝐫}_2)$$
(11)
the collective fixed-surface octupole response function is given by

$$_3(s)=\frac{_3^0(s)}{1\kappa _3_3^0(s)}.$$
(12)
The parameter $`\kappa _3`$ specifies the strength of the residual interaction, its value can be estimated in a self-consistent way, giving (, p.557),

$$\kappa _{BM}=\frac{4\pi }{3}\frac{m\omega _0^2}{AR^4}\mathrm{1.\hspace{0.33em}10}^5\mathrm{MeV}/\mathrm{fm}^6$$
(13)
for the octupole case. The parameter $`\omega _0`$ is given by $`\omega _041A^{\frac{1}{3}}\mathrm{MeV}`$ . Since this estimate is based on a harmonic oscillator mean field and we are assuming a square-well potential instead, we expect some differences. Hence we shall determine the parameter $`\kappa _3`$ phenomenologically, by requiring that the peak of the high-energy octupole resonance agrees with experiment. This requirement implies $`\kappa _32\kappa _{BM}`$, which is in agreement with the prescription obtained in the quadrupole case . In Fig.2 we report the collective octupole strength function given by Eq. (12), with this value of $`\kappa _3`$. We can clearly see the effects of collectivity that result in a shift and concentration of the strength into two sharp peaks around 20 Mev and 6-7 Mev. The experimentally observed  concentration of isoscalar octupole strength in the two regions usually denoted by HEOR (high energy octupole resonance) and LEOR (low energy octupole resonance) is qualitatively reproduced, however the considerable strength experimentally observed at lower energy (low-lying collective states) is absent from our fixed-surface response function. Like for the quadrupole response , we need to consider a different solution of the linearized Vlasov equation (in which the nuclear surface is allowed to vibrate ) in order to account for this feature of the response.

The poles of the collective response function are determined by the vanishing of the denominator in Eq. (12). A solution of this equation at zero or purely imaginary frequency could be interpreted as an instability of the spherical shape against an octupole-type deformation of the ground state, however we can see that this could happen only for interactions stronger than

$$\kappa _{3cr}=\frac{1}{_3^0(0)}\mathrm{3.14\hspace{0.17em}10}^5\mathrm{MeV}/\mathrm{fm}^6.$$
(14)
Our value of $`\kappa _3`$ is smaller (in absolute value) than this, so in the present model the spherical shape is stable. This is valid for any value of $`A`$, as long as we take $`\kappa _3=2\kappa _{BM}`$, since the product $`\kappa _{BM}_3^0(0)`$ is $`A`$-independent.

The triangular orbits do not enter directly into the discussion of this collective response, but the fact that they do not generate pathologies at the level of the single-particle response $`_3^0(s)`$ simplifies the discussion of the stability in the collective response also.

## IV Collective response (moving surface)

The collective fixed-surface response function (12) does reproduce two important features of the octupole response: the HEOR and the LEOR, however it misses the experimentally observed low-lying collective $`3^{}`$ states. In order to account for these low-lying states, a different solution of the linearized Vlasov equation has been proposed in Ref. . In that approach the nuclear surface is allowed to vibrate according to the usual liquid-drop model expression

$$R(\theta ,\phi ,t)=R+\underset{LM}{}\delta R_{LM}(t)Y_{LM}(\theta ,\phi ).$$
(15)
A self-consistency condition involving the surface tension is then used to determine the time-dependence of the additional collective variables $`\delta R_{LM}`$(t). Then, always for a separable residual interaction of the kind (11), the collective fixed-surface response function (12) is replaced by the following moving-surface response function

$$\stackrel{~}{}_3(s)=_3(s)+𝒮_3(s),$$
(16)
with $`_3(s)`$ still given by Eq. (12), while $`𝒮_3(s)`$ gives the moving-surface contribution. For the separable interaction (11) the function $`𝒮_3(s)`$ can be evaluated explicitly as

$$𝒮_3(s)=\frac{R^6}{1\kappa _3_3^0(s)}\frac{[\chi _3^0(s)+\kappa _3\varrho _0R^3_3^0(s)]^2}{[C_3\chi _3(s)][1\kappa _3_3^0(s)]+\kappa _3R^6[\chi _3^0(s)+\varrho _0R^3]^2},$$
(17)
with $`C_3=10\sigma R^2+(C_3)_{coul}`$ ($`\sigma 1\mathrm{M}\mathrm{e}\mathrm{V}\mathrm{fm}^2`$ is the surface tension parameter obtained from the mass formula, $`(C_3)_{coul}`$ gives the Coulomb contribution to the restoring force) and $`\varrho _0=A/\frac{4\pi }{3}R^3`$ the equilibrium density.

The functions $`\chi _3^0(s)`$ and $`\chi _3(s)`$ are defined as

$$\chi _3^0(s)=\frac{9A}{4\pi }\underset{nN}{}C_{3N}^2_0^1𝑑xx^2s_{nN}(x)\frac{()^n(Q_{nN}^{(3)}(x)/R^3)}{s+i\epsilon s_{nN}(x)},$$
(18)
and

$$\chi _3(s)=\frac{9A}{2\pi }ϵ_F(s+i\epsilon )\underset{nN}{}C_{3N}^2_0^1𝑑xx^2\frac{1}{s+i\epsilon s_{nN}(x)},$$
(19)
their structure is similar to that of the zero-order propagator (4).

In Fig.3 we report the octupole strength function given by the moving-surface response function (16) and compare it to the collective fixed-surface response given by Eq. (12). The most relevant change induced by the moving surface is the large double hump appearing at low energy. This feature is in qualitative agreement both with experiment  and with the result of RPA-type calculations , moreover is rather similar to that found in Ref.  within the same model, but with different excitation operator and residual interaction. We interpret this low-energy double hump as a superposition of surface vibrations and LEOR. In order to explain why we do not obtain one or more sharp $`3^{}`$ states at low energy, we have to analyze our moving-surface response function in some detail. Since the explicit expression (17) looks rather involved, we consider the limit of non-interacting nucleons. If we let $`\kappa _30`$, the function (17) becomes

$$𝒮_3^0(s)=R^6\frac{[\chi _3^0(s)]^2}{C_3\chi _3(s)},$$
(20)
and the full moving-surface response function (16)

$$\stackrel{~}{}_3^0(s)=_3^0(s)+𝒮_3^0(s).$$
(21)
In Fig.4 we show the octupole response function evaluated both for $`\kappa _3=2\kappa _{BM}`$ (solid curve) and for $`\kappa _3=0`$ (dashed curve). We can see from this figure that, while the residual interaction changes drastically the response in the giant resonance region, at low energy the octupole response is affected only very slightly by the interaction (11). Thus the low-energy octupole response can be analyzed by using the simpler formula (20), rather than Eq. (17).

The eigenfrequencies of low-energy collective modes are approximately determined by the vanishing of the denominator in Eq. (20), moreover at low frequency the function $`\chi _3(\omega )`$ can be expanded as

$$\chi _3(\omega )=i\omega \gamma _3+D_3\omega ^2+\mathrm{},$$
(22)
implying that the parameters $`\delta R_{3M}(t)`$, that describe octupole surface vibrations in Eq. (15), approximately satisfy an equation of motion of the damped oscillator kind:

$$D_3\delta \ddot{R}_{3M}(t)+\gamma _3\delta \dot{R}_{3M}(t)+C_3\delta R_{3M}(t)=0.$$
(23)
It has already been pointed out in  that for $`A=208`$ the numerical values of parameters in this equation are such that this oscillator is actually overdamped, moreover, always in , the coefficients $`\gamma _L`$ have been evaluated analyticallty in the low-frequency limit, giving (for a generic $`L`$)

$$\gamma _L=\gamma _{wf}2\frac{(4\pi )^2}{2L+1}\underset{N=1}{\overset{L}{}}\frac{1}{N}Y_{LN}(\frac{\pi }{2},\frac{\pi }{2})^2\underset{n=1}{\overset{+\mathrm{}}{}}\mathrm{cos}\alpha _{nN}\mathrm{sin}^3\alpha _{nN}\mathrm{\Theta }(\frac{\pi }{2}\alpha _{nN}),$$
(24)
with $`\gamma _{wf}=\frac{3}{4}\varrho _0p_FR^4`$ and $`\alpha _{nN}=\frac{n}{N}\pi `$. The angles $`\alpha _{nN}`$ are related to the nucleon trajectories as discussed in Sect. II. In our case the coefficient $`\gamma _{L=3}`$ gets a contribution only from the term with $`n=1`$ and $`N=3`$, thus we see that only nucleons moving along closed triangular trajectories can contribute to the damping of octupole surface vibrations. In order to check that indeed the closed triangular orbits are the main source of Landau damping in the low-energy octupole response, we study also the response that is obtained when the contribution of these orbits is excluded. This exclusion can be made simply by avoiding a small interval of length $`2\mathrm{\Delta }`$ about the value $`x=\frac{\sqrt{3}}{2}`$ in the $`x`$-integration in Eqs. (4), (18) and(19). In Fig.5 we show the response obtained with a value of $`\mathrm{\Delta }=0.02`$. The shape of the low-energy hunp changes dramatically because of the lack of damping due to the missing triangular trajectories. A very sharp peak is now developed at low energy. The position of this peak is influenced by the value of the mass parameter $`D_3`$, omitting the triangular orbits changes also the numerical value of this parameter. Once again we see that the closed triangular orbits, rather than generating a shape instability in the octupole channel, are the main source of Landau damping in this response function. Of course quantum effects can play an essential role at this level since in quantum mechanics the angular momentum is quantized and there could be no value of angular momentum corresponding to triangular trajectories.

We would also like to point out that our moving-surface response function can display a shape instability in the octupole channel, which can arise if the restoring force parameter $`C_3`$ vanishes. This could happen if the repulsive Coulomb term $`(C_3)_{coul}`$ exactly balances the attractive surface-tension part of the restoring force. In this case our moving-surface response function developes a pole at the origin, as shown in Fig.6. For odd multipolarities the instability condition in our model is the same as for the liquid-drop model .

Before concluding, in Fig.7 we compare our moving surface response function with the response function of the overdamped oscillator (23). This figure supports our interpretation of the lower-energy hump as due to overdamped surface vibrations.

## V Conclusions

We have analyzed the octupole response function of a hypothetical heavy nucleus containing $`A=208`$ nucleons in a square-well potential by employing a semiclassical model that relates the response to features of the classical trajectories of nucleons. Triangular trajectories are particularly relevant for the octupole response because, at the single-particle level, a simple combination of their characteristic frequencies, appearing in this response function, can vanish. This fact had already been noticed long ago and some speculations had been made on the possible connection between these classical orbits and the onset of a shape instability in heavy spherical nuclei against octupole-type deformations . Our detailed calculations performed in a semiclassical approach show that the vanishing of the eigenfrequencies associated with triangular trajectories has no consequences on the stability of the spherical shape in heavy nuclei, rather, triangular orbits are essential in providing a damping of low-energy octupole excitations.


text:
# Tailoring the photonic bandgap of porous silicon dielectric mirror

## Abstract

A systematic method to fabricate porous silicon one dimensional photonic crystals has been engineered to have a photonic bandwidth up to 2000nm. The observation of the tailorability of the photonic bandgap (PBG) underscores the requirement of the large refractive index contrast for making broad PBG structures. In this letter, we present the fabrication and characteristics of such structures that may be promising structures for a large variety of applications.

Dielectric mirrors can reflect light for any direction of propagation in specific wavelength regionsjohn \- john1 . This property will allow the creation of all-optical integrated circuits as it can be used to confine, manipulate and guide the photons. The advantage of a dielectric mirror over metallic ones is due to the dispersive and absorbing regions in the visible and infrared spectra of metals handbook . This kind of dielectric photonic band gap (PBG) structures have already been used in vertical-cavity surface-emitting lasers breiland , dielectric interference filters Hetch , sensors sensor  and other devices.

Alternating layers of material with different dielectric constant are the simplest possible photonic crystals known so far. These dielectric multilayer films act as a perfect mirror for light with a frequency within a sharply-defined gap. For a fixed number of periods, this PBG increases, as the ratio of the refractive indices of the two layers ($`n_H/n_L`$, high, H, and low, L, refractive index $`n`$) increases john ,Pavesi , and for a fixed ratio of refractive indices, upto a certain limit, the PBG increases with the number of periods Pavesi . In this letter we propose an idea to tailor the PBG structures, ranging from few nanometers to 2000 nm or more, for two different values of refractive indices, ($`n_H,n_L`$) without the prerequisite of a huge contrast. The idea has been successfully implemented experimentally on porous silicon.

Porous Silicon (pSi) is already considered as a promising material for photonic applications selena . This material can be fabricated with electrochemical etching of silicon with HF. In this material, the porosity is a linear function of the current density for a specific HF concentration and anodization time Pavesi  and is known to show a refractive index contrast for different porosites theiss . Therefore, a periodic pulse alternating between two different current densities has been a convinient procedure to fabricate multilayer pSi films Pavesi ,selena , Canham . The change of the current density does not affect the previously formed pSi layer because silicon dissolution occurs at the silicon electrolyte interface Theib .

The multilayer pSi one dimensional-PBG structures known till now had the PBG of less than 500nm Pavesi  selena . What we require is an easy method to modify the PBG of the structure according to our requirements. For the first time, through this letter, a technique to tailor the PBG structures, is successfully materialized on pSi multilayer mirrors. Stacks of $`\lambda /4`$ ($`n_Hd_H=n_Ld_L=\lambda /4`$, with $`d_i`$ being the real thickness) mirrors have been prepared for different wavelengths in the visible and near infrared region (one submirror after the other submirror). The number of submirrors fabricated in a single structure varies from 2 to 80. Each mirror is designed for a different wavelength and can consist of 2.5 to 20 periods according to the requirements and the region we want to work with. Thus, the number of pSi layers can go up to 700. The exact selection of the $`\lambda `$-values depends on the application which will be discussed in the latter part of the paper.

Our method has been inspired by the ladder structure simulated by Kavokin et al. kavokin  in order to observe photonic Bloch oscilations. The structures proposed are coupled microcavities, but here we present the structures without cavities, considering only dielectric Bragg mirrors (DBM). Of course, ladder structures can be prepared by different techniques, i.e., Molecular Beam Epitaxy, Metal Oxide chemical Vapor Deposition, etc. However in this letter, we use a cheap, easy and fast techinque to produce the same.

We have fabricated pSi multilayers by wet electrochemical etching viv  of highly boron doped substrate (p+, 0.001-0.005 $`\mathrm{\Omega }cm`$). In order to have better interfaces, an aqueous HF/ethanol/glycerol electrolyte with 15% HF, 75% ethanol and 10% glycerol servidori  concentration is used to anodize the silicon substrate. In addition, in order to maintain a constant HF concentration over the interface between Si and pSi under chemical attack, during the etching process a peristatical pump is used to circulate the electrolyte within the Teflon cell. Anodization begins when a constant current is applied between the silicon wafer and the electrolyte by means of an electronic circuit controlling the anodization process. To stabilize the pSi multilayer structure, all the samples are thermally oxidized in an oxygen ambient at 900 0C for 10 min fauchet . Oxidation induces a blue shift in the peak reflectivity due to the decrease in the refractive indices of the layers. In our structures, the multilayers have alternating porosities of the order of $`3555\%(n_H)`$ and $`6085\%(n_L)`$. The time taken in the etching process varies from few minutes to 3 hours. Hence to make the surface stable, the first layer is always a low porosity layer and made 10-20 nm thicker than usual for thick structures (more than 20 microns). The refractive indices of the pSi layers have been estimated using reflectivity spectra of 2 $`\mu `$m thick single layers at 1500nm. Scanning electron Microscopy (SEM) was used to examine the structural features of the films. The reflectivity spectra of the samples are taken with a Shimadzu UV3010 ultraviolet-visible-near infrared spectrophotometer at 5 degree incidence.

The schematic of the structures can be seen in Fig.1. The idea is to make one submirror after another submirror like a continuum array of submirrors. The specific range of $`\lambda `$-values depends directly on the application. For instance, a mirror reflecting the complete near infrared (NIR) range (800 to 2500 nm) only; for this we need to select the different wavelengths in this specific range in a way that the PBGs of respective submirrors overlap and give us a continuous PBG reflecting the full NIR range. For designing the mirrors for the latter part of the visible region and the NIR region, we have empirically selected the $`\lambda `$ values as follows. First select the starting $`\lambda _1`$-value and then the subsequent $`\lambda `$-values obey a simple relation $`\lambda _{i+1}\lambda _i=2+i`$, where $`i`$ represents the submirror number. This selection has been performed and shown in the latter part of the letter.

We have made various structures composed of single mirror with 10 periods (one period consists of one low porosity and one high porosity layer) up to a mirror composed of a set of 80 submirrors. The total number of layers varies from 20 to 650 respectively. Taking into account that in the wide PBG structures the anodization time is of the order of hours, to make the structure more stable the first mirror always starts with a high refractive index layer, consisting of 2.5 to 4.5 periods. The subsequent submirrors consists of 3 to 5 periods reducing the etching time without affecting the flatness of the photonic band gap.

Fig. 2(a) shows cross sectional SEM images of 9 out of the 54 submirrors in the complete structure designed for the partial visible and complete near infrared region. This image shows a continuum array of submirrors prepared using $`\lambda _{i+1}\lambda _i=2+i`$, designed by taking $`\lambda _1=700`$ with a refractive index contrast given by $`n_L/n_H=1.4/2.2`$. The reflectance spectrum of this structure is shown in Fig.3(b).

The absolute reflectance spectra of some of the structures have been shown in Fig. 3. The absolute reflectivity of commercially available Aluminium (Al) mirrors with a SiO overcoat (as bare silicon is prone to oxidation and loses its reflectivity) is shown in Fig. 3(a) oriel . The figure clearly shows an absorption in the visible region. From 400nm to 2500 nm the reflectivity varies between 85-90% with a dip to 75% at around 825nm. Taking this mirror as a reference (available with the Shimadzu UV3101 UV-vis-NIR spectrophotometer) the reflectivity spectra of our structures have been measured. Hence in the visible region, the reflectance shown by most of our structures was higher than that for the Al mirror. Thus, the data have been corrected for the Al mirror response. Fig.3(b)-(d) show the absolute reflectance of few examples demonstrating the successful implementation of the idea. In fig 3(b) we present the response of the mirror shown in Fig.2. It has the PBG for the partial visible and complete NIR region i.e. from 615nm to 2450nm. The absolute reflectance can be clearly seen within 85-100% which is clearly 2-20 % more than that for the Al mirror from 615-1200nm. In fig 3(c) we plot the absolute reflectance of a mirror designed for the complete visible range. The absolute reflectance of this mirror is higher than the Al mirror from 390 to 650 nm range, as can be seen, from the inset which shows the raw data coming from Shimadzu UV3101 UV-vis-NIR spectrophotometer having Al mirror as reference. This mirror is composed by 77 submirrors, with the starting wavelength 450nm and the subsequent DBMs are determined by the relation $`\lambda _{i+1}\lambda _i=3`$ having 4 periods in each submirror. The refractive index contrast used was $`1.4/2.2`$.

The reflectivity of a mirror reflecting two different wavelength ranges has been shown in Fig. 3(d). This particular mirror is composed by two submirrors, one prepared for 600 nm ($`n_L/n_H=1.3/2.3`$) and another for 1700 nm ($`n_L/n_H=1.6/2`$) with 10 periods each. The idea is to show how we can change the refractive index contrast from one submirror to another submirror if required. Theoretical simulations have been done for all the tailored PBG structures and seen to be in agreement with the experimental reflectivity spectra va

Here we can clearly see the enhancement of the PBG width with the increase in the number of submirrors. We can very well understand that this broad PBG is due to the overlapping of all the PBGs of the submirrors designed for different wavelengths. The advantage of these mirrors over metallic mirrors (Al or silver) is their selective nature and higher reflectivity. For example, designing a mirror for the complete NIR range without significant reflection from the visible range or viceversa will be excellent for different applications.

It is important to stress that Pavesi group too has developed photonic structures with wide photonic band gap using aperiodic porous silicon layers oton .

Apart from the structures shown in Fig. 3 we have tailored the PBG for many other mirrors, for example in the spectral ranges of 356-650nm, 400-600nm, 600-1200nm, 600-2000nm, etc. The idea is to design a PBG of desired width and at required wavelengths. Here we present just few examples with a lot of scope for improvement. It will also stimulate the interest of making better filters based on single microcavity microcavities , because with the method presented here we are able to increase the high reflecting range on both sides of the cavity.

For the first time, it has been made possible to tailor the PBG of the dielectric mirrors, for the partial UV and complete vis-NIR range (from 356nm to 2500nm) to have desired width at required wavelenths. These mirrors are easy to prepare, selective and more reflective than the commercially available Al mirrors.

We gratefully acknowledge the useful discussions with Dr. Miguel Robles and help provided by Jose Campos, Dr. Oswaldo Flores (CCF-UNAM).