Generalized Langevin equation with colored noise description of the stochastic oscillations of accretion disks

We consider a description of the stochastic oscillations of the general relativistic accretion disks around compact astrophysical objects interacting with their external medium based on a generalized Langevin equation with colored noise, which accounts for the general memory and retarded effects of the frictional force, and on the fluctuation-dissipation theorem. The presence of the memory effects influences the response of the disk to external random interactions, and modifies the dynamical behavior of the disk, as well as the energy dissipation processes. The generalized Langevin equation of the motion of the disk in the vertical direction is studied numerically, and the vertical displacements, velocities and luminosities of the stochastically perturbed disks are explicitly obtained for both the Schwarzschild and the Kerr cases. The Power Spectral Distribution (PSD) of the disk luminosity is also obtained. As a possible astrophysical application of the formalism we investigate the possibility that the Intra Day Variability (IDV) of the Active Galactic Nuclei (AGN) may be due to the stochastic disk instabilities. The perturbations due to colored/nontrivially correlated noise induce a complicated disk dynamics, which could explain some astrophysical observational features related to disk variability.


I. INTRODUCTION
In an astrophysical environment, compact general relativistic objects like black holes or neutron stars are often surrounded by an accretion disk. The physics of accretion disks around compact objects is thought to explain a number of observations related to X-ray binaries or Active Galactic Nuclei (AGN). The disks are usually considered as being composed from massive test particles that move in the gravitational field of the central compact object. Numerical simulations have shown that the accretion induced collapse of a white dwarf may lead to a neutron star surrounded by a disk with mass up to 0.8M ⊙ [1]. Waves and normal-mode oscillations in geometrically thin and thick disks around compact objects have been studied extensively both within Newtonian gravity, and within a relativistic framework (see [2] for a review of accretion disk properties). It was also suggested that the Intra Day Variability (IDV) of AGNs can also be related to the vertical oscillations of the accretion disks around black holes [3,4].
Stochastic processes and methods play a fundamental role in our understanding of many physical/astrophysical phenomena [5][6][7]. The study of the general relativistic oscillations of thin accretion disks around compact astrophysical objects interacting with the surrounding medium through non-gravitational forces was considered * Electronic address: t.harko@ucl.ac.uk † Electronic address: csinleung@polyu.edu.hk ‡ Electronic address: gabriela.mocanu@ubbcluj.ro in [8]. The interaction of the accretion disk with the external medium was modelled via a friction force and a random force, respectively. The general equations describing the stochastically perturbed disks were derived by considering the perturbations of trajectories of the test particles in equatorial orbits, assumed to move along the geodesic lines. By taking into account the presence of a viscous dissipation and of a stochastic force the dynamics of the stochastically perturbed disks can be formulated in terms of a general relativistic Langevin equation. The stochastic energy transport equation was also obtained. The vertical oscillations of the disks in the Schwarzschild and Kerr geometries were considered in detail, and the vertical displacements, velocities and luminosities of the stochastically perturbed disks were explicitly obtained for both the Schwarzschild and the Kerr cases by numerically integrating the corresponding Langevin equations.
The Langevin equation, used to describe the stochastic oscillations of the accretion disks in [8], provides a correct phenomenological and statistical description of the Brownian motion only in the long time limit, for times long as compared to the characteristic relaxation time of the velocity autocorrelation function [9,10]. In order to describe the dynamics of a homogeneous system without restriction on a time scale, one should generalize the Langevin equation by introducing, instead of the simple friction term, a systematic force term with an integral kernel [9,10]. The convolution term expresses the memory, or retardation effect, which can make the accretion disk response to an external perturbing stochastic force quite complicated. As shown in [8], a simple Brownian motion model for the disk oscillations produces a Power Spectrum Distribution (PSD) P (f ) ∼ f −α , with α = −2. Completely uncorrelated interaction of the disk with the external medium would produce white noise, with a PSD P (f ) = f 0 = const. However, the observational data for a large percentage of AGN show that their optical IDV has PSDs for which α is neither 0 nor −2, as shown by discrete Fourier transform analysis [11], structure function analysis [12] and fractal dimension analysis [13]. Therefore, if this type of variability is produced and has its origin in the accretion disks rotating around compact astrophysical objects, a consistent mathematical description of the stochastic phenomena in accretion disks would require to go beyond the simple description of the complex accretion disk -cosmic environment interaction in terms of a simple friction and a stochastic force, respectively. The Langevin type equation with a damping term and stochastic force was used in [14] to describe the stochastic oscillations on the vertical direction of the accretion disks around a black hole, and to calculate the luminosity and PSD for the oscillating disk. The stochastic resonance phenomenon in PSD curves for different parameter values of viscosity coefficient, accretion rate, mass of black hole and outer radius of the disk was also studied. The results show that the simulated PSD curves of luminosity for disk oscillation have the same profile as the observed PSD of black hole X-ray binaries, and the stochastic resonance of accretion disk oscillation may be an alternative interpretation of the persistent low-frequency quasi-periodic oscillations.
A description of the stochastic oscillations of the general relativistic accretion disks around compact astrophysical objects based on the generalized Langevin equation, which accounts for the general retarded effects of the frictional force, and on the fluctuation-dissipation theorems, was introduced in [15]. The vertical displacements, velocities and luminosities of the stochastically perturbed disks were explicitly obtained for both the Schwarzschild and the Kerr cases. The PSD of the simulated light curves was determined and it was found that the spectral slope α has values that correspond with observations. The theoretical predictions of the model were compared with the observational data for the luminosity time variation of the BL Lac S5 0716+714 object. The influences of the friction parameter, spin parameter and mass M of the central compact object on the stochastic resonance in PSD curves was discussed, by using the generalized Langevin equation, in [16]. The results show that a large spin parameter can enhance the stochastic resonance phenomenon, but the larger the friction coefficient or the central mass is, the weaker the stochastic resonance phenomenon becomes. The simulated PSD curves of the output luminosity of stochastically oscillating disks have the same profile as the observed PSDs of X-ray binaries. Hence the resonance peak in the PSD curve can explain the quasi-periodic oscillations.
It is the purpose of the present paper to consider the memory/retardation effects in the description of the stochastic oscillations of the particles in the accretion disks around compact objects under the influence of the external environment, generating some random interactions between the disk and the cosmic medium. More specifically, we consider two astrophysical accretion disk -external medium interaction models that can be described by means of the Langevin type stochastic differential equations. If the black hole -accretion disk system is located in a dense cosmic environment, like, for example, a stellar cluster, then the gravitational force exerted by the neighbouring stars can be described as an effective stochastic, time dependent force [5,6]. The stochastic gravitational force appears due to the local fluctuations of the numbers of stars in the clusters, as shown in [5,6], where these processes have been studied in much detail. Therefore for a massive central object -accretion disk system located in a dense stellar cluster, the accretion disk is subjected to two gravitational forces, one produced by the black hole, and a second one representing an effective random force due to the neighboring stars. In this case the dynamical evolution of the accretion disk can be described by means of the Langevin equation.
As a second astrophysical process that can be modelled by a Langevin type stochastic differential equation we consider the dynamical interaction between the accretion disk and the remnant matter from a star tidally disrupted by the central black hole. In this case, by assuming that the mass of the disk is much smaller than the mass of the inflowing gas, one can model the system as a single particle (the disk) located into the hydrodynamical flow of the remnant gas. Under these simplifying assumptions the dynamics of the system can be described again by a stochastic Langevin type equation.
Moreover, in the present analysis we assume that the stochastic perturbations of the disk (the particle), assumed to randomly interact with an external medium (fluctuating gravitational force and gas remnants of a tidally disrupted star), are described by the generalized Langevin equation, introduced in [17].
Mathematically, the stochastic equation of motion of the disk becomes a differential -integral equation, whose integral kernel describes the memory effects. We use this equation to provide a detailed description of the vertical stochastic oscillations of an accretion disk around a Schwarzschild black hole and also of an accretion disk around a Kerr black hole. The solution of the generalized Langevin equation is obtained numerically, and the basic physical parameters of the disk (vertical displacement, velocity and luminosity) are obtained in a numerical form. As a possible application of the present formalism we investigate the possibility that the IDV of the AGNs has its origin in some random/stochastic effects at the accretion disk level. In the present study we will focus only on supermassive black holes, both static and rotating, which can be described by the Schwarzschild and Kerr geometries, and which are the best candidates to explain the electromagnetic radiation emission properties of the AGNs, assumed to be powered by accretion of mass onto massive black holes, with masses of the order of 10 6 M ⊙ to 10 10 M ⊙ .
The present paper is organized as follows. The generalized Langevin equation, describing the stochastic oscillations of the accretion disks with memory is presented in Section II. Two specific examples of astrophysical interest of disk -external medium interaction that can be described by a Langevin type equation are also considered. The dimensionless form of the disk equation of motion, as well as the numerical procedure for solving the generalized Langevin equation are presented in Section III. The stochastic oscillations of the accretion disks in the presence of colored noise are studied for both Schwarzschild and Kerr geometries in Section IV. The application of the present formalism for the study of the IDV's is considered in Section V. We conclude and discuss our results in Section VI.

II. THE GENERALIZED LANGEVIN EQUATION FOR STOCHASTIC OSCILLATIONS OF ACCRETION DISKS
In the present Section we consider the equation of motion of the stochastically oscillating accretion disks interacting with their external cosmic environment. We consider two specific cases of astrophysical interest, namely, the stochastic perturbations of the accretion disks, located in dense stellar clusters, by the gravitational fields of their neighboring stars, and the interaction between the accretion disk with the gaseous hydrodynamic flow generated by the disruption of a star by the tidal forces of the central supermassive black hole. By assuming that the particles in the disk move on geodesic lines, we show first that the effects of a stochastic perturbation can be described, in the Newtonian approximation, by a Langevin type equation. The Langevin equation with memory is introduced, and the connection between the generalized Langevin equation and the fluctuating hydrodynamics [9,18] is also established.

A. Equation of motion of massive test particles in an accretion disk
We model the accretion disk as a single thermodynamic system, consisting of a large number of particles, with overall energy density ǫ and pressure p, respectively. The energy-momentum tensor T µν of the disk is given by where the four-velocity u µ of the particles with mass m satisfies the normalization condition u µ u µ = 1, and the relation u µ ∇ ν u µ = 0, respectively, where ∇ ν represents the covariant derivative with respect to the metric. We also introduce the projection operator h µλ = g µλ − u µ u λ , which satisfies the relation h µλ u µ = 0. The equations of motion of the particles in the disk follow from the equation of conservation of the energy-momentum tensor, ∇ µ T µ ν = 0, and, after contraction with h µλ u µ , leads to the equations Finally, contraction with g αλ gives the equation of motion of the fluid element of the disk as [19] where Γ µ αβ are the Christoffel symbols associated to the metric, and a comma denotes the partial derivative with respect to the coordinate x ν .
We consider that the particles in the disk are in contact with the external environment (stochastic gravitational field of neighboring compact astrophysical objects, remnant gas from tidally disrupted stars flowing towards the central black hole etc.), whose action on the disk can be described generally by means of a stochastic force ξ µ . If the disk is perturbed, a point particle in the nearby position x ′µ must satisfy the generalized equation of motion where ξ µ is the external, random force, acting on the particle. The coordinates x ′µ are given by x ′µ = x µ + δx µ , where δx µ is a small quantity. Thus, in the first approximation, we obtain for the Christoffel symbols where Γ µ αβ,λ = ∂Γ µ αβ /∂x λ . Moreover, and respectively. In the following we introduce the fourvelocity of the perturbed motion as V µ = dδx µ /ds. V µ satisfies the condition u µ V µ = 0. Therefore the equation of motion for δx µ becomes In obtaining Eq. (9) we have assumed that the matter content of the disk can be modelled as a perfect fluid, and we have neglected the possibility of energy dissipation, either by viscous processes, or by radiation emission, considered to have a negligible effect on the disk structure and dynamics. Moreover, we have adopted a full general relativistic approach, which takes into effect the modifications of the space-time geometry in the disk due to the presence of the central supermassive black hole. For realistic astrophysical systems, the study of Eq. (9) can be done by using numerical methods only. In the following Sections we will obtain some approximate, but still realistic, forms of Eq. (9), which allow a simpler numerical treatment of the problem.

B. The generalized Langevin equation
Since the self-gravity of the disk is negligibly small as compared to the gravitational field of the central supermassive black hole, the infinitesimal transformations x ′µ = x µ + δx µ do not change the metric. Therefore the coordinate transformations satisfy the Killing equations Hence it follows that in this approximation the equation of motion of the disk particles perturbed by an external force is given by In the case of a disk consisting of particles rotating with azimuthal angular velocity Ω = dφ/dt, the components of the four velocity are given by u µ = u t (1, 0, Ω/c, 0). If the disk consists of dust particles, with negligible pressure, moving along the geodesic lines, the term containing the pressure gradient can be neglected, and we obtain the equation of motion considered in [20]. In the adopted approximation of non-zero pressure, the deviation equation contains some new terms, which are proportional to the ratio of the second derivative of the pressure and the sum of the energy density and pressure, respectively, and the product of the pressure gradient and the deviation velocity V µ , respectively.
In the non-relativistic case the friction force is given by f i f r = −νmv i , where ν is the friction coefficient, and v i are the components of the non-relativistic velocity [10]. The relativistic generalization of the friction force requires the introduction of the friction tensor ν µ α , given in [21][22][23] as Therefore the friction force can be expressed as [21][22][23] In the following we will consider only the vertical oscillations of the disk, which for an arbitrary axisymmetric metric are described by the equation where is the frequency of the gravitational oscillations of the disk [8].
In the following we assume that the gravitational effects, which enter in Eq. (9) via the derivatives of the Christoffel symbols Γ µ αβ,λ are much stronger than the hydrodynamic pressure effects, and therefore ω 2 By considering the non-relativistic limit of Eq. (10), the perturbed stochastic geodesic equation describing the vertical oscillations of the disk is equivalent with the Langevin type equation [8] Eq. (15) practically represents a linear oscillator, with a viscosity term, and a random perturbing force. In obtaining Eq. (15) we have assumed that the external perturbative force is a function of time only, and we have neglected its possible spatial coordinate dependence. A fundamental assumption of non-relativistic Galilean physics is the existence of a universal time t. Therefore, within the non-relativistic Langevin theory, it is quite natural to identify this universal time t with the time parameter of the stochastic driving process [23]. In general, the external stochastic noise ξ z can be represented as a function of a t-parameterized two dimensional stochastic process, ξ z = ξ z [x µ (t), p µ (t)], where p µ is the fourmomentum of the external source. By assuming that the motion is non-relativistic, with v/c << 1, and | p| << mc, where | p| is the absolute value of the momentum, we have E ≈ mc 2 ≈ constant, and consequently the dependence on the four-momentum of the stochastic force can be neglected. For a spatially inhomogeneous external perturbing cosmic environment, the function ξ z would also depend on z. However, if the geometric dimension of the external stochastic source is much larger than the size of the disk, the inhomogeneities in the generated stochastic force can be neglected. Hence, in the following the external noise ξ z is modeled by the 1-dimensional standard homogeneous Wiener process, with any possible momentum or space-like coordinate dependence neglected.
The Langevin equation Eq. (15) provides a correct phenomenological and statistical description of the Brownian motion only in the long time limit, for times long as compared to the characteristic relaxation time of the velocity autocorrelation function [9,10]. In order to describe the dynamics of a homogeneous system without restriction on a time scale, one should generalize the Langevin equation by introducing, instead of the simple friction term, a systematic force term with an integral kernel [9,10]. The convolution term expresses the memory, or the retardation effects, and, in the specific case of the stochastic oscillations of the accretion disks, the effects of the non-zero disk pressure and energy density.
In the case of the vertical oscillations of a thin accretion disk, the generalized Langevin equation, introduced in [17], which takes into account the retarded effects of the frictional force, is given by and respectively. In Eq. (17) τ d represents a characteristic disk memory time, while α 1 is interpreted as a disk damp- From a physical point of view the coefficient α 1 describes the strength of energy dissipation in the disk. The autocorrelation function of the stochastic force satisfies the condition where β is a constant. It is worthwhile to summarize now the physical assumptions, implicitly assumed by Eqs. (16), (17), and (19): a) The external cosmic environment interacting with the accretion disk around a supermassive black hole is spatially homogeneous and stationary; i.e., relaxation processes within the external environment occur on time scales much shorter than the relevant dynamical time scales associated with the motion of the accretion disk. Interaction with a spatially inhomogeneous nonstationary cosmic environment can be modelled by considering friction and noise amplitude functions depending on the spatial coordinates, and disk momenta. b) On a macroscopic level, the interaction between the disk particles and the external cosmic environment is sufficiently well described by the generalized friction coefficient with memory γ(|t−τ |), and the stochastic Langevin force ξ z (t). c) Stochastic collisions between the disk particles and the constituents of the external environment occur virtually uncorrelated.
In the next Section we will use Eq. (16) to study the vertical oscillations of the accretion disks around supermassive black holes in both Schwarzschild and Kerr geometries.
C. Equation of motion of a thin disk in the presence of a stochastic gravitational field As a first example of the stochastic interaction between accretion disks around supermassive black holes and the cosmic environment we consider the situation in which the disk, and its central object, is a member of a larger stellar system, in which all components interact gravitationally. We assume that at the center of the disk we have a compact object of mass M . The force acting on the disk can be divided into two components, the gravitational force of the central object, which we denote by F G , and the force F R generated by the influence of the immediate local neighborhood, which may consist of ordinary stars and black holes. The force F G is a smoothly varying function of position, while F R can be subject to relatively rapid fluctuations [5,6]. The fluctuation in the external force is a direct consequence of the changing of the local stellar distribution that makes the influence of the near neighbors on the disk time variable. It can be shown [5,6] that the fluctuations in the force acting on the central object disk system, due to the changing local stellar distribution, do occur with extreme rapidity as compared to the rate at which any of the other physical parameters change. Let p be the pressure of the matter in the disk. Therefore we can write for the force acting on the disk the expression where F p is the force due to the pressure variation in the disk, F visc is the force generated by the viscous dissipative effects in the disk, and the fluctuating random force F R , due to the gravitational interaction of the disk with the nearby stellar systems, can be considered as a function of time only [5,6].
In order to obtain the equation of motion of the disk in the presence of the gravitational force of the central object, and of the external stochastic forces, we consider that the disk is thin in the sense discussed in [24,25]. We approximate the distribution of the surface density Σ by the expressions [26] where R in is the innermost radius of the disk, R adj is an adjustment radius in the disk, and R out is the outer radius of the disk. The index χ of the surface density can be either 3/5 or 3/4 [26]. The mass of the disk M d from R in to R out is given by where we have used the assumptions R in << R out and R adj ≈ R out , respectively. We assume that the disk as a whole is perturbed, and that it deviates from the equatorial plane by a small distance δz. The restoring gravitational force F G caused by the attraction of the central object is [26] F G (δz) = GM δz Rout Rin 2πrΣ (r) After integration, using the density distribution given by Eqs. (20), we obtain We define the mean vertical density of the disk as Therefore the vertical oscillations of the disk as a whole can be described by the equation of motion As for the acceleration due to the vertical pressure distribution we can approximate it as a press = − (1/ρ d ) ∂p/∂z ≈ −κ z (r), where the vertical epicyclic frequency κ z (r) depends on the horizontal radius r only.
The simplest way to obtain such a relation is by considering a power series expansion of the pressure so that The constant term in the acceleration can be eliminated by a rescaling of the acceleration, and the epicyclic frequency can be obtained as κ z (r) = ∂ 2 p(r, z)/∂z 2 z=0 /ρ d . Therefore for the equation of motion of the disk we obtain where If some dissipative forces, like, for example, dissipation due to radiation emission, or viscosity friction effects in the disk matter, are also acting on the disk, one should include the dissipative forces in the equation of motion. The simplest phenomenological approach would be to assume that the viscous/heat dissipation is proportional to the velocity of the vertical oscillations, F visc ∝ δz/dt. However, a better phenomenological description of the viscous friction can be obtained by assuming an integral kernel description of the dissipative processes. Therefore the equation of motion of the disk becomes In this model the disk as a whole undergoes stochastic oscillations determined by the interaction between the disk, and the fluctuating stellar/black hole neighborhood. The equation of motion is given by the generalized Langevin Eq. (28).
In order to obtain the oscillations frequency of the disk we assume that the innermost radius of the disk R in is located at R in = 6GM/c 2 , corresponding to the innermost stable orbits in a general relativistic Schwarzschild geometry. We represent the adjustment radius of the disk R adj as R adj = f 1 R in and the outer radius of the disk as R out = g 1 R in , where f 1 > 1 and g 1 >> 1 are constants. Then for the oscillation frequency of the disk we obtain or The period of the free disk oscillations is given by In the present paper we will consider only the case ω 2 0 >> κ 2 z , and therefore our study will be limited to the consideration of the gravitational effects in the disk.
D. Stochastic equation of motion of accretion disks in a fluctuating environment due to stellar capture Close encounters between stars from a nuclear cluster and a massive black hole lead to the disruption of stellar bodies due to the gravitational tidal forces of the black holes, and the subsequent accretion of the remnant gas [27,28]. These astrophysical processes are among the likely mechanisms for feeding black holes that are embedded in a dense stellar cluster. Stars on highly elongated orbits are susceptible to tidal disruption and hence they provide a natural source of material to replenish the inner disk. Such events can produce debris and lead to recurring episodes of enhanced accretion activity, with the infalling matter interacting both gravitationally and mechanically with the accretion disk. Intermediate-mass black holes have a high chance of capturing stars through tidal energy dissipation within a few core relaxation times [27]. The interaction between the accretion disk and the remnant gas falling towards the central black hole can be considered as a random process, in which the particles of the disk interact with the particles of the tidally disrupted star. In order to model the interaction between the disk and the remnant matter from the disrupted star we assume that the disk is thin, and its mass M d is much smaller than the mass M F falling towards the black hole, M d << M F , condition valid at least for finite time interval 0 < t < t F , where t F is the in-falling time of M F . Hence, if this condition is satisfied, one can obtain a simplified description of the disk -remnant stellar matter interaction by considering the disk as a single physical system randomly moving in a fluid represented by the matter accreted by the black hole. This simple physical description allows the description of the dynamical evolution of the disk in terms of a Langevin type equation, obtained from the fluctuating hydrodynamic equations of motion of the interacting particle -fluid system. Therefore in the following we consider the accretion disk as a macroscopic system of mass M d in an incompressible fluctuating fluid (the debris from stellar disruption). The position of a point on the surface of the accretion disk at time t is described through the position vector r 1 (t). We denote by r the position vector of a point in the fluid and by R the position vector of the center of the accretion disk. The motion of the fluid (the gaseous remnant of the star) in the presence of an exterior (gravitational) force F ext is described by the linearized stochastic Landau-Lifshitz equations of motion [18,[29][30][31], where ↔ P ( r, t) is the total fluid stress tensor, with components where p is the pressure in the fluid, η is the shear viscosity coefficient of the matter falling towards the black hole, and ↔ Π is the random stress tensor, with components Π ij , which has the following stochastic properties [18,29] where the square brackets denote averages over an equilibrium ensemble, k B is the Boltzmann constant and Θ is the temperature of the external cosmic environment. We have also denoted ∇· ↔ P ( r, t) = ∂P ij /∂x j . We assume that the accreted gas is incompressible, and therefore its velocity satisfies the condition The motion of the accretion disk, modelled as a single physical system, immersed in the fluid accreted by the central black hole in the presence of an external gravitational force F ext ( r, t) is governed by the equation [9,18,[29][30][31] where ↔ P ( r, t)· n = P ij n j , and n is the unit vector perpendicular to the thin accretion disk, and S(t) is the surface area of the accretion disk. Eqs. (32) and (37) simply express the conservation of the momentum. In the fully linearized scheme one may neglect the time dependence of R(t) and S(t) in these equations, and take the origin of the coordinates in the center of the accretion disk, so that R(t) = 0. The stochastic boundary value problem associated with the motion of the accretion disk of radius r 1 (t) in the fluid accreted by the black hole is specified by the boundary condition where r 1 (t) is the position vector of a point on the surface of the accretion disk at time t having the velocity V (t). r 1 (t) is given by r 1 (t) = r 1 (0)+ t 0 V (t ′ ) dt ′ . The components of the force acting on the surface of the accretion disk due to the fluid falling towards the black hole are given by [9,18,[29][30][31] By estimating the perturbations of the velocities in the hydrodynamic flow of the gas remnant due to the presence of the accretion disk it follows that the equation of motion of the perturbed accretion disk can be written in the form of a generalized Langevin equation of the form [31] where the function γ (t − t ′ ) represents a retarded effect of the frictional force and ξ(t) is the random force generated by the stochastic fluctuations in the fluid accreted by the black hole. Assuming that the flow is unsteady, one should also take into account the so-called retarded viscous force effect, which is due to an additional term to the Stokes drag, related to the history of the particle acceleration. This additional drag force, nowadays referred to as the Basset history force, has the form [9] where r p is the particle radius. However, in the present paper we ignore the effects of the Basset history force on the accretion disk dynamics.

III. THE NUMERICAL APPROACH TO THE GENERALIZED LANGEVIN EQUATION
In the present Section we write down the equation of motion of the stochastically oscillating accretion disks in a dimensionless form. The numerical procedure to solve the Langevin equation in the presence of a colored noise is also described. To characterize the luminosity of the stochastically oscillating disks from a global point of view we will use the Power Spectral Distribution (PSD), whose main properties are briefly reviewed.

A. Dimensionless form of the generalized Langevin equation
By introducing a set of dimensionless parameters (θ, σ, q) defined as where M D is the mass of the disk, the equation of motion of the stochastically oscillating disk given by Eq. (16) can be written as where we have denoted and respectively, where M is the mass of the central compact object.
For the luminosity calculation, the energy is the sum of kinetic plus potential energy and the output luminosity is the time variation of the energy The luminosity, written in dimensionless form as has an evolution given by B. The numerical procedure To solve the dimensionless generalized Langevin Eq. (43), the algorithms already developed in the literature, see e.g. [32], were used to obtain a set of three differential equations as followṡ where the noise Υ has the property and the initial value of z is drawn from a distribution with second moment The dimensionless form of the disk luminosity is given by The set of the stochastic linear differential equations (50)-(52) is solved by implementing the fourth order Runge-Kutta integrator developed in [32]. Time is discretized with timestep h. At each timestep j + 1 the values of the variables in the set T j+1 = {q, p, z} are calculated as i.e., as the sum of a deterministic part and a random part, where both parts depend on the full values of all the parameters in the set at the previous timestep. The deterministic part follows by evolving the "normal" Runge-Kutta algorithm in time and has the form where D is a matrix with components and The random part is calculated as detailed in Section II.B of [32], and for the case studied here has the form where the Z i 's are linear combinations of four independent Gaussian variables {a 1 , a 2 , a 3 , a 4 } ∈ N (0, h) (71)

C. The Power Spectral Distribution
As a parameter allowing a global characterization of the stochastic behavior of the disk we use the slope of the Power Spectral Distribution (PSD) of the physical parameters of the disk. The numerical values of the slopes of the PSD provide insight to the nature of the mechanism leading to the observed variability. In statistical analysis, if X is some fluctuating stationary quantity, with mean µ X and variance σ 2 X , then an autocorrelation function for X is defined as [47,48] where X s is the values of X measured at time s and denotes averaging over all values s. The PSD is defined based on the correlation function as [47,48] and it is straightforward to see its importance in terms of the "memory" of a given process. For example, if X is the B band (optical) magnitude of the disk, the slope of the PSD of a time series of X provides insight to the degree of correlation the underlying physical process has with itself. The system needs additional energy to fluctuate and this mechanism is historically best explained for Brownian motion, in which case the energy is thermal. As an example of the connection between the process generating the light curve and the observed spectral slope, a Brownian motion process would produce a PSD of the shape P (f ) ∼ f −2 while a completely uncorrelated evolution of a system would produce white noise, i.e., a PSD of the shape P (f ) = f 0 = const. From this it is important to remember that the slope of the PSD can be regarded as an indicator of the type of underlying process generating the light curves. As already stated in Section I, the observed IDV light curves have spectral slopes that do not fit in the simple models of Brownian motion or uncorrelated evolution. However, in the case of experimental curves, one cannot know before-hand if the process is purely stationary. The experimental and simulated light curves are analysed in this paper with the software .R which produces spectral slopes in accordance with the algorithm described in [36].
In order to trust an interpretation that the spectral slope represents a sign of correlation in the signal, the investigated signal should be stationary. By definition, a stationary stochastic process is a process whose properties are constant in time. More precisely, the mean value and the variance of such a signal are time independent and the autocorrelation of the signal at different times depends only on the time difference and not on the actual time at which the correlation is evaluated.
In our model, the perturbing noise is a stationary stochastic process, as can be easily seen because of its prescribed property Eq. (56). To show that the resulting light curveL = −pz is stationary, one should look at Equations (56)-(71). The stochastic part of the signals involved are given by linear combinations of Gaussian random variables (the Z i s), so their properties are stationary. We thus expect that our interpretation of the spectral slope is somewhat accurate.

IV. COLORED STOCHASTIC OSCILLATIONS OF THE GENERAL RELATIVISTIC DISKS
In the present Section we will obtain the basic physical parameters (displacement, velocity, and luminosity), of the vertically oscillating accretion disks in both static Schwarzschild geometry, and in the rotating Kerr geometry, respectively.
(74) For the Schwarzschild solution the functions U and ζ are given by where M is the mass of the central compact object in geometrical units [33]. The radial frequency at infinity Ω 2 ⊥ of the particles in the disk, and the proper free oscillation frequency ω 2 ⊥ of the disk can be obtained as [8,34,35] and respectively. In the case of the static Schwarzschild geometry, the free vertical oscillation frequency of general relativistic disks can be written as [8,35] where M is expressed in natural units as M = GM/c 2 , and ρ is the distance from the black hole center in the radial plan of the disk. If we express ρ in terms of the mass of the central object as ρ = nM , then For M = 10 8 M ⊙ , n = 10, ω 2 ⊥ = 0.001 and fixed values of the integration step h = 1, dimensionless initial displacement q(0) = 0, dimensionless velocity p(0) = 0.01 and b = 4, we analysed two cases. All trajectories are obtained by mediating over 1000 stochastic trajectories (for which the variable z 0 is redrawn from its distribution).
First for fixed A = 0.0015 we varied a and obtained the behaviour of the dimensionless displacement, velocity and luminosity. Then for fixed a = 0.03 we varied A and obtained the behaviour of the same physical parameters. The results of the numerical simulations of the displacements, velocities and luminosities for the stochastic oscillations of the accretion disks are presented in Figs. 1, 2, and 3, respectively.
The PSD of luminosity for the cases a = 0.03, A = 0.0015 and a = 0.03 and A = 1.5 are shown in Fig. 4, respectively. The PSD curves were obtained by using the .R software [36].

B. Stochastic disk oscillations in the presence of colored noise in the Kerr geometry
The line element of a stationary axisymmetric spacetime is given by the Weyl-Lewis-Papapetrou metric [37] where ψ i , i = 1, 2, 3 are functions of z and ρ only. The metric functions generating the Kerr black hole solution in axisymmetric form are given by By defining ρ = nM as the radial distance from the origin in the plane of the disk, and the angular momentum of the compact object as J = kM , and denoting δ = M √ n 2 + 1 − k 2 , with n > 0 and k ∈ [0, 1], the frequency of the free vertical oscillations of the general relativistic disks in the Kerr geometry can be represented as [8,38] where ω 2 1 (n, k) = 3n 4 + n 2 √ 1 + δ × 2k (5 + 3δ) + 3 √ 1 + δ 6 + 3δ + δ 2 − and respectively. By adopting for the mass of the central object a value of M = 10 8 M ⊙ , and by fixing n = 10, k = 0.9, and ω 2 ⊥ = 0.0012, respectively, and considering for the step and the initial conditions the values h = 1, q(0) = 0, p(0) = 0.01, b = 4, we have analysed two cases. All trajectories are obtained by mediating over 1000 stochastic trajectories (for which the variable z 0 is redrawn from its distribution).
For fixed A = 0.0015 we varied a and obtained the behaviour of the dimensionless displacement, dimensionless velocity, and dimensionless luminosity. Then, for The PSDs of the luminosity for the cases a = 0.03, A = 0.0015 and a = 0.03, A = 1.5, respectively, are shown in Fig. 8. Values for the spectral slope obtained for a = 0.03, A = 0.0015 and different configurations of the system producing the simulated light curves are given in Table I. The fourth column in Table I is the Bayesian probability that the source actually behaves like in the attempted fit, where p B close to 1 suggests that the fit is correct. The definition and calculation of p B is quite lengthy, but we sketch it here and refer the reader to [36] for full detail. In cases such as investigating data from an AGN, when there is only one light curve available (call it x obs ), where reproducibility of the experiment is not an option, the .R software with the Bayes script offers the possibility to asses whether or not the recorded data was produced by a process behaving like a null hypothesis, H 0 . For this purpose, the script simulates numerous sets of data based on the statistical properties of x obs , stores them in x rep and treats these new data sets as they would be multiple measurements of the same process. Afterwards it defines Ξ( x), called a test statistics, which is a real valued function of the data chosen such that extreme values are unlikely when the null hypothesis H 0 is true. With these parameters, the Bayesian probability is defined as where Pr = {x|y} means probability of event x given that event y has occurred.

V. COMPARISON WITH AGN IDV DATA
Extensive observational and theoretical efforts have been made in order to explain IntraDay Variability (IDV) in some classes of Active Galactic Nuclei (AGN) [3,[39][40][41][42]. IDV manifests as the fast change (in less than one day) in the luminosity output of an object, and for supermassive black holes IDV occurs in the optical domain. The PSD of these light curves are found to be nontrivial [11][12][13].
In the present paper we assume that the IDV emission occurs from the disk of a super-massive black hole -accretion disk system [3,4], and it is due to the oscillations of the disk perturbed by a sudden interaction with the external environment. Such a transient interaction can be due to the stochastic gravitational force generated by the fluctuations of the star number in the near neighborhood of the black hole -accretion disk system, if the system is located in a dense stellar cluster, or to the interaction of the disk with the debris produced by the tidal disruption of a star captured by the central black hole.
There are indications that the relationship between the root-mean-square-deviation (rms) and the flux provides more information regarding the source-process of fast variability, at least in X-Ray Binaries [45]. For these objects, for which IDV is exhibited in the X-Ray domain, it was found that this relationship is linear. However, a test performed on optical IDV for BL Lac S5 0716+714 on observational data from different epochs showed that the rms-flux relation was not linear [46].
In order to test if the optical IDV flux data of BL Lac S5 0716+714 are log-normally distributed, in Fig. 9 we compare the observational data with the normal distribution. The bins represent the logarithm of the observational data, and the smooth line is a superimposed Gaussian with the same mean and variance as the data. If the two would have coincided (i.e. if the flux would have been log-normally distributed), the rms-flux relation would have been linear. We find that the simulated data, obtained from the generalized Langevin equation description of the luminosity of the stochastically oscillating disks, also exhibit this feature, i.e., that the rmsflux relation is not linear, for neither of the Schwarzschild or the Kerr case (an example is shown in Figs. 10).
The model presented in this work is thus consistent with observational data, from the point of the statistics exhibited by the light curves. This is very encouraging and supports the idea that perturbations due to colored/nontrivially correlated noise consistently explain the statistical properties of the observed IDV light curves.  [41]. If the bins (the distribution of the data) can be fitted by the continuous curve (a normal distribution with the same mean and variance as the data), the distribution of the data is log-normal. The data was recorded in the B filter, on JD 2454824.
accretion disks based on the generalized Langevin equation with colored noise and memory effects. From an astrophysical point of view this approach may be relevant for modeling of accretion disks around compact general relativistic objects randomly interacting with the cosmic environment surrounding the disk. Events like the encounter of a star, or a massive compact object, by the central black hole, or the interaction between the disk and the cosmic dust resulted from the tidal disruption of a star captured by the central black hole, may trigger such interactions, which are essentially random in their nature. In our description we have assumed that the oscillation velocities of the massive test particles composing the disk are small, and therefore they can be described by the Newtonian equation of motion. Hence in our approach we have assumed that the motion of the disk is non-relativistic in the sense that the velocity δz/dt of the particles is much smaller than the speed of light, δz/dt << c. However, since the gravitational field in the vicinity of the black hole is strong, we have adopted for the oscillation frequency of the disk in the gravitational field the general relativistic values. The cases of the static and rotating black holes were considered independently, and the effect of the rotation on the disk luminosity was explicitly considered. The analysis of the PSD of the luminosity curves has shown that in the presence of memory effects and colored noise non-trivial values can be obtained for the spectral slope. Therefore the presence of the colored noise/memory effects could allow a better description of the astrophysical processes involving random/stochastic factors.
In the present paper we have considered that accretion disks are composed from massive test particles forming a fluid system, and moving on non -geodesic lines around the central black hole. By considering the equation of motions of the disk perturbations, and assuming that the gravitational effects dominates the hydrodynamic ones, in the first order of approximation the effects of the nonzero pressure and energy density give in the equation of motion a correction term proportional to the magnitude of the vertical perturbations δz. In the non-relativistic limit the disk pressure effects and the external gravitational perturbation terms can be combined in a single, effective coefficient, proportional to the vertical displacement, and which determines the frequency of the disk oscillations. In order to construct realistic physical models one should assume an equation of state for the disk matter components, and adopt some particular disk structures. Since there is no general consensus on the physical properties of the accretion disks, there is no guarantee that adopting a particular disk model would produce observationally relevant results. On the other hand, in the presence of a strong gravitational field, the gravitational effects dominate the pressure effects, which can be neglected with a very good approximation. In this case the motion of the particles in the disk can be considered as geodesic. This approximation is certainly valid for particles located in the inner region of the disk. As for the viscous dissipative effects, in the present paper we have assumed that they can be described by a single integral kernel term, whose functional form has been fixed phenomenologically, so that the fluctuation -dissipation theorems can be satisfied. However, some essential features of the accretion disks dynamics can be obtained, and analyzed, even by the present simplified theoretical model.
As a possible astrophysical application of our model we have compared a set of IDV observational light curve data with the results of the simulations of the vertical stochastic oscillations of the accretion disks, described by the generalized Langevin equation. Some common features of the theoretical model and observations, namely, the possibility of obtaining non-trivial PSD values, and the deviation in both cases of data from a log-normal distribution, were identified. Hence the generalized Langevin equation with colored noise could provide a powerful tool in the description of the physical processes in the vicinity of compact general relativistic objects.