© ESO, 2015

1. Introduction

For bright galaxies, Dressler (1980) showed that the population fraction of elliptical galaxies in clusters increases with increasing local environmental galaxy density, while at the same time the fraction of spiral galaxies decreases. An analogous morphology-density relation is observed for dwarf galaxy types (Binggeli et al., 1987). Dwarf irregulars and low-mass spirals dominate outside of dense environments, while early-type dwarf galaxies (dE, subsuming objects termed dwarf lenticulars) are the most common type of galaxy in clusters. However, the origin of this large population of dEs is still a topic of ongoing debate, partly because they were found to follow a number of continuous relations in structure and colour with bright elliptical and/or lenticular galaxies (e.g. Côté et al., 2007; Chen et al., 2010; Kormendy & Bender, 2012), but also since their identification and detailed analysis are not yet possible at higher shifts (see Barazza et al. 2009 for dEs in a z = 0.165 cluster). Clues to the origin of dEs may be given by the presence of weak disk features in a significant fraction of them (Jerjen et al., 2000; De Rijcke et al., 2003; Lisker et al., 2006b), the presence of residual star formation activity in dEs outside of the densest environments (Lisker et al., 2006a; Tully & Trentham, 2008; Pak et al., 2014) and the finding that these have significantly flatter shapes than dEs with low relative velocities in the Virgo cluster core (Lisker et al., 2009). Taken together, this may be regarded as supporting the hypothesis that many dEs were formed from disk galaxies that were affected and eventually transformed by a dense environment (for further observational aspects related to this hypothesis, see also Kormendy, 1985; Binggeli et al., 1990; De Rijcke et al., 2005; Beasley et al., 2006; Smith et al., 2012; Ryś et al., 2014). In this scenario, dEs with disk features (dEdis) could be galaxies in an incomplete transformation phase (Graham et al., 2003; Mastropietro et al., 2005).

Different processes have been discussed to explain the morphological transformation. In general, they can be divided into two main types based on the underlying physics: hydrodynamical scenarios and tidal scenarios. In reality, a mix of both types will act on an infalling galaxy, but it is still an open and relevant question how strong the impact is for which kind of interaction. One of the most popular tidal scenarios is galaxy harassment Moore et al. (1996). It is characterised as the combined effect of tidal interactions between an infalling galaxy and the global cluster potential and of the tidal interaction with close high speed encounters with other cluster members. Moore et al. (1996) used a set of N-body simulations to show that these tidal interactions are suitable for causing a morphological transformation of spiral galaxies. This outcome is generally supported by studies of various authors, even though debating the actual efficiency (Moore et al., 1998, 1999; Mastropietro et al., 2005; Aguerri & González-García, 2009; Smith et al., 2010).

The efficiency of the harassment scenario is expected to depend on the orbit of an infalling galaxy within a galaxy cluster (Smith et al., 2010). The deeper an orbit is placed inside a cluster, the stronger the tidal force of the background mass of the cluster and the higher the probability of close encounters with other cluster members. Mastropietro et al. (2005) show in their simulations of bulgeless disk galaxies with masses of 7 × 10¹⁰M_⊙ that combined tidal interactions are efficient at harassing galaxies in the innermost regions of a cluster on orbits with pericentres very close to the cluster centre (between 40 kpc and 100 kpc) and less efficient on orbits in the outer regions of a cluster. Smith et al. (2010) show that late-type dwarf irregular galaxies with masses of 10¹⁰M_⊙ that enter a cluster from its virial radius are only slightly affected and not harassed, but they find that dwarfs on orbits very close to the cluster centre comparable to those of Mastropietro et al. (2005) could be harassed efficiently. Both studies have shown what happens to late-type galaxies of different mass ranges on orbits very close to the cluster centre and far away from it. But especially on more moderate orbits, the question is still open whether galaxies could be efficiently harassed or not. With the study presented here, we intend to add a piece to that puzzle.

Complementary to the orbital parameters, Moore et al. (1999) show that the properties of the progenitor galaxies are relevant for the efficiency of transformation by tidal interactions. They show that compact high surface brightness galaxies (HSB) are hard to affect, and more extended low surface brightness galaxies (LSB) can be affected efficiently. A comparative result was shown by Chang et al. (2013) by probing the efficiency of tidal stripping for different morphologies of satellite galaxies in host galaxy of Milky Way-like mass. In their simulations of disk galaxies orbiting inside dark matter haloes with group masses, Villalobos et al. (2012) have additionally found that the inclination between the orbital plane and the disk of the infalling galaxy is an important parameter for the effect of tidal stripping. However, tidal stripping is only the most extreme effect of tidal interactions, and one has to expect that an inclination of the infalling galaxy may affect the reshaping of the galaxy.

To sum up, the parameter space for the efficiency of the tidally induced transformation in a galaxy cluster spans the orbital parameters, the type of progenitor, and the orientation of the galactic disk. Furthermore, the properties of the cluster itself, such as its mass, its galaxy density, and the velocity dispersion of its galaxies, are expected to be relevant. Some parameters are partly degenerated; e.g., a variation in the orbit and the cluster properties could end up with a similar tidal strength. To determine the contribution of harassed galaxies to the population of dEs in clusters, it is necessary to understand the efficiency of harassment depending on these parameters.

In this study, we investigate the efficiency of tidally induced transformations in a Virgo-like cluster by varying the orbital parameters. We adopt three different orbits. On the innermost orbit, we use three different orientations of the disk with respect to the orbital plane. Furthermore, on the innermost orbit, we probe how the efficiency changes if the progenitor galaxy is less extended. We use a simple stellar evolution approach to compare our galaxy remnants to the observed dE population of the Virgo cluster. The aim of this study is to show the relevance of the parameters, as well as their general trends with respect to transformation efficiency. In Sect. 2 we outline the numerical method and the different setups. Results of the different simulations are presented in Sect. 3. In Sect. 4 we compare our simulations with observations. In Sect. 5 the results are discussed and interpreted, leading to our conclusions in Sect. 6.

2. Numerical methods

There are many parameters that affect the tidal interaction of an infalling galaxy inside a galaxy cluster: the orbital parameters, the structural parameters of the infalling galaxy, and the structural parameters of the cluster. For this study we decided to keep the cluster fixed and concentrate mostly on the orbital parameters, but including a variation of the structural parameters of the infalling galaxy. We used a combination of N-body simulations and a Monte-Carlo method to simulate the harassment scenario. For these simulations the direct N-body integrator NBODY6 was taken (Aarseth, 1999, 2003), with GPU acceleration by Nitadori & Aarseth (2012) and the massively parallel version using MPI (Spurzem, 1999; Spurzem et al., 2008). This code, called NBODY6++GPU (Wang et al., in prep.), was modified to treat strong and time-varying tidal fields (see Appendix) and to simulate encounters with massive point-like particles for the harassment scenario.

While the current version of the code can simulate up to 10⁶ particles now on a few dozen GPU nodes (order 10⁵ GPU cores for Kepler GPUs), we used only much more moderate particle numbers; any unwanted two-body relaxation occurring due to too small a particle number has been corrected for (see Sect. 2.3). The code contained the single-stellar evolution (SSE) package of Hurley et al. (2000) and has been combined with the GALEV spectral synthesis code (Kotulla et al. 2009; Pang et al., in prep.). Although NBODY6 and its offsprings were originally designed for star-by-star cluster simulations, the code turned out to be very useful in this project, also thanks to the GPU acceleration. We have not resolved our galaxies star-by-star, so a single body in the simulation typically represents many stars.

2.1. Cluster model

The gravitational force of the cluster was modelled as an analytical force field. An Navarro-Frenk-White profile (Navarro et al., 1996) was adopted to describe the density profile of the galaxy cluster, and the fit parameters ρ₀ = 3.2 × 10^-4M_⊙/pc^-3 and R_s = 560 kpc from McLaughlin (1999) were used to reproduce the force field of a Virgo-like cluster with M₂₀₀ = 4.2 × 10¹⁴M_⊙ and R₂₀₀ = 1.55 Mpc. Nbody6++GPU is using Hill’s approximation to take external force-fields into account. In this approximation it is assumed that the galaxy is moving on an circular orbit around the cluster centre, which allows the Euler forces to be neglected. We have added some code to Nbody6++GPU to include the Euler forces, so that we are able to use eccentric orbits for the galaxies around the cluster centre. An explanation of the equations of motions that are used is given in Appendix A. The orbit of the infalling galaxy around the cluster centre was calculated using a leap-frog algorithm. Free parameters of the orbit were the apocentre and pericentre; we kept the pericentre fixed at 0.3 Mpc and used three different apocentres (0.5 Mpc, 1.0 Mpc, and 1.5 Mpc).

2.2. Encounters with other cluster members

Flyby events of cluster members were simulated by using a Monte-Carlo method. A random number test was used to check periodically if an encounter happened or not. In the case of an encounter, an additional particle was injected into the simulation, which represents the flying-by cluster member. After that we waited half of the typical time between two encounter events and checked again whether another encounter had happened.

The typical time between two encounters depends on the local galaxy density and the local velocity dispersion of the cluster. In this way the encounters with other cluster members are more frequent near the cluster centre. The local galaxy density was derived from the radial galaxy distribution of those galaxy clusters of the Millennium II simulation (Boylan-Kolchin et al., 2009) that have masses between 2.4−4 × 10¹⁴M_⊙. We made use of the corresponding semi-analytic model of Guo et al. (2011) to select cluster member galaxies with total (dark+stellar) masses of more than 0.1 M_gal, where M_gal is the initial total mass of our simulated infalling galaxy. The local velocity dispersion was calculated numerically as a solution of the spherical Jeans-equation for the cluster potential.

The masses and velocities of the injected particles were chosen randomly. The masses following the mass function of the aforementioned galaxy clusters of the Millennium II simulation (Boylan-Kolchin et al., 2009). The minimum mass of a flying-by galaxy that was assumed to be relevant strong enough to disturb the infalling galaxy, hence is worth considering, was 0.1 M_gal. The velocities follow a Maxwell-Boltzmann distribution, and the velocity orientation is chosen randomly.

The impact parameters of the encounters were chosen by following the geometrical cross section. Two different galaxy models with different extents were used in this study. The minimum impact parameter for simulations with the more extended galaxy model was set to 8 kpc and for simulations with the smaller galaxy model to 4 kpc, the maximal impact parameter was set to 60 kpc. Owing to gravitational focusing, the distance of minimal approach will be closer than the impact parameter at infinity. To take this into account and to correct the impact parameter at infinity, we enlarged it by using the relation of the closest approach and impact parameter of Spurzem et al. (2009). A more detailed description on how the parameters of encounters were chosen is given in Appendix B.

2.3. Galaxy models

The galaxy models were three component models with a stellar disk and bulge and a dark matter halo created by the GalactICS code of Kuijken & Dubinski (1995). This code uses a King model (King, 1962) for the bulge, an exponential Shu model (Shu, 1969) for the disk, and an Evans model (Evans, 1993) for the dark matter halo.

Two different galaxy models were used, one with a disk scale length r₀ of 3 kpc and one with a disk scale length of 1.2 kpc. The 3 kpc model was our main model, which was used in 20 of 24 simulations and is referred to as the “large model”. The 1.2 kpc model is referred to as the “small model”. The large model is based on the typical values of dust-corrected k-band observations of Sc/Sd galaxies (Graham & Worley, 2008), which we assume to be suitable late-type progenitor galaxies. The mass-to-light ratios for the conversion of the k-band values were 0.46 for the disk and 0.82 for the bulge component. We determine these conversion factors with the GALEV evolutionary synthesis models from Kotulla et al. (2009) with Sandage-functions (Sandage, 1986) for the star formation histories and delay times of τ = 3 Gyr for the bulge and τ = 9 Gyr for the disk.

The small model is constructed to be a dwarf version of an Sc/Sd galaxy. The properties of both are listed in Table 1. With a thickness of r₀/z₀ = 9.2 for the small model and r₀/z₀ = 6.7 for the large model, both models have a comparatively thin disk component (Bizyaev & Mitronova, 2009). The models were represented by 10⁵ particles, with 3 × 10⁴ disk particles, 5 × 10³ bulge particles, and 6.5 × 10⁴ halo particles. Because Nbody6++GPU is a collisional N-Body code, one needs to correct the simulation results with respect to two-body relaxation. To determine the evolution under two-body relaxation, a control simulation was made for each of both models. In this control simulation, the evolution of the galaxy models without any tidal field or galaxy encounters was simulated, yielding the evolution of the model galaxy parameters that is caused by two-body relaxation alone. Following the idea of Velazquez & White (1999), one can correct the measured parameters from the simulations by subtracting the reference values of the parameter change caused by two-body relaxation, as determined by the control simulation. In that way one gets the relaxation-free change of the model parameters. We have cross-checked this assumption, by simulating two of our simulations again, doubling the total number of particles, and it turned out that the correction works fine. The differences in the relaxation-free values were in the range 5%−10%.

2.4. Setups

Three different orbits were used in this study: their apocentres were 1.5 Mpc, 1.0 Mpc, and 0.5 Mpc, and all of them have the same pericentre at 0.3 Mpc. The different orbits are shown in Fig. 1. All simulations started at their apocentre and spanned a time of 5 Gyr. The time span of 5 Gyr is similar to previous studies (Moore et al., 1996, 1998, 1999; Mastropietro et al., 2005), and it was chosen to follow the evolution of galaxies inside a cluster environment since z ≃ 0.5.

To probe the importance of the initial structural parameters of a galaxy, we ran one simulation setup with a smaller galaxy model. It was shown by Villalobos et al. (2012) that the inclination of the disk to the orbital plane could change the evolution of a galaxy in a group or cluster environment, in the way that galaxies with an inclination could resist the tidal forces for a longer time. To investigate how strong this impact is, we took three different inclinations for the galaxies on the innermost orbit: the disk plane parallel to the orbital plane 0°, the disk plane perpendicular to the orbital plane 90°, and one orientation between 45°.

We made three simulations per setup, each with a different random number seed to get different realisations of the random processes. In addition, one simulation per setup was made without encounters with other cluster members, to study the influence of the tidal field alone and to investigate the stochasticity of the harassment scenario. Our setups are listed in Table 2.

Fig. 1 Orbits with a pericentre of 0.3 Mpc and an apocentre of 0.5 Mpc (red), 1.0 Mpc (green dashed), and 1.5 Mpc (blue dotted).

2.5. Bound mass

The following algorithm was used to determine the bound particles. We calculated the kinetic and potential energy of each particle and checked whether the absolute value of the potential energy was higher than its kinetic energy. In that case the particle was bound energetically, but it remained unclear whether it was bound to the galaxy or only to the cluster. Particles that were outside the Lagrangian points L1 and L2 of the galaxy are not bound to the gravitational centre of the galaxy but to the cluster and will leave the galaxy following the attraction of the cluster. To distinguish between both cases, we used the tidal radius R_T, which is an approximation of the distance from the galaxy centre to the Lagrangian points L1 and L2, in order to clarify whether a particle lies beyond the Lagrangian points or not. Because the tidal radius itself is only an approximation of this problem, we used the criteria for the discrimination that a particle has to be inside of 1.5 R_T around the density centre of the bulge, in order to be bound to the galaxy.

All simulations started at the apocentre of their orbits. At this point the initial tidal radius was mostly already smaller than the extent of the dark matter halo and therefore only part of the whole mass of the galaxy was bound initially. The last column of Table 2 shows the initially bound mass fractions.

2.6. Virtual photometry

The photometric properties of different galaxy types are not randomly distributed in phase space, but the different galaxy types populate specific regions of the phase space. If the tidal transformation of late-type galaxies by the harassment scenario is one of the main production channels of dE galaxies in galaxy clusters, the remnants of harassed galaxies have to match the photometric properties of the dE galaxies. To investigate whether the simulated remnants fulfil this condition, the photometric properties of the remnants were calculated and compared to those of dE galaxies in the Virgo galaxy cluster.

To compute the brightness and the half-light radius of the simulated galaxies, the stellar evolution of the bulge and disk component were simulated by using the single stellar evolution code of Hurley et al. (2000). To do so, we adopted a Sandage-function (Sandage, 1986) for the star formation history, following Gavazzi et al. (2002), with a delay time τ = 3 Gyr for the bulge and τ = 9 Gyr for the disk. As initial mass function, the Kroupa IMF (Kroupa, 2001) in a mass range of 0.08 M_⊙−8 M_⊙ was used. The metallicities of the stars were calculated by using the chemically consistent evolution of the GALEV population synthesis models (Kotulla et al., 2009). Each stellar particle was assigned an individual mass, age, and metallicity, and the evolution of the star was traced. Finally parts of the GALEV population-synthesis code (Kotulla et al., 2009) were used to convert the stellar properties into the r-band brightness of the stars by using the BaSeL-library (Lejeune et al., 1997, 1998).

Furthermore, it was assumed that the star formation ended abruptly when the galaxy entered the cluster. While it is expected that a low-mass galaxy loses most of its gas relatively quickly by ram pressure (Gunn & Gott, 1972; Boselli et al., 2008) – even outside the virial radius when entering along filaments (Bahé et al., 2013) – this is certainly a simplified approach. Models employing an instantaneous removal of the extended gas reservoir of galaxies that become satellites led to an overabundance of passive galaxies in clusters (Weinmann et al., 2010); a recent study combining observations and simulations instead points to a picture in which “satellite quenching is delayed-then-rapid” (Wetzel et al., 2013). On the other hand, the relevance of pre-processing (Villalobos et al., 2012) may imply an earlier decline of the star formation rate. Since we are primarily interested in simulating the effects on the stellar structure, we consider the star formation rate to be zero from the starting point of our simulation, 5 Gyr ago.

To compute the total brightness and the half-light radius of the galaxies, projections of the galaxies were created under different projection angles Θ ∈ [ 0°,22.5°,45°,67.5°,90° ] along the x and y-axes. We extracted radial profiles by fitting ellipses to the isophotes following the formula of Bertin & Arnouts (1996). As total brightness of the galaxy, the brightness inside of an elliptical aperture of two Petrosian semi-major axes was used (Petrosian, 1976; Lisker et al., 2007). The major axis of the elliptical profile which encloses half of the total luminosity was adopted as the half-light radius of the galaxy. From the total brightness and the half-light radius one can compute the effective surface brightness.

For the comparison with observations, data from the Sloan Digital Sky Survey (SDSS; Adelman-McCarthy et al. 2006) were used. The sample consists of 864 members and possible members of the Virgo cluster, classified as dE, dEdi, S0, Sc, Sd, Sm, and Im galaxies based on the Virgo Cluster Catalogue (Binggeli et al., 1985) down to the completeness limit of a B-band apparent magnitude 18 and with velocities below 3500 km s^-1. The measurements of r-band brightness, surface brightness, and half-light radius were taken from Lisker et al. (2007), Janz & Lisker (2008, 2009), and Meyer et al. (2014). The sample comprises 449 dE , 41 dEdi, 73 S0, 114 Sc/Sd, 26 Sm, and 161 Im galaxies.

3. Results

In the following section, we want to analyse the impact of orbital parameters, the structural parameters of a galaxy, and the inclination of a galaxy. In Table 3, the initial properties of the galaxy models, as well as the properties of the galaxy remnants of the different simulations, are listed. The values of Table 3 are already relaxation-corrected as explained in Sect. 2.3.

3.1. Effect of orbital parameters and encounters

As expected, the orbital parameters were essential for the evolution of a galaxy inside of a galaxy cluster environment. The evolution of the bound stellar mass fractions in the simulations with different apocentres (0.5 Mpc, 1.0 Mpc, and 1.5 Mpc) is shown in Fig. 2. On the farthest orbit, the stellar mass loss was between 4% and 6%, whereas on the innermost orbit, the stellar mass loss was 24% to 50%. On the intermediate orbit, it reached 8% to 13%. The scatter of the mass loss in the different setups is a result of the stochastical nature of the harassment scenario.

Figures 3 and 4 show the evolution of the galaxies in the first and second runs of the setup with the closest orbit by a series of snapshots with intervals of 625 Myr. In the first run of this setup, no encounter with another cluster member resulted in additional mass loss. As opposed to that, in the second run at a very early stage of the simulation, an encounter pulled nearly 5% of all stellar particles out of the disk and disturbed the system strongly. This resulted in a different mass loss history and different shape of the remnant. For a better understanding of the mass loss histories and evolution of the galaxies, it is reasonable to separate the effects of the global tidal field and close encounters, as we do in the following.

3.1.1. Mass loss on different orbits without encounters

Using the simulations without encounters we could follow the evolution of the infalling galaxy introduced only by the global tidal field of the cluster. This evolution is a kind of minimal effect that each galaxy experiences if it falls on a given orbit into a cluster. That effect could be amplified by encounters, but not weakened.

The stellar mass loss by the global tidal field (dotted lines in Fig. 2) started when the galaxy passed the pericentre the first time. At the pericentre the tidal forces reached their maximum strength and therefore the tidal radius reached a local minimum. The local minimum of the tidal radius is equivalent to a local minimum in the bound mass fraction. The small tidal radius at the pericentre enabled the growth of tidal tails. Along these tidal tails, stars also left the galaxy more than one Gyr after the galaxy passed the pericentre. The smaller the apocentre of an orbit, the more frequently the galaxy passes the pericentre and the smaller the tidal radius averaged over the simulation time. The tidal force becomes more efficient in removing material the smaller the averaged tidal radius gets relative to the averaged half-mass radius of the stellar component of the galaxy. Figure 5 shows that for the galaxies without an inclination, the ratio of the averaged tidal radius along the orbit and the averaged half-mass radius of the stellar component is an adequate way to quantify the efficiency of tidal forces. The mass loss caused by the persistent global tidal forces of the cluster’s background mass represents the minimum mass loss by the harassment scenario.

Fig. 2 Time evolution of the bound stellar mass fraction in simulations with different orbital parameters of the large galaxy model (coloured) and with the small model (grey). The dotted lines represent simulations without encounters with other cluster members.

Fig. 3 Snapshots from the first simulation run of the setup with the closest orbit, i.e. an apocentre of 0.5 Mpc and pericentre of 0.3 Mpc. The timestep between snapshots is 625 Myr; they start at the upper left at 0 Gyr and reach the lower right after 5 Gyr. Each snapshot image shows the galaxy face-on and from two perpendicular edge-on viewpoints. In the upper right corner of each image, the actual position on the orbit around the cluster centre is marked. The red dots are bound disk particles, the green dots are bound bulge particles, and the black dots are unbound stellar particles.

Fig. 4 Like Fig. 3, but for the second simulation run of the setup with the closest orbit.

A huge part of the dark matter halo was also stripped at the first pericentre passage. For the large model, the loss of initially bound total mass was close to or more than 50%, for the small model this loss was in most cases only 6%.

3.1.2. Mass loss on different orbits with encounters

Within one setup, the runs show a large scatter in their stellar mass loss. This scatter results from different encounters with cluster members. Not every encounter has the same impact, it depends on the minimal distance, the mass ratio and the relative velocity of both galaxies. To describe the impact of encounters in a similar way to how we described the impact of the global tidal field before, we calculated the tidal radius of the galaxy with respect to the flying-by objects at the moment of closest approach. This tidal radius is expressed in units of the half-mass radius of the stellar component of the galaxy model. This quantity shows whether a flying-by cluster member could have any influence on the galaxy. However, one also needs to take the relative velocity of both galaxies into account, because the slower a galaxy is as it fies by, the higher the energy and angular momentum transfer from this object. Figure 6 shows the distribution of the above properties for all flying-by galaxies in our simulations. The symbol colour indicates the approximate energy that the flying-by cluster members inject into the stellar component of the galaxy, in units of the total energy of the stellar component. To calculate this deposited energy, we used the approximation of Spitzer (1958) for impulsive high-speed encounters. This is a simplification of the real galaxy-galaxy encounter, but it is accurate enough to give an impression of the strength of the encounter. The figure illustrates that only a small number of all encounters with other cluster members fulfilled both criteria: to have a strong influence and be slow in relative velocity. Consequently, only a few encounters result in significant energy injection and additional mass loss. Therefore several encounters with cluster members are needed to have a strong impact on the infalling galaxy. This small fraction of particularly strong encounters with cluster members is also why the mass loss curves of several simulations are close to the mass loss curves of simulations without encounters. The different evolution of a galaxy with and without strong encounter is illustrated in Figs. 3 and 4 for the same simulation snapshots. The total number of encounters and the number of strong encounters (with ΔE ≥ 10^-2) in each simulation run is given in Cols. (6) and (7) of Table 3. The total number of encounters in simulations with the small-galaxy model is two times higher than for simulations with the large model on the same orbit, because for the encounter probability, we take cluster members with masses above 0.1 M_gal into account (M_gal is the initial mass of the galaxy model), therefore more galaxy encounters occur to the smaller and lighter model. But it is also harder to affect the small model by an encounter. It has to come closer to the galaxy to strip the disk stars of the smaller model away, and the probability of such close encounters declines quadratically with the impact parameter.

The scatter of mass loss within a given setup declines with the increasing apocentre of the orbit. The reason is that, in setups with larger apocentre, the galaxies have spent a smaller amount of time in regions of the cluster with a high galaxy density, and therefore fewer encounters occurred.

3.1.3. Transformation of shapes and kinematics

The top panel of Fig. 7 shows the thickness of the galaxy remnant against the ratio of the average tidal radius to the average half-mass radius. The thickness of the galaxies is measured as the ratio of the square root of the quadratic mean height of stars and the square root of the quadratic mean radius of stars ($\mathit{c}\mathit{/}\mathit{a}\mathrm{=}\sqrt{⟨{\mathit{z}}^{\mathrm{2}}⟩}\mathit{/}\sqrt{⟨{\mathit{r}}^{\mathrm{2}}⟩}$). When averaging over the runs of each setup, we find that the stronger the tidal force, the stronger the average thickening of the galaxies. At a fixed ratio of averaged tidal radius to averaged half-mass radius, the thickening was increased by the individual encounters with cluster members. Within a given setup, the stronger the stellar mass loss induced by such an encounter, the stronger the induced thickening. But also without additional mass loss, encounters were able to disturb the disk and induce a small thickening of the galaxy. The bottom panel of Fig. 7 shows the ratio of maximal rotation velocity and the mean velocity dispersion of the stellar component. The galaxy model started with a rotation-dominated stellar component, but the rotational support decreased with time by mass stripping and disturbances. Even if only small or no additional stellar mass loss happened – as in the case of the simulations with an apocentre of 1.5 Mpc – the ratio of peak velocity and mean dispersion could be decreased by perturbations by ~10% or more. In the two simulations where the stellar mass loss was above 30%, the ratio of the peak velocity and the mean velocity dispersion drops below one and the systems become dominated by pressure.

3.2. Effect of initial structural parameters

Changing the initial extent and the mass of a galaxy model is another way to change the ratio of the averaged tidal radius and the averaged half-mass radius of the stellar component while keeping the orbital parameters fixed. For a better comparison we plot the results of the simulations with the small model, together with those of Sect. 3.1, where we used the large model in Figs. 2, 5, and 7. It illustrates how important the initial structural parameters of a galaxy model are. They are as relevant as the orbital parameters.

In the simulation without encounters, the small model lost no stellar mass, experienced no thickening, its ratio of peak velocity and mean dispersion decreased only by ~6% and only ~6% of its dark matter mass got lost. The small model also had a smaller cross section for encounters with other cluster members, which made it hard to be influenced by such encounters. Only one encounter in one of the runs was strong enough to disturb the galaxy, and it took only 2% of the stellar mass away. However, this single encounter increased the thickness of the disk by a factor of two and decreased the ratio of maximal rotation velocity and mean dispersion by an additional 8%. These results fit well with the finding of Sects. 3.1 and 3.2 that the extent of a galaxy has to be close to its tidal radius and the galaxy has to be in a region with high galaxy density to get tidally transformed efficiently.

3.3. Effect of inclination

Fig. 5 Bound stellar mass fraction versus the averaged tidal radius in units of the averaged half-mass radius of the model galaxy. Open symbols represent simulations without encounters with other cluster members.

Fig. 6 Relative velocity of the flying-by cluster members versus the tidal radius with respect to the flying-by cluster members in units of the half-mass radius of the model galaxy at the moment of closest approach. The colour scale denotes the approximate energy that the flying-by cluster members inject into the stellar component of the galaxy, in units of the total energy of the stellar component.

Fig. 7 Top: axis ratio of the stellar component of the galaxy remnants versus the averaged tidal radius in units of the averaged half-mass radius. Bottom: ratio of the rotational peak velocity and the mean velocity dispersion of the galaxy remnants versus the averaged tidal radius in units of the averaged half-mass radius of the stellar component. Both panels only contain simulations without any inclination between orbital plane and galactic disk. Open symbols represent simulations without encounters with other cluster members.

Fig. 8 Bound stellar mass fraction in the simulations with different orientations of the galactic disk. The dotted lines represent simulations without encounters with other cluster members.

Fig. 9 As in Fig. 7, but for simulations with an inclination between orbital plane and galactic disk.

The inclination of the galactic disk relative to the orbital plane is an additional parameter that influences the efficiency of the tidal transformation of galaxies. We ran simulations for galaxies with different inclinations (0°, 45°, and 90°) on those of our orbits that are farthest inside. An inclination is expected to decrease the efficiency of mass stripping by the global tidal field (Villalobos et al., 2012), but in general it is not expected to change the effect of close encounters. Figure 8 shows the evolution of the bound stellar mass fraction in the simulations with different inclinations and the closest orbit. As expected the mass loss induced by the global tidal field alone (dotted lines in Fig. 8) decreases with increasing inclination. In the case of a perpendicular orientation of the disk (90°), the stellar mass loss induced by the global tidal field was reduced by half compared to the case where the disk is parallel to the orbital plane (0°). The additional mass loss caused by encounters with other cluster members still shows a large scatter.

Figure 9 shows the thickness and the ratio of peak velocity to the mean velocity dispersion against the bound stellar mass fraction. Galaxies with the same stellar mass loss show stronger thickening if they do not move parallel to their disk. Under the influence of the global tidal field and without encounters, the galaxies with an inclination thicken more. Especially with an inclination of 45°, the galaxy thickens by a factor of five, whereas in the case of 0°, the galaxy thickens only by a factor of two. The simulations without encounters show that the galaxy with a perpendicular orientation of the disk to the orbital plane could preserve its rotational support much better than the galaxies with other inclinations. Nevertheless, under the influence of encounters, a galaxy with perpendicular orientation of the disk also became pressure-dominated. This highlights the stochastical nature of the harassment scenario again.

One has to stress that the effects an inclination has on the interaction with the global tidal field of a cluster are diverse: it could reduce the mass loss, amplify the thickening, and reduce the kinematical transformation. However, it seems that the influence on each of these quantities differs between the different orientations. Therefore, to describe the evolution of galaxies in a cluster environment, one has to pay attention to the different possible orientations of the disk.

4. Comparison with observations

Figure 10 shows the distribution of different galaxy types of the Virgo cluster (Lisker et al., 2007; Janz & Lisker, 2008, 2009; Meyer et al., 2014) and the distribution of our initial galaxy models (green), as well as their remnants (red) in the ⟨μ⟩_eff − r_eff (top panel) and in the M_r − ⟨μ⟩_eff (bottom panel) plane. The points for the initial and remnant galaxies are the mean values averaged over different projection angles (edge-on = 90°, 67.5°, 45°, 22.5°, 0° = face-on) along the x-axis and along the y-axis; the error bars span the range from the minimum to the maximum value of these projections of the galaxies. The small scatter of the initial models originates in the different random-number seeds that are used to assign masses and ages to the stellar particles. The large error bars in the surface brightness of the initial models stems from the fact that the surface brightness of the edge-on projection of the model is nearly 2 mag / arcsec² brighter than the face-on projection. The point clouds of initial and remnant galaxies that belong to simulations with the large model can be identified by their clear separation in effective radius and magnitudes in comparison to the point clouds of the smaller model.

Fig. 10 Top: distribution of Virgo galaxies, as well as the initial galaxy models and remnants in the ⟨μ⟩eff − reff plane. Bottom: distribution of Virgo galaxies, as well as the initial galaxy models and remnants in the Mr − ⟨μ⟩eff plane. The points for the initial and remnant galaxies are the mean values averaged over different projection angles (edge on = 90°, 67.5°, 45°, 22.5°, face on = 0°) along the x-axis and along the y-axis. The error bars span over the range from the minimum to the maximum value of the projections.

Compared to the distribution of the Sc/Sd galaxies in the Virgo cluster as given by Lisker et al. (2007), Janz & Lisker (2008, 2009), and Meyer et al. (2014), the initial large models that are based on the Sc/Sd parameters of Graham & Worley (2008) tend to have larger radii. The initial small models were smaller than the typical extent of the Virgo Sc/Sd galaxies, which have a mean effective radius of $\overline{{\mathit{r}}_{\mathrm{eff}}}\mathrm{=}\mathrm{2.63}$ kpc. The point cloud of remnants of the large model overlaps partially with the more extended dEs and dEdis. But on average the remnants of the large model tend to have radii that are too large to fit with the distribution of dEs well. Their magnitudes are close to –18 mag in r, which is often used as the threshold for classifying bright dwarf galaxies, but in general they are too bright for ordinary dEs or dEdis. This trend indicates that the stripping of stellar material, hence the truncation of the disk and the dimming of light in simulations with the large model, was not efficient enough to produce typical dE or dEdi galaxies. However, the direction of evolution of the large model implies that harassment of galaxies could potentially produce dE and dEdis if the effects were amplified in some way.

The remnants of the small model are all placed inside the distribution of bright dE and dEdi galaxies. But the dimming of light of these remnants was only caused by the truncation of the star formation at the beginning of the simulation, not by stripping of stellar matter. We have to stress that it has been shown before that these galaxies were nearly unaffected by the environment, so these remnants are still disk galaxies with a small thickness (c/a = 0.09–0.19) and are strongly dominated by rotation (v/σ = 1.99–2.14). They could not be used to explain the distribution of dE galaxies because they did not undergo a morphology transformation. Nonetheless, we can imagine that if the effect of the tidal interaction was amplified, then the direction of evolution of these galaxies would be similar to that of the large model, and they would undergo a morphological transformation and end up as dE galaxies.

This comparison shows that there is already a partial overlap of harassed galaxies with dE/dEdi galaxies after an evolution of 5 Gyr in a cluster environment. From these results one can expect that if the environmental effects were amplified, the harassment scenario could be one production channel of dE and dEdi galaxies.

5. Discussion

5.1. Orbital parameters

By varying the apocentre of the orbit of an infalling galaxy between 0.5 Mpc and 1.5 Mpc at a fixed pericentre of 0.3 Mpc, it turned out that the tidally induced morphological transformation of low-mass late-type galaxies to bright dE galaxies became more efficient, the closer the orbit of the galaxies was. With decreasing apocentre, the mass loss, the thickness, and the pressure support increased (also see Mastropietro et al., 2005; Smith et al., 2010). However, the comparison of the galaxy remnants with observations of dE galaxies shows that even on the innermost orbit, the environmental effects were still not efficient enough to transform the infalling galaxies into typical dE galaxies within an evolution time of 5 Gyr. In general, the remnants are still too extended and too bright, which indicates that a stronger environment would be necessary to strip more material and complete the transformation, or the progenitor galaxies would have to be late-type disks of lower luminosity and size than our model (see Graham et al., 2003).

If one takes dynamical friction into account, galaxies that enter a cluster will sink deeper into the cluster core. We neglected this effect in our simulations and probed what happened to galaxies that did not end up near the cluster centre. This is a relevant question because one observes dE galaxies not only in cluster cores but also in the outskirts of clusters, so a main production process of dE galaxies would also have to work in these regions and not just in the very cluster centre.

On our more extended orbit (apocentre 1.5 Mpc) in the outskirts of the cluster, the tidal forces of the cluster were not strong enough and the galaxy encounters not frequent enough to harass the galaxies. On this orbit, the galaxies passed the pericentre only once in 5 Gyr, and the mass loss, as well as the thickening of the galaxies, was marginal. These results agree with the study of Smith et al. (2010), who found in their simulations that low-mass dwarf irregular galaxies that enter a galaxy cluster from its present-day virial radius were mostly unaffected by the tidal forces of the cluster or by interactions with outer cluster members. Only on those orbits that are very close to the cluster centre did Smith et al. find the tidal interactions to be efficient to transform dwarf irregular galaxies. They conclude that galaxies that are transformed in a harassment scenario have to have entered a galaxy cluster a long time ago to reach the innermost regions of a galaxy cluster and get transformed.

De Lucia et al. (2012) have shown by using the data of the Millennium simulation (Springel et al., 2005) that galaxies with halo masses below 10¹¹M_⊙ that are in a region closer than 0.5 R_vir today became satellites of bigger systems 9 Gyr ago (median value). Nine billion years are nearly double the time that we simulated. However, we could expect from our simulation results that a galaxy that spends such a long time in an environment like a cluster core could complete the transformation to a dE galaxy.

One point that needs to be taken into account when simulating such a long time period is that the cluster itself will evolve on this time scale. A galaxy that enters a cluster nine billion years ago enters a protocluster environment with a mass that is more like present-day groups. Villalobos et al. (2012) have done simulations of the evolution of galaxies in a group-like tidal field to study the relevance of preprocessing of galaxies in groups. They find that the group tidal field can induce a morphological transformation (also see Mayer et al., 2001). Taking the preprocessing in galaxy groups into account could thus also be an answer for how harassed galaxies could populate the outskirts of clusters, because not all protocluster members will end up in the cluster core.

Lisker et al. (2007) mention that the different infall times and the associated cluster centric distances could be an explanation of the origin of different dE subpopulations like dEdis and could give a natural answer to the question of why the different subpopulations of dEs follow a morphology-density relation. In this sense one could interpret dEdis not only as galaxies in an uncompleted transformation phase (Mastropietro et al., 2005) but also as a separate subpopulation that evolved largely in a different environment than ordinary dEs. In the same sense, all the different subpopulations of dEs may have evolved in parallel in different environments and would then be the present-day product of this parallel evolution (Lisker et al., 2013).

The next step in clarifying the relevance of the harassment scenario as a production channel of dE galaxies would have to be studies that simulate a longer period in time and take the evolution of the environment itself into account.

5.2. Structural parameters

We have probed the influence of structural parameters of the progenitor galaxy by using two different models: one with a small and one with a large disk component, which differ by a factor 2.5 in effective radius. While the size of the large model becomes smaller by 26% on average, the size of the small model remains unchanged (see Table 3). The small model was sitting so deep inside its tidal radius that it remained nearly untouched by the global tidal field. The slight increase in the half-light radius (Fig. 10) is due to the older stellar population of the bulge component.

Small galaxies are much harder to affect by the tidal force of the cluster or by encounters with other cluster members, in agreement with the findings of Moore et al. (1999). They studied the efficiency of the harassment scenario for low and high surface brightness galaxies and found that high surface brightness galaxies are very resistant to environmental effects.

A galaxy like the small model has to be on an orbit closer to the cluster centre to undergo a similar evolution as the large model. To reach the same ratio of mean tidal radius and half-mass radius, the small model would have to be about two times closer to the cluster centre than the closest orbit we used for the large model if we assume that the enclosed cluster mass stays the same. If one takes the radial mass profile of the cluster into account, one can in fact estimate that the small model would have to be about seven times closer to the cluster centre. This means it would have to be on an orbit with an apocentre closer than 100 kpc to the cluster centre. This emphasizes again that the initial structural parameters of the progenitor galaxy are as relevant as the orbit of the galaxy.

5.3. Inclination

The inclination of the disk of an infalling galaxy has been found to be an important parameter for the influence of the cluster’s tidal field. The induced stellar mass loss was reduced by half when the disk was perpendicular to the orbital plane as compared to a parallel orientation on the same orbit. On the other hand, the thickening of galaxies with an inclination was stronger than for parallel oriented galaxies. Especially with an inclination of 45°, the thickening induced by tidal forces was more than twice as large as the parallel oriented galaxy. Additionally, it turned out that galaxies with a perpendicular orientation could shield themselves better against kinematical transformations by the tidal field.

Villalobos et al. (2012) have found in their simulations of the evolution of galaxies in a group-like tidal field that galaxies that are initially face-on or retrograde-oriented lost less mass and preserved their disk structure and kinematic for a longer time than prograde galaxies¹. They conclude that this could be an explanation for disk galaxies in the inner regions of groups and clusters. We agree with Villalobos et al. (2012) that galaxies with an inclination experienced less mass loss and kinematical transformations.

We speculate that the effect of inclination could be an origin for dEdi galaxies. dEdis appear to be less round than ordinary dEs (Lisker et al., 2006b), but with a thickness of 0.35–0.4, they are somewhat thicker than typical late-type disk galaxies at similar stellar mass (Sánchez-Janssen et al., 2010). The rotation parameter v_max/σ of dEdi galaxies is smaller than for disk galaxies but higher than for ordinary dE galaxies (Toloba et al., 2011, 2015). We found that our simulated galaxies with an inclination of 45° or 90° thicken more and reach the regime of dEdi galaxies. Galaxies with a disk perpendicular to their orbital plane (90°) could shield their kinematics much better against the influence of the tidal field. Both observations together lead to the speculation that disk galaxies that enter a cluster at a high inclination 45°–90° have a good chance of ending up with a thick disk with high v_max/σ values. In that case, harassed high-inclination disk galaxies may be an origin for dEdi galaxies.

6. Conclusions

Harassment of galaxies is a complex process. Its efficiency depends on a combination of many parameters, such as the orbit, the inclination, and the structural parameters of the progenitor galaxy, as well as on the properties of the galaxy cluster itself. It is likely that harassment is one of the most important engines of galaxy evolution inside galaxy clusters, and is often referred to in past studies. However, this process would have had to start at a very early epoch in protocluster environments to complete the transformation by today. To fully understand the overall importance of harassment, more simulations that reproduce the complicated interaction between galaxies and their environment are needed.

Online material

Appendix A: Equations of motions

The reference frame that was used to describe the equations of motion of a particle in our simulation was centred on the density peak of the bulge. The galaxies themselves moved on an eccentric orbit around the cluster centre. Additionally, we rotated our coordinate system in such a way that the x-axis of the reference frame always pointed away from the centre of the cluster. This means that our coordinate system is not an inertial system, but an accelerated and rotated one.

The equation of motion of a particle is $\begin{array}{ccc}\begin{array}{c}{¨}\\ {r}\end{array}\mathrm{=}{{f}}_{\mathrm{gal}}\mathrm{+}{{f}}_{\mathrm{cl}}\mathrm{+}{{a}}_{\mathrm{0}}\mathrm{+}{{f}}_{\mathrm{ce}}\mathrm{+}{{f}}_{\mathrm{co}}\mathrm{+}{{f}}_{\mathrm{eu}}\mathit{.}& & \end{array}$(A.1)The first two terms are the gravitational forces of the galaxy f_gal and the cluster f_cl, which acts on a particle. The gravitational force of the galaxy was calculated by direct summation of the two-body interaction of the particles, and the gravitational force of the cluster was calculated by an analytical force field as described in Sect. 2.1: $\begin{array}{ccc}{{f}}_{\mathrm{gal}}& \mathrm{=}& \mathrm{-}{\nabla }{\mathrm{\Phi }}_{\mathrm{gal}}\\ {{f}}_{\mathrm{cl}}& \mathrm{=}& \mathrm{-}{\nabla }{\mathrm{\Phi }}_{\mathrm{cl}}\mathit{.}\end{array}$The third term of Eq. (A.1) is the force of inertia a₀. It has the same strength but inverse direction as the force that accelerates the galaxy on its orbit around the cluster centre. The distance of the cluster centre and the origin of ordinates is r₀. The direction of the force of inertia points away from the cluster centre, which means that in our rotating reference frame along the direction of the x-axis, that $\begin{array}{ccc}{{a}}_{\mathrm{0}}\mathrm{=}{\left[{{e}}_{\mathit{x}}\frac{\mathrm{d}{\mathrm{\Phi }}_{\mathrm{cl}}}{\mathrm{d}{\mathit{r}}_{\mathrm{cl}}}\right]}_{{\mathit{r}}_{\mathrm{0}}}\mathrm{·}& & \end{array}$(A.4)Since it is a rotating reference frame, we have to take the centrifugal f_ce, Coriolis f_co and Euler force f_eu into account. The orbital plane of the galaxy within the cluster is chosen so that the angular velocity of the motion of the galaxy around the cluster centre is pointing in the z direction. This angular velocity is the same as the angular velocity that is used to rotate the reference frame: $\begin{array}{ccc}{{f}}_{\mathrm{ce}}& \mathrm{=}& \mathrm{-}{\omega }\mathrm{\times }\mathrm{(}{\omega }\mathrm{\times }{{r}}^{\mathrm{)}}\mathrm{=}{\mathit{\omega }}^{\mathrm{2}}\mathrm{(}\mathit{x}{{e}}_{\mathit{x}}\mathrm{+}\mathit{y}{{e}}_{\mathit{y}}\mathrm{)}\\ {{f}}_{\mathrm{co}}& \mathrm{=}& \mathrm{-}\mathrm{2}{\omega }\mathrm{\times }\dot{{r}}\mathrm{=}\mathrm{2}\mathit{\omega }\mathrm{(}\mathit{ẏ}{{e}}_{\mathit{x}}\mathrm{-}\mathit{ẋ}{{e}}_{\mathit{y}}\mathrm{)}\\ {{f}}_{\mathrm{eu}}& \mathrm{=}& \mathrm{-}\dot{{\omega }}\mathrm{\times }{r}\mathrm{=}\mathit{\omega ̇}\mathrm{(}\mathit{y}{{e}}_{\mathit{x}}\mathrm{-}\mathit{x}{{e}}_{\mathit{y}}\mathrm{)}\mathit{.}\end{array}$After an expansion of the cluster potential in a Taylor series to the first order around the origin of ordinates, the equations of motions are $\begin{array}{ccc}\begin{array}{c}\mathrm{¨}\\ \mathit{x}\end{array}& \mathrm{\approx }& \mathrm{-}\frac{\mathit{\partial }}{\mathit{\partial x}}{\mathrm{\Phi }}_{\mathrm{gal}}\mathrm{+}\left({\mathit{\omega }}^{\mathrm{2}}\mathrm{-}{\left[\frac{{\mathit{\partial }}^{\mathrm{2}}{\mathrm{\Phi }}_{\mathrm{cl}}}{\mathit{\partial }{\mathit{r}}_{\mathrm{cl}}^{\mathrm{2}}}\right]}_{{\mathit{r}}_{\mathrm{0}}}\right)\mathit{x}\mathrm{+}\mathrm{2}\mathit{\omega }\mathit{ẏ}\mathrm{+}\mathit{\omega ̇}\mathit{y}\\ \begin{array}{c}\mathrm{¨}\\ \mathit{y}\end{array}& \mathrm{\approx }& \mathrm{-}\frac{\mathit{\partial }}{\mathit{\partial y}}{\mathrm{\Phi }}_{\mathrm{gal}}\mathrm{+}\left({\mathit{\omega }}^{\mathrm{2}}\mathrm{-}{\left[\frac{\mathit{\partial }{\mathrm{\Phi }}_{\mathrm{cl}}}{\mathit{\partial }{\mathit{r}}_{\mathrm{cl}}}\right]}_{{\mathit{r}}_{\mathrm{0}}}\frac{\mathrm{1}}{{\mathit{r}}_{\mathrm{0}}}\right)\mathit{y}\mathrm{-}\mathrm{2}\mathit{\omega }\mathit{ẋ}\mathrm{-}\mathit{\omega ̇}\mathit{x}\\ \begin{array}{c}\mathrm{¨}\\ \mathit{z}\end{array}& \mathrm{\approx }& \mathrm{-}\frac{\mathit{\partial }}{\mathit{\partial y}}{\mathrm{\Phi }}_{\mathrm{cl}}\mathrm{-}{\left[\frac{\mathit{\partial }{\mathrm{\Phi }}_{\mathrm{cl}}}{\mathit{\partial }{\mathit{r}}_{\mathrm{cl}}}\right]}_{{\mathit{r}}_{\mathrm{0}}}\frac{\mathit{z}}{{\mathit{r}}_{\mathrm{0}}}\mathrm{·}\end{array}$

If one assumes a circular orbit for the galaxy , these equations of motion simplify to the Hill’s approximation.

Appendix B: Properties of encounter galaxies

The typical time between two encounters depends on the local galaxy density n(r) and the local velocity dispersion σ(r) of the cluster, which are functions of the cluster centric distance r. We take all those encounters into account that pass the infalling galaxy within a distance of p_max = 60 kpc: $\begin{array}{ccc}\mathit{\tau }\mathrm{=}{\left(\mathrm{\Sigma }\mathit{n}\mathrm{\left(}\mathit{r}\mathrm{\right)}\mathit{\sigma }\mathrm{\left(}\mathit{r}\mathrm{\right)}\right)}^{-1}\mathit{,}& & \end{array}$(B.1)where Σ is the geometrical cross section of such a galaxy-galaxy encounter $\mathrm{\Sigma }\mathrm{=}\mathit{\pi }{\mathit{p}}_{\max }^{\mathrm{2}}$. The local galaxy density was derived by assuming a galaxy distribution that follows a power law: $\begin{array}{ccc}\frac{\mathrm{d}\mathit{N}}{\mathrm{d}\mathit{r}}\mathrm{=}{\mathit{N}}_{\mathrm{0}}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathit{r}}^{\mathit{\beta }}\mathrm{·}& & \end{array}$(B.2)The parameters N₀ = 154 and β = 0.6 were determined by a fit to the radial galaxy distribution of those galaxy clusters of the Millennium II simulation (Boylan-Kolchin et al., 2009) that have a mass between 2.4−4 × 10¹⁴M_⊙.

Cluster members with total masses down to 0.1 M_gal were taken into account, where M_gal was the mass of the infalling galaxy. The total masses of the galaxies were taken from the corresponding semi-analytic model of (Guo et al., 2011). The local velocity dispersion was calculated numerically as a solution of the spherical Jeans-equation for the cluster potential. The occurrence of an encounter was simulated by randomly placing the particle that represents the perturber galaxy on a sphere with a radius of 200 kpc around the infalling galaxy. The velocities of the flying-by galaxies followed a Maxwell-Boltzmann distribution. The velocity vector of the flying-by galaxy pointed to the centre of the sphere, but offset by the impact parameter of the encounter. The impact parameter was chosen randomly based on the geometrical cross section, so that the probability of an impact parameter p ∈ [ p_min,p_max ] is given by $\begin{array}{ccc}\mathit{f}\mathrm{\left(}\mathit{p}\mathrm{\right)}\mathrm{=}\mathrm{2}\frac{\mathit{p}\mathrm{d}\mathit{p}}{\mathrm{(}{\mathit{p}}_{\max }^{\mathrm{2}}\mathrm{-}{{\mathit{p}}_{\min }^{\mathrm{2}}}^{\mathrm{)}}}\mathrm{·}& & \end{array}$(B.3)One also has to consider that the flying-by galaxy will by focused on the infalling galaxy by gravitation. We corrected the offset to take this gravitational focusing into account by following Spurzem et al. (2009).

The masses of the perturber galaxies were chosen randomly by following the mass distribution of the aforementioned galaxy clusters of the Millennium II simulation (Boylan-Kolchin et al., 2009). To determine the mass distribution, a power-law ansatz was used: $\begin{array}{ccc}\frac{\mathrm{d}\mathit{N}}{\mathrm{d}\mathit{m}}\mathrm{=}{\mathit{N}}_{\mathrm{0}}\hspace{0.17em}\hspace{0.17em}{\mathit{m}}^{\mathrm{-}\mathit{\alpha }}\mathrm{·}& & \end{array}$(B.4)The fit yielded α = 2.0.

Acknowledgments

We thank Isabel Franco for her efforts with a previous simulation study that eventually led to the project presented here. D.B. acknowledges funding by the IMPRS for Astronomy & Cosmic Physics at the University of Heidelberg. D.B. and R.S. acknowledge support by Chinese Academy of Sciences through the Silk Road Project at NAOC, through the Chinese Academy of Sciences Visiting Professorship for Senior International Scientists, grant number 2009S1−5 (R.S.), and through the “Qianren” special foreign experts program of China. T.L. was supported within the framework of the Excellence Initiative by the German Research Foundation (DFG) through the Heidelberg Graduate School of Fundamental Physics (grant number GSC 129/1). C.O. appreciates funding by the German Research Foundation (DFG) grant OL 350/1-1 and SFB 881. The work of C.O. and R.S. was supported in part by the National Science Foundation under grant PHYS-1066293 and by the hospitality of the Aspen Center for Physics. R.K. gratefully acknowledges financial support from STScI theory grant HST-AR-12840.01-A. Support for Program number HST-AR-12840.01-A was provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract NAS5-26555. We used computing resources at the Center of Information and Computing at National Astronomical Observatories, Chinese Academy of Sciences, funded by the Ministry of Finance of the People’s Republic of China under grant ZDYZ2008 − 2, and at the ZAH, Heidelberg University, GPU clusters hydra and kepler, funded under grants I/80041-043 and I/84678/84680 of the Volkswagen Foundation.

References

All Tables

All Figures

	Fig. 1 Orbits with a pericentre of 0.3 Mpc and an apocentre of 0.5 Mpc (red), 1.0 Mpc (green dashed), and 1.5 Mpc (blue dotted).
In the text

	Fig. 2 Time evolution of the bound stellar mass fraction in simulations with different orbital parameters of the large galaxy model (coloured) and with the small model (grey). The dotted lines represent simulations without encounters with other cluster members.
In the text

Fig. 3 Snapshots from the first simulation run of the setup with the closest orbit, i.e. an apocentre of 0.5 Mpc and pericentre of 0.3 Mpc. The timestep between snapshots is 625 Myr; they start at the upper left at 0 Gyr and reach the lower right after 5 Gyr. Each snapshot image shows the galaxy face-on and from two perpendicular edge-on viewpoints. In the upper right corner of each image, the actual position on the orbit around the cluster centre is marked. The red dots are bound disk particles, the green dots are bound bulge particles, and the black dots are unbound stellar particles.

In the text

	Fig. 4 Like Fig. 3, but for the second simulation run of the setup with the closest orbit.
In the text

	Fig. 5 Bound stellar mass fraction versus the averaged tidal radius in units of the averaged half-mass radius of the model galaxy. Open symbols represent simulations without encounters with other cluster members.
In the text

Fig. 6 Relative velocity of the flying-by cluster members versus the tidal radius with respect to the flying-by cluster members in units of the half-mass radius of the model galaxy at the moment of closest approach. The colour scale denotes the approximate energy that the flying-by cluster members inject into the stellar component of the galaxy, in units of the total energy of the stellar component.

In the text

Fig. 7 Top: axis ratio of the stellar component of the galaxy remnants versus the averaged tidal radius in units of the averaged half-mass radius. Bottom: ratio of the rotational peak velocity and the mean velocity dispersion of the galaxy remnants versus the averaged tidal radius in units of the averaged half-mass radius of the stellar component. Both panels only contain simulations without any inclination between orbital plane and galactic disk. Open symbols represent simulations without encounters with other cluster members.

In the text

	Fig. 8 Bound stellar mass fraction in the simulations with different orientations of the galactic disk. The dotted lines represent simulations without encounters with other cluster members.
In the text

	Fig. 9 As in Fig. 7, but for simulations with an inclination between orbital plane and galactic disk.
In the text

Fig. 10 Top: distribution of Virgo galaxies, as well as the initial galaxy models and remnants in the ⟨μ⟩eff − reff plane. Bottom: distribution of Virgo galaxies, as well as the initial galaxy models and remnants in the Mr − ⟨μ⟩eff plane. The points for the initial and remnant galaxies are the mean values averaged over different projection angles (edge on = 90°, 67.5°, 45°, 22.5°, face on = 0°) along the x-axis and along the y-axis. The error bars span over the range from the minimum to the maximum value of the projections.

In the text