text:
# Quantum optical versus quantum Brownian motion master equation in terms of covariance and equilibrium properties

## I INTRODUCTION

A theory of quantum dissipation, even restricted to the Markovian case, is a subject of major interest for many different scientific communities, ranging from mathematicians to physicists and chemists, according to various perspectives. Among these in first line experimental applications in phenomena where spontaneous emission, decoherence and dissipation play an important role, but also theoretical studies regarding the connection between quantum and classical description of dynamics, since, thanks to the lack of simple quantization recipes, such as the correspondence principle, dissipative systems become a fruitful working area where typical quantum structures may emerge. This interest has led to a huge number of proposals of Markovian master equations for the description of such dissipative phenomena, based on microphysical, phenomenological or purely mathematical approaches (see references in ), not always accompanied by clear statements with regard to obeyed physical and mathematical properties, thus often leading to amendments of these models in view of some missing desired feature. This is in particular true with respect to the property of positivity or complete positivity , proper distinction between Hamiltonian and dissipative part , translational invariance , and decoherence effects .

As a result some efforts have been made to compare and characterize the different proposals in view of relevant mathematical and physical properties: preservation of the positivity of the statistical operator, existence of a suitable canonical stationary state, and translational invariance . The main starting point for this research work was the result of Lindblad for the most general structure of bounded generator of a completely positive quantum dynamical semigroup , together with his paper on quantum Brownian motion , in which some of these issues were already considered. Complete positivity, actually equivalent to positivity for the considered Markovian quasi-free systems , ensures that the statistical operator preserves positivity during the time evolution, and has emerged as a typical quantum feature, corresponding to the requirement that positivity of the time evolution is preserved under entanglement. It is in fact by now an essential property in the realm of quantum communication , even though its origins lie in the theory of quantum measurement . The result of Lindblad rigorously holds for bounded operators, though it is usually exploited as a starting point also for unbounded operators, leaving in this case the task open, to show that the considered structure is a proper generator of a completely positive quantum dynamical semigroup.

In this paper we will try to further clarify the situation, showing that if besides complete positivity and equipartition, i.e., existence of a suitable canonical stationary solution, proper covariance properties of the generator of the dynamics, reflecting the relevant symmetries of the reservoir (and therefore not necessarily only translational invariance), are taken into account, a suitable characterization of two different classes of master equations can be given in a very neat way. In particular with the two one-dimensional Lie groups $`\mathrm{U}(1)`$ and $`\mathrm{R}`$ two distinct type of master equations are associated, describing respectively the damped harmonic oscillator (the so-called quantum optical master equation) and the motion of a Brownian particle (the so-called quantum Brownian motion master equation), which, despite the fact that their physical realm of validity is essentially well-understood , are sometimes mixed up , their characterization in connection to underlying symmetry of the reservoir being usually neglected . Recent work on completely positive quantum dynamical semigroups has shown that, asking for suitable covariance properties, a characterization of generators of such semigroups can be given also in the case of unbounded operators, where very few results are available , so that one can check whether a proposed formal Lindblad structure is indeed a proper generator of a quantum dynamical semigroup. The room left by these mathematical and physical requirements should be covered by microphysical approaches, determining their relevance and predictive power. It turns out that the quantum optical master equation is in the essence fixed by these requirements, while in the quantum Brownian motion case there is a non trivial freedom left, thus explaining the huge, sometimes contradictory literature devoted to the quantum Brownian motion master equation. In this connection a recently obtained  microphysical model for the quantum description of the motion of a Brownian particle is presented, derived from a quantum version of the linear Boltzmann equation , extending previous phenomenological models  where dissipation effects leading to the correct stationary solution could not be accounted for. For further work relying on symmetry properties in the case of a fermionic oscillator see .

The paper is organized as follows: in Sec. II we recall the most general Lindblad structure corresponding to a quantum Fokker-Planck equation, further showing the expressions that come out asking for shift-covariance or translation-covariance; in Sec. III we outline the results of a microphysical model for the description of quantum Brownian motion, obtained from a quantum linear Boltzmann equation expressed in terms of the operator-valued dynamic structure factor of the reservoir; in Sec. IV we briefly comment on our results, indicating possible future extensions.

## II GENERAL STRUCTURE OF QUANTUM FOKKER-PLANCK EQUATION WITH A LINDBLAD STRUCTURE

In order to give a quantum description of dissipative phenomena which at the classical level are described by second-order Fokker-Planck equations, we first concentrate on generalizations of the quantum Liouville equation which preserve trace and positivity of the statistical operator, also possibly accounting for friction effects. This can be done considering the most general expression of generator of a completely positive quantum dynamical semigroup (complete positivity being for these Markovian quasi-free systems necessary in order to preserve positivity ), given by a Lindblad structure in which position and momentum operator for the microsystem (respectively $`\widehat{𝗑}`$ and $`\widehat{𝗉}`$) appear restricted to bilinear expressions, according to the classical approximations leading to a friction term proportional to velocity. The general expression of a Markovian quantum Fokker-Planck equation preserving the positivity of the statistical operator is thus given, in the one-dimensional case to which we will restrict for simplicity, by

$$_{\widehat{𝗑}\widehat{𝗉}}[\widehat{\rho }]=\frac{i}{\mathrm{}}[H_0(\widehat{𝗑},\widehat{𝗉})+\frac{\mu }{2}\{\widehat{𝗑},\widehat{𝗉}\},\widehat{\rho }]+\underset{i=1}{\overset{2}{}}[\widehat{𝖵}_i,\widehat{\rho }\widehat{𝖵}_i^{}\frac{1}{2}\{\widehat{𝖵}_i^{}\widehat{𝖵}_i,\widehat{\rho }\}]$$
(1)
$$\widehat{𝖵}_i=\alpha _i\widehat{𝗉}+\beta _i\widehat{𝗑},\alpha _i,\beta _i\mathrm{C},\mu \mathrm{R},$$

where $`H_0`$ is a self-adjoint operator given by a quadratic expression in $`\widehat{𝗑}`$ and $`\widehat{𝗉}`$ describing the free system, and the added Hamiltonian term proportional to $`\mu `$ has been introduced for later convenience. For the sake of comparison with classical Fokker-Planck equations, in order to make the intuitive physical meaning of the different contributions clear, (1) is usually conveniently written in the following form using nested commutators and anticommutators :

$`_{\widehat{𝗑}\widehat{𝗉}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝗑},\widehat{𝗉}),\widehat{\rho }]{\displaystyle \frac{i}{\mathrm{}}}{\displaystyle \frac{(\mu \gamma )}{2}}[\{\widehat{𝗑},\widehat{𝗉}\},\widehat{\rho }]{\displaystyle \frac{i}{\mathrm{}}}\gamma [\widehat{𝗑},\{\widehat{𝗉},\widehat{\rho }\}]`$ (2)
$`{\displaystyle \frac{D_{pp}}{\mathrm{}^2}}[\widehat{𝗑},[\widehat{𝗑},\widehat{\rho }]]{\displaystyle \frac{D_{xx}}{\mathrm{}^2}}[\widehat{𝗉},[\widehat{𝗉},\widehat{\rho }]]+{\displaystyle \frac{D_{px}}{\mathrm{}^2}}[\widehat{𝗑},[\widehat{𝗉},\widehat{\rho }]]+{\displaystyle \frac{D_{xp}}{\mathrm{}^2}}[\widehat{𝗉},[\widehat{𝗑},\widehat{\rho }]],`$ (3)
where due to $`[\widehat{𝗉},[\widehat{𝗑},]]=[\widehat{𝗑},[\widehat{𝗉},]]`$ actually $`D_{xp}=D_{px}`$ and the new coefficients are related to $`\alpha _i`$ and $`\beta _i`$ through the equations

$$D_{xx}=\frac{\mathrm{}}{2}\underset{i=1}{\overset{2}{}}|\alpha _i|^2D_{pp}=\frac{\mathrm{}}{2}\underset{i=1}{\overset{2}{}}|\beta _i|^2D_{px}=\frac{\mathrm{}}{2}\mathrm{}\underset{i=1}{\overset{2}{}}\alpha _i^{}\beta _i\gamma =\frac{\mathrm{}}{2}\mathrm{}\underset{i=1}{\overset{2}{}}\alpha _i^{}\beta _i,$$

so that the following inequalities hold

$$D_{xx}0D_{pp}0D_{xx}D_{pp}D_{px}^2\frac{\gamma ^2\mathrm{}^2}{4},$$
(4)
which are necessary and sufficient conditions for an expression of the form (2) to be cast in Lindblad form, corresponding to the requirement that the matrix of coefficients

$$𝐃=\left(\begin{array}{cc}D_{xx}& D_{px}+i\frac{\mathrm{}}{2}\gamma \\ D_{px}i\frac{\mathrm{}}{2}\gamma & D_{pp}\end{array}\right)$$
(5)
has a nonnegative determinant.

An alternative but equivalent expression for (2) can be given introducing, with the aid of a length $`l`$, whose physical meaning and expression will depend on the system to be described, creation and annihilation operators $`\widehat{𝖺}`$ and $`\widehat{𝖺}^{}`$:

$$\widehat{𝖺}=\frac{1}{l\sqrt{2}}(\widehat{𝗑}+\frac{i}{\mathrm{}}l^2\widehat{𝗉})\widehat{𝖺}^{}=\frac{1}{l\sqrt{2}}(\widehat{𝗑}\frac{i}{\mathrm{}}l^2\widehat{𝗉}),$$
(6)
satisfying the commutation relation $`[\widehat{𝖺},\widehat{𝖺}^{}]=1`$. One thus obtains the expression:

$`_{\widehat{𝖺}\widehat{𝖺}^{}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝖺},\widehat{𝖺}^{}),\widehat{\rho }]{\displaystyle \frac{(\mu \gamma )}{2}}[\widehat{𝖺}^2\widehat{𝖺}^{}{}_{}{}^{2},\widehat{\rho }]`$
$`{\displaystyle \frac{\gamma }{2}}\left[[\widehat{𝖺},\{\widehat{𝖺},\widehat{\rho }\}][\widehat{𝖺},\{\widehat{𝖺}^{},\widehat{\rho }\}]\right]+\mathrm{h}.\mathrm{c}.`$
$`{\displaystyle \frac{1}{2}}\left({\displaystyle \frac{D_{xx}}{l^2}}+{\displaystyle \frac{D_{pp}l^2}{\mathrm{}^2}}\right)[\widehat{𝖺}^{},[\widehat{𝖺},\widehat{\rho }]]+{\displaystyle \frac{1}{2}}\left({\displaystyle \frac{D_{xx}}{l^2}}{\displaystyle \frac{D_{pp}l^2}{\mathrm{}^2}}2i{\displaystyle \frac{D_{px}}{\mathrm{}}}\right)[\widehat{𝖺},[\widehat{𝖺},\widehat{\rho }]]+\mathrm{h}.\mathrm{c}.,`$
or equivalently, collecting terms as in (1)

$`_{\widehat{𝖺}\widehat{𝖺}^{}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝖺},\widehat{𝖺}^{}),\widehat{\rho }]{\displaystyle \frac{\mu }{2}}[\widehat{𝖺}^2\widehat{𝖺}^{}{}_{}{}^{2},\widehat{\rho }]`$ (7)
$`+\left({\displaystyle \frac{D_{xx}}{l^2}}+{\displaystyle \frac{D_{pp}l^2}{\mathrm{}^2}}+\gamma \right)\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}^{}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^{}\widehat{𝖺}\}\right]`$ (8)
$`+\left({\displaystyle \frac{D_{xx}}{l^2}}+{\displaystyle \frac{D_{pp}l^2}{\mathrm{}^2}}\gamma \right)\left[\widehat{𝖺}^{}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}\widehat{𝖺}^{}\}\right]`$ (9)
$`\left({\displaystyle \frac{D_{xx}}{l^2}}{\displaystyle \frac{D_{pp}l^2}{\mathrm{}^2}}2i{\displaystyle \frac{D_{px}}{\mathrm{}}}\right)\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^2\}\right]+\mathrm{h}.\mathrm{c}.`$ (10)
The recalled expressions essentially give the possible Lindblad structures at most bilinear in the operators $`\widehat{𝗑}`$ and $`\widehat{𝗉}`$ or $`\widehat{𝖺}`$ and $`\widehat{𝖺}^{}`$.

### A Shift-covariance

We now analyze the behavior of the considered expressions with respect to suitable symmetry transformations. Consider a locally-compact group $`G`$ and a unitary representation $`\widehat{U}(g)`$, with $`gG`$, on the Hilbert space of the system: following  we say that a mapping $``$ in the Schrödinger picture is $`G`$-covariant if it commutes with the mapping $`𝒰_g[]=\widehat{U}(g)\widehat{U}^{}(g)`$ for all $`gG`$

$$[𝒰_g[]]=𝒰_g[[]].$$
(11)
Let us now consider the following unitary representation of the group $`\mathrm{U}(1)`$

$$\widehat{U}_\varphi =e^{i\varphi \widehat{𝖭}}$$

where $`\widehat{𝖭}=\widehat{𝖺}^{}\widehat{𝖺}`$ is the number operator and $`\varphi [0,2\pi ]`$. If we now ask for the general expression (7) invariance under the action of the group $`\mathrm{U}(1)`$, i.e., shift-covariance according to

$$_{\widehat{𝖺}\widehat{𝖺}^{}}[e^{i\varphi \widehat{𝖭}}e^{i\varphi \widehat{𝖭}}]=e^{i\varphi \widehat{𝖭}}_{\widehat{𝖺}\widehat{𝖺}^{}}[]e^{i\varphi \widehat{𝖭}},$$
(12)
then the Hamiltonian has to be a function of the generator of the transformation $`\widehat{𝖭}`$, and the following stringent requirements appear:

$$D_{xx}=D_{pp}\frac{l^4}{\mathrm{}^2}D_{px}=0\mu =0.$$

Note that the condition $`\mu =0`$ appears here as a necessary condition for shift-covariance or invariance under the relevant symmetry group and not as a natural or most simple choice as sometimes advocated . From a physical point of view the master equation is expected to reflect the relevant symmetry of the reservoir the microscopic system is interacting with. Considering for example a single mode (harmonic oscillator) interacting with the electromagnetic field, one has a $`\mathrm{U}(1)`$ symmetry and condition (12) is actually equivalent to the rotating wave approximation, essentially saying that the master equation is invariant under the transformation

$$\widehat{𝖺}\widehat{𝖺}e^{i\theta }\widehat{𝖺}^{}\widehat{𝖺}^{}e^{+i\theta }$$

or equivalently in terms of $`\widehat{𝗑}`$ and $`\widehat{𝗉}`$

$$\widehat{𝗑}\widehat{𝗑}\mathrm{cos}\theta +\widehat{𝗉}\mathrm{sin}\theta \widehat{𝗉}\widehat{𝗑}\mathrm{sin}\theta +\widehat{𝗉}\mathrm{cos}\theta .$$
(13)
The rotating wave approximation is therefore strictly linked to a $`\mathrm{U}(1)`$ symmetry and one cannot expect or try to obtain translational invariance in this case .

If now one further makes the choice $`H_0(\widehat{𝖭})=\mathrm{}\omega (\widehat{𝖭}+\frac{1}{2})`$ corresponding to a single mode (harmonic oscillator) and asks that an operator with the canonical structure $`\widehat{\rho }_0=e^{\beta H_0(\widehat{𝖭})}`$ be a stationary solution, i.e., $`_{\widehat{𝖺}\widehat{𝖺}^{}}[\widehat{\rho }_0]=0`$, a further connection between the coefficients of the master equation appears

$$2\frac{D_{pp}l^2}{\mathrm{}^2}=\gamma \mathrm{coth}\left(\frac{1}{2}\beta \mathrm{}\omega \right),$$

where $`\beta `$ is the inverse temperature characterizing the thermal electromagnetic field the mode is interacting with. The requirements of complete positivity, shift-covariance and existence of the expected canonical stationary solution then fix the quantum Fokker-Planck equation to be of the form

$`_{\widehat{𝖺}\widehat{𝖺}^{}}^{\mathrm{QO}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝖭}),\widehat{\rho }]`$ (14)
$`+\gamma \left[\mathrm{coth}\left({\displaystyle \frac{1}{2}}\beta \mathrm{}\omega \right)+1\right]\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}^{}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^{}\widehat{𝖺}\}\right]`$ (15)
$`+\gamma \left[\mathrm{coth}\left({\displaystyle \frac{1}{2}}\beta \mathrm{}\omega \right)1\right]\left[\widehat{𝖺}^{}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}\widehat{𝖺}^{}\}\right],`$ (16)
or in terms of the average of the number operator over a thermal distribution

$$N_\beta (\omega )=\frac{1}{e^{\beta \mathrm{}\omega }1}=\frac{1}{2}\left[\mathrm{coth}\left(\frac{1}{2}\beta \mathrm{}\omega \right)1\right]$$

setting $`\eta =2\gamma `$

$$_{\widehat{𝖺}\widehat{𝖺}^{}}^{\mathrm{QO}}[\widehat{\rho }]=\frac{i}{\mathrm{}}[H_0(\widehat{𝖭}),\widehat{\rho }]+\eta \left(N_\beta (\omega )+1\right)\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}^{}\frac{1}{2}\{\widehat{\rho },\widehat{𝖺}^{}\widehat{𝖺}\}\right]+\eta N_\beta (\omega )\left[\widehat{𝖺}^{}\widehat{\rho }\widehat{𝖺}\frac{1}{2}\{\widehat{\rho },\widehat{𝖺}\widehat{𝖺}^{}\}\right],$$

i.e., the well-known quantum optical master equation for the description of a damped harmonic oscillator (for a recent review see ), where the only free parameter is the decay rate $`\eta `$ and a further freedom appears in the commutator term, where a function of $`\widehat{𝖭}`$ corresponding to a frequency shift may be considered. The quantum optical master equation is therefore essentially fixed by formal requirements, well in accordance with its stability with respect to microphysical derivations, which are in fact predictive in so far as they give explicit expressions for $`\eta `$ and the energy shift. As a last step, using as natural length of the problem $`l=\sqrt{\frac{\mathrm{}}{M\omega }}`$, (14) may be written in terms of $`\widehat{𝗑}`$ and $`\widehat{𝗉}`$ as

$`_{\widehat{𝗑}\widehat{𝗉}}^{\mathrm{QO}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[{\displaystyle \frac{\widehat{𝗉}^2}{2M}}+{\displaystyle \frac{1}{2}}M\omega ^2\widehat{𝗑}^2,\widehat{\rho }]`$
$`{\displaystyle \frac{i}{\mathrm{}}}{\displaystyle \frac{\gamma }{2}}\left([\widehat{𝗑},\{\widehat{𝗉},\widehat{\rho }\}][\widehat{𝗉},\{\widehat{𝗑},\widehat{\rho }\}]\right)`$
$`{\displaystyle \frac{1}{\mathrm{}}}{\displaystyle \frac{\gamma }{2}}\mathrm{coth}\left({\displaystyle \frac{1}{2}}\beta \mathrm{}\omega \right)\left(M\omega [\widehat{𝗑},[\widehat{𝗑},\widehat{\rho }]]+{\displaystyle \frac{1}{M\omega }}[\widehat{𝗉},[\widehat{𝗉},\widehat{\rho }]]\right),`$
where invariance under (13) can be easily checked.

That $`\mathrm{U}(1)`$ symmetry or shift-covariance may lead under suitable restrictions to the quantum optical master equation can also be seen considering the recently obtained most general structure of a proper generator of a shift-covariant completely positive quantum dynamical semigroup given in , where also the unboundedness of the operators appearing in the formal Lindblad structure has been taken in due account, using the notion of form-generator. The formal operator expression associated to the form-generator is:

$`[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H(\widehat{𝖭}),\widehat{\rho }]+\left[A_0(\widehat{𝖭})\widehat{\rho }A_0^{}(\widehat{𝖭}){\displaystyle \frac{1}{2}}\{\widehat{\rho },A_0^{}(\widehat{𝖭})A_0(\widehat{𝖭})\}\right]`$ (17)
$`+{\displaystyle \underset{m=1}{\overset{\mathrm{}}{}}}\left[\widehat{𝖶}^{}{}_{}{}^{m}A_{m}^{}(\widehat{𝖭})\widehat{\rho }A_m^{}(\widehat{𝖭})\widehat{𝖶}^m{\displaystyle \frac{1}{2}}\{\widehat{\rho },A_m^{}(\widehat{𝖭})\widehat{𝖶}^m\widehat{𝖶}^{}{}_{}{}^{m}A_{m}^{}(\widehat{𝖭})\}\right]`$ (18)
$`+{\displaystyle \underset{m=1}{\overset{\mathrm{}}{}}}\left[\widehat{𝖶}^mA_m(\widehat{𝖭})\widehat{\rho }A_m^{}(\widehat{𝖭})\widehat{𝖶}^{}{}_{}{}^{m}{\displaystyle \frac{1}{2}}\{\widehat{\rho },A_m^{}(\widehat{𝖭})A_m(\widehat{𝖭})\}\right],`$ (19)
where $`H()=H^{}()`$, $`A_i()`$ are functions of the number operator, and $`\widehat{𝖶}=_{n=0}^{\mathrm{}}|n+1n|`$ it the so-called shift-operator . Recalling the polar decompositions $`\widehat{𝖺}^{}=\widehat{𝖶}\sqrt{\widehat{𝖭}+1}`$ and $`\widehat{𝖺}=\widehat{𝖶}^{}\sqrt{\widehat{𝖭}}`$, for the following simple choice of functions:

$`H(n)`$ $`=`$ $`H_0(n)=\mathrm{}\omega \left(n+{\displaystyle \frac{1}{2}}\right)`$
$`A_m(n)`$ $`=`$ $`0m=0,|m|>1`$
$`A_1(n)`$ $`=`$ $`\sqrt{\gamma (\beta )}\sqrt{n+1}`$
$`A_1(n)`$ $`=`$ $`e^{\frac{\beta }{2}}\sqrt{\gamma (\beta )}\sqrt{n},`$
where

$$\gamma (\beta )=\gamma \left[\mathrm{coth}\left(\frac{1}{2}\beta \mathrm{}\omega \right)1\right],$$

one recovers from (17) the quantum optical master equation.

It is interesting to observe that complete positivity, shift-covariance and the requirement of a canonical stationary solution also allow as a proper generator of a quantum dynamical semigroup the following expression

$`[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝖭}),\widehat{\rho }]\gamma _0[\widehat{𝖭},[\widehat{𝖭},\widehat{\rho }]]`$ (20)
$`+{\displaystyle \underset{m=1}{\overset{\mathrm{}}{}}}\gamma _m\{[\mathrm{coth}\left({\displaystyle \frac{1}{2}}\beta \mathrm{}\omega \right)+1]^m\left[\widehat{𝖺}^m\widehat{\rho }\widehat{𝖺}^{}{}_{}{}^{m}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^{}{}_{}{}^{m}\widehat{𝖺}_{}^{m}\}\right]`$ (21)
$`+[\mathrm{coth}\left({\displaystyle \frac{1}{2}}\beta \mathrm{}\omega \right)1]^m[\widehat{𝖺}^{}{}_{}{}^{m}\widehat{\rho }\widehat{𝖺}^m{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^m\widehat{𝖺}^{}{}_{}{}^{m}\}]\},`$ (22)
corresponding to the functions:

$`A_0(n)`$ $`=`$ $`\sqrt{\gamma _0}n`$
$`A_m(n)`$ $`=`$ $`\sqrt{\gamma _m(\beta )}\sqrt{{\displaystyle \frac{(n+m)!}{n!}}}`$
$`A_m(n)`$ $`=`$ $`e^{m\frac{\beta }{2}}\sqrt{\gamma _m(\beta )}\sqrt{{\displaystyle \frac{n!}{(nm)!}}},`$
where

$$\gamma _m(\beta )=\gamma _m\left[\mathrm{coth}\left(\frac{1}{2}\beta \mathrm{}\omega \right)1\right]^m.$$

Eq. (20) provides a generalization of the quantum optical master equation in which both phase-diffusion, related to the coefficient $`\gamma _0`$, and $`m`$-photon processes with decay rates $`\eta _m=2^m\gamma _m`$, which should quickly approach zero in order to allow for non-explosion of the associated Markov process, can be considered. Thus far we have dealt with the case of shift-covariance, corresponding to $`\mathrm{U}(1)`$ symmetry of the system with many degrees of freedom acting as reservoir and determining the non-Hamiltonian dynamics of the microsystem.

### B Translation-covariance

We now consider the case in which the relevant symmetry is invariance under translations, corresponding to a homogeneous reservoir. Given the unitary representation of the translation group

$$\widehat{U}(b)=e^{\frac{i}{\mathrm{}}b\widehat{𝗉}},$$

where $`\widehat{𝗉}`$ is the momentum operator of the microsystem and $`b\mathrm{R}`$, translation-covariance according to  amounts to the requirement

$$_{\widehat{𝗑}\widehat{𝗉}}[e^{\frac{i}{\mathrm{}}b\widehat{𝗉}}e^{+\frac{i}{\mathrm{}}b\widehat{𝗉}}]=e^{\frac{i}{\mathrm{}}b\widehat{𝗉}}_{\widehat{𝗑}\widehat{𝗉}}[]e^{+\frac{i}{\mathrm{}}b\widehat{𝗉}}.$$
(23)
Invariance under the group $`\mathrm{R}`$ of translations thus implies for the structure of the quantum Fokker-Planck equation that the Hamiltonian has to be a function of the generator of the transformation $`\widehat{𝗉}`$ and the following simple requirement in the coefficients appearing in (2)

$$\mu =\gamma .$$

Considering for example a free particle interacting with a homogeneous reservoir, one has $`\mathrm{R}`$ symmetry corresponding to invariance under translations, which is reflected by the fact that according to (23) the master equation is invariant under the transformation

$$\widehat{𝗑}\widehat{𝗑}+b\widehat{𝗉}\widehat{𝗉}$$

or equivalently in terms of $`\widehat{𝖺}`$ and $`\widehat{𝖺}^{}`$

$$\widehat{𝖺}\widehat{𝖺}+\frac{1}{\sqrt{2}}\frac{b}{l}\widehat{𝖺}^{}\widehat{𝖺}^{}+\frac{1}{\sqrt{2}}\frac{b}{l}.$$

Further making the obvious choice $`H_0(\widehat{𝗉})=\frac{\widehat{𝗉}^2}{2M}`$, corresponding to a free particle of mass $`M`$ and asking that an operator with the canonical structure $`\widehat{\rho }_0=e^{\beta H_0(\widehat{𝗉})}`$ be a stationary solution, i.e., $`_{\widehat{𝗑}\widehat{𝗉}}[\widehat{\rho }_0]=0`$, one has the condition

$$D_{pp}=\gamma \frac{2M}{\beta },$$

with $`\beta `$ the inverse temperature of the homogeneous reservoir. The requirement of complete positivity, translation-covariance and existence of the expected stationary solution thus constrain the quantum Fokker-Planck equation to be of the form

$`_{\widehat{𝗑}\widehat{𝗉}}^{\mathrm{QBM}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝗉}),\widehat{\rho }]{\displaystyle \frac{i}{\mathrm{}}}\gamma [\widehat{𝗑},\{\widehat{𝗉},\widehat{\rho }\}]\gamma {\displaystyle \frac{2M}{\beta \mathrm{}^2}}[\widehat{𝗑},[\widehat{𝗑},\widehat{\rho }]]`$ (24)
$`{\displaystyle \frac{D_{xx}}{\mathrm{}^2}}[\widehat{𝗉},[\widehat{𝗉},\widehat{\rho }]]+2{\displaystyle \frac{D_{px}}{\mathrm{}^2}}[\widehat{𝗑},[\widehat{𝗉},\widehat{\rho }]],`$ (25)
i.e., the typical structure one finds in the extensive physical literature aiming at the description of quantum Brownian motion (see  for a review). With respect to the case of the quantum optical master equation (14) the remaining freedom in the structure is much bigger: apart from the friction coefficient $`\gamma `$ and a real function of $`\widehat{𝗉}`$ correcting the Hamiltonian, the coefficients $`D_{xx}`$ and $`D_{px}`$ are undetermined except for the relation

$$\gamma \frac{2M}{\beta }D_{xx}D_{px}^2\frac{\gamma ^2\mathrm{}^2}{4}$$

stemming from (4). This is reflected by the fact that a much wider literature has been devoted to the subject, looking for microphysical derivations of the quantum Brownian motion master equation, in order to obtain expressions for the undetermined parameters.

We now compare the result (24) with the general structure of a proper generator of a translation-covariant quantum dynamical semigroup, given in , where also the case of an unbounded generator has been considered, introducing the notion of form-generator and specifying a suitable domain for the mapping. Restricting to the continuous component of the generator, corresponding to a quantum Fokker-Planck equation describing friction and diffusion, one has for the formal operator expression associated to the form-generator the following result

$$[\widehat{\rho }]=\frac{i}{\mathrm{}}[\beta \widehat{𝗑}+H(\widehat{𝗉}),\widehat{\rho }]+\alpha \widehat{𝖵}\widehat{\rho }\widehat{𝖵}^{}\widehat{𝖪}\widehat{\rho }\widehat{\rho }\widehat{𝖪}^{}\beta \mathrm{R},\alpha 0$$
(26)
$$\widehat{𝖵}=\widehat{𝗑}+L(\widehat{𝗉})\widehat{𝖪}=\frac{\alpha }{2}\left[\widehat{𝗑}^2+2\widehat{𝗑}L(\widehat{𝗉})+L^{}(\widehat{𝗉})L(\widehat{𝗉})\right],$$

$`H()=H^{}()`$, $`L()`$ being functions of the momentum operator $`\widehat{𝗉}`$, and $`\beta 0`$ implying e.g., a constant gravitational or electric field. According to (26) one has a single operator of the form $`\widehat{𝖵}=\widehat{𝗑}+L(\widehat{𝗉})`$ (or one for each Cartesian coordinate considering higher dimensions) instead of two as considered in (1). Expressing (26) in terms of nested commutators and anticommutators as in (2), restricting to the case in which $`L()`$ is a linear function, according to the fact that we are considering friction effects at most linear in the velocity, the inequalities in (4), corresponding to the fact that the determinant of the matrix given in (5) be zero or positive, now become more restrictive:

$$D_{xx}0D_{pp}0D_{xx}D_{pp}D_{px}^2=\frac{\gamma ^2\mathrm{}^2}{4},$$

corresponding to $`det𝐃=0`$. Coming back to (24) this implies the further restriction

$$D_{xx}=\gamma \frac{\beta \mathrm{}^2}{8M}+\frac{\beta }{2\gamma M}D_{px}^2,$$

so that apart from the overall multiplying coefficient $`\gamma `$ only another coefficient $`D_{px}`$ is left free, and one has

$`_{\widehat{𝗑}\widehat{𝗉}}^{\mathrm{QBM}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝗉}),\widehat{\rho }]{\displaystyle \frac{i}{\mathrm{}}}\gamma [\widehat{𝗑},\{\widehat{𝗉},\widehat{\rho }\}]\gamma {\displaystyle \frac{2M}{\beta \mathrm{}^2}}[\widehat{𝗑},[\widehat{𝗑},\widehat{\rho }]]\gamma {\displaystyle \frac{\beta }{8M}}[\widehat{𝗉},[\widehat{𝗉},\widehat{\rho }]]`$ (27)
$`{\displaystyle \frac{\beta }{2\gamma M}}{\displaystyle \frac{D_{px}^2}{\mathrm{}^2}}[\widehat{𝗉},[\widehat{𝗉},\widehat{\rho }]]+2{\displaystyle \frac{D_{px}}{\mathrm{}^2}}[\widehat{𝗉},[\widehat{𝗑},\widehat{\rho }]].`$ (28)
Therefore a predictive microphysical model of quantum Brownian motion essentially has to indicate an explicit expression for the coefficients $`\gamma `$ and $`D_{px}`$.

It is interesting to express (27) in terms of the creation and annihilation operators given in (6). Setting $`\frac{\beta \mathrm{}^2}{4M}=\lambda _{\mathrm{th}}^2`$, the square of the thermal wavelength of the microsystem undergoing Brownian motion, one has:

$`_{\widehat{𝖺}\widehat{𝖺}^{}}^{\mathrm{QBM}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝖺},\widehat{𝖺}^{}),\widehat{\rho }]{\displaystyle \frac{\gamma }{2}}[\widehat{𝖺}^2\widehat{𝖺}^{}{}_{}{}^{2},\widehat{\rho }]`$ (29)
$`+\left[{\displaystyle \frac{\gamma }{2}}\left({\displaystyle \frac{\lambda _{\mathrm{th}}^2}{l^2}}+{\displaystyle \frac{l^2}{\lambda _{\mathrm{th}}^2}}+2\right)+{\displaystyle \frac{1}{\gamma }}{\displaystyle \frac{2}{\mathrm{}^2}}{\displaystyle \frac{\lambda _{\mathrm{th}}^2}{l^2}}D_{px}^2\right]\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}^{}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^{}\widehat{𝖺}\}\right]`$ (30)
$`+\left[{\displaystyle \frac{\gamma }{2}}\left({\displaystyle \frac{\lambda _{\mathrm{th}}^2}{l^2}}+{\displaystyle \frac{l^2}{\lambda _{\mathrm{th}}^2}}2\right)+{\displaystyle \frac{1}{\gamma }}{\displaystyle \frac{2}{\mathrm{}^2}}{\displaystyle \frac{\lambda _{\mathrm{th}}^2}{l^2}}D_{px}^2\right]\left[\widehat{𝖺}^{}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}\widehat{𝖺}^{}\}\right]`$ (31)
$`\left[{\displaystyle \frac{\gamma }{2}}\left({\displaystyle \frac{\lambda _{\mathrm{th}}^2}{l^2}}{\displaystyle \frac{l^2}{\lambda _{\mathrm{th}}^2}}\right)+2{\displaystyle \frac{D_{px}}{\mathrm{}}}\left({\displaystyle \frac{1}{\gamma }}{\displaystyle \frac{\lambda _{\mathrm{th}}^2}{l^2}}{\displaystyle \frac{D_{px}}{\mathrm{}}}i\right)\right]\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^2\}\right]+\mathrm{h}.\mathrm{c}.`$ (32)
Eq. (29) strongly simplifies if one takes for the length $`l`$, used in order to introduce the operators $`\widehat{𝖺}`$ and $`\widehat{𝖺}^{}`$ in terms of the operator position and momentum of the particle, the value $`l=\lambda _{\mathrm{th}}=\sqrt{\frac{\beta \mathrm{}^2}{4M}}`$, naturally suggested by the underlying physics. With $`l`$ the thermal de Broglie wavelength (29) reduces to

$`_{\widehat{𝖺}\widehat{𝖺}^{}}^{\mathrm{QBM}}[\widehat{\rho }]`$ $`=`$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝖺},\widehat{𝖺}^{}),\widehat{\rho }]{\displaystyle \frac{\gamma }{2}}[\widehat{𝖺}^2\widehat{𝖺}^{}{}_{}{}^{2},\widehat{\rho }]+2\gamma [\widehat{𝖺}\widehat{\rho }\widehat{𝖺}^{}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^{}\widehat{𝖺}\}]`$ (33)
$`+{\displaystyle \frac{2}{\gamma }}{\displaystyle \frac{D_{px}^2}{\mathrm{}^2}}\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}^{}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^{}\widehat{𝖺}\}+\widehat{𝖺}^{}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}\widehat{𝖺}^{}\}\right]`$ (34)
$`2{\displaystyle \frac{D_{px}}{\mathrm{}}}\left({\displaystyle \frac{1}{\gamma }}{\displaystyle \frac{D_{px}}{\mathrm{}}}i\right)\left[\widehat{𝖺}\widehat{\rho }\widehat{𝖺}{\displaystyle \frac{1}{2}}\{\widehat{\rho },\widehat{𝖺}^2\}\right]+\mathrm{h}.\mathrm{c}.,`$ (35)
where the last three contributions can only vanish if the real coefficient $`D_{px}`$ is equal to zero, corresponding to $`L()=L^{}()`$ in (26). Eq. (27) or equivalently (33) express the general structure of a quantum Fokker-Planck equation which is invariant under translation, warrants the existence of the expected canonical expression as a stationary solution, thus recovering equipartition, and is furthermore a proper generator of a completely positive quantum dynamical semigroup.

### C Covariance and uniqueness of the stationary solution

In recent work  aiming at comparing and clarifying different approaches to quantum dissipation, which relies on  but neglects the more recent and thorough results of , the statement can be found that no Markovian theory can satisfy all three criteria of positivity, translational invariance, and asymptotic approach to the canonical equilibrium state $`e^{\beta H_0}`$, except in special cases. This statement is always correct in view of the last observation, and these simple exceptional cases, which are usually not spelled out, can be read in (17): the microsystem has to be a free particle apart from an effective correction to the Hamiltonian, given by a real function of $`\widehat{𝗉}`$ describing for example an effective mass, and a potential term depending linearly on position (such as, e.g., a constant gravitational or electric field). These cases are often neglected, having in mind that translational invariance, mainly seen as an abstract property rather than the expression of homogeneity of the reservoir, should be always asked for. In this spirit translational invariance, i.e., $`\mathrm{R}`$ symmetry, is asked for also for the damped harmonic oscillator.

The three physical requirements one can reasonably ask together for Markovian systems in the weak-coupling limit are complete positivity, existence of the stationary solution predicted by equipartition and invariance under the relevant symmetry, not necessarily translational invariance. Having translational invariance apart from the potential term is not physically significant since the potential term actually breaks this invariance, and would furthermore lead to high non-uniqueness of the stationary state, as argued below. The dynamics of the microsystem is driven by both the potential term (which could also arise as a mean field effect) and the contributions describing decoherence and dissipation, so that a physically relevant symmetry should pertain to the system as a whole. Let us in fact consider a mapping $``$ covariant with respect to a given symmetry group $`G`$ according to (11), which admits a stationary solution $`\widehat{\rho }_0`$, i.e., $`[\widehat{\rho }_0]=0`$. If the operator $`\widehat{\rho }_0`$ is not invariant under the unitary representation $`𝒰_g`$ of $`G`$, so that

$$\widehat{\rho }_g=𝒰_g[\widehat{\rho }_0]$$

is linearly independent from $`\widehat{\rho }_0`$ at least for some $`g`$ in $`G`$, due to the $`G`$-covariance of $``$, $`\widehat{\rho }_g`$ still is a stationary solution

$$[\widehat{\rho }_g]=[𝒰_g[\widehat{\rho }_0]]=𝒰_g[[\widehat{\rho }_0]]=0,$$

so that one cannot have the expected uniqueness of the solution. In the case of the one-dimensional Lie groups considered in Sec. II for example a stationary solution can be unique only if it commutes with the generator of the group. Note that this simple argument is independent on whether the mapping $``$ ensures complete positivity of the time evolution or not, so that the clash between the requirement of translation invariance and the existence of the correct stationary state corresponding to equipartition is not due to the requirement of complete positivity of the time evolution mapping.

## III QUANTUM FOKKER-PLANCK EQUATION FOR THE MOTION OF A BROWNIAN PARTICLE

As we have shown in Sec. II A, in the case of an underlying $`\mathrm{U}(1)`$ symmetry formal requirements are enough to essentially fix structure and coefficients of the quantum Fokker-Planck equation, apart from an energy shift and an overall coefficient. The same is not true in the case of $`\mathrm{R}`$ symmetry describing translational invariance. In this paragraph we therefore quickly recall a recently proposed quantum Fokker-Planck equation for the description of the motion of a heavy Brownian particle interacting through collisions with a homogeneous fluid made up of much lighter particles (see  for details). This result has been obtained within a kinetic approach, where the dynamics is driven by single events described as collisions in which one generally has momentum and energy transfer, alternative to the Zwanzig Caldeira Leggett approach  (recently criticized in ), where one describes the reservoir as a collection of harmonic oscillators coupled to the microsystem through its position, usually performing calculations in terms of path-integral techniques.

In the aforementioned approach one obtains in the first instance a kinetic equation for the statistical operator analogous to the classical linear Boltzmann equation given by the simple expression

$`[\widehat{\rho }]=`$ $``$ $`{\displaystyle \frac{i}{\mathrm{}}}[H_0(\widehat{𝗉}),\widehat{\rho }]`$ (36)
$`+`$ $`{\displaystyle \frac{2\pi }{\mathrm{}}}(2\pi \mathrm{})^3n{\displaystyle _{𝐑^3}}d^3𝐪|\stackrel{~}{t}(q)|^2\left[e^{\frac{i}{\mathrm{}}𝐪\widehat{𝗑}}\sqrt{S(𝐪,\widehat{𝗉})}\widehat{\rho }\sqrt{S(𝐪,\widehat{𝗉})}e^{\frac{i}{\mathrm{}}𝐪\widehat{𝗑}}{\displaystyle \frac{1}{2}}\{S(𝐪,\widehat{𝗉}),\widehat{\rho }\}\right],`$ (37)
where $`\stackrel{~}{t}(q)`$ is the Fourier transform of the T-matrix describing the microphysical collisions and the function $`S(𝐪,𝐩)`$ appearing operator-valued in (36) is a positive two-point correlation function known in the physical community as dynamic structure factor , usually expressed as a function of momentum and energy transfer, $`𝐪`$ and $`E`$, according to

$$S(𝐪,𝐩)S(𝐪,E)E(𝐪,𝐩)=\frac{(𝐩+𝐪)^2}{2M}\frac{𝐩^2}{2M},$$

with $`M`$ the mass of the Brownian particle. The dynamic structure factor $`S(𝐪,E)`$ is the Fourier transform of the two-point time dependent density autocorrelation function of the fluid

$$S(𝐪,E)=\frac{1}{2\pi \mathrm{}}_𝐑𝑑t_{𝐑^3}d^3𝐱e^{\frac{i}{\mathrm{}}[E(𝐪,𝐩)t𝐪𝐱]}\frac{1}{N}_{𝐑^3}d^3𝐲N(𝐲)N(𝐲+𝐱,t),$$
(38)
and is always positive since it is proportional to the energy dependent scattering cross-section of a microscopic probe off a macroscopic sample , giving the spectrum of its spontaneous fluctuations. Equation (36) also has the three properties of complete positivity, translational invariance and canonical stationary solution, and it actually gives a physical example of the general structure of generator of a translation-covariant quantum dynamical semigroup , going beyond the diffusive case considered in (17). In order to recover from the general integral kinetic equation (36) a quantum Fokker-Planck equation for the description of the Brownian motion of a heavy particle one has to consider the limit of small momentum transfer and small energy transfer (corresponding to the Brownian limit in which the mass of the test particle is much heavier than the particles of the fluid). Considering a gas of free particles obeying Maxwell-Boltzmann statistics, due to

$$S_{\mathrm{MB}}(𝒒,𝒑)=\frac{1}{(2\pi \mathrm{})^3}\frac{2\pi m^2}{n\beta q}z\mathrm{exp}\left[\frac{\beta }{8m}\frac{(2mE(𝒒,𝒑)+q^2)^2}{q^2}\right]$$
(39)
one has

$$[\widehat{\rho }]=\frac{i}{\mathrm{}}[H_0(\widehat{𝗉}),\widehat{\rho }]\frac{i}{\mathrm{}}\gamma \underset{i=1}{\overset{3}{}}[\widehat{𝗑}_i,\{\widehat{𝗉}_i,\widehat{\rho }\}]\frac{D_{pp}}{\mathrm{}^2}\underset{i=1}{\overset{3}{}}[\widehat{𝗑}_i,[\widehat{𝗑}_i,\widehat{\rho }]]\frac{D_{xx}}{\mathrm{}^2}\underset{i=1}{\overset{3}{}}[\widehat{𝗉}_i,[\widehat{𝗉}_i,\widehat{\rho }]].$$
(40)
In this kinetic case the free parameters in (27) are now determined as:

$$D_{px}=0\gamma =\frac{1}{3}z\frac{\pi ^2m^2}{\beta \mathrm{}}_{𝐑^3}d^3𝐪|\stackrel{~}{t}(q)|^2qe^{\frac{\beta }{8m}q^2}$$

with $`z=e^{\beta \mu }`$ the fugacity of the gas , while $`D_{xx}`$ and $`D_{pp}`$ are expressed in terms of $`\gamma `$ as can be read in (27)

$$D_{xx}=\frac{\beta \mathrm{}^2}{8M}\gamma D_{pp}=\frac{2M}{\beta }\gamma .$$

Also the expression of (40) in terms of the operators $`\widehat{𝖺}_i`$ and $`\widehat{𝖺}_i^{}`$ according to (33) takes in this case a particularly simple form

$$[\widehat{\rho }]=\frac{i}{\mathrm{}}[H_0(\widehat{𝖺}_i,\widehat{𝖺}_i^{}),\widehat{\rho }]\frac{\gamma }{2}\underset{i=1}{\overset{3}{}}[\widehat{𝖺}_i^2\widehat{𝖺}_i^{}{}_{}{}^{2},\widehat{\rho }]+2\gamma \underset{i=1}{\overset{3}{}}[\widehat{𝖺}_i\widehat{\rho }\widehat{𝖺}_{}^{}{}_{i}{}^{}\frac{1}{2}\{\widehat{\rho },\widehat{𝖺}_{}^{}{}_{i}{}^{}\widehat{𝖺}_i\}],$$
(41)
so that one has a single $`\widehat{𝖺}_i=\sqrt{\frac{2M}{\beta \mathrm{}^2}}(\widehat{𝗑}_i+i\frac{\beta \mathrm{}^2}{4M}\widehat{𝗉}_i)`$ operator for each Cartesian coordinate.

## IV CONCLUSIONS AND OUTLOOK

The main scope of this paper was to show how relevant symmetries can be in the determination of structures of quantum Fokker-Planck equation. While in the literature only translational invariance is considered and asked for also in the case of the damped harmonic oscillator, leading to no go statements regarding the possibility of having quantum Fokker-Planck equations with all the physically relevant features (complete positivity, covariance and equipartition), we here considered both one-dimensional Lie groups $`\mathrm{U}(1)`$ and $`\mathrm{R}`$, showing that they are connected to two distinct classes of quantum Fokker-Planck equations: the quantum optical master equation associated to shift-covariance and the quantum Brownian motion master equation associated to translation-covariance. That these two classes of models actually correspond to different physics can also be seen in connection with recent studies on their properties with respect to decoherence . In particular, independently of complete positivity of the time evolution, it has been shown that covariance properties put severe restrictions on the structure of the stationary solution, provided uniqueness is asked for.

In the case of shift-covariance a generalization of the quantum optical master equation has been proposed, in which also $`m`$-photon processes can be considered. Moreover recent mathematical results by Holevo have been considered, concerning the structure of generators of completely positive quantum dynamical semigroups and leading to further restrictions in the quantum Brownian motion case. For the Brownian motion of a test particle in a fluid, where formal requirements are not enough to essentially fix the structure of the quantum Fokker-Planck equation describing the phenomenon, a recent microphysical approach has been briefly recalled, based on scattering theory, where the quantum Fokker-Planck equation is obtained as the small momentum and energy transfer limit of a quantum generalization of the classical linear Boltzmann equation. Covariance properties with respect to some physically relevant group, typically reflecting a symmetry under certain transformations of the given reservoir, can therefore be a most useful requirement in the determination of structures of quantum Fokker-Planck equations or more generally of linear kinetic equations.

An interesting extension of this work would entail the study of generators of the dynamics of systems in which one has an important correlation between internal and translational degrees of freedom, both coupled through different interactions to some reservoir, a problem recently considered in , where the translational degrees of freedom are treated in a classical way, assuming decoherence is strong enough. Such models could be of interest for the implementation of quantum computing, where indeed some experimental scheme actually relies on a coupling between internal and center of mass degrees of freedom .

## ACKNOWLEDGMENTS

The author would like to thank Prof. L. Lanz for many useful suggestions during the whole work and Prof. A. Barchielli for useful discussions. He also thanks Dr. F. Belgiorno for careful reading of the manuscript. This work was partially supported by MIUR under Cofinanziamento and Progetto Giovani.


format: text

text:
# Universality of the Crossing Probability for the Potts Model for 𝑞={1,2,3,4}

## I Introduction

The concept of the universality and scaling relations St  is general concept of the modern phase transition theory. The main point of the scaling theory is that in the vicinity of the critical point for a system of linear size $`L`$, the critical behavior of thermodynamical quantities can be expressed as a universal function of two variables: reduced temperature $`\tau =\frac{TT_c}{T_c}`$ and external field $`h`$. The finite-size scaling of thermodynamical functions of spin models was studied theoretically and numerically FB1 ; FB2 ; Landau1 ; Landau2 . In Ref. PF , Privman and Fisher proposed the idea of the universal finite-size scaling with nonuniversal metric factors. For example, for the free-energy density of system size $`L`$ and dimension $`d`$

$$f(T,h;L)=L^dY(C_1\tau L^{\frac{1}{\nu }},C_2hL^{(\beta +\gamma )/\nu })$$
(1)
where $`\nu `$ is the scaling index for correlation length, $`\beta `$ is the scaling index of magnetization and $`\gamma `$ is the scaling index of magnetic susceptibility. For systems of different boundary conditions, aspect ratios and geometries (square, honeycomb, triangle), the scaling function $`Y(x,y)`$ is universal and only metric factors $`C_1,C_2`$ are system dependent. Scaling properties of thermodynamical functions of the Potts model were investigated in Refs. KL ; dQ .

Langlands, Pichet, Pouliot and Saint-Aubin Lang1  show that for site and bond percolation on square, honeycomb and triangle lattices with the aspect ratios $`a`$, $`a\sqrt{3}`$ and $`a\sqrt{3}/2`$ respectively, the crossing probability $`\pi _h`$ of a system to percolate in the horizontal direction is the universal function of $`a`$. Hu, Lin show that by choosing a very small number of nonuniversal metric factors, all scaled data for percolation functions and the number of percolating clusters on square, honeycomb and triangle lattices may fall on the same universal scaling functions Hu1 ; Hu2 . Their scaling argument was $`x=(pp_c)L^{\frac{1}{\nu }}`$ where $`p_c`$ is the critical point, $`L`$ is the lattice size, $`\nu `$ is the correlation length index. The scaling of crossing probabilities for the three-dimensional percolation was investigated in Ref. Hu3 .

The continuum limit of the crossing probability $`\pi _h(p_c)`$ was investigated by J. Cardy by conformal field methods Cardy0 ; Cardy1 ; Cardy2 . The analogous formula for the crossing probability $`\pi _{hv}=\pi _h\pi _{hs}`$ was found by Watts Watts . The works of Smirnov Sm1 ; Sm2  analytically proved that the crossing formula holds for the continuum limit of site percolation on the triangle lattice Sm1 ; Sm2 .

The $`q`$-state Potts model can be represented as the correlated site-bond percolation in terms of Fortuin-Kasteleyn (FK) clusters FK . At the critical point of the second order phase transition Potts model, FK-clusters exhibit the percolation transition. So there is an intrinsic relationship between critical properties of the Potts model and percolation properties of FK clusters. The universality of the crossing probability for the Ising model on rectangle lattices of square, honeycomb and triangle geometries was investigated by Langlands et all. Lang2 ; Lang3 . Hu, Chen and Lin show the universality of the crossing probability $`\pi _h`$ and a number of percolation clusters for the correlated site-bond percolation $`q=2`$ Hu4 .

In this paper the the crossing probabilities of FK clusters is studied numerically. We investigate the universality of the crossing probability with respect to a number of spin states $`q`$ of the Potts model. We show numerically that the probability of a system to percolate only in horizontal direction $`\pi _{hs}`$ is an universal function of the scaling variable $`z=\left[b(q)L^{\frac{1}{\nu }}(pp_c)\right]^{\zeta (q)}`$ for $`q=1,2,3,4`$ where $`p=1\mathrm{exp}(\beta )`$ is the probability of bond to be closed, $`\beta =1/Tk_B`$ is an inverse temperature, $`p_c=1\mathrm{exp}(\beta _c)`$ is the location of critical point in the $`p`$ scale, $`b(q)`$ is a nonuniversal metric factor and $`\zeta (q)`$ is a scaling index, described dependence of the form of the crossing probability on $`q`$.

We show that for each value $`q=1,2,3,4`$ on the square, lattice the index $`\zeta `$ is practically does not depend on the lattice size. Therefore, by introducing this index $`\zeta (q)`$, we can lay on the same curve the points for different $`q`$ in the critical region.

## II  Crossing probabilities for the Potts model

In the Potts model, each spin $`\sigma _i`$ can take one of the $`q`$ different values $`1,\mathrm{},q`$ and the Hamiltonian is written as $`=J\underset{(i,j)}{}\delta (\sigma _i,\sigma _j)`$, where $`J`$ is the ferromagnetic coupling constant, which we set it equal 1. The partition function of the Potts model Baxter  is

$$Z=\underset{\sigma }{}\mathrm{exp}(\beta (\omega ))=\underset{\sigma }{}\underset{(i,j)}{}\left[(1p)+p\delta (\sigma _i,\sigma _j)\right]$$
(2)
where $`\beta =1/Tk_B`$ is the inverse temperature, $`p=1\mathrm{exp}(\beta )`$ is the probability of bond to be closed and $`1p=\mathrm{exp}(\beta )`$ is the probability of bond to be open; summation is performed over all spin configurations $`\sigma `$. Sometimes the term in square brackets is expressed via $`v=\mathrm{exp}(\beta )1`$: $`\left[(1p)+p\delta (\sigma _i,\sigma _j)\right]=\mathrm{exp}(\beta )\left[1+v\delta (\sigma _i,\sigma _j)\right]`$, but we write Eq (2) via $`p`$ to emphasize the fact that in correlated site-bond percolation FK , the probability of bond to be closed is $`p`$ and to be open is $`1p`$.

For the square lattice of linear size $`L`$ with periodic boundary conditions the total number of bonds is $`N=2L^2`$. Let us define by $``$ the graph of all the edges on the lattice. The product over all bonds $`(i,j)`$ consists of $`N`$ terms, so the product may be expanded as sum of $`2^N`$ terms. Let us associate to each of these $`2^N`$ terms subgraph $`G`$ of graph $``$, by following rule. Each of $`2^N`$ terms can be considered as the product of $`N`$ factors. Each factor $`[(1p)+p\delta (\sigma _i,\sigma _j)]`$ corresponds to some edge $`(i,j)`$ of the graph $``$. To construct $`G`$, we perform the following procedure. If this factor for edge $`(i,j)`$ is equal to $`1p`$, we delete edge $`(i,j)`$ from subgraph. If this factor for edge $`(i,j)`$ is equal to $`p`$, we leave this edge in the subgraph. So we associate to each term in sum (2) subgraph $`G`$. Each subgraph $`G`$ consists of $`b(G)`$ edges (closed bonds) and $`C(G)`$ connected components. The term in (2) corresponding to $`G`$ contains factor $`(1p)^{Nb(G)}p^{b(G)}`$, and delta-functions guarantee equivalence of spins in each connected component. As a result of the summation over all possible spin configurations, the contribution of this subgraph $`G`$ into the partition function is $`q^{C(G)}(1p)^{Nb(G)}p^{b(G)}`$. Therefore, we can replace the summation over all spin configurations in Eq. (2) by summation over all possible subgraph $`G`$ on $``$

$$Z=\underset{G}{}q^{C(G)}(1p)^{Nb(G)}p^{b(G)}$$
(3)
We shall keep in mind that this definition is valid even for the non-integer value of $`q`$. The partition function of the Potts model for $`q1`$ corresponds to the percolation, where $`p`$ is a probability of bond to be occupied.

Now we introduce the crossing probability $`\pi _{hs}`$, the probability of a system to percolate only in the horizontal direction while the percolation in the vertical direction is absent. We must distinguish it from the probability of a system to percolate in the horizontal direction irrespective to percolation in the vertical direction $`\pi _h`$. We define the indicator function $`I_{hs}(G)`$ for each subgraph $`G`$ in accordance with the rule

$$I_{hs}(G)=\{\begin{array}{cc}1\hfill & \mathrm{if}G\mathrm{percolates}\mathrm{only}\mathrm{in}\mathrm{horizontal}\mathrm{direction}\hfill \\ 0\hfill & \mathrm{in}\mathrm{all}\mathrm{others}\mathrm{cases}\hfill \end{array}$$
(4)
We mean that $`G`$ percolates only in the horizontal direction if it contains at least one connected component touching left and right sides of the lattice, and there are no components joining top and bottom sides. Therefore, the crossing probability $`\pi _{hs}(\beta ;L,q)`$ of a system of size $`L`$ with $`q`$ possible spin states to percolate only in the horizontal direction at $`\beta `$ can be written as

$$\pi _{hs}(\beta ;L,q)=\frac{1}{Z}\underset{G}{}I_{hs}(G)q^{C(G)}(1p)^{Nb(G)}p^{b(G)}$$
(5)
Here, $`p=1\mathrm{exp}(\beta )`$. We can introduce the crossing probability for vertical direction $`\pi _{vs}`$ by the same way. But on the square lattice $`\pi _{hs}=\pi _{vs}`$ and later on, investigate only $`\pi _{hs}`$.

Let us notice that $`\pi _{hs}(\beta ;L,q)`$ has maximum near the critical point $`\beta _c(L,q)`$ because in the ordering phase $`\beta >\beta _c(L,q)`$, the most probable configurations contain large clusters touching top and bottom sides of the lattice, and do not contribute into $`\pi _{hs}`$.

## III Numerical results

We use the Wolff cluster algorithm Wolff  to generate the different spin configurations. For each spin configuration we generate a bond configuration in accordance with FK . Then, we break the lattice into independent clusters of connected sites by using the Hoshen-Kopelman HK  algorithm. Then, we analyze crossing properties $`I_{hs}(G)`$ of this configuration. Between checking the crossing we skip 5 Monte-Carlo steps. For each value of $`\beta `$, the averaging is performed over $`10`$ series each of $`10^5`$ configurations. The total number of configurations is $`10^6`$. Sets of configurations are used for the estimation of the numerical inaccuracy. Quantity $`I_{hs}`$ is an indicator function. It means that it takes values 0 or 1 for each configuration. Therefore, the resolution of our computations for $`\pi _{hs}`$ is $`10^6`$. We compute data for $`q=1(\mathrm{percolation}),2,3,4`$ and lattice sizes $`L=16,32,48,64,80,96,112,128`$. For each pair $`(q,L)`$, we perform computation for 200 values of $`\beta `$ (or $`p`$ for percolation) in the critical region.

There is a very important question, how the scaling variable should be chosen?

For the percolation $`q=1`$ the scaling variable is defined as the deviation of the bond concentration from critical point $`r=pp_c`$. For the Potts model, we can take

1. $`t=\frac{TT_c}{T_c}`$
2. $`\tau =\frac{\beta \beta _c}{\beta _c}`$
3. $`r=(1\mathrm{exp}(\beta ))(1\mathrm{exp}(\beta _c)=pp_c`$

The first variant $`t=\frac{TT_c}{T_c}`$ is usual for the investigation of thermodynamic quantities St ; Baxter , the second one is widely used for the approximation of magnetic susceptibility $`\chi `$ of the Ising model in the critical region Wu , and the third variant $`r=pp_c`$ may be argued by the fact that for the site-bond correlated percolation, the probability of bond to be occupied is $`p=1\mathrm{exp}(\beta )`$ – see equations (2), (3).

Let us evaluate the symmetry of the crossing probability as a function of all these variables by comparing third moments of the function $`\pi _{hs}`$. We assume that $`\pi _{hs}`$ must be the most symmetric function as a function of ”right” variable. We calculate four moments of $`\pi _{hs}(x)`$ as a functions of $`x=t,\tau ,r`$ by numerical integration in accordance with formula $`\mu _0=\pi _{hs}(x)𝑑x`$, $`\mu _1=\frac{1}{\mu _0}x\pi _{hs}(x)𝑑x`$, $`\mu _k=\frac{1}{\mu _0}(x\mu _1)^k\pi _{hs}(x)𝑑x`$. Results of computations $`\mu _3`$ for the lattice size $`L=128`$ are given in Table 1.

We see that for all the cases $`\mu _3`$ is smaller for scaling variable $`r`$, with one exception $`q=4`$, $`\tau `$, when numerical error is approximately thirty times greater than value $`\mu _3(\tau )`$ and six times greater than $`\mu _3(r)`$. For all other lattice sizes $`L=16,\mathrm{},112`$ the third moment $`\mu _3`$ calculated by using variable $`r`$, is smaller than for $`t`$ or for $`\tau `$. Therefore, we will work with the scaling variable $`p=1\mathrm{exp}(\beta )`$, the probability of bond to be closed for $`q=2,3,4`$. The critical point of the $`q`$-state Potts model on the infinite lattice is Baxter  $`p_c(q)=\frac{\sqrt{q}}{\sqrt{q}+1}`$. The bond percolation critical point $`p_c(q=1)=\frac{1}{2}`$ can be obtained from this formula. Variable $`p=1\mathrm{exp}(\beta )`$ naturally provides crossover from the percolation to the Potts model. We can plot the crossing probability $`\pi _{hs}`$ as a function of $`p`$ for the percolation and the Potts model on the same graph.

In Fig.1, data for $`\pi _{hs}(p;L,q)`$ are plotted for $`q=1,2,3,4`$ and $`L=32,64,128`$. We can see shift of the critical point $`p_c(q,L)`$ for finite lattice sizes $`L`$. We also see the change of the width of the function due to size scaling.

Let us try to identify the shape of $`\pi _{hs}`$ by analyzing the ratio $`\frac{\mu _4}{\mu _2^2}`$. We know that at least for the percolation, function $`\pi _{hs}`$ looks like the Gaussian function vas  $`A\mathrm{exp}(x^2)`$, and for the Gaussian function we expect $`\frac{\mu _4}{\mu _2^2}(Gaussian)=3`$. But, in reality, crossing probabilities are not Gaussian. The moments ratio for function $`\pi _h^{}(p)`$ (the derivative of crossing probability $`\pi _h(p)`$ with respect to $`p`$), for the percolation model was found by Ziff $`\frac{\mu _4}{\mu _2^2}(q=1)3.20(5)`$ Ziff1  and more precisely $`\frac{\mu _4}{\mu _2^2}(q=1)3.176(4)`$ Ziff2  and by Hovy and Aharony $`\frac{\mu _4}{\mu _2^2}(q=1)3.174(25)`$ HA . The explanation of the non-Gaussian crossing probability shape was first given by Berlyard and Wehr BW  and was also discussed in Ref. Wester . Newman and Ziff verify that the tails of the crossing probability falls of as $`\mathrm{exp}(c(pp_c)^{4/3})`$, therefore the tail of the distribution is characterized by the correlation length exponent NZ1 ; NZ2 . However, as the first approximation we introduce an additional scaling index $`\zeta `$ and check moment ratios for function $`g(x;\zeta )=A\mathrm{exp}(x^\zeta )`$. For this function $`g(x;\zeta )`$ the moment ratios are

$$\frac{\mu _4}{\mu _2^2}(\zeta )=\frac{\mathrm{\Gamma }(\frac{1}{\zeta })\mathrm{\Gamma }(\frac{5}{\zeta })}{\mathrm{\Gamma }(\frac{3}{\zeta })^2}$$
(6)
We calculate moments for several values of $`\zeta `$ and put these data into the first ($`\zeta `$) and second ($`\frac{\mu _4}{\mu _2^2}(\zeta )`$) rows of the Table 2. The choice of values $`\zeta `$ will be argued later. Results of computation of the moments ratio $`\frac{\mu _4}{\mu _2^2}(q)`$ for the crossing probability $`\pi _{hs}(p;L,q)`$ for $`L=32,128`$ is placed in the fourth and fifth row of the Table 2. We see that the moments ratio practically does not depend upon the lattice size $`L`$. We check this fact for others lattice sizes.

Moment ratios for $`g(x;\zeta )`$ slightly differ from moment ratios for $`\pi _{hs}`$ (this fact will be explained later), but we can try to approximate the crossing probability by the function $`g(x,\zeta )`$ and then compare results of approximation with numerical data.

Below, we describe the fitting procedure. As we can see from Fig.1, there are many nonuniversal scaling factors for $`\pi _{hs}`$: the amplitude of the crossing probability $`A(L,q)`$, the finite size critical point $`p_c(L,q)`$, which differs from $`p_c(L=\mathrm{},q)=\frac{\sqrt{q}}{\sqrt{q}+1}`$, the nonuniversal scaling factor $`B(L,q)`$, which provide the finite size scaling of the function, and the additional scaling index $`\zeta (L,q)`$.

We perform the four-parametric fit of $`\pi _{hs}(p;L,q)`$ by the function (7)

$$F(p;L,q)=A(q,L)\mathrm{exp}(\left[B(L,q)(pp_c(L,q))\right]^{\zeta (L,q)})$$
(7)
We use points $`\mathrm{log}(\pi _{hs})>9`$ for this fit and the log scale for the ordinate axis. As a result of the approximation we get a set of scaling amplitudes $`A(q,L)`$, nonuniversal metric $`B(q,L)`$, critical points $`p_c(q,L)`$ and scaling indices $`\zeta (q,L)`$. Then, we use this scaling factors to adjust our numerical data onto one line for $`q=1,2,3,4`$ and $`L=32,128`$. We plot $`f=|\mathrm{log}(\pi _{hs}(p;L,q))\mathrm{log}(A(L,q))|`$ as a function of the new scaling variable $`z=\left[B(L,q)\left(pp_c(L,q)\right)\right]^{\zeta (L,q)}`$ in Fig.2.

The scaling straight lines $`|z|`$ shown in Fig.2, correspond to the fitting function $`F(z)=\mathrm{exp}(|z|)`$. We see that all points lie on the one curve, and this curve is very close to $`|z|`$ in the range $`0.1<z<3`$. So, we might expect that our fitting procedure is correct. However, on the graph, the tail points deviate from the line $`|z|`$, and this deviation explains the fact that moment ratios for $`\pi _{hs}(p;L,q)`$ in the Table 2 differ from analytical values for $`g(x,\zeta )`$. Results of approximation for the crossing amplitude $`A(L,q)`$ are placed in Table 3. For each value of $`q=1,2,3`$ amplitudes $`A(L,q)`$ are equal within our numerical accuracy of the approximation. Therefore, we can conclude that the scaling amplitude depends upon $`q`$ and depends weakly upon $`L`$.

Results of the approximation for the scaling index $`\zeta (L,q)`$ are represented in the Table 4. We can see that this scaling index practically does not depend upon lattice size.

We assume that for each value $`q`$ we have the additional scaling index $`\zeta (q)`$, which is nondependent upon the lattice size. Test values of $`\zeta `$ in the first row of Table 2 are chosen to be close to $`\zeta (q)`$. It seems, that $`\zeta (3)=4/3`$, $`\zeta (4)=6/5`$.

We know that in accordance with the scaling theory for each fixed value $`q`$ the crossing probability must be a universal function of variable $`x=L^{\frac{1}{\nu (q)}}(pp_c(L,q))`$. Therefore, we approximate the numerical data for $`B(L,q)`$ (see Fig.3) by function $`b(q)L^{y(q)}`$ and represent results in Table 5. Results of the approximation $`b(q)L^{y(q)}`$ are also shown in Fig.3 by lines.

We can see that the thermal scaling index $`y(q)`$ is very close to the analytical value of the inverse correlation length index $`\frac{1}{\nu (q)}`$ Baxter , which is represented in the Table 5 for the comparison. This fact once more confirms our approximation procedure. The exception is the case $`q=4`$, for which there is a difference between $`y(q=4)`$ and $`\frac{1}{\nu (q=4)}`$.

Many critical quantities in the Potts model $`q=4`$ exhibit logarithmic corrections CNS ; SS ; AA ; CZ . These logarithmic corrections explain the difference between analytical value $`\frac{1}{\nu (q=4)}=1.5`$ and numerical approximation for the scaling index $`y=1.39(6)`$.

## IV Conclusions

To summarize results, we can conclude that the crossing probability $`\pi _{hs}(p;L,q)`$ is a universal function $`Q(z)`$

$$\pi _{hs}(p;L,q)=A(q)Q\left(\left[b(q)L^{\frac{1}{\nu (q)}}(pp_c(L,q))\right]^{\zeta (q)}\right)$$
(8)
of the scaling variable $`z=\left[b(q)L^{\frac{1}{\nu (q)}}(pp_c(L,q))\right]^{\zeta (q)}`$. As we can conclude from Fig.2, the function $`Q(z)`$ looks like the exponent $`Q(z)\mathrm{exp}(|z|)`$ on the interval $`0.1<z<3`$, but deviates from it in the vicinity of 0 and on the tails $`|z|>3`$.

Let us pay attention to the new details of this work. We consider the universality of the crossing probability for different values of $`q`$ by adding the scaling index $`\zeta (q)`$. We work in the scale $`p`$, where $`p=1\mathrm{exp}(\beta )`$ is the probability of a bond to be closed instead of the usual scale $`t=(TT_c)/T_c`$ to make the crossing probability symmetric. We find numerically the scaling index $`\zeta (q)`$. The universal function looks like $`Q(z)\mathrm{exp}(|z|)`$ on the interval $`0.1<z<3`$.

The author would like to thank Robert M. Ziff for helpful remarks, comments and discussion. The author is grateful to S. Nechaev for revision of the manuscript and to the Joint SuperComputer Center RAS (www.jscc.ru) for providing computational resources.


text:
# Interstellar Scintillation Velocities of a Relativistic Binary PSR B1534+12 and Three Other Millisecond Pulsars

## 1 INTRODUCTION

Interstellar scintillation (ISS) of pulsars is caused by diffraction and refraction of radio waves from electron density irregularities in the interstellar plasma (Rickett, 1990). Measurements of the ISS related modifications of the pulsar signal, such as the angular source broadening, temporal broadening of pulses and their intensity variations in time and frequency, have been used to study the small-scale structure of the interstellar plasma (e.g. Cordes, Weisberg, & Boriakoff, 1985; Lambert & Rickett, 2000).

Pulsars are fast moving objects with mean velocities of $``$460 km s-1 (Lyne & Lorimer, 1994) and $``$85 km s-1 (Toscano et al., 1999) for the normal and the millisecond pulsar populations, respectively. Therefore, measurements of a net speed of the interstellar diffraction pattern derived from decorrelation scales of pulse intensity variations in time and frequency are dominated by the pulsar motion and can be used to estimate a transverse component of the pulsar velocity in space. In fact, the ISS technique has been routinely used as one of the means to measure pulsar velocities (e.g. Lyne & Smith, 1982; Cordes, 1986; Gupta, Rickett, & Lyne, 1994), in addition to the more direct methods of pulse timing (e.g. Toscano et al., 1999) and radio interferometry (e.g. Fomalont et al., 1997). As the latter two techniques are limited to low timing noise, bright, nearby objects, the ISS observations provide a very helpful, additional way to obtain velocity estimates for a relatively large number of sources. Such measurements, in combination with simulation studies (Tauris & Bailes, 1996; Cordes & Chernoff, 1997) are important for global characterization of the galactic population of neutron stars.

The ISS velocity measurements, when combined with timing or interferometric astrometry data, provide useful constraints on the distribution of scattering material along the pulsar line of sight and provide a way to verify a reliability of the ISS method itself (Cordes & Rickett, 1998; Deshpande & Ramachandran, 1998). Recent analyses of a growing body of such measurements have shown that, despite their obvious model dependence, the ISS velocities are remarkably accurate (Gupta, 1995; Johnston, Nicastro, & Koribalski, 1998; Nicastro et al., 2001).

In principle, the ISS observations of binary pulsars can be used to track their orbital velocities projected onto the plane of the sky and to deduce some orbital parameters that would be impossible to measure in any other way (Lyne, 1984; Dewey, Cordes, & Wolszczan, 1988). This is particularly true for relativistic binaries, in which case the orbital velocity variations are exceptionally large and may dominate the ISS measurements.

In this paper, we present ISS velocity measurements of four millisecond pulsars obtained from long-term monitoring observations with the 305-m Arecibo radio telescope. PSR J1640+2224 (Foster et al., 1995) and PSR J1713+0747 (Foster, Wolszczan, & Camilo, 1993) are binary millisecond pulsars in long period, nearly circular orbits with white dwarf companions. PSR B1257+12 is a fast moving solitary millisecond pulsar with a system of at least three planets around it (Wolszczan & Frail, 1992; Wolszczan, 1994). Finally, PSR B1534+12 is a relativistic 10-hour binary pulsar in orbit with another neutron star (Wolszczan, 1991; Stairs et al., 1998). Here, we focus our attention on PSR B1534+12 and the measurements of changes in its ISS velocity caused by orbital motion. We also verify the earlier ISS speed measurements for PSR B1257+12 and PSR J1713+0747 (Johnston et al., 1998; Gothoskar & Gupta, 2000) and provide the first such measurement for PSR J1640+2224.

## 2 OBSERVATIONS AND ANALYSIS OF THE DYNAMIC SPECTRA

The four pulsars discussed in this paper have been regularly monitored with the Arecibo radio telescope to measure their pulse arrival times. The most recent results of timing modeling of these objects will be presented elsewhere. We have reprocessed the timing measurements made at multiple epochs between 1994 and 2001 at 430 MHz to obtain the dynamic spectra of pulsar scintillations and analyzed them to estimate the ISS velocities of the pulsars themselves.

Typically, PSR B1257+12, PSR J1640+2224 and PSR J1713+0747 were observed for 45-120 min. in a series of 3 min. integrations carried out synchronously with the Doppler-corrected pulsar period. In the case of PSR B1534+12, the dynamic spectra were analyzed in 45-60 min. sections to avoid excessive averaging over changes in the frequency and time decorrelation scales induced by the pulsar’s orbital motion. Similarly, to eliminate biases in the decorrelation scales introduced by refractive scintillation (Bhat, Rao, & Gupta, 1999) the dynamic spectra strongly affected by drifts of spectral features and quasiperiodic fringing were excluded from the analysis.

The left- and right-hand circularly polarized pulsar signals were fed into the Penn State Pulsar Machine (PSPM), which is a 2 $`\times `$ 128 $`\times `$ 60 kHz channel, computer controlled, on-line pulsar processor. After detection, the dual-channel signals were added together, smoothed, 4-bit quantized and pulse phase averaged by the analysis computer in each of the 128 frequency channels. In the off-line processing, the 128 $`\times `$ 60 kHz point spectra of the pulsar signal were formed from the ON-pulse and OFF-pulse level measurements in each frequency channel by calculating the corresponding (ON-OFF)/OFF spectral intensities. A dynamic spectrum of the pulsar intensity variations in time and frequency was then constructed from a set of 3-minute averaged, contiguous spectra for each object and each observing session.

Examples of the dynamic spectra of the four pulsars are shown in Figures 1 and 2. Also shown are the two-dimensional autocorrelation functions of these spectra, which were used to measure the decorrelation scales of pulsar intensity variations, $`\mathrm{\Delta }t_d`$ and $`\mathrm{\Delta }\nu _d`$. Following the established convention, $`\mathrm{\Delta }t_d`$ and $`\mathrm{\Delta }\nu _d`$ were measured as the respective $`1/e`$ widths and half-widths of elliptical Gaussians least-squares fitted to the autocorrelation functions. The decorrelation widths for PSR B1257+12, PSR B1534+12, PSR J1640+2224 and PSR J1713+0747, averaged over all observing epochs are listed in Table 1.

The errors of $`\mathrm{\Delta }t_d`$ and $`\mathrm{\Delta }\nu _d`$ measurements consist of statistical uncertainty due to a finite number of “scintles” in the dynamic spectra, the signal strength, and the long-term changes of these parameters caused by refractive effects (Bhat et al., 1999). We found that the refraction induced variations in $`\mathrm{\Delta }t_d`$ and $`\mathrm{\Delta }\nu _d`$ occurring on the time scales from months to years were a dominant source of error for all the four pulsars discussed here. Consequently, the errors quoted in Table 1 were calculated from the scatter in the estimates of the decorrelation widths over the entire time span of observations.

## 3 THE ISS VELOCITY ANALYSIS

The ISS velocities of pulsars can be calculated from their dynamic spectra, because the two observables, the characteristic decorrelation scales of pulse intensity in time, $`\mathrm{\Delta }t_d`$, and in frequency, $`\mathrm{\Delta }\nu _d`$ are directly related to the net speed of a diffractive scintillation pattern (Lyne & Smith, 1982). For $`\mathrm{\Delta }t_d`$ and $`\mathrm{\Delta }\nu _d`$ measured in seconds and megahertz, respectively, the velocity $`V_{ISS}`$ (km s-1) is derived from the equation (Gupta et al., 1994):

$$V_{ISS}=3.85\times 10^4\frac{\sqrt{\mathrm{\Delta }\nu _ddx}}{\nu \mathrm{\Delta }t_d},$$
(1)
where $`x`$ is the ratio of the distances of the scattering screen to the observer and to the pulsar, $`d`$ is the pulsar distance (in kpc), and $`\nu `$ is the observing frequency (in GHz). The numerical constant is valid for the assumed homogeneous turbulent medium characterized by a Kolmogorov turbulence spectrum.

In principle, the observed ISS velocity of a pulsar represents the net effect of the pulsar’s motion, the Earth’s orbital motion around the Sun, and the motion of the scattering medium. Thus, $`V_{ISS}`$ can be written as the magnitude of the vector sum of all the velocities contributing to it (Gupta et al., 1994):

$$V_{ISS}=|x𝑽_𝒐+x𝑽_{𝒑𝒎}+𝑽_{\mathbf{}}(1+x)𝑽_{𝒊𝒔𝒎}|_{}$$
(2)
where, for binary pulsars, $`𝑽_𝒐`$ is the pulsar’s orbital velocity, $`𝑽_{𝒑𝒎}`$ is its proper motion velocity, $`𝑽_{\mathbf{}}`$ is the velocity of the Earth, and $`𝑽_{𝒊𝒔𝒎}`$ is the velocity of the scattering screen. The subscript $``$ indicates that the only contributions to the ISS velocity come from the components of the velocity vectors perpendicular to the line of sight. Choosing the origin of a coordinate system to coincide with the center of mass of a binary with its x-y plane tangent to the sky, with the x-axis parallel to the orbit’s line of nodes and the z-axis pointing from the solar system barycenter to that of the binary system, $`V_{ISS}`$ can be expressed in terms of the projections of the component velocity vectors onto the plane of the sky:

$`V_{ISS}^2=`$ $`(xV_{o,x}+xV_{pm,x}+V_{,x})^2`$ (3)
$`+(xV_{o,y}+xV_{pm,y}+V_{,y})^2`$
In this expression we have assumed that the contribution of $`V_{ism}`$ to $`V_{ISS}`$ is negligible (Gupta, 1995). The component velocities in Eqn. (3) are defined as follows. The orbital velocity components are:

$`V_{o,x}=\xi `$ $`\{[e\mathrm{sin}\varphi \mathrm{cos}\varphi (1+e\mathrm{cos}\varphi )\mathrm{sin}\varphi ]\mathrm{cos}\omega `$ (4)
$`[e\mathrm{sin}\varphi \mathrm{sin}\varphi +(1+e\mathrm{cos}\varphi )\mathrm{cos}\varphi ]\mathrm{sin}\omega \}`$
$`V_{o,y}=\xi `$ $`\{[e\mathrm{sin}\varphi \mathrm{cos}\varphi (1+e\mathrm{cos}\varphi )\mathrm{sin}\varphi ]\mathrm{sin}\omega `$ (5)
$`+[e\mathrm{sin}\varphi \mathrm{sin}\varphi +(1+e\mathrm{cos}\varphi )\mathrm{cos}\varphi ]\mathrm{cos}\omega \}(\pm \mathrm{cos}i)`$
where $`\xi =(2\pi a\mathrm{sin}i)/(P_b\sqrt{1e^2}\mathrm{sin}i)`$, $`a`$ is the semi-major axis, $`e`$ is the eccentricity, $`P_b`$ is the orbital period, $`\varphi `$ is the orbital phase, $`\omega `$ is the longitude of periastron, and $`i`$ is the orbital inclination. The two possible signs of the $`\mathrm{cos}i`$ term arise from the fact that we have no way of determining which of the two nodes of the orbit is the ascending node so both cases need to be considered. The proper motion components are defined by the familiar relationships:

$$V_{pm,x}=4.74d(\mu _\alpha \mathrm{sin}\mathrm{\Omega }+\mu _\delta \mathrm{cos}\mathrm{\Omega })$$
(6)
$$V_{pm,y}=4.74d(\mu _\alpha \mathrm{cos}\mathrm{\Omega }+\mu _\delta \mathrm{sin}\mathrm{\Omega })$$
(7)
where $`\mathrm{\Omega }`$ is the longitude of the ascending node of the orbit, defined as the position angle (reckoned from north through east) of the line of nodes measured in the plane tangent to the sky, $`\mu _\alpha `$ and $`\mu _\beta `$ are the proper motions in right ascension and declination, respectively, in units of mas yr-1, and the distance $`d`$ is given in kpc. Finally, the Earth velocity components are given by:

$`V_{,x}=`$ $`1.73\times 10^6`$ $`[\dot{X}(\mathrm{cos}\mathrm{\Omega }\mathrm{sin}\delta \mathrm{cos}\alpha \mathrm{sin}\mathrm{\Omega }\mathrm{sin}\alpha )`$
$`+\dot{Y}(\mathrm{cos}\mathrm{\Omega }\mathrm{sin}\delta \mathrm{sin}\alpha +\mathrm{sin}\mathrm{\Omega }\mathrm{cos}\alpha )`$
$`+\dot{Z}(\mathrm{cos}\mathrm{\Omega }\mathrm{cos}\delta )]`$
$`V_{,y}=`$ $`1.73\times 10^6`$ $`[\dot{X}(\mathrm{sin}\mathrm{\Omega }\mathrm{sin}\delta \mathrm{cos}\alpha +\mathrm{cos}\mathrm{\Omega }\mathrm{sin}\alpha )`$
$`+\dot{Y}(\mathrm{sin}\mathrm{\Omega }\mathrm{sin}\delta \mathrm{sin}\alpha \mathrm{cos}\mathrm{\Omega }\mathrm{cos}\alpha )`$
$`+\dot{Z}(\mathrm{sin}\mathrm{\Omega }\mathrm{cos}\delta )]`$
where $`\alpha `$ and $`\delta `$ are the pulsar coordinates and the parameters $`\dot{X}`$, $`\dot{Y}`$, and $`\dot{Z}`$ are the components of Earth’s velocity in the Solar System Barycenter coordinate system expressed in units of AU d-1.

Expression (3) conveniently reduces to one for a solitary pulsar by simply neglecting the orbital motion terms and setting $`\mathrm{\Omega }`$ equal to zero:

$`V_{ISS}^2=`$ $`\{x(4.74d\mu _\delta )`$ $`+1.73\times 10^6[\dot{X}(\mathrm{sin}\delta \mathrm{cos}\alpha )+\dot{Y}(\mathrm{sin}\delta \mathrm{sin}\alpha )+\dot{Z}\mathrm{cos}\delta ]\}^2`$ (10)
$`+\{x(4.74d\mu _\alpha )`$ $`+1.73\times 10^6[\dot{X}\mathrm{sin}\alpha \dot{Y}\mathrm{cos}\alpha ]\}^2`$
In what follows, Eq. (1) is used to derive the observed ISS velocities from ISS measurements, whereas Eq. (3) or Eq. (10) give the expected ISS speeds from the known pulsar proper motion and Earth velocities.

## 4 ISS VELOCITY AND ORBITAL MOTION OF PSR B1534+12

Johnston et al. (1998), Gothoskar & Gupta (2000), and Lommen (2001) have measured the ISS velocity for this pulsar to be 190 km s-1, 191 km s-1, and 149 km s-1, respectively. In principle, these estimates can be significantly biased by a varying transverse component of the pulsar’s orbital velocity, whose semi-amplitude exceeds $``$100 km s-1. Here we have performed a rigorous analysis of this effect by taking into account all the motions contributing significantly to the observed ISS velocity of PSR B1534+12, as well as the dependence of the ISS velocity on the pulsar’s orbital phase.

Compared to the results of similar analyses carried out for the binary pulsars PSR B0655+64 by Lyne (1984) and PSR B1855+09 and PSR B1913+16 by Dewey et al. (1988), the ISS velocities for PSR B1534+12 shown in Fig. 3 exhibit a much more pronounced dependence on orbital phase and allow an unambiguous modeling of this effect with the aid of Eqs. (3-9). In fact, the problem can be conveniently parametrized in terms of only two unknowns, the longitude of the ascending node $`\mathrm{\Omega }`$ and the screen location $`x`$, with the Keplerian orbital parameters and proper motion of the pulsar already determined from long-term timing measurements (Stairs et al., 1998, and references therein). For the pulsar distance we have adopted a value of $`d=1.1\pm 0.2`$ kpc derived by Stairs et al. (1998) on the assumption that General Relativity (GR) is the correct theory of gravity.

Although the longitude of periastron $`\omega `$ of the pulsar’s orbit advances at a rate of 1.75 deg yr-1, it does not change fast enough to have a significant effect on our analysis. Consequently, we have used the value of $`\omega `$ at the mid-point of the actual observing period. Finally, because the transverse component of the Earth’s velocity varies and its direction relative to the pulsar’s orbital motion cannot be determined without knowing $`\mathrm{\Omega }`$, we have ignored the $`V_{}`$ terms in Eq. (3) and included them in the uncertainties of the $`V_{ISS}`$ measurements.

Since both the model (Eq. 3) and Eq. (1) include a dependence on the parameter $`x`$, the best fit values of $`\mathrm{\Omega }`$ and $`x`$ were obtained by manually searching the $`\chi ^2`$-space for the minimum. The fit yielded $`\mathrm{\Omega }_1=70\pm 20`$ and $`\mathrm{\Omega }_2=290\pm 20`$ degrees for the longitude of the ascending node and $`x=1.3\pm 0.4`$, which corresponds to a screen located at a distance of $`d_s=630\pm 200`$ pc from the Earth. Fig. 3 shows the ISS velocity measurements of PSR B1534+12 and the best fitting curve given by Eq. (3) plotted against orbital phase.

¿From our ISS velocity measurements of PSR B1534+12 it is also possible to give an estimate of the pulsar’s transverse systemic velocity. Since the transverse component of the pulsar’s orbital velocity never vanishes the best we can do is to consider only the measurements made around orbital phases 0.2 and 0.7 at which it reaches the minimum value. Using these measurements we obtain a value of 192 km s-1, which is an overestimate of the actual proper motion due to contributions from the pulsar’s orbital motion as well as the location of a scattering screen. Nonetheless, this value is remarkably close to those published by Johnston et al. (1998) and Gothoskar & Gupta (2000).

## 5 ISS VELOCITIES OF PSR B1257+12, PSR J1640+2224, AND PSR J1713+0747

Scintillation velocities of PSR B1257+12 were calculated using the distance based on the Taylor & Cordes (1993) model of galactic electron distribution (the TC model hereafter). For the pulsar distance of 0.62 kpc, we obtained an average value of $`V_{ISS}=197\pm 57`$ km s-1. The ISS velocity estimate of 225 km s-1 using the same distance has been derived by Gothoskar & Gupta (2000). We find these two estimates to be remarkably similar given the large, long-term variations of the ISS velocity observed for this pulsar. The measured velocity variations ranging from 117 km s-1 to 290 km s-1 on a time scale of $``$2 years are consistent with modifications of the scintillation decorrelation scales by the refractive scintillation. We have found the ISS velocity of this pulsar to be consistently lower than the effective transverse velocity of $``$278 km s-1, computed using the proper motion values given by Wolszczan et al. (2000), with the average ratio $`V_{ISS}/V_{eff}=0.67`$. This difference can be plausibly accounted for by assuming that the scattering screen is located closer to the observer than to the pulsar ($`x<1`$). It is also likely that the distance used is an underestimate of the true distance to the pulsar.

Our measurement of $`V_{ISS}`$ of 38 $`\pm `$ 8 km s-1 for PSR J1640+2224 represents the first ISS velocity estimate for this binary millisecond pulsar. In the calculations we adopted the TC model distance of 1.18 kpc. In comparing this velocity with the value of 75 km s-1 obtained from Eq. (10), computed with proper motion measurements provided by Wolszczan et al. (2000), we have ignored the transverse component of the pulsar’s orbital velocity, because it is unlikely to be greater than $`10`$ km s-1. On the other hand, the transverse component of the Earth’s orbital velocity accounts for as much as 40% of the effective ISS velocity of the pulsar and it must be included in the analysis. Since the ISS velocity for PSR J1640+2224 is lower than the effective transverse velocity with the ratio $`V_{ISS}/V_{eff}=0.52`$, it is likely that the effective location of the scattering screen is closer to the observer than to the pulsar. However, this discrepancy may also mean that the distance of 1.18 kpc we used may be an underestimate of the actual pulsar distance.

ISS velocity estimates of 29 km s-1 and 118 $`\pm `$ 40 km s-1 for PSR J1713+0747 have already been published by Johnston et al. (1998), using a distance of 0.89 kpc, and by Lommen (2001), respectively. Even if corrected for the TC distance of 1.1 kpc used here, the value given by Johnston et al. is still significantly lower than our velocity estimate of 82 $`\pm `$ 16 km s-1 and that of Lommen (2001). In absence of observations made at sufficiently similar epochs, the source of large discrepancies between the three values remains uncertain, although the refraction induced $`V_{ISS}`$ variability provides a possible explanation. Since the transverse component of the orbital velocity of PSR J1713+0747 does not contribute more than $`15`$ km s-1 to the pulsar’s ISS velocity, we included it in the effective uncertainty of the measurement and used Eq. (10) instead of Eq. (3). With a net effective velocity of 52 km s-1 for this pulsar, calculated using the proper motion values given by Camilo, Foster, & Wolszczan (1994), our measurements of the ISS velocity are consistently larger than the effective transverse velocity. In this case, the ratio $`V_{ISS}/V_{eff}`$ = 1.71 suggests either the presence of an asymmetrically located scattering screen which is closer to the pulsar than to the observer (i.e. $`x>1`$) or that the pulsar distance is overestimated.

## 6 DISCUSSION

We have presented measurements of interstellar scintillation velocities for four pulsars together with a complete set of relationships that should be used to compare these velocities with the proper motion measurements. Our $`V_{ISS}`$ estimates are generally consistent with the corresponding proper motion based values derived from timing measurements of these objects. As shown by the results of similar analyses conducted in the past (e.g Johnston et al., 1998; Gothoskar & Gupta, 2000), the discrepancies between these two kinds of the transverse pulsar velocity measurements can be accounted for by either an asymmetric location of the effective scatterer along the line of sight, or by an inaccurate estimate of the pulsar distance. In principle, these two factors can be decoupled by the analysis that involves the angular and temporal broadening data, in addition to the scintillation and proper motion measurements of a pulsar (Cordes & Rickett, 1998; Deshpande & Ramachandran, 1998). It is likely, that these “hybrid” methods will offer an important means to study a line-of-sight distribution of the scattering material, as timing and astrometric proper motion measurements for many more pulsars become available.

Our analysis of the relativistic binary pulsar PSR B1534+12 provides the best demonstration so far that the changes of a transverse component of the pulsar orbital velocity can be unambigously detected in the ISS velocity data. The most important practical consequence of this detection is the measurement of the longitude of the ascending node $`\mathrm{\Omega }`$ of the binary orbit. This parameter cannot be determined from pulsar timing, except in very special circumstances (Kopeikin, 1996; van Straten et al., 2001). One of the instances where the knowledge of $`\mathrm{\Omega }`$ is important is the effect of aberration caused by the relativistic orbital motion of a pulsar (Damour & Taylor, 1992). In this case, $`\mathrm{\Omega }`$ has to be taken into account in the determination of aberration-induced variations of the orbital phase dependent offset of the linear polarization angle of the pulse profile (Damour & Taylor, 1992, Eqs. 2.16-2.24). The value of $`\mathrm{\Omega }`$ is also needed for a precise measurement of a proper motion related rate of change of the projected size of the pulsar orbit (Kopeikin, 1996, and references therein):

$$\delta \dot{x}=1.54\times 10^{16}x\mathrm{cot}i(\mu _\alpha \mathrm{sin}\mathrm{\Omega }+\mu _\delta \mathrm{cos}\mathrm{\Omega })$$
(11)
where, $`x`$ is the size of the projected semi-major axis expressed in seconds. From this equation we obtain values of $`\delta \dot{x}_1=9.43\times 10^{16}`$ s s-1 and $`\delta \dot{x}_2=1.25\times 10^{15}`$ s s-1, corresponding to $`\mathrm{\Omega }_1`$ and $`\mathrm{\Omega }_2`$, respectively. With the currently available timing data for PSR B1534+12, this effect is not yet measurable.

Aside from the ambiguity in $`\mathrm{\Omega }`$ resulting from the unknown direction of the pulsar’s orbital motion, the two values of this parameter obtained from our analysis are very robust in the sense that small changes of $`\mathrm{\Omega }`$ lead to the correspondingly large changes of the shape of the $`V_{ISS}`$ curve (Fig. 3). On the other hand, the screen location is less well determined, because of the uncertain estimate of the pulsar distance. However, if the GR derived distance to PSR B1534+12 is reliable, it is probably safe to assume that the screen is indeed located closer to the pulsar than to the observer.

We would like to thank Ingrid Stairs for her continuing help with the PSPM observations of PSR B1534+12 and Maciej Konacki for insightful discussions. This work was supported by a National Science Foundation grant AST-9988217 (SB, AW) and a Polish Committee for Scientific Research grant 2-P03D-006-16 (WL, MP, AW). The Arecibo Observatory is part of the National Astronomy and Ionosphere Center, which is operated by Cornell University under a cooperative agreement with the NSF.