Advances in Astronomy

Volume 2015 (2015), Article ID 196304, 19 pages

http://dx.doi.org/10.1155/2015/196304

## Hydrodynamic Modeling of the Interaction of Winds within a Collapsing Turbulent Gas Cloud

Departamento de Investigación en Física, Universidad de Sonora, 83000 Hermosillo, SON, Mexico

Received 6 January 2015; Accepted 20 April 2015

Academic Editor: William Reach

Copyright © 2015 Guillermo Arreaga-García and Julio Saucedo-Morales. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

By using the particle-based code Gadget2, we follow the evolution of a gas giant molecular cloud, in which a set of gas particles representing the wind are created by a Monte Carlo scheme and suddenly move outwards from the cloud’s center. The particles representing the gas cloud initially have
a velocity according to a turbulent spectrum built in a Fourier space of 64^{3} grid elements. The level of turbulence and the temperature of the cloud are both adjusted so that a gravitational collapse of the cloud is initially induced. All the winds are activated in a very early stage of evolution of the cloud. We consider only two kinds of winds, namely, one with spherical symmetry and the second one of a bipolar collimated jet. In order to assess the dynamical change in the cloud due to interactions with the winds, we show isovelocity and isodensity plots for all our simulations. We also report on the accretion centers detected at the last simulation time available for each model.

#### 1. Introduction

Stars are born in large gas structures made of molecular hydrogen. According to [1] a cloud is a gas structure with mass and size around and 2–15 pc, respectively.

The physical process by which the molecular gas is transformed from a gas structure into stars is mainly gravitational collapse, whose main effect on the cloud is a reduction in size with the corresponding increase in density at various points of the cloud. These small overdensities within the larger gas structure are defined as cores, which have a typical mass and size around and pc, respectively, [1].

The basic idea of* prompt fragmentation* is that during gravitational collapse a molecular cloud may spontaneously break into several cores. At some point during the transformation process from clouds to cores and then to stars, the densest gas structures settle down in more or less dynamically stable gas objects called protostars. If the cores do not undergo further subfragmentation upon collapsing to higher densities, they will become protostars.

To achieve a better understanding of the huge scale changes due to this transformation process, we may recall that the number density of a typical star is around molecules per . Meanwhile, the number density of a typical typical cloud structure ranges around molecules per and that of a typical protostar ranges around molecules per .

Many other phenomena occurring in the interstellar medium can have a great influence on the evolution of the cloud, among others: (i) the highly ionized gas ejected by supernovae explosions; (ii) the bipolar collimated winds ejected by massive protostars; (iii) rapidly expanding regions which hit the slower and less dense gas structures.

For instance, in [2] a set of SPH simulations were conducted to study star formation triggered by an expanding region within a spherical gas structure of uniform density that was not globally collapsing. The expanding shock plays the role of a snowplow, which sweeps out the gas of the surrounding medium. The density of the swept out gas increases as a consequence of this agglomeration process, and a gravitational collapse may then be locally initiated. Small gas overdensities can be formed in this way, which may achieve the protostar stage. Because these two processes are complementary in forming protostars, this star formation scenario is known as* the collect and collapse model*. The one-dimensional model of freely expanding wind into an uniform ambient medium was analitically and numerically solved some time ago; see [3] and references there in.

In this paper, we investigate the change in the dynamical configuration of a typical turbulent cloud when a wind of particles is outwardly ejected from the central region of the cloud. The most important differences between the present paper and [2] are the turbulent nature of our cloud and its induced global tendency to gravitational collapse. Turbulence makes a big change in the spatial distribution as the cloud becomes filamentary and flocculent. The tendency to collapse favors the condensation of gas towards the central region of the cloud. Thus, a gas wind flowing through this collapsing and turbulent gas cloud is a very difficult hydrodynamical phenomenon that can only be solved by numerical simulations. It is shown that the presence of winds (i) keeps the collapsing nature of the turbulent cloud, but they cause a delay in the runway collapse; (ii) changes the number and location of accretion centers detected in the densest central region of the cloud.

Furthermore, in order to study the effects of protostellar outflows on the turbulence of a star forming region, in [4, 5], magneto-hydrodynamics simulations (MHD) were conducted with a mesh-based code which implements the adaptive mesh refinement technique (AMR). As we use here the fully parallelized Gadget2 code which implements the SPH technique, this will make it possible to carry out an interesting comparison between the two computational techniques.

#### 2. The Physical System

In this section we briefly describe the physics of the cloud and the winds which will be considered in the following sections.

##### 2.1. The Cloud

We consider here a typical spherical cloud with a radius pc and mass . Initially, it has a radially uniform density distribution with an average density given by gr , which is equivalent to a number density molecules for molecular hydrogen with molecular mass gr/mol. The size and mass of this cloud are chosen here to be typical in statistical sense, in accordance with [1].

The free fall time is defined as the time needed for an external particle to reach the center of the cloud when gravity is the only force pulling the particle. In this idealized gravitational collapse, we havewhere is Newton’s universal gravitational constant. For our cloud, we have years.

For a spherical cloud, the approximate total gravitational potential energy is . The average total thermal energy (kinetic plus potential interaction terms of the molecules) is , where is the Boltzmann constant, is the equilibrium temperature, is the total number of molecules in the gas, and is the speed of sound; see Section 3.5 for a more precise definition.

The kinetic energy can be estimated by , where is the average translational velocity of the cloud. In order to have both energies of the same order of magnitude, , the gas elements of the cloud must attain average velocities within the rangeor km/sec, for a speed of sound given byso that the corresponding temperature associated with the cloud is K.

It is possible to define the crossing time by means ofwhich sets a time scale for a sound wave to travel across the cloud. To make the crossing time comparable in magnitude to the free fall time of (1), the front wave must have velocities around or km/sec, which are velocities a little bit slower than the ones estimated above; see (2). In this paper, we will treat propagation velocities of gas particles ranging around Mach.

##### 2.2. The Wind

The dynamical characteristics of the wind strongly depend on its type of source. All stars eject winds of the first kind, which are driven by stellar radiation. For instance, in cool stars, like the ones observed in the AGB (asymtotic giant branch) of the Galaxy, the winds cause a mass loss in the range yr whereas the terminal wind velocities are around km/sec. In OB stars, the mass loss ranges over /yr and the terminal wind velocities can go up to thousands of km/sec.

Supernovae dump around joules of thermal and kinetic energy into the interstellar medium. But there are several types of supernovas, so that the mass losses and terminal velocities are very different; see [6]. For example, for a supernova whose progenitor was a He star, the mass loss and terminal velocities are within the ranges /yr and km/sec, respectively. When the progenitor was a RSG star, their values range over /yr and km/sec, respectively.

It seems that all protostars eject highly collimated jets of gas during their formation process by gravitational collapse. The origin of these jets is still unclear, but it is very likely that the accretion disk and magnetic field around protostars play a crucial role in determining the velocities and degree of collimation of the jets. For the molecular winds associated with protostarts of Class 0 and Class 1, the characteristic velocities are around 20 km/sec. However, for optical jets of highly ionized gas, the typical jet velocities are a few hundred km/sec; see [7] and the references therein.

#### 3. The Computational Methods and the Models

In this section we briefly describe the way we set up the physical system outlined in Section 2 in computational terms as well as the models we consider in this paper.

##### 3.1. The Initial Configuration of Particles

We set million SPH particles for representing the gas cloud. By means of a rectangular mesh we make the partition of the simulation volume in small elements each with a volume ; at the center of each volume we place a particle (the th, say), with a mass determined by its location according to the density profile being considered, that is, with . Next, we displace each particle from its location by a distance on the order of in a random spatial direction.

As stated earlier, in this paper we only consider a uniform density cloud, for which , for all the simulations (see Section 2.1). Therefore, all the cloud particles have the same mass.

##### 3.2. The Initial Turbulent Velocity of Particles

To generate the turbulent velocity spectrum for the cloud particles* we follow a procedure based on the papers [8, 9]*. We set a second mesh with the size of each element given by , , and . In Fourier space, the partition is , , and . Each Fourier mode has the components , , and , where the indices , , and take values in , , and , respectively. The wave number magnitude is , and so and . The Fourier wave can be equally described by a wave length , and then we see that and .

The components of the particle velocity arewhere the spectral index was fixed at and thus we have . The vector whose components are denoted by takes values obeying a Rayleigh distribution. The wave phase vector, , given by takes random values on the interval . The components of the vector are calculated by means of , where is a random number in . is a fixed parameter with value .

##### 3.3. The Set-Up of the Particle Wind

We generate the initial wind particles by a traditional Monte Carlo scheme, in which particles are initially located randomly in the central volume of the cloud. Let , , and be random uniform variables taking real values within the interval . According to the fundamental probability conservation law for a system with spherical symmetry, we must have . By means of integration we obtain that the spherical coordinates of the particles are related to the uniform random variables by the following equations:where is the total mass contained in the central cloud region of radius and its mass density. It must be noticed that the first relation should be integrated numerically to obtain once has taken an allowed uniform random value.

In order to assign a velocity vector to each wind particle, let us consider the equation of mass conservation for a set of particles moving radially outwards; that is,In order to determine the wind velocities according to (7), where the radius ranges in the interval , the mass loss and the mass density will be considered as free parameters in the simulations; see Section 3.8. As the velocity magnitude diverges for particles around , we must then set up a cut velocity value, . Therefore, the velocity vector for the wind particle is given by , where is the speed of the particle obtained by (7) and is the unit vector pointing radially outwards.

Of course, there are other possibilities, which will be considered elsewhere: one is to fix the radial density and/or the velocity profile in order to obtain the mass loss as a result. For example, for modeling an expanding region, the authors of [2] proposed a complicated velocity function, which gives a constant expansion velocity at the last stages of time evolution, so that the average velocity of the shocked shell considered by [2] is or km/sec.

##### 3.4. Resolution and Thermodynamical Considerations

Following [10, 11], in order to avoid artificial fragmentation, the code must fulfil certain resolution criteria, imposed on the Jeans wavelength , which is given bywhere is the instantaneous speed of sound and is the local density. To obtain a more useful form for a particle based code, the Jeans wavelength is transformed into the Jeans mass given byIn this paper, the values of the density and speed of sound are updated according to the following equation of state:as proposed by [12], where , and for the critical density we assume the value .

For the turbulent cloud under consideration, we have g, where we take .

In this paper, the mass of an SPH particle is g, so that and therefore the Jeans resolution requirement is very easily satisfied.

##### 3.5. Initial Energies

Following [13], the dynamical properties of the initial distribution of gas are usually characterized by , the ratio of the thermal energy to the gravitational energy, and , the ratio of the rotational energy to the gravitational energy.

In a particle based code, we approximate the thermal energy of the cloud by calculating the sum over all the particles described in Section 3.1; that is,where and are the pressure and density associated with particle by means of the equation of state given in (10). In a similar way, the approximate potential energy iswhere is the gravitational potential of particle .* For all the simulations considered in this paper*, the values of the speed of sound (see (3)) and the level of turbulence are chosen so that the energy ratios have the numerical values

##### 3.6. The Evolution Code

We carry out the time evolution of the initial distribution of particles with the fully parallel code, which is described in detail by [14]. is based on the method for computing the gravitational forces and on the standard SPH method for solving the Euler equations of hydrodynamics. incorporates the following standard features: (i) each particle has its own smoothing length ; (ii) the particles are also allowed to have individual gravitational softening lengths , whose values are adjusted so that for every time step is of order unity. fixes the value of for each time-step using the minimum value of the smoothing length of all particles; that is, if for , then .

The code has an implementation of a Monaghan-Balsara form for the artificial viscosity; see [15, 16]. The strength of the viscosity is regulated by setting the parameters and ; see in [14]. Here we fix the Courant factor to .

##### 3.7. The Detection of Accretion Centers

Let us now describe the code implemented to detect accretion centers in the models. Any gas particle with density higher than is a candidate to be an accretion center. We localize all the candidate particles for a given snapshot at time .

We then test the separation between candidate particles: if there is one candidate with no other candidates closer than , then this particle is identified as an accretion center at time . We define as the neighbor radius for an accretion center, given by , where is the minimum value of the smoothing length of all particles; that is, if with being the smoothing length of each particle , for in this way determines a set of particles which are within the sphere having this radius and whose center is the accretion center itself. All those particles will contribute to the mass and momentum of the accretion center.

We carry out the detection of accretion centers for all the models, using a fixed value of given by , as was done by [5].

The code described here is a first step towards the full sink-particle technique implementation already in use worldwide. Our code is still incomplete as we have not yet included all the particle tests suggested in [17, 18].

##### 3.8. The Models

The initial configuration of particles has many parameters and their values determine the models to consider as follows. To sum up the meaning of these parameters in order of importance, we make a list in Table 1, without any reference to some specific physical system. In the first line, determines the time when the wind appears within the turbulent cloud. In the second and third lines, and determine the velocity by (7). Next, the radius determines the initial spatial extension of the wind particles. () is the total number (the total mass) of wind particles within . is the peak wind velocity allowed at the initial simulation time to avoid divergences in the velocity given by (7). Finally, we need to define two angles in order to characterize the collimated jets: the polar angular aperture and the azimuthal angular aperture , centered around the -axis with , and the equator with , respectively.