Change Light Theme Back to Top

P464 - MHD Project

• Implemented by Adhilsha A (2011006)
• Last updated on 28-04-2024

This website follows the updates on my project on MHD. The project is about the study of the Galactic Dynamos. Each task in the project gets one step closer to the final goal. The project is divided into 3 tasks.

---

Table of Contents

---

• P464 - MHD Project
• TASK 1: Solving the diffusion equation
  • 1.1 Introduction to Diffusion equation
  • 1.2 Equation forms
  • 1.3 Boundary conditions
  • 1.4 Implementation
    • 1.4.1 Spatial derivative
    • 1.4.2 RK4 time stepping
  • 1.5 Constants and parameters for simulation
    • 1.5.1 Study of evolution of Magnetic field strength
    • 1.5.2 Study of evolution of pitch angle.
    • 1.5.3 Study of Magnetic field evolution with different seed fields and different boundary conditions
• TASK 2: (TBD)
  • 2.1 Introduction and equations
  • 2.3 Implementation
  • 2.4 Results
    • 2.4.1 Sample implementation
    • 2.4.2 Decay rate at different spatial steps
    • 2.4.3 Pitch angle evolution
    • 2.4.4 Magnetic field evolution with different seed fields
    • 2.4.5 Magnetic field evolution with different boundary conditions
    • 2.4.6 Finding critical dynamo number
• Main Project
  • Abstract
  • 1. Introduction
  • 2. Methods
    • 2.1. Finite difference
    • 2.2. Time-stepping using RK4
    • 2.3. Implementation of BCs
    • 2.4. Calculation of Global growth rate
    • 2.5. Intuition behind Velocity profiles
    • 2.6. Addition of Disk-flaring
  • 3. Simulation and results
    • 3.1. Effect of non-zero vertical outflow
    • 3.2. Effect of different profiles of vertical outflow
    • 3.3. Disk flaring and vertical outflow
  • Conclusions
• References

---

TASK 1: Solving the diffusion equation

1.1 Introduction to Diffusion equation

The first task is to solve the diffusion equation in 1D. The induction equation is given by:

$$\frac{\partial \overline B_r}{\partial t} = - \frac{\overline V}{r} \frac{\partial}{\partial r}(r\overline B_r) - \frac{\partial}{\partial z}(\overline V_z \overline B_r) - \frac{\partial}{\partial z}(\alpha \overline B_{\phi}) + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_r)\right) + \frac{\partial ^2 \overline B_r}{\partial z^2} \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = - q\Omega \overline B_{r} - \frac{\partial}{\partial r}(\overline V_r \overline B_{\phi}) - \frac{\partial}{\partial z}(\overline V_z \overline B_{\phi}) + \frac{\partial}{\partial z}(\alpha \overline B_{r}) + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_{\phi})\right) + \frac{\partial ^2 \overline B_{\phi}}{\partial z^2} \right]$$

where $\overline B_r$, $\overline B_{\phi}$ and $\overline B_z$ are the components of the magnetic field, $\overline V_r$, $\overline V_{\phi}$ and $\overline V_z$ are the components of the velocity field.

For this task, we solve the diffusion equation omitting the $\nabla \times (\overline V \times \overline B)$ and $\alpha$ terms. The equation we have to solve becomes:

$$\frac{\partial \overline B_r}{\partial t} = \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_r)\right) + \frac{\partial ^2 \overline B_r}{\partial z^2} \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_{\phi})\right) + \frac{\partial ^2 \overline B_{\phi}}{\partial z^2} \right]$$

Since both $\overline B_r$ and $\overline B_{\phi}$ have same form of equation, we generalize the form of the equation to:

$$\frac{\partial \overline B}{\partial t} = \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B)\right) + \frac{\partial ^2 \overline B}{\partial z^2} \right]$$

1.2 Equation forms

In case of $\overline B_r(z)$ and $\overline B_{\phi}(z)$, the equations become:

$$\frac{\partial \overline B}{\partial t} = \eta_t \frac{\partial ^2 \overline B}{\partial z^2}$$

In case of $\overline B_r(r)$ and $\overline B_{\phi}(r)$ under no-$z$ appproximation, the equations become:

$$\frac{\partial \overline B}{\partial t} = \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B)\right) - \frac{\pi ^2}{4 h^2}\overline B \right]$$

On expanding the equation, we get:

$$\frac{\partial \overline B}{\partial t} = \eta_t \left[ \frac{\partial ^2 \overline B}{\partial r^2} + \frac{1}{r} \frac{\partial \overline B}{\partial r} - \frac{\overline B}{r^2} - \frac{\pi ^2}{4 h^2}\overline B \right]$$

These are the equations to be solved in this task.

1.3 Boundary conditions

We will be using the simplest boundary conditions for this task. The boundary conditions are:

$$\overline B_r = \overline B_{\phi} = 0 \quad \text{at} \quad (z = -h, h)\text{ or }(r = 0,r_{max})$$

We will be taking dimensionless units with $h = 1$, $\eta_t = 1$.

In $r$ implementation, we will be starting with $r$ = 0.01 since the value $r = 0$ is can possibly cause division by zero errors.

1.4 Implementation

For our implementation we will be using finite difference method for defining the spatial derivatives and RK4 for the time stepping / evolution. The code can be observed at task1_code.

1.4.1 Spatial derivative

The spatial derivatives are calculated using the finite difference method. here, we are using 6th order finite difference method for the spatial derivatives. The 6th order finite difference method is given by:

$$\frac{d f}{d x} = \frac{1}{60\delta x}\left( - f_{i-3} + 9f_{i-2} - 45f_{i-1} + 45f_{i+1} - 9f_{i+2} + f_{i+3} \right)$$

$$\frac{d^2 f}{d x^2} = \frac{1}{180\delta x^2}\left( 2f_{i-3} - 27f_{i-2} + 270f_{i-1} - 490f_{i} + 270f_{i+1} - 27f_{i+2} + 2f_{i+3} \right)$$

Other orders can be used and is available for use in the code.

For the sake of obtaining smooth derivatives as well for the sake of maintaining the boundary conditions, we use ghost zones. By default, 'relative anti-symmetric' ghost zones are used where the ghost cell at $i$ distance away from boundary $b$ is given by:

$$f_{b-i} = 2f_{b} - f_{b+i} \quad \text{at the start of the spatial domain}$$

$$f_{b+i} = 2f_{b} - f_{b-i} \quad \text{at the end of the spatial domain}$$

where $i = 1, 2, 3, ..., n$ ($n$ is the half the order of the finite difference method).

Other ghost zones types, such as 'symmetric' and 'anti-symmetric' can also be used and is available for use in the code.

Using the 6th order finite difference method and 'relative anti-symmetric' ghost zones, a following magnetic field and the spatial derivatives are obtained (test case):

As you can see the spatial derivatives are smooth and the second derivative is zero at the boundaries as expected from the 'relative anti-symmetric' ghost zones. This works exceptionally well for the $z$ implementation. For the $r$ implementation, the ghost zones cannot help much in keeping the boundary conditions constant as the equation contains first and second derivatives as well as the function itself. Hence, the boundary conditions are not maintained well in the $r$ implementation. This is a known issue and is being worked on.

1.4.2 RK4 time stepping

The time stepping is done using the RK4 method. Our evolution function can be given of form:

$$\frac{\partial \overline B}{\partial t} = f(\overline B, t)$$

Then, the RK4 method evolution for a single time step can be given by:

$$\kappa_1 = \delta t f(\overline B, t)$$

$$\kappa_2 = \delta t f(\overline B + \frac{\kappa_1}{2}, t + \frac{\delta t}{2})$$

$$\kappa_3 = \delta t f(\overline B + \frac{\kappa_2}{2}, t + \frac{\delta t}{2})$$

$$\kappa_4 = \delta t f(\overline B + k_3, t + \delta t)$$

$$\overline B(t + \delta t) = \overline B(t) + \frac{1}{6}(\kappa_1 + 2\kappa_2 + 2\kappa_3 + \kappa_4)$$

The thing to note here about the generalized form, with regards to the RK4 method, is that:

• The spatial derivatives do not change values with the change ($\overline B \rightarrow \overline B + \kappa/2$) in the RK4 method. This is due to the symmetry of the spatial derivative equations. Hence, the spatial derivatives are calculated only once at the start of the time step. Thus, it is not mentioned in the equation form.

Our method, as mentioned above, is to find the spatial derivatives using finite difference method at a certain time step and then use the RK4 method to evolve the system to the next time step. This is done repetitively to get the solution at all time steps.

1.5 Constants and parameters for simulation

For $z$ implementation, we will be using the following parameters:

• $z\ \epsilon\ (-1, 1)$ # range of z
• $t\ \epsilon\ [0, 1]$ # range of time
• $N_z = 80$ # number of spatial steps
• $N_t = 5000$ # number of time steps
• $order = 6$ # order of the finite difference method
• ghostzone_type = 'relative anti-symmetric' # type of ghostzone

For $r$ implementation, we will be using the following parameters:

• $r\ \epsilon\ (0.01, 8.01)$ # range of r
• $t\ \epsilon\ [0, 1]$ # range of time
• $N_r = 100$ # number of spatial steps
• $N_t = 4000$ # number of time steps
• $order = 6$ # order of the finite difference method
• ghostzone_type = 'relative anti-symmetric' # type of ghostzone

The seed field used for the z-implementation is:

$$B_{\phi}(z) = 3 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 1 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 2 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

whereas the seed field used for the r-implementation is:

$$B_{\phi}(r) = 3 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 1 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 2 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

The results are:

1.5.1 Study of evolution of Magnetic field strength

For this part, I am setting $\overline B_r = 0$ and evolving $\overline B_{\phi}$ over time. The magnetic field strength, $|B_{\phi}|$, is given at different points in the spatial domain as a function of time. The results are (N.B. The decay rate is given in the legend off the plots):

from these results, taking the middle point in the spatial domain, the magnetic field strength, $|B_{\phi}|$, is given as a function of time in log-scale. Then it can be fitted to find the decay rate as well. The results are:

1.5.2 Study of evolution of pitch angle.

For this part, I am setting the seed $\overline B_r$ and $\overline B_{\phi}$ to non-zero values and evolving them over time. The pitch angle is given by:

$$p = \tan^{-1}\left(\frac{\overline B_{r}}{\overline B_{\phi}}\right)$$

The seed fields taken for this simulation are, for the z-implementation:

$$\text{{B}}_{\phi}(z) = 3 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 1 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 2 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

$$\text{{B}}_{r}(z) = 3 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 2.5 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 1 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

and for the r-implementation:

$$\text{{B}}_{\phi}(r) = 3 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 1 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 2 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

$$\text{{B}}_{r}(r) = 3 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 2.5 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 1 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

The time evolution of the fields are given below for reference:

The results of the pitch angle evolution at certain specific spatial steps were obtained as well. The results:

To get another perspective, the pitch angle variation over spatial domain at different time steps is also found. The results are:

1.5.3 Study of Magnetic field evolution with different seed fields and different boundary conditions

Here, I did two trials, first one with zero boundary conditions. Then, the second one with non-zero boundary conditions and different seed fields from the previous trial.

The trial 1 seed field with zero boundary conditions is, for the z-implementation:

$$B_{\phi}( z) = 3 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 1 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 2 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

$$B_{\phi}( z) = 3 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 2.5 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 1 \sin \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

and for the r-implementation:

$$B_{\phi}( r) = 3 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 1 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 2 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

$$B_{r}( z) = 3 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 2.5 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 1 \sin \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

The results are:

Here, with zero boundary conditions, the boundaries are maintained well in the z-implementation. The reason for this is 'relative anti-symmetric' ghost zones which keeps the $\frac{\partial^2 \overline B}{\partial z^@}$ zero at the boundaries.

However, the boundaries are not maintained well in the r-implementation, even though the it is hardly noticeable in the plots. This is because the magnetic field is decaying and the boundary condition is already zero. Hence, the boundary condition is 'seen' to be maintained well. But, the difference is noticeable in the next trial.

With the different seed fields used in the implementations, all magnetic field are seen to decay with time, as expected.

The trial 2 seed field with non-zero boundary conditions is, for the z-implementation:

$$B_{\phi}( z) = 3 \cos \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 1 \cos \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 2 \cos \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

$$B_{r}(z) = 3 \cos \left( \frac{{z - z_i}}{{z_f - z_i}} \pi \right) + 2.5 \cos \left( \frac{{z - z_i}}{{z_f - z_i}} 2\pi \right) + 1 \cos \left( \frac{{z - z_i}}{{z_f - z_i}} 5\pi \right)$$

and for the r-implementation:

$$B_{\phi}( r) = 3 \cos \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 1 \cos \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 2 \cos \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

$$B_{r}( r) = 3 \cos \left( \frac{{r - r_i}}{{r_f - r_i}} \pi \right) + 2.5 \cos \left( \frac{{r - r_i}}{{r_f - r_i}} 2\pi \right) + 1 \cos \left( \frac{{r - r_i}}{{r_f - r_i}} 5\pi \right)$$

The results are:

Here in $z$ implementation, even the non-zero boundaries are maintained well as expected.

As you can see in the $r$ implementation, the non-zero boundary conditions are not maintained at all. This is because the ghost zones cannot help much in keeping the boundary conditions constant as the equation contains first and second derivatives as well as the function itself. The 'relative anti-symmetric' ghost zones can only keep the second derivative zero at the boundaries.

The decay with of magnetic field with time is evident here, even with different seed fields used in the implementations.

---

TASK 2

2.1 Introduction and equations

From the original induction equation, we omit the terms involving $\overline V$ from both equations and the $\alpha$ term from the equation for $B_{\phi}$, in order to get a $\alpha\omega$ dynamo. Task 2 will be focused on solving these equations. THe equations are:

$$\frac{\partial \overline B_r}{\partial t} = - \frac{\partial}{\partial z}(\alpha \overline B_{\phi}) + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_r)\right) + \frac{\partial ^2 \overline B_r}{\partial z^2} \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = - q\Omega \overline B_{r} + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_{\phi})\right) + \frac{\partial ^2 \overline B_{\phi}}{\partial z^2} \right]$$

By applying no-$z$ approximation and expanding, equations become:

$$\frac{\partial \overline B_r}{\partial t} = - \frac{2}{\pi h}(\alpha \overline B_{\phi}) + \eta_t \left[ \frac{\partial ^2 \overline B}{\partial r^2} + \frac{1}{r} \frac{\partial \overline B}{\partial r} - \frac{\overline B}{r^2} + \frac{\pi^2 }{4h^2}\overline B_r \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = - q\Omega \overline B_{r} + \eta_t \left[ \frac{\partial ^2 \overline B}{\partial r^2} + \frac{1}{r} \frac{\partial \overline B}{\partial r} - \frac{\overline B}{r^2} + \frac{\pi^2 }{4h^2}\overline B_{\phi} \right]$$

2.2 Implementation

As done before, we will be using the dimensionless form of the equations. Hence it becomes:

$$\frac{\partial \tilde B}{\partial \tilde t} = - \frac{2}{\pi}(\alpha \tilde B_{\phi}) + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2 }{4}\tilde B \right]$$

$$\frac{\partial \tilde B}{\partial \tilde t} = - q\Omega \tilde B_{r} + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2}{4}\tilde B \right]$$

We will use the following definitions in our implementation:

$$\alpha = \alpha_0(r)\hat\alpha(z);\ \ \ \hat \alpha = sin\left(\pi \tilde z\right);\ \ \ q=-\frac{r}{\Omega}\frac{\partial \Omega}{\partial r}; \ \ \ \Omega = \frac{\Omega_0}{\sqrt{1 + (r/r_\omega)^2}}$$

$$R_\alpha = \frac{\alpha_0h}{\eta_t};\ \ \ R_\omega = \frac{-q\Omega h^2}{\eta_t}$$

We will be keeping $\alpha$ a constant to properly study the dynamo, instead of zero at midplane.

The final equations are:

$$\frac{\partial \tilde B}{\partial \tilde t} = - \frac{2}{\pi}(R_\alpha\hat\alpha \tilde B_{\phi}) + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2 }{4}\tilde B \right]$$

$$\frac{\partial \tilde B}{\partial \tilde t} = R_\omega B_{r} + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2}{4}\tilde B \right]$$

For solving these equations, we will use the same spatial derivative function using finite difference as given before in Section 1.4.1. For time-stepping, we implement coupled RK4 which solved the coupled differential equations including every single spatial position.

Find the code at task2_code.

2.3 Results

2.3.1 Sample implementation

To test the solver as well as obtain preliminary results, we ran the following simulation. THe results shows amplification as expected.

2.3.2 Decay rate at different spatial steps

To find the decay rate, we take the, we first plotted the magnetic field strength evolution over time in a log-plot, which shows straight lines as expected. Some local variation (as seen as one deviant curve in the figure below) could be due to the choice of the seed field. as it is not seen in every trial done.

To find the decay rate, we take the slope of the log-plot of the magnetic field strength at the middle of the spatial domain, assuming exponential decay. The results are:

2.3.3 Pitch angle evolution

The evolution of pitch angle was also tracked and is given below. The choice of seed field play a major role here since these simulations were of short time periods.

2.3.4 Magnetic field evolution with different seed fields

For studying the evolution with different seed fields, we ran two different simulations as shown below. The magnetic strength evolution of the different seed fields are:

2.3.5 Magnetic field evolution with different boundary conditions

To study the effect of boundary conditions, we ran the simulation with symmetric, anti-symmetric and relative anti-symmetric ghost zones, each of which ensures a different boundary condition. The symmetric boundary condition keeps the first derivative zero, the anti-symmetric boundary condition keeps the boundary at zero itself, whereas the relative anti-symmetric keeps the second derivative zero. We found that the anti-symmetric boundary conditions suits the current task and will be followed in the next section.

The evolutions of magnetic field strength for each of the implementation above is:

2.3.6 Finding critical dynamo number

To find the critical dynamo number at some time, we plot the local growth rate at that time and check at which spatial position does it cross zero (move from decay to growth). THe values used in this implementation is inspired from the numerical simulation in (Chamandy et.al. 2014). It was used such that those $R_omega$ and $D_C$ values are present in our value range across $r$.

The variation of $R_\omega$ and Dynamo number $D$ across $r$ was obtained as: The magnetic field strength evolution of this implementation turned out to be: The local growth rate and corresponding Critical Dynamo Number at unit time is marked in the figure below.

The above method was done, drawing inspiration from SS21 Figure 11.1 (Shukurov & Subramanian 2021) for our experimental setting. Reading into another paper (Chamandy et.al. 2014), I believe it might have been better to run the simulation with a range of values for global dynamo number, each simulation extending to the stable regime. Then, the critical dynamo number can be found by checking the value of fynamo number at which the steady decay changes into steady growth. This required a large amount of time and thus could not be explored right now.

For doing a similar comparison, our simulation will be done for a longer time period till a global decay rate is obtained and will be compared with the theoretical growth rate predicted by (Chamandy et.al. 2014). The values there were given by:

$$D_C = - \frac{\pi^5}{32};\ \ \ \gamma = \frac{2}{\pi} t_d^{-1} (\sqrt{-D} - \sqrt{-D_C})$$

With the theoretical values at $t_d = 1$ (unit diffusion time, since time is scaled) and the global decay rate obtained at a longer simulation, the below figure was plotted.

---

Main Project

---

Abstract

Studying galactic mean-field dynamos involves exploring the intricate dynamics governing cosmic magnetic fields, which provide invaluable insights into the formation and evolution of galaxies. Our research focuses on understanding the impact of vertical outflows within these systems, as they greatly influence the distribution of gas and magnetic fields within galaxies and are crucial to the star formation rates and matter-energy transport processes. In this study, our objective was to investigate how vertical outflows, characterized by various velocity profiles, influence the evolution of galactic magnetic fields.

Utilizing a combination of finite differencing and RK4 methods, we simulated the kinematic regime to observe the behavior of magnetic fields over time. Our findings demonstrate that the incorporation of velocity outflows, particularly stronger ones, hinders the growth of magnetic fields, and in some cases, even leads to decay. This highlights the significant role that vertical outflows play in shaping the evolution of magnetic fields in galactic dynamos, particularly emphasizing their interaction with radial terms.

These insights can not only deepen our understanding of the dynamics of galactic magnetic fields but also contribute to broader discussions surrounding galaxy formation. By studying the intricate interplay between various physical processes within cosmic structures, this project sheds light on aspects of galactic evolution that are closer to our comprehension of the universe's complex mechanisms.

1. Introduction

Galactic magnetic fields are one of the intriguing phenomena that is spread across the vast expanse of space within galaxies, stretching across thousands of light-years. Unlike magnetic fields on smaller scales, such as those around individual stars or planets, these fields display coherent patterns that span vast regions. Often organized in spiral or toroidal configurations, they reflect the overall shape of the galaxy itself.

Various mechanisms have been proposed to explain the generation and sustenance of these large-scale magnetic fields. A leading theory, the Dynamo theory, posits that magnetic fields arise from the movement of charged particles within galaxies, particularly the rotation and turbulence of interstellar gas and cosmic rays. This dynamo mechanism has the capability to amplify weak initial magnetic fields to the strengths observed in galaxies today.

Dynamo theory explains how magnetic fields are continuously generated and develop over immense time spans in the cosmos. It originated from foundational studies by Parker and Vainshtein and Ruzmaikin, who introduced mean-field dynamo models to account for the emergence of magnetic fields in spiral galaxies. Dynamos are vital in shaping the magnetic environment of various astrophysical systems, including galaxies, stars, and accretion disks. Understanding the mechanisms behind dynamo processes is essential for deciphering cosmic magnetism, star formation dynamics, and galactic evolution.

In Galactic mean-field theory, the evolution of magnetic fields is thoroughly described by the induction equation, which is crucial for explaining the dynamics of these fields within cosmic structures. This equation, fundamental to magnetohydrodynamics (MHD), clarifies how magnetic fields evolve under the influence of various physical processes.

From the galactic mean field theory, we have magnetic field evolution equations that are given by the induction equation. The induction equation (Chamandy 2024) is given by:

$$\frac{\partial \overline B_r}{\partial t} = - \frac{\overline V_r}{r} \frac{\partial}{\partial r}(r\overline B_r) - \frac{\partial}{\partial z}(\overline V_z \overline B_r) - \frac{\partial}{\partial z}(\alpha \overline B_{\phi}) + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_r)\right) + \frac{\partial ^2 \overline B_r}{\partial z^2} \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = - q\Omega \overline B_{r} - \frac{\partial}{\partial r}(\overline V_r \overline B_{\phi}) - \frac{\partial}{\partial z}(\overline V_z \overline B_{\phi}) + \frac{\partial}{\partial z}(\alpha \overline B_{r}) + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_{\phi})\right) + \frac{\partial ^2 \overline B_{\phi}}{\partial z^2} \right]$$

where $\overline B_r$, $\overline B_{\phi}$ and $\overline B_z$ are the components of the magnetic field, $\overline V_r$, $\overline V_{\phi}$ and $\overline V_z$ are the components of the velocity field, $\alpha$ is the dynamo alpha effect, $\Omega$ is the angular velocity of the galaxy, $q$ is the shear parameter, and $\eta_t$ is the turbulent magnetic diffusivity.

In this project, we will be focusing on the impact of vertical outflows on the evolution of galactic magnetic fields compared to task 2 where we we were studying the impact of $\alpha$-$\omega$ dynamo effects. The velocity field is given by $\overline V_z$ and terms involving $\overline V_r$ as well as $\alpha$ term from the equation for $B_{\phi}$ are omitted from the equations. Task 3 will be focused on solving these equations. THe equations are:

$$\frac{\partial \overline B_r}{\partial t} = - \frac{\partial}{\partial z}(\overline V_z \overline B_{r}) - \frac{\partial}{\partial z}(\alpha \overline B_{\phi}) + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_r)\right) + \frac{\partial ^2 \overline B_r}{\partial z^2} \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = - \frac{\partial}{\partial z}(\overline V_z \overline B_{\phi}) - q\Omega \overline B_{r} + \eta_t \left[ \frac{\partial}{\partial r} \left( \frac{1}{r} \frac{\partial}{\partial r}(r \overline B_{\phi})\right) + \frac{\partial ^2 \overline B_{\phi}}{\partial z^2} \right]$$

By applying no-$z$ approximation and expanding, equations become:

$$\frac{\partial \overline B_r}{\partial t} = - \frac{2}{\pi h}(\overline V_z \overline B_{r}) - \frac{2}{\pi h}(\alpha \overline B_{\phi}) + \eta_t \left[ \frac{\partial ^2 \overline B}{\partial r^2} + \frac{1}{r} \frac{\partial \overline B}{\partial r} - \frac{\overline B}{r^2} + \frac{\pi^2 }{4h^2}\overline B_r \right]$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = - \frac{2}{\pi h}(\overline V_z \overline B_{\phi}) - q\Omega \overline B_{r} + \eta_t \left[ \frac{\partial ^2 \overline B}{\partial r^2} + \frac{1}{r} \frac{\partial \overline B}{\partial r} - \frac{\overline B}{r^2} + \frac{\pi^2 }{4h^2}\overline B_{\phi} \right]$$

We can convert these equations into dimensionless form and the same will be used in the simulations. Hence it becomes:

$$\frac{\partial \tilde B_r}{\partial\tilde t} = - \frac{2}{\pi\tilde h}(\tilde V_z \tilde B_{r}) - \frac{2}{\pi\tilde h}(\tilde\alpha \tilde B_{\phi}) + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2 }{4\tilde h^2}\tilde B_r \right]$$

$$\frac{\partial \tilde B_{\phi}}{\partial\tilde t} = - \frac{2}{\pi\tilde h}(\tilde V_z \tilde B_{\phi}) - q\tilde\Omega \tilde B_{r} + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2 }{4\tilde h^2}\tilde B_{\phi} \right]$$

To make it dimensionless, we use the following definitions in our implementation (Chamandy et.al. 2014):

$$\eta_t = \frac{1}{3}lu;\quad t_d = h^2/\eta_t;$$

$$\Omega = \frac{\Omega_0}{\sqrt{1 + (r/r_\omega)^2}};\quad q=-\frac{r}{\Omega}\frac{\partial \Omega}{\partial r};\quad \tilde \Omega = \frac{\Omega h^2}{\eta_t}$$

$$\alpha = \alpha_0(r)\hat\alpha(z);\quad \hat \alpha = sin\left(\pi \tilde z\right); \quad \alpha_0(r) = \frac{\Omega l^2}{h};\quad \tilde \alpha_0 = \frac{\alpha_0 h}{\eta_t}$$

$$R_\alpha = \tilde\alpha_0;\quad R_\omega = -q\tilde \Omega; \quad D = R_\alpha R_\omega$$

$$\tilde h = 1;\quad \tilde t_d = 1; \quad\tilde \eta_t = 1;\quad \tilde r = \frac{r}{h};\quad \tilde h = 1;\quad \tilde t_d = 1; \quad\tilde B = \frac{B}{B_0}; \quad \tilde V_z = \frac{V_z h}{\eta_t}$$

Here, $l$ is the correlation length of the turbulence, $u$ is the turbulent velocity, $h$ is the scale height of the disk, $r_\omega$ is the radial scale of the angular velocity profile, and $B_0$ is the initial magnetic field strength. The values of the parameter used ((Chamandy et.al. 2014)) will be given in Section 3. We will be keeping $\hat\alpha = 1$ a constant instead of zero at midplane, for our study. For disk flaring, we will use $\tilde h = h(r) / h$ and the same will affect the $\alpha_0$ term as well. This won't affect the conversion to dimensionless form, as we are using the same $h$ in the conversion.

The final equations are:

$$\frac{\partial \tilde B_r}{\partial\tilde t} = - \frac{2}{\pi\tilde h}(\tilde V_z \tilde B_{r}) - \frac{2}{\pi\tilde h}(R_\alpha \tilde B_{\phi}) + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2 }{4\tilde h^2}\tilde B_r \right]$$

$$\frac{\partial \tilde B_{\phi}}{\partial\tilde t} = - \frac{2}{\pi\tilde h}(\tilde V_z \tilde B_{\phi}) + R_\omega \tilde B_{r} + \left[ \frac{\partial ^2 \tilde B}{\partial\tilde r^2} + \frac{1}{\tilde r} \frac{\partial \tilde B}{\partial\tilde r} - \frac{\tilde B}{\tilde r^2} + \frac{\pi^2 }{4\tilde h^2}\tilde B_{\phi} \right]$$

For solving these equations, we will use the same spatial derivative function using finite difference as given before in Section 2.1. For time-stepping, we implement coupled RK4 which solved the coupled differential equations including every single spatial position.

Similar equations have been solved in the past by others (Shukurov et.al. 2016) where though significantly different in terms of implementation and settings, the resutls can be compared in the kinematic regime. A result from the same is shown in Figure 1 (Shukurov et.al. 2016).

Figure 1: The inset shows different levels of vertical outflow: absence (solid), weak (dashed), moderate (dotted) and strong (dot-dashed)

Here, in the inset of Figure 1, the magnetic field strength evolution in the kinematic regime can be seen for different vertical outflows. Stronger outflows are seen to hinder the growth of magnetic fields, and in some cases, even lead to decay. This is the main focus of our study.

Find the code at task3_code.

2. Methods

2.1. Finite difference

The spatial derivatives are calculated using the finite difference method. here, we are using 6th order finite difference method for the spatial derivatives. The 6th order finite difference method is given by (Bradenburg 2019) as follows:

$$\frac{d f}{d x} = \frac{1}{60\delta x}\left( - f_{i-3} + 9f_{i-2} - 45f_{i-1} + 45f_{i+1} - 9f_{i+2} + f_{i+3} \right)$$

$$\frac{d^2 f}{d x^2} = \frac{1}{180\delta x^2}\left( 2f_{i-3} - 27f_{i-2} + 270f_{i-1} - 490f_{i} + 270f_{i+1} - 27f_{i+2} + 2f_{i+3} \right)$$

Other orders can be used and is available for use in the code.

By using finite difference, we can calculate the spatial derivatives at every spatial point. This approximation is used to convert the partial differential equations into ordinary differential equations, which can be solved using time-stepping methods (RK4 in our case). This approximation method is called method of lines.

2.2. Time-stepping using RK4

The time stepping is done using the RK4 method. Our evolution function can be given of form:

$$\frac{\partial \overline B}{\partial t} = f(\overline B, t)$$

But, since our $\overline B$ has two components, we have two equations to solve.

$$\frac{\partial \overline B_r}{\partial t} = f_r(\overline B_r, \overline B_{\phi}, t)$$

$$\frac{\partial \overline B_{\phi}}{\partial t} = f_\phi(\overline B_r, \overline B_{\phi}, t)$$

Then, the RK4 method evolution for a single time step $\delta t$ from the current $\overline B(t)$ can be given by the coupled RK4 equations:

$$\kappa_{1r} = \delta t f_r(\overline B_r(t), \overline B_{\phi}(t), t)$$

$$\kappa_{1\phi} = \delta t f_\phi(\overline B_r(t), \overline B_{\phi}(t), t)$$

$$\kappa_{2r} = \delta t f_r(\overline B_r(t) + \frac{\kappa_{1r}}{2}, \overline B_{\phi}(t) + \frac{\kappa_{1\phi}}{2}, t + \frac{\delta t}{2})$$

$$\kappa_{2\phi} = \delta t f_\phi(\overline B_r(t) + \frac{\kappa_{1r}}{2}, \overline B_{\phi}(t) + \frac{\kappa_{1\phi}}{2}, t + \frac{\delta t}{2})$$

$$\kappa_{3r} = \delta t f_r(\overline B_r(t) + \frac{\kappa_{2r}}{2}, \overline B_{\phi}(t) + \frac{\kappa_{2\phi}}{2}, t + \frac{\delta t}{2})$$

$$\kappa_{3\phi} = \delta t f_\phi(\overline B_r(t) + \frac{\kappa_{2r}}{2}, \overline B_{\phi}(t) + \frac{\kappa_{2\phi}}{2}, t + \frac{\delta t}{2})$$

$$\kappa_{4r} = \delta t f_r(\overline B_r(t) + \kappa_{3r}, \overline B_{\phi}(t) + \kappa_{3\phi}, t + \delta t)$$

$$\kappa_{4\phi} = \delta t f_\phi(\overline B_r(t) + \kappa_{3r}, \overline B_{\phi}(t) + \kappa_{3\phi}, t + \delta t)$$

$$\overline B_r(t + \delta t) = \overline B_r(t) + \frac{1}{6}(\kappa_{1r} + 2\kappa_{2r} + 2\kappa_{3r} + \kappa_{4r})$$

$$\overline B_{\phi}(t + \delta t) = \overline B_{\phi}(t) + \frac{1}{6}(\kappa_{1\phi} + 2\kappa_{2\phi} + 2\kappa_{3\phi} + \kappa_{4\phi})$$

The thing to note here about the generalized form, with regards to the RK4 method, is that:

• The spatial derivatives do not change values with the change ($\overline B \rightarrow \overline B + \kappa/2$) in the RK4 method. This is due to the symmetry of the spatial derivative equations. Hence, the spatial derivatives are calculated only once at the start of the time step. Thus, it is not mentioned in the equation form.
• The order of calculation of $\kappa$ values is important. The order is $\kappa_{1}$, $\kappa_{2}$, $\kappa_{3}$, $\kappa_{4}$. Both components of $\overline B$ are calculated at each step before moving to the next step. This is important to maintain the integrity of the coupled equations.

2.3. Implementation of BCs

For the sake of obtaining smooth derivatives as well for the sake of maintaining the boundary conditions, we use ghost zones. To implement ghost zones, we do the following:

• At the start and end of the spatial domain, we add ghost cells outside the domain of interest. We are free to choose the number of ghost cells to add and to set the values of the ghost cells.
• The number of ghost cells to add is half the order of the finite difference method. Hence, for 6th order finite difference method, we add 3 ghost cells at the start and end of the spatial domain. THis will be useful for calculating the derivatives at the boundaries smoothly.
• The values of the ghost cells are set according to the type of ghost zones used and the boundary conditions.

By default, 'anti-symmetric' ghost zones are used where the ghost cell at $i$ distance away from boundary $b$ is given by:

$$f_{b-i} = -f_{b+i} \quad \text{at the start of the spatial domain}$$

$$f_{b+i} = -f_{b-i} \quad \text{at the end of the spatial domain}$$

This is used to implement the boundary conditions $B_r = 0$ and $B_{\phi} = 0$ at the boundaries.

Symmetric ghost zones are used where the ghost cell at $i$ distance away from boundary $b$ is given by:

$$f_{b-i} = f_{b+i} \quad \text{at the start of the spatial domain}$$

$$f_{b+i} = f_{b-i} \quad \text{at the end of the spatial domain}$$

This is used to implement the boundary conditions $\frac{\partial B_r}{\partial r} = 0$ and $\frac{\partial B_{\phi}}{\partial r} = 0$ at the boundaries.

Relative anti-symmetric ghost zones are used where the ghost cell at $i$ distance away from boundary $b$ is given by:

$$f_{b-i} = 2f_{b} - f_{b+i} \quad \text{at the start of the spatial domain}$$

$$f_{b+i} = 2f_{b} - f_{b-i} \quad \text{at the end of the spatial domain}$$

This is used to implement the boundary conditions $\frac{\partial^2 B_r}{\partial r^2} = 0$ and $\frac{\partial^2 B_{\phi}}{\partial r^2} = 0$ at the boundaries.

where $i = 1, 2, 3, ..., n$ ($n$ is the half the order of the finite difference method). These ghost zones are used to maintain the boundary conditions as well as to obtain smooth derivatives at the boundaries and are free to be changed in the code. By default, 'anti-symmetric' ghost zones are used to implement the boundary conditions $B_r = 0$ and $B_{\phi} = 0$ at the boundaries.

2.4. Calculation of Global growth rate

In the calculation of the global growth rate and the iterative process for running simulations until the desired growth rate is achieved, several systematic steps are followed:

1. Initial Evolution Period: During the first 10 diffusion times, the magnetic field evolves without monitoring growth rates.

2. Collection of Growth Rates: At the end of each subsequent diffusion time, the growth rate is determined at every spatial point. This involves analyzing the evolution of magnetic field strength ($|B|$) over the last few hundred time steps at each spatial point and calculating the slope of the logarithm of $|B|$. This process is repeated for all spatial points, generating a distribution of growth rates ($\gamma$) across the radius ($r$).

3. Mean Growth Rate Calculation: The mean growth rate ($\overline{\gamma}$) across all radii is calculated from the obtained distribution of growth rates.

4. Tolerance Check: The absolute difference between the mean growth rate ($\overline{\gamma}$) and the individual growth rates ($\gamma$) is compared against a predefined tolerance threshold.

5. Simulation Continuation: The simulation continues until the absolute difference falls below the specified tolerance. Furthermore, this condition must be maintained for an additional 3 diffusion times to ensure stability.

By meticulously following these steps, we can ensure accurate determination of global growth rates in galactic dynamo simulations.

2.5. Intuition behind Velocity profiles

Typical choice of the Velocity profile for vertical outflow as given in (Chamandy et.al. 2014) is:

$$V_z = \frac{V_zz}{h}$$

Hence, at a certain $z$, the Vertical outflow will be a constant.

In addition to exploring different values of constant profiles, I was also thinking to explore radial variations of the Vertical outflow due to reasons discussed before in Section 1.

Since the star formations are more in smaller radius, a radially decreasing profile seemed correct. In addition to this, if we consider the effect of higher gravity in smaller radius, we can also try a velocity profile which decreases till some radius $r_v$ and then increase due to the lower effect of gravity (henceforth mentioned as U-profile). To explore both using similar equation, I defined it in the following way:

$$V_z = V_0 \times \left(1 + 0.5 \times \left(\left(\frac{r - r_{v}}{r_{v}}\right)^2 - 1\right)\right)$$

The last two profiles mentioned above are obtained by putting $r_v$ = $r_f$ (final radius) and $r_v$ = 10 $kpc$ respectively. The profiles are shown below in Figures 2 and 3.

Figure 2: Velocity profile that decreases till the end of spatial domain.
Figure 3: Velocity U-profile that decreases till 10 kpc, then increases
Note that this is just two sample profiles mimicking the shapes intuitive to the reasoning above.

2.6. Addition of Disk-flaring

To include disk flaring, we changed applied the radial profile $h(r)$ on the dimensionless form and propagated it through the other affected parameters ($R_\alpha, D, etc/$) and evolution. The form of disk flaring used is:

$$h(r) = h \times \sqrt{1 + \left(\frac{r}{r_h}\right)^2}$$

Here, $r_h$ controls the disk flaring and is set to $10 kpc$ as per (Chamandy et.al. 2014) The form of $h(r)$ is given below in Figure 4.

Figure 4: Visualization of disk flaring

Note that the scaling and making the parameters are kept unaffected by this change as we use $h$ and the corresponding $t_d$ for that purpose. See the implementation in task_3-code.

3. Simulation and results

All simulation paramters, trial settings, and the results are available in Github and task_3-code.

We used the following parameters and the forms discussed in Section 1 and (Chamandy et al. 2014) for the implementation.

The profiles of $R_\omega$, $R_\alpha$ and $D$ across $r$ are given below in Figures 5 to 7 for reference.

Figure 5: Variation of R_omega with respect to r
Figure 6: Variation of R_alpha with respect to r
Figure 7: Variation of Dynamo number with respect to r

Note that the disk flaring is not used for studies in Sections 3.1 and 3.2 for the sake of simplicity and to focus on the effect of vertical outflow. The disk flaring is used in the study in Section 3.3.

3.1. Effect of non-zero vertical outflow

To study the effect of vertical outflow compared to the previous results from Task 2, we ran the simulations with $V_z$ = 0, 2, 4, 6 , 8, 10 $km/s$ for the modified equations and implementation in Section 1. We will first focus on results from $V_z$ = 0 and $V_z$ = 4 $km/s$ (representing non-zero $V_z$) for comparison of different aspects of magnetic field evolution.

The results from $V_z$ = 0 and $V_z$ = 4 $km/s$ are shown below. Below, in Figures 8 to 11, are the final eigen modes and the magnetic field evolution of both simulations.

Figure 8: Final eigen mode of magnetic field in the absence of vertical outflow
Figure 9: Final eigen mode of magnetic field in the presence of vertical outflow (V_z = 4 km/s)
Figure 10: Evolution of magnetic field strength in the absence of vertical outflow
Figure 11: Evolution of magnetic field strength in the presence of vertical outflow (V_z = 4 km/s)

As you can see the final eigen modes are the same and is not affected by the vertical outflow here.

We can also look at the evolution of growth rate across $r$ for both simulations given below in Figures 12 and 13.

Figure 12: Evolution of growth rate in the absence of vertical outflow
Figure 13: Evolution of growth rate in the presence of vertical outflow (V_z = 4 km/s)

The global growth rate $\Gamma$ achieved by the simulation with Vertical outflow is lower than that without the vertical outflow. This is expected as the vertical outflow removes the magnetic field from the system that was growing. The pitch angle evolution is not severely affected by the vertical outflow as seen below in Figures 14 and 15.

Figure 14: Final form of pitch angle in the absence of vertical outflow
Figure 15: Final form of pitch angle in the presence of vertical outflow (V_z = 4 km/s)

Looking at the general trend across all the above range of values for $V_z$, the important one is the change in the global growth rate $\Gamma$ obtained. The variation of $\Gamma$ across various values of $V_z$ is shown below in Figure 16 through animation and in Figure 17 through a plot of magnetic field evolution.

Figure 16: Variation of global growth rate with respect to V_z
Figure 17: Evolution of magnetic field strength for various values of V_z

As you can see, the global growth rate is negatively affected by the outflow, and the same is expected the outflow removes the magnetic field from the system that was growing. As $V_z$ further increases, the $\Gamma$ crosses zero and becomes negative, i.e. decaying magnetic field. This is also expected as the vertical outflow can remove magnetic faster than it evolves, causing decay.

Since, all other parameters are the same, we can see how the vertical outflow affected the growth of the magnetic field. The results can be compared to the result mentioned in Figure 1 (Shukurov et.al. 2006) and discussed in Section 1.

3.2. Effect of different profiles of vertical outflow

Apart from constant velocity profile across $r$, we will also explore the decreasing and U-shaped profile effect on the magnetic field evolution.

The magnetic field evolution, pitch angle evolution and final eigen functions are given below for both simulations in Figures 18 to 23.

Figure 18: Final eigen mode of magnetic field in the decreasing vertical outflow profile
Figure 19: Final eigen mode of magnetic field in the U-shaped vertical outflow profile
Figure 20: Evolution of magnetic field strength in the decreasing vertical outflow profile
Figure 21: Evolution of magnetic field strength in the U-shaped vertical outflow profile
Figure 22: Final form of pitch angle in the decreasing vertical outflow profile
Figure 23: Final form of pitch angle in the U-shaped vertical outflow profile

Looking at their evolution of $\gamma$ given below, we can see that the decreasing profile has increased the global growth rate which can be explained looking at how the growth rates evolve and interact radially. See the animation below in Figure 24.

Figure 24: Evolution of growth rate in different Velocity profiles; constant (blue), decreasing (green) and U-shaped (red)

As the growth at the radius 2 $kpc$ (maximum Dynamo number as shown in the beginning of Section 3) is influencing and spreading to other radial positions, the decreasing profile has the same effect as a decreased $V_z$ value which implies higher growth rate. This also implies how a different profile (actual profile as well) would interact different with magnetic field evolution.

3.3. Disk flaring and vertical outflow

With the implementation of disk flaring as mentioned in Section 2.6, we can see the interaction of the disk flaring with the vertical outflow. For this, we ran the simulation with disk-flaring and $V_z$ ranging from 0 to 10 $km/s$ in steps of 2 $km/s$. Looking at their evolution of growth rates given below, we can see how disk flaring affected the impact of vertical outflows in the system. See the animation below in Figure 25.

Figure 25: Evolution of growth rate in different Velocity profiles with disk flaring

For higher vertical outflows, the growth rates are spread slowly in radius as the eariler interactions are affected by the disk flaring. We can compare these to the results from simulations without disk flaring to see the effect of disk flaring on the magnetic field evolution. In comparison to results shown in Figure 16, the disk flaring has affected the growth rates and the spread of growth rates across the radius. See the animation below in Figure 26.

Figure 26: Comparison of evolution of growth rate in different Velocity values with (denoted by DF in animation) and without disk flaring

As you can see in Figure 26, the radial interactions are affected by disk flaring and affecting the spread of growth. This shows that with additional terms dependent on $r$, the vertical outflows impact in various rates or in even more complicated manner which is further not explored here.

Conclusions

In our comprehensive exploration of galactic mean-field dynamos, we have meticulously observed and drawn significant conclusions from the simulations and results obtained. Across the span of our investigations, several key observations have emerged, shedding light on the intricate dynamics of cosmic magnetic fields within galactic systems.

• Our analysis indicates that the presence of vertical outflows exerts a notable influence on the growth of magnetic fields, distinctly contrasting with scenarios where such outflows are absent. Notably, we have observed a linear hindrance to magnetic field growth in the presence of vertical outflows, with stronger outflows removing magnetic fields at a rate surpassing their growth, eventually leading to decay. This fundamental insight, clearly seen in Figure 17, underscores the crucial role of vertical outflows in modulating magnetic field evolution within galaxies.

• Different profiles of vertical outflows, as discussed in Section 3.2, have revealed intriguing correlations between the evolution of growth rates and Dynamo number profiles (See Figure 24). Notably, we have observed that decreasing profiles exhibit effects akin to weaker velocity outflows, underscoring the intricate interplay between velocity profiles and magnetic field evolution. This highlights the nature of interactions within galactic dynamo systems, where variations in velocity profiles can show diverse responses in magnetic field evolution.

• Our inclusion of radially changing systems, as detailed in Section 3.3 and seen in Figure 26, has shown the impact of radial variations on the influence of outflows, affecting convergence rates and growth dynamics. This underscores the complexity inherent in galactic dynamo systems, where radial variations can significantly modulate the effects of vertical outflows on magnetic field evolution.

The significance of our findings extends beyond these systems, as they hold profound implications for understanding real-world galactic systems. Traditionally, velocity outflow terms have been neglected in favor of simplifying system dynamics. However, our research underscores the critical relevance of considering these terms, particularly in light of their impact on shaping the magnetic field evolution and their relation to star formation rates and matter transport within and from the system.

Moving forward, our study suggests promising aspects to be explored for future research. Specifically, the formulation and exploration of various profiles of velocity outflows, accounting for effects such as star formation rates and disk flaring at different radii, hold considerable potential. Additionally, the possibility of coupling these considerations with the relaxation of various other assumptions and exploring system dynamics beyond the kinematic regime presents exciting opportunities for advancing our understanding of galactic dynamo theory.

In conclusion, our comprehensive investigation into the interactions of vertical outflows within galactic mean-field dynamos has yielded profound insights into the complex dynamics governing magnetic field evolution. Through meticulous analysis and simulation, we have illuminated the intricate interplay between velocity outflows and magnetic field dynamics, paving the way for further advancements in galactic dynamo theory and our understanding of cosmic structures.

References

1. Brandenburg, A. 2019. Computational aspects of astrophysical MHD and turbulence. In Advances in nonlinear dynamos (CRC Press) [link]
2. Chamandy, L., Shukurov, A., Subramanian, K. & Stoker, K. 2014. MNRAS 443. [link]
3. Chamandy, L. 2024. P464 Class Notes. (NISER)
4. Shukurov, A., Sokoloff, D., Subramanian, K., & Brandenburg, A. 2006. Galactic dynamo and helicity losses through fountain flow. Astronomy & Astrophysics, 448(2), L33-L36. [link]
5. Shukurov, A., & Subramanian, K. 2021. Astrophysical Magnetic Fields: From Galaxies to the Early Universe. (Cambridge: Cambridge University Press)[link]