Magnetically driven accretion disc winds: the role of gas thermodynamics and comparison to ultra-fast outflows

Winds are commonly observed in luminous active galactic nuclei (AGNs). A plausible model of those winds is magnetohydrodynamic (MHD) disc winds. In the case of disc winds from a thin accretion disc, isothermal or adiabatic assumption is usually adopted in such MHD models. In this work we perform two-dimensional MHD simulations implementing different thermal treatments (isothermal, adiabatic and radiative) to study their effects on winds from a thin accretion disc. We find that both the isothermal model and the adiabatic model overestimate the temperature, underestimate the power of disc winds, and cannot predict the local structure of the winds, compared to the results obtained by solving the energy equation with radiative cooling and heating. Based on the model with radiative cooling and heating, the ionization parameter, the column density and the velocity of the disc winds have been compared to the observed ultrafast outflows (UFOs). We find that in our simulations the UFOs can only be produced inside hundreds of Schwarzschild radius. At much larger radii, no UFOs are found. Thus, the pure MHD winds cannot interpret all the observed UFOs.

In the aspect of theories, the presence of wind can affect the dynamics and spectrum of the accretion flow. For example, the angular momentum of an accretion disc can be removed by winds, which drives gas fall onto the central object (e.g. Konigl 1989;Stone & Norman 1994;Bai & Stone 2013;Li & Begelman 2014;Zhu & Stone 2018;Li & Cao 2019). Winds can also influence the density and temperature profiles of an accretion disc, and subsequently change the emitted spectrum (e.g. Yuan 2003). Moreover, it is generally believed that winds from AGNs can affect the evolution of its host galaxy via the "AGN feedback" process (e.g. ★ E-mail:wxwang@shao.ac.cn (WW); dfbu@shao.ac.cn (DB); fyuan@shao.ac.cn (FY) Fabian 2012; King & Pounds 2015;Chen et al. 2022;He et al. 2022). Thus, the wind production and propagation in and from AGNs are of great interest.
It is well known for a long time that magnetic field plays a prominent role in driving winds. When the poloidal magnetic field at the disc surface is inclined at an adequately large angle from the rotation axis, winds can be driven by magneto-centrifugal force (e.g. Blandford & Payne 1982;Cao & Spruit 1994;Contopoulos & Lovelace 1994). If its toroidal component is dominant, winds can be driven by the magnetic pressure gradient force (e.g. Lynden-Bell 1996; Kato et al. 2004;Yuan et al. 2015). Alternatively, radiatively driven mechanism (e.g. Proga et al. 2000;Proga & Kallman 2004) and thermally driven mechanism (e.g. Begelman et al. 1983;Dyda et al. 2017;Waters & Proga 2018) are plausible for launching disc winds. In this paper, we only focus on the magnetically driven mechanism.
From the view of thermodynamics of disc winds, early analytical literature of MHD disc winds had assumed an isothermal/adiabatic process or a polytropic process with a parametrized index (e.g. Contopoulos & Lovelace 1994;Cao & Spruit 1994;Kudoh & Shibata 1997). It results in a quite simple description of the thermal properties. They also assumed self-similar solutions to make it feasible to analytically solve the set of MHD equations. Even in semi-analytical work the treatments for the thermal properties are limited (e.g. Li & Cao 2021). When an energy equation including radiative cooling and heating is considered, the MHD equations cannot be normalized to a dimensionless form; self-similar solutions are no longer valid. Therefore, numerical studies are invoked to self-consistently obtain the wind solutions when considering an appropriate energy equation.
In the case of winds launched from a cold thin disc, technical difficulties in numerically simulating a geometrically thin accretion flow have limited studies and understandings of the disc winds. Lots of efforts have been made to self-consistently solve the accretion process by global simulations (e.g. Zhu & Stone 2018;Mishra et al. 2020). However, the computational domain in these simulations is too narrow to obtain complete wind solutions. In addition, the scale height / of a typical standard thin disc with characteristic parameters 1 in a luminous AGN is in the order of 10 −3 , which has hitherto never been approached in global simulations. In the case of minimal scale height, the disc can be regarded as a boundary without being resolved if it is assumed that the presence of winds do not influence the disc structure.
In this work, we will show the difference of wind properties between the isothermal model, the adiabatic model and the radiative model (as our fiducial model). The purpose is to investigate whether the isothermal or adiabatic assumption on the disc winds are valid and adequate.
Apart from the thermodynamics of MHD disc winds, their observable properties are of great interest. The MHD wind model has been employed to explain the observed UFOs. Fukumura et al. (2015Fukumura et al. ( , 2018 utilize the self-similar MHD wind solution to interpret the UFOs in nearby quasars PG 1211+143 and PDS 456. The authors adjust wind model parameters to fit the observed spectra. They find that for proper parameters, the self-similar MHD winds can well produce the observed UFOs. Therefore, they believe that the UFOs are MHD driven accretion disc winds. Nevertheless, the self-similar solutions oversimplify the wind model. Whether the MHD winds are radially self-similar is unclear. In this case, numerical studies of the disc winds considering radiations are still necessary for the realistic case especially where self-similarity is violated.
According to these considerations, the first aim of the present paper is to study the effects of different thermodynamics treatment on the properties of winds from a thin disc. The second aim is to study whether UFOs can be launched by magnetic field. The paper is organized as follows. Section 2 provides basic equations and assumptions of the model. Section 3 describes numerical setup. The results are presented in Section 4. We summarize our work in Section 6.

Governing Equations
Cold accretion disc in luminous AGNs is geometrically thin. Technically, it is arduous to resolve the thin accretion disc in global simulations. In this paper, we simulate winds launching from the surface of a thin disc. Same as the way that Ustyugova et al. (1999), Proga et al. (2000), Nomura et al. (2016) and Nomura & Ohsuga (2017) implement the disc, the main body of the thin disc is not resolved and not included in the computational domain; rather, the surface of the accretion disc is treated as the boundary of our simulation. Spherical coordinates ( , , ) are employed and we only simulate the region above the midplane. The plane of = 90 • corresponds to the disc surface. Following Nomura et al. (2016) and Nomura & Ohsuga (2017), we assume that the black hole is located at the location = 180 • , cos = −3.1 s , where is the Eddington ratio of the AGN. Observations imply that a compact hot corona should exist sandwiching the disc in the innermost region. We assume that the corona is a point source and is located at the origin ( = 0). 1 Black hole mass BH = 10 8 , the prescription of viscosity = 0.1 and mass accretion rate BH = 0.1 Edd .
The dynamics of the winds are governed by the following equations in general, where is mass density, v is the velocity and is the gas specific internal energy, is the gas pressure, f B = ( B,r , B, , B, ) is the Lorentz force per unit volume, g = ( , ) is the gravity whosecomponent arises because the black hole is not located at the origin and L is the net heating/cooling function. A Newtonian potential is utilized. The equation of state for perfect gas is adopted: Without loss of generality, the energy equation is represented in the form of Equation 5. When L = 0, it describes adiabatic flow; if Equation 5 is not solved but replaced by ∝ , it displays the case of isothermal process.

Radiative heating and cooling
The X-ray radiation from the corona heats, cools and photoionizes gas in winds through radiative processes of Compton scattering and photoionization. Additionally, winds lose energy via bremsstrahlung, line cooling and recombination cooling. Here we use exactly the same cooling function as in Blondin (1994). In Equation (5) the net cooling rate is written as where Γ Compton is the heating or cooling rate of Compton scattering Γ X is the photoionization heating/recombination cooling, and Λ is the cooling rate of bremsstrahlung and line cooling In the above equations, X is the characteristic temperature of the X-ray radiation and is a parameter to represent optically thin or thick atmosphere. X is set to 10 8 K corresponding to 10 keV bremsstrahlung radiation from the central corona. Assuming the wind is optically thin to its own radiation, we have = 1.

Accretion disc treatment
In this paper, we simulate winds launching from the surface of a thin disc. The surface of the disc is implemented in the first layer of grids of the computational domain above the equatorial plane (with = 0). The black hole is located below the surface of accretion disc and following Nomura & Ohsuga (2017) we put the black hole at the location ( = 0 = −3.1 s , = ). 0 corresponds to the scale height of the accretion disc at the inner radial boundary, which is 30 s . According to the standard thin disc model (Shakura & Sunyaev 1973), the azimuthal velocity equals to the Keplerian velocity at every radius. We assume that the density of the surface of the disc does not vary with time. In this sense, the mass taken away by wind will be immediately replenished. This setting is similar as that in Nomura et al. (2016). In that paper, the authors assume that the density of the surface of the disc at the midplane is constant with both time and radius. In our present paper, we assume that the density of the surface of the accretion disc varies with radius, but not with time. In order to set a reasonable non-evolving accretion disc surface density, we employ the accretion disc theory. We set it according to the standard thin accretion disc value given by (e.g. Nomura & Ohsuga 2017;Kato et al. 2008 where crit = 18( BH / ) 2/21 ( / ) 16/21 s is the critical radius. Inside crit , radiation pressure dominates gas pressure; outside crit , gas pressure dominates radiation pressure. According to the parameters adopted in the present paper, crit ∼ 154 s . is the radiative efficiency of the accretion disc. We note that the radial velocity at the midplane is not calculated according to the thin disc theory. Initially, we set it to be 0. It can evolve with time. Therefore, the accretion rate of the midplane disc is not a constant with time. In our fiducial model, we set = 0.1. This does not mean that the accretion disc has an accretion rate of 0.1 times Eddington rate. We just use Equation (11) to set a reasonable density of the surface of a thin disc. The temperature profile is, where in is the inner radius of the disc and is set to be 3 s , in is the effective temperature at in , which is set to satisfy D = ∫ out in 2 4 d , where is the Stefan-Boltzmann coefficient. Under the assumption made in Section 2.3, the outer radius out can be derived from in ( out / in ) −3/4 = 3 × 10 3 K. The evolution of disc is neglected. The density and the temperature of the disc are fixed, but we allow the velocities ( , , ) and fields ( , , ) to float.

Initial and boundary conditions
The gas in the computational domain is initially set to be vertically in hydro-static equilibrium without radiative force. Then the initial density distribution above the surface of the thin disc reads, Here is the gravitational constant and is the sound speed. Since the density in the vertical direction is exponentially decreasing, the ratio of density at equatorial plane to density at high latitude can increase with height and approaches infinity. In numerical simulations if the density discrepancy in the computational domain is too large, the simulation is usually quite hard to run. In order to overcome this problem, we set a density floor to be floor = 1.49 × 10 −20 g cm −3 . Note that the value of density floor is 11 orders of magnitude smaller than the density at the equatorial plane, and therefore, it has negligibly small effects on our results. Initially, the gas has Keplerian rotational velocity. The initial radial and directions velocities are set to be zero. So the initial velocity is 0 = (0, 0, sin √︃ 2 0 + 2 +2 0 cos ). The initial temperature is vertically isothermal so that 0 ( , ) = = /2 ( sin ). The initial internal energy density can be derived by using 0 ( , ) and 0 ( , ). Additionally, initial ionization parameter is needed to determine the opacity at the first time step of simulations, and the initial opacity is needed to determine the ionization X-ray flux at the fist time step. Due to the low density of background, the gas is assumed to be initially photoionized by X-ray without attenuation.
The initial magnetic field is defined by using the vector potential below (Zanni et al. 2007), The parameter 0 determines the strength of the initial magnetic field. The inclination angle of the initial magnetic field lines with respect to the equatorial plane is determined by the parameter . In this paper, we set = 0.4. The initial plasma , defined as the ratio of gas pressure to magnetic pressure, is shown in Fig. 1. We impose outflow boundary conditions at the inner and outer radial boundaries. Axi-symmetry boundary condition are set at = 0. At = /2, we use reflecting boundary conditions.

Models and parameters
We assume the mass of a non-rotating black hole BH = 10 8 and the radiative efficiency of the accretion disc is = 0.06. The disc luminosity is set to be 0.1 times the Eddington luminosity. In this paper, we do not simulate the accretion process. The effects of winds on the accretion disc are neglected. In reality, in an accretion system, the presence of winds can affect the accretion rate of the disc. In the hydrodynamic simulations, Nomura et al. (2020) have examined the effects of mass loss on the disc. The authors utilized iterative method to adjust the accretion rate. The iteration which can reach convergence is a plausible method. However, in the MHD case, the situation becomes complicated, since the magnetic field evolves with the fluid. It is difficult to determine a self-consistent magnetic field when one does iteration.In the present MHD simulations, we do not use the iterative method to adjust the accretion rate.
In order to compare the different thermodynamics influencing the disc winds, we consider adiabatic and isothermal processes as well. For the isothermal model, the internal energy does not evolve so that  the cooling function L is set to zero. Without radiation, the ionization parameter can not be defined. The set of equations remain closed. Models with referred parameters are listed in Table 1.

RESULTS
In this section, we present simulation results and comparison with observations. It is organized as follows. Section 4.1 describes the time-averaged quantities; real winds are identified in Section 4.2; in Section 4.3 we show the detailed properties of real winds and demonstrate the dependence on different thermodynamics; Section 4.4 describes the mechanisms of driving disc winds for the fiducial model FID; Section 4.5 compares the observable features of disc winds with the UFOs.

Time-averaged quantities
The simulations are executed in several tens of free-fall time-scale ff,o at the outer radial boundary o = 1500 s . The mass loss rate defined as for each model is shown in Fig. 2. The reason for the integral range extending to 80 • instead of /2 is that the flow near the equatorial plane varies dramatically as turbulence. As the flow fluctuates, 80 • is Figure 2. Evolution of mass loss rate at the outer radial boundary in unit of Eddington limit Edd = Edd / / 2 , where Edd is the Eddington luminosity for BH = 10 8 . a rough but tolerant approximation. It reflects more realistic outflow rate outside of the turbulent region. The mass loss rate fluctuates within one order of magnitude, which implies that the system evolves to a quasi-steady state. Thus, it is reasonable to analyse the timeaveraged quantities. Fig. 4 displays time-averaged quantities of each model. As the density and streamlines in Fig. 4 show, the basic morphology in the simulated region is similar to each other in different thermodynamic models. The whole region is divided into two main parts: laminar wind close to the pole and turbulent flow close to the equator. The turbulent region is gas puffed up from the disc surface and the laminar wind is launched from the turbulent region. The fraction of each part differs in different models.
The temperature in isothermal models (ISO6 and ISO7) and adiabatic model (ADB) is much higher than that in the fiducial model (FID) globally. The turbulent region in the former three models is more like hot corona. The compression work accounts for the high temperature of model ADB, which is a natural result of adiabatic process. However, in the fiducial model the turbulent gas resembles cold flow. It indicates that the radiative cooling is rather efficient in dense region.
The magnetic fields are ordered where the gas is sparse, while they are entangled where the gas is dense. The ordered large-scale field lines corresponding to outflow streamlines reflect the ideal-MHD scheme. As no physical diffusion of magnetic field is considered in the present paper, the entangled field might be adjusted by the numerical diffusion. The gas internal energy is much smaller than magnetic energy. At first glance, one would expect that the magnetic field can not be deformed. However, for a thin disc, the gas temperature is much smaller than the Virial temperature. The sound speed is much smaller than the Keplerian velocity. Even though the gas internal energy is much low, the kinetic energy of gas is much high. We compare the gas kinetic energy to the magnetic energy in Fig. 5. It can been seen in our models, in most region of the computational domain, the kinetic energy is much higher than the magnetic energy. Therefore, the magnetic field can be deformed by the motions of the gas. Generally, the essential morphologies of the magnetic fields as well as the fluids are similar in models ISO6, ISO7 and FID. Model ADB has hotter and more stochastic flow than others, which leads to larger spacial fraction of turbulent region.
The vertical magnetic field can subtract angular momentum from the disc. The timescale for the angular momentum subtraction is, where denotes the disc scale height. In Fig. 3, we plot the ratio of the angular momentum subtraction timescale to the disc orbital timescale. The result shows that the angular momentum subtraction timescale is at least more than 2 orders of magnitude times larger than the orbital timescale. For most of our computational domain, the simulation time is much shorter than the angular momentum subtraction timescale, which indicates that the specific angular momentum of disc is not subtracted too much within the evolution time.

Real wind identification
The streamlines overlaying on density in Fig. 4 suggest that not all materials from the disc surface would continuously flow outwards. Only those materials that with positive Bernoulli constant will flow out to infinity. In this sense, to identify the real winds, we plot Bernoulli parameter calculated by where ℎ ≡ + is the specific enthalpy. For isothermal case, the specific enthalpy is determined by the sound speed solely ℎ = −1 2 s . As Fig. 6 shows, the Bernoulli parameter close to the rotational axis is nearly constant with radius. It coincides with the radially laminar streamlines of winds shown in Fig. 4, indicating the winds are in quasi-steady state. Positive values (reddish color) represents potential outflows that have sufficient energy to surpass gravitational potential energy. The combination of the positive and laminar streamlines distinguishes the real wind from the turbulent gas which might fall back rather than flow out to the infinity. The density contour of = 10 −15 g cm 3 in Fig. 6 roughly separates the whole region into the turbulent flow and the laminar outflow. Thus, the position of this density contour is regarded as the wind base. For the model FID, to avoid the bulk of flow with negative at around (900 s , 600 s ), we constrain the density contour by substituting the polar angle outside 600 s with the angle at 600 s (see the right-bottom panel in Fig. 6). Fig. 7 shows radial profiles of outflow mass flux (mass loss rate), momentum flux and energy flux of real winds. They are calculated as follows,

Wind properties dependence on the thermodynamics
where and are time-averaged density and radial velocity, and max ( ) is the angle of the interface between the laminar flow and the turbulent region at each radius. Quantities at each radius are integrated from the polar axis ( = 0) to the wind base ( max ( )) which is roughly identified by the density contours demonstrated in Fig. 6 in Section 4.2.
The outflow mass flux rapidly increases within about 200 s and remains nearly constant at outer radii. The model FID with cooling function gives the largest value of each flux than other models, which suggests more efficient wind production in model FID. The fluxes in model ISO6 and ISO7 are several times lower than those in model FID. The deviation of model ADB can go up to several orders of magnitude. The large deviation is consistent with the morphologies of magnetic fields and flows that implied in Fig. 4. Here note that the ratio of outflow rate and accretion rate exceeds unity in models ISO6, ISO7, and FID. It is not contradictory with mass conservation. In the simulation, the accretion rate is not self-consistently calculated from an evolving disc but fixed to a constant value corresponding to the assumption of a steady disc. The assumed constant density profile of the disc surface only works as gas supply. The fixed density profile can produce a corresponding accretion rate much larger than 0.1 Edd , which leads to a mass loss rate greater than 0.1 Edd .
Initially, our magnetic field has vertical component. Such a configuration is expected when the strong disc wind is emanated. The upper panel in Fig. 7 shows the radial dependence of the wind mass flux. It is shown that the outflow rate keeps increasing outside 100rs. In the fiducial model for instance, the mass flux from 100 s to 200 s increased from around 0.1 to 0.5 Eddington rate. Inside 100 s , the contribution of wind mass flux is smaller than 0.1 Eddington rate. In other words, most of the wind mass flux comes from the region outside 100 s . For the adiabatic model, the mass flux of wind keeps increasing from the inner boundary to the outer boundary. In this model, winds can be generated at every radii. Therefore, the assumption that the initial magnetic field has a vertical component is reasonable. Fig. 8 demonstrates the angular distribution of winds at the outer radial boundary. The winds in model FID can be accelerated up to 0.3 ; the maximum velocity of winds in other models are also in the order of 0.3 . The winds with sufficiently large velocity have low density. They are close to the polar axis rather than the disc plane, which is a distinguished feature compared to the line-driven winds Nomura et al. (2016); Nomura & Ohsuga (2017).
The adiabatic model gives the lowest velocity, mass flux and spatial fraction of winds and overestimates the temperature. It deviates the most from the fiducial model. The essential reason would be that the overheated flow suppresses the effects of magnetic fields, leading to a more stochastic flow and more tangled field lines (see the third column in Fig. 4). However, large-scale ordered field lines are substantial to the acceleration of winds via the Blandford & Payne   mechanism. Less powerful wind is attributed to the lack of large-scale ordered field lines.
According to the first law of thermodynamics, it is apparent that the increase in the internal energy is attributed to the compression work for an adiabatic model; for isothermal models, constant internal energy constrains an artificial cooling/heating to exactly balance the compression/expansion work. To further examine the thermodynamics, we calculate the divergence of the fluid velocity to trace the compression and expansion. The results are shown in Fig. 9. Fig. 9 illustrates that most compression occur in the turbulent region, especially at the wind base. For model ADB, the compression operates in almost the whole region. It confirms that compression plays a dominant role during the evolution in model ADB leading   Fig. 4. Compared to that in models ADB, ISO6 and ISO7, the gas within >∼ 60 • in model FID is much cooler even there is compression work. It indicates the radiative cooling in the turbulent region is dominant in model FID compared to the compression work and radiative heating.
In the wind region ( <∼ 60 • ) in all models, fluids are expanding. The expansion does negative work and cools the flow. In the model FID, compared to the turbulent region, the relatively high temperature of winds indicates that radiative heating is operative in the wind region. It also coincides with the region of high ionization parameter (see the colormap in Fig. 10a). This is a natural consequence of X-ray heating. Another eminent feature in the model FID is that the temperature of winds is radially stratified, which cannot be described by the adiabatic and isothermal models.
To summarise, cooling function plays an essential role in determining the temperature of disc winds. Adiabatic assumption incorrectly overestimates the temperature. Isothermal assumption erases the local structure and underestimate the power of winds. Although the isothermal model with an appropriate temperature (10 7 K in this work) produces similar wind dynamics to the fiducial model, this unique value of temperature cannot be provided as a priori parameter. In the following sections, we represent further detailed results of disc winds for model FID in which a cooling function is employed to mimic a more realistic case.

Mechanisms of driving disc winds in model FID
To quantitatively study the winds accelerations in model FID, we have calculated various forces per unit mass exerted along the two streamlines "A" and "B" (refer to Fig. 10a; the latter is closer to the surface of turbulent region). The results are shown in Fig. 10b. The Lorentz force is The first term on the right hand side of this equation is the magnetic tension force and the second term is the magnetic pressure gradient force. The magnetic tension force is negligibly small, so only the Lorentz force and the magnetic pressure gradient force are presented. Along both streamlines, the winds are accelerated by both the centrifugal and magnetic pressure gradient forces. The centrifugal force is slightly stronger than the magnetic pressure gradient force. The gas pressure gradient force is negligible as expected. It is a natural result for the cold accretion disc, where gas temperature is significantly smaller than the Virial temperature. The acceleration along the streamline B reminds us the wind acceleration model in Blandford & Payne (1982). We find that the magnetic field lines do have an appropriate inclination angle with the rotational axis (recall the magnetic field lines in Fig. 4), as required by the Blandford & Payne mechanism to work.

Comparison to observations: UFOs
The disc winds produced in our simulations have inhomogeneous properties and extensive range of the ionization parameter. Not all winds meet the conditions of the observed UFOs. To compare the observable variables, we exploit the constraints of the observed features from Gofford et al. (2015) on the disc winds. Generally, UFOs have velocity higher than 10 4 km s −1 , ionization parameter in the range 2.5 < log < 5.5 and column density in the range 22 < logN H < 24. The angular distribution of these UFOs is shown in Fig. 11. Similarly, the UFOs locate close to the pole rather than the equator. They cover the region between = 4.7 • and 45.6 • . The corresponding solid angle Ω ≈ 3.85; the fraction of Ω to a whole sphere, Ω/4 ≈ 30.6 per cent. It is close to the lower limit of the fraction of UFO in AGNs (Gofford et al. 2015).
To obtain the mass flux, momentum flux and kinetic power of UFOs, we integrate the Equation (18)-(20) over the flows which meet the above conditions. The comparison of the UFOs found in model FID to observations is shown in Fig. 12. The cut-off at 600 s implies that in model FID no wind meets the condition of UFOs outside 600 s . The properties of UFOs within that radius are well in the range given by observations and close to the averaged values.
Combining the ionization parameter (colormap in Fig. 10a) and the temperature (second column in Fig. 4) of model FID, the reason of the missing UFOs at large radii is the low ionization state. However, it has been reported that the UFOs can exist within 10 2 − 10 4 s (Gofford et al. 2015). Thus, some physics might be missing in MHD disc winds models. A plausible conjecture for the missing physics would be line force. As the winds are sparse and cool, line force might operate on the winds to change their properties.

Dependence on midplane gas density and X-ray flux
In our fiducial model, the gas density at the midplane is set corresponding to that of a thin disc which has Eddington ratio of 0.1. We do a test simulation EPS-2. In this simulation, the midplane density is also obtained by Equation (11) by setting to 0.01. Also, in the test simulation, the X-ray luminosity is 10 times lower than that in the fiducial model. In Fig. 13, we plot the midplane density in our fiducial model (blue line) and test simulation EPS-2 (yellow line). It can be seen that inside 200 , where the UFOs are launched, the midplane density in model EPS-2 is much higher than that in our fiducial model. Outside 200 , the density in model EPS-2 is much lower. Fig. 14 (black lines) shows the UFOs in model EPS-2. The UFOs are present only inside 100 . The reason of the absence of UFOs outside 100 is as follows. We find that the wind density in the region < 45 • in model EPS-2 is comparable to that in the fiducial model.  . Column density (black solid line), column density of gas with 2.5 < log < 5.5 (black dashed line), and maximum velocity of the gas with 2.5 < log < 5.5 (red line). The horizontal axis is viewing angle (or angle). The vertical axis on the left-hand side shows the value of velocity (red line). The vertical axis on the right-hand side shows the value of column density (black solid and black dashed lines). The horizontal solid line marks the velocity of 10000km/s. The horizontal dashed lines mark column density of 10 22 cm −2 and 10 24 cm −2 . The vertical dotted lines mark = 4.7 • and 45.6 • where the disc winds satisfy the UFO conditions in between.
However, the X-ray flux in model EPS-2 is 10 times lower. Thus, the winds at large radii in model EPS-2 have too low ionization parameter to satisfy the conditions of the observed UFOs. Both the mass flux and kinetic power of UFOs in model EPS-2 are several times lower than those in the fiducial model. The reason is the same that in model EPS-2, the ionization parameter is relatively low. Therefore, the winds which satisfy the conditions of UFOs are reduced.
Observations show that the UFOs are located at the region of 100 -10000 (Gofford et al. 2015). In model EPS-2, the UFOs are located inside 100 . Therefore, it is not consistent with observations. We also find that the maximum velocity of UFOs in model EPS-2 is quite similar as that in the fiducial model. This is not consistent with observation which show that the UFOs velocity is faster for higher luminosity.

Dependence of UFOs on magnetic field strength
In our fiducial model, when the quasi-steady state is achieved, at the midplane, the plasma beta ∼ 10 −4 . We perform a test simulation to study the dependence of UFOs on the magnetic field strength. In our test simulation, when the quasi-steady state is achieved, the midplane plasma beta ∼ 10 −2 . The magnetic field strength in the test simulation is much lower. The mass flux, momentum flux and kinetic power of UFOs are plotted in Fig. 14 (blue curve). Comparing figures 12 with 14, we can see that, in the test simulations, the mass flux of UFOs is one order of magnitude smaller than that in the fiducial model. The kinetic power of UFOs in the test simulation is more than 1 order of magnitude smaller than that in the fiducial model. In future, further complete parameter survey is necessary to draw a quantitatively solid conclusion of the dependence of UFOs on magnetic field strength.  Gofford et al. (2015) are located inside the shade region. In other words, the upper and lower boundaries of the shade region correspond to the maximum and minimum values of properties of the observed UFOs. The horizontal dotted line is the averaged values from observations (Gofford et al. 2015). The vertical dotted line denotes the innermost radii of oberved UFOs (Gofford et al. 2015), i.e. about 100 s .

SUMMARY
We have studied the difference of wind properties by considering different thermodynamics: isothermal, adiabatic and radiative. The wind properties strongly depend on the treatment of thermodynamics. The isothermal and adiabatic models underestimates the wind power and overestimates the wind temperature. Additionally, the isothermal and adiabatic assumptions erase the local structure of disc winds. We suggest future work consider energy equation with radiative cooling and heating to obtain accurate wind properties. In this way the wind model can be more appropriate to directly compared to the observations. The fiducial model with radiative cooling and heating can produce powerful radially stratified winds. The winds inside 600 s are consistent with the observations of UFOs, while no UFOs are found outside 600 s . The lack of UFOs with appropriate ionization pa-  rameter, column density and velocity at large radii implies that extra physics need to be considered. In future, seeking these extra physics would be meaningful to the interpretation of the observed winds from a thin accretion disc.