User:Lothar.brendel/belt-notes: Difference between revisions
(new) |
(→SEDNA: Ionization → Rate Equations: +unstable y=1) |
||
| (4 intermediate revisions by the same user not shown) | |||
| Line 1: | Line 1: | ||
= SEDNA: Ionization → Rate Equations = | = SEDNA: Ionization → Rate Equations = | ||
The rate equations (44) and (45) in [https://ui.adsabs.harvard.edu/abs/2020ApJS..250...13K/abstract Kuiper etal 2020] are non-linear, namely quadratic in \(x\) and \(y\), respectively. Solving them with an implicit Euler scheme due to stability reasons (the rate equations seemingly are not considerd in the \(\Delta t\) control), yields quadratic equations, (46) and (53), for the values at \(t+\Delta t\). For each, the correct solution must be chosen, namely the one for which \(x^{n+1}-x^n=\mathcal O(\Delta t)\). | |||
For \(B>0\), the function <code>Sedna1.1/RateEquations.c:QuadraticEquation(A,B,C,&x1,&x2)</code> assigns the solution with larger absolute value, which happens to be negative, to \(x_2\). For \(B<0\), the solution with larger absolute value is positive and is assigned to \(x_1\). The sign of the other solution, with smaller absolute value, depends on the sign of \(C\). | For \(B/A>0\), the function <code>Sedna1.1/RateEquations.c:QuadraticEquation(A,B,C,&x1,&x2)</code> assigns the solution with larger absolute value, which happens to be negative, to \(x_2\). For \(B/A<0\), the solution with larger absolute value is positive and is assigned to \(x_1\). The sign of the other solution, with smaller absolute value, depends on the sign of \(C\). | ||
== \(x(t)\) == | == \(x(t)\) == | ||
* \(0<b+c=A\propto\Delta t\) | * \(0<b+c=A\propto\Delta t\) | ||
* \(B=1+a-b>0\) for not too large \(\Delta t\) | * \(B=1+a-b>0\) (for not too large \(\Delta t\)!) | ||
* \(C=-(x^n+a)<0\) | * \(C=-(x^n+a)<0\) | ||
| Line 14: | Line 14: | ||
x^{n+1}=\frac{-B+\sqrt{\dots}}{2A}~, | x^{n+1}=\frac{-B+\sqrt{\dots}}{2A}~, | ||
\] | \] | ||
which is the solution with smaller absolute value. <code>QuadraticEquation()</code> yields \(x_2<0 | which is the solution with smaller absolute value. <code>QuadraticEquation()</code> yields \(x_2<0\le x_1<\vert x_2\vert~\Rightarrow~x^{n+1}=x_1\). | ||
== \(y(t)\) == | == \(y(t)\) == | ||
| Line 27: | Line 25: | ||
y^{n+1}=\frac{-B-\sqrt{\dots}}{2A}~, | y^{n+1}=\frac{-B-\sqrt{\dots}}{2A}~, | ||
\] | \] | ||
which again is the solution | which again is the solution of smaller absolute value. <code>QuadraticEquation()</code> yields \(0\le x_2<x_1~\Rightarrow~y^{n+1}=x_2\). | ||
== Caveats == | |||
Without radiation (\(a=0\)), \(1-x=y=1\) is a stationary solution, but it's unstable. Due to round off errors, \(y^n=1\Rightarrow y^{n+1}=1-\epsilon\), triggering the instability. On the other hand, \(x^n=1\Rightarrow x^{n+1}=1\) due to \(C=0\). Hence, in opaque regions, an \(x/y\)-inconsistency will we reported. | |||
Latest revision as of 19:45, 18 May 2026
SEDNA: Ionization → Rate Equations
The rate equations (44) and (45) in Kuiper etal 2020 are non-linear, namely quadratic in \(x\) and \(y\), respectively. Solving them with an implicit Euler scheme due to stability reasons (the rate equations seemingly are not considerd in the \(\Delta t\) control), yields quadratic equations, (46) and (53), for the values at \(t+\Delta t\). For each, the correct solution must be chosen, namely the one for which \(x^{n+1}-x^n=\mathcal O(\Delta t)\).
For \(B/A>0\), the function Sedna1.1/RateEquations.c:QuadraticEquation(A,B,C,&x1,&x2) assigns the solution with larger absolute value, which happens to be negative, to \(x_2\). For \(B/A<0\), the solution with larger absolute value is positive and is assigned to \(x_1\). The sign of the other solution, with smaller absolute value, depends on the sign of \(C\).
\(x(t)\)
- \(0<b+c=A\propto\Delta t\)
- \(B=1+a-b>0\) (for not too large \(\Delta t\)!)
- \(C=-(x^n+a)<0\)
\[
x^{n+1}=\frac{-B+\sqrt{\dots}}{2A}~,
\]
which is the solution with smaller absolute value. QuadraticEquation() yields \(x_2<0\le x_1<\vert x_2\vert~\Rightarrow~x^{n+1}=x_1\).
\(y(t)\)
- \(0<b+c=A\propto\Delta t\)
- \(B=-1-a-b-2c<0\)
- \(C=y^n+a>0\)
\[
y^{n+1}=\frac{-B-\sqrt{\dots}}{2A}~,
\]
which again is the solution of smaller absolute value. QuadraticEquation() yields \(0\le x_2<x_1~\Rightarrow~y^{n+1}=x_2\).
Caveats
Without radiation (\(a=0\)), \(1-x=y=1\) is a stationary solution, but it's unstable. Due to round off errors, \(y^n=1\Rightarrow y^{n+1}=1-\epsilon\), triggering the instability. On the other hand, \(x^n=1\Rightarrow x^{n+1}=1\) due to \(C=0\). Hence, in opaque regions, an \(x/y\)-inconsistency will we reported.