New Heating Model
A More Generalized Arber Formula
Last week I read the Hockney papers to understand where the heating time formulas come from which are used by Arber for his formula. I realized that there is no particular reason why $(\lambda_D / \Delta x)$ has to be squared, so we should allow the possibility that it may not be squared. Maybe, this relationship breaks down when the simulation gets sufficiently under-resolved. Regardless, I decided to run some simulations that would allow me to explore how effective such a model to be. Let us assume, from the Arber Paper, that the heating time is only dependent on the dimensionless debye ratio and the inverse of the plasma frequency.
\begin{equation} \tau_H \sim \frac{n_{ppc}}{\omega_{pe}} ( \lambda_D / \Delta x)^\alpha \end{equation}
The two other dimensionless parameters are $\omega_{pe} \Delta t$ and $v_{th} \Delta t / \Delta x$ which we are going to ignore for now, but if necessary could be included in a more advanced model. From my previous sims, I am convinced that the dependence on $n_{ppc}$ is correct, so I will be leaving that as a constant. This would correspond to an Arber-like formula of
\begin{equation} \frac{dT_{eV}}{dt_{ps}} = C_H T_{0,eV}^{1 - \alpha/2} \Delta x_{nm}^\alpha n_{23}^{(\alpha + 1)/2} \end{equation}
where for a spline shape function, $C_H = \alpha_{H} / n_{ppc} = 2/25 = 0.08$. Here, I assumed that $n_{ppc}$ = 25 which is what I used for the sims in the next section. I want to test this model for a variety of different ratioss.
Simulation Setup
From my previous sims, it seems like the ratio is the most important factor in determining how much heating there is. As a result, I decided to fix the ratio at certain values and adjust $\Delta x$ , $n_{23}$ and $T_0$ in such a way to keep the ratio the same. For example, using $T_0 = 100$ eV, $n_{23} = 1$, $\Delta x = 3 $ nm (and $n_{ppc} = 25$ for all the sims), I did the following 13 simulations:
$({T_0, n_{23}, \Delta x})$ $({4 T_0, 4 n_{23}, \Delta x})$ $({9 T_0, 9 n_{23}, \Delta x})$ $({\frac{1}{4} T_0, \frac{1}{4} n_{23}, \Delta x})$ $({\frac{1}{9} T_0, \frac{1}{9} n_{23}, \Delta x})$ $({4 T_0, n_{23}, 2 \Delta x})$ $({9 T_0, n_{23}, 3 \Delta x})$ $({\frac{1}{4} T_0, n_{23}, \frac{1}{2} \Delta x})$ $({\frac{1}{9} T_0, n_{23}, \frac{1}{3} \Delta x})$ $({T_0, 4 n_{23}, \frac{1}{2} \Delta x})$ $({T_0, 9 n_{23}, \frac{1}{3} \Delta x})$ $({T_0, \frac{1}{4} n_{23}, 2 \Delta x})$ $({T_0, \frac{1}{9} n_{23}, 3 \Delta x})$
this corresponded to a ratio of $0.0783$. I repeated this process, but just varied $\Delta x$ from $\Delta x = {1/4, 1/2, 1, 3, 6, 9, 12}$ nm. In setting up the simulations this way, I can make several assessment of the data in order to confirm that an $\alpha$ model is correct.
also, I did all the sims for 1 ps. As a result, for shorthand, I will represent $\frac{dT}{dt} \rightarrow \Delta T$ which is assumed to be the change in temperature in electron volts in 1 ps.
Assume $\Delta x$ is constant
If we assume $\Delta x$ is constant, if we are at a constant ratio and increase the initial temperature by some factor $\beta$, we will get
$T_0’ = \beta T_0$ leads to $n’ = \beta n$.
So, if $\Delta T = C_H T_0^{1 - \alpha/2} \Delta x^\alpha n^{(\alpha+1)/2}$, then we get
\begin{equation} \Delta T’ = \beta^{3/2} \Delta T \end{equation}
which is indepedent of $\alpha$.
Assume $n$ is constant
If $n$ is constant, $\Delta x’ = \beta \Delta x$ leads to $T_0’ = \beta^2 T_0$ so that
\begin{equation} \Delta T’ = \beta^2 \Delta T \end{equation}
Assume $T_0$ is constant
If $T_0$ is constant, $\Delta x’ = \beta \Delta x$ leads to $n’ = \frac{1}{\beta^2} n$ so that
\begin{equation} \Delta T’ = \beta^{-1} \Delta T \end{equation}
Results for $\Delta x = 3$ nm
Keeping $\Delta x$ fixed, I got the expected $\beta^{3/2}$ relationship
Keeping $T_0$ fixed, I got the expected $\beta^{-1}$ relationship
Keeping $n$ fixed, I got the expected $\beta^2$ relationship
So, in combination, these three plots show that the $\alpha$-Model should be a very good fit to the data from my sims. Then, I iterated over a range of $\alpha \in ({0, 4})$ and $C_H \in ({0.001, 100})$ and picked them so that mean squared error between the $\alpha$-Model and PIC sims was the lowest. In this case, I got $\alpha = 2.5$ and $C_H = 0.735$. This is in contrast to Arber’s prediction of $\alpha = 2$ and $C_H = 0.08$. Below, I have a plot that shows how the Arber Model compares to my Model in terms of the difference in the amount of heating $\Delta T$ compared to the PIC simulations.
To me, it is clear that this model performs much better. However, it is not as clear for the other models whether the improvement is due to a good choice for hyper-parameters $\alpha$ and $C_H$ or possibly just to overfitting to my data. I say this because the predictions for $C_H$ and $\alpha$ has no clear pattern when I changed the cell size.
Summary of Results for the other simulations
For an example of the smaller simulation size simulations, for $\Delta x = \frac{1}{4}$, we get $\beta$ model graphs that look like this:
When we keep $\Delta x$ and $T_0$ fixed we do not get the expected relationships. To me, this looks like the case because the small grid of 32x32 combined with only using 25 ppc may make the simulations too noisy. Additionally, Since the Debye-Length isn’t too badly under-resolved, there is not enough heating to make a good comparison. On the other end of the spectrum, when $\Delta x$ is large (for example $\Delta x = 12$ nm), the $\beta$ models do match up moderately well. However, when we try to fit the alpha model, we get that the best prediction has $\alpha \sim 0$ which means that all simulations would just exhibit a constant amount of heating $C_H$ in one ps. This seems related to the observation from previous weeks that after some value of $\Delta x$, the amount of numerical heating does not continue to go up. I wonder if these issues could be solved by using more particles per cell and more grid cells.
I have summarized my findings in the Table below
Things to Do
- Explore ‘Optimal Path’ Formula from Hockney Paper
- Work on Modifying Model and a better selection criteria for the $C_H$ and $\alpha$ parameters
- Determine if Hockney paper was specific to 2D or if it generalizes to 1D and 3D.
- Calculate $\omega_{pe} \Delta t$ for my sims.