| 74 | |
| 75 | Corner Transport Upwind attempts to modify predictor interface states with transverse flux corrections. |
| 76 | |
| 77 | 1. Spatial reconstruction and averaging of solution over domain $\lambda \Delta t$ to get interface states |
| 78 | 2. Application of source term to interface states over $\Delta t/2$ (or estimate source terms using cell centers) |
| 79 | 3. Calculate Fluxes from interface states |
| 80 | 4. Update interface states with transverse fluxes |
| 81 | 5. Calculate new Fluxes |
| 82 | 6. Update Cells using those fluxes |
| 83 | 7. Calculate source term using time averaged cell center |
| 84 | 8. Update final states using source term for $\Delta t$ |
| 85 | |
| 86 | |
| 87 | === Variant for source terms involving force fields (ie gradients) === |
| 88 | If the source terms can be directionally split, then the above algorithm can be modified as follows |
| 89 | |
| 90 | 1. Spatial reconstruction and averaging of solution over domain $\lambda \Delta t$ to get interface states |
| 91 | 2. Application of parallel source term to interface states over $\Delta t/2$ (or estimate source terms using cell center) |
| 92 | 3. Calculate Fluxes from interface states |
| 93 | 4. Update interface states with transverse fluxes |
| 94 | 5. Apply perpendicular source terms to interface states over $\Delta t/2$ (or estimate source terms using cell center) |
| 95 | 6. Calculate new Fluxes |
| 96 | 7. Update Cells using those fluxes |
| 97 | 8. Calculate source term using time averaged cell center |
| 98 | 9. Update final states using source term for $\Delta t$ |
| 99 | |
| 100 | === Variant for Self Gravity === |
| 101 | Same first 7 steps |
| 102 | |
| 103 | 8. Use conservative formalism for momentum transport (see https://astrobear.pas.rochester.edu/trac/blog/johannjc10222013) to calculate momentum fluxes associated with self-gravity |
| 104 | 9. Apply momentum fluxes from s-g to cells |
| 105 | 10. Use mass fluxes to calculate energy flux |
| 106 | 11. Make it second order by updating fluxes with new potential following elliptic solve. |
| 107 | |
| 108 | === Variant for Cylindrical Geometry === |
| 109 | |
| 110 | Most geometric source terms can be removed by adjusting the underlying algorithm and fluxes appropriately. See [http://iopscience.iop.org/article/10.1088/0067-0049/188/1/290/pdf Skinner & Ostriker 2010] - This has yet to be done for AstroBEAR |
| 111 | |
90 | | |
91 | | |
92 | | |
93 | | |
94 | | We can use the midpoint rule to estimate |
95 | | |
96 | | $Q_{L_{i+1/2}} = \displaystyle \frac{\int_{0}^{\Delta t} q_L(t) dt}{\Delta t} \approx q_L \left ( \frac{\Delta t}{2} \right ) + O(\Delta t^3)$ |
97 | | |
98 | | Usually some form of spatial reconstruction is used to calculate $q_L(0)$ and then the evolution equation for q is used to approximate time derivatives using spatial derivatives. |
99 | | |
100 | | $q_L(\Delta t/2) = q(x_L,0) + \displaystyle \int_0^{t/2} \frac{\partial q(x_L,t')}{\partial t} $ |
101 | | |
102 | | which we then use the evolution equation to replace time derivatives with spatial ones |
103 | | |
104 | | $q_L(\Delta t/2) = q_L(x_L,0) + \displaystyle \int_0^{t/2} \left [ -\frac{\partial f(x_L,t'))}{\partial x} + s(q(x_L,t')) \right ]$ |
105 | | |
106 | | and then we can linearize the flux function about $q$ |
107 | | |
108 | | $\frac{\partial f(q(x,t))}{\partial x} = \frac{\partial f}{\partial q} \frac{\partial q(x,t)}{\partial x} = \lambda \frac{\partial q(x,t)}{\partial x} $ |
109 | | |
110 | | |
111 | | $q_L(\Delta t/2) = q_L(0) + \displaystyle \int_0^{t/2} \lambda \frac{\partial (q(x_L,t'))}{\partial x} dt'+ \int_0^{t/2} s(q(x_L,t')) dt' $ |
112 | | |
113 | | Now if we ignore the source term, we can directly solve for |
114 | | |
115 | | $q(x,t)=q(x-\lambda t, 0)$ |
116 | | |
117 | | and under the change of variables $x'=x_L-\lambda t'$ we can rewrite the integrals as |
118 | | |
119 | | $q_L(\Delta t/2) = q_L(0) + \displaystyle \int_{x_L}^{x_L-\lambda t/2} \frac{\partial (q(x',0))}{\partial x} dx' + \int_{x_L}^{x_L-\lambda t/2} s(q(x',0)) dx'$ |
120 | | |
121 | | which just gives us |
122 | | |
123 | | $q_L(\Delta t/2) = q_L(0) - (q(x_L,0) - q(x_L-\lambda \Delta t/2)) + s(\bar{q})\Delta t/2$ |
124 | | |
125 | | where $\bar{q} = \displaystyle \int_{x_L-\lambda \Delta t/2}^{x_L}q(x',0) dx'$ |
126 | | |
127 | | |
128 | | |
129 | | $\frac{\partial f(q(x,t))}{\partial x} = \frac{\partial f}{\partial q} \frac{\partial q(x,t)}{\partial x} = A \frac{\partial q(x,t)}{\partial x} = R \Lambda L \frac{\partial q(x,t)}{\partial x} = \displaystyle \sum_i R_i \lambda_i L_i \frac{\partial q(x,t)}{\partial x} \approx \displaystyle \sum_i R_i \lambda_i L_i \frac{\partial q(x-\lambda_i t,0)}{\partial x} $ |
130 | | |
131 | | |
132 | | and use the characteristics (eigenvectors of the Jacobian of the flux) to linearize the flux |
133 | | |
134 | | $q_L(\Delta t/2) = q_L(0) + \displaystyle \int_0^{t/2} \left [ R \Lambda L \frac{\partial (q(x_L,t'))}{\partial x} + s(q(x_L,t')) \right ] $ |
135 | | |
136 | | $q_L(\Delta t/2) = q_L(0) + \displaystyle \int_0^{t/2} \left [ \sum R_i \lambda_i L_i \frac{\partial (q(x_L,t'))}{\partial x} + s(q(x_L,t')) \right ] $ |
137 | | |
138 | | $q_L(\Delta t/2) = q_L(0) + \displaystyle \int_0^{t/2} \left [ \sum R_i \lambda_i \delta l (L_i \frac{\partial (q(x_L,t'))}{\partial x} + s(q(x_L,t')) \right ] $ |
139 | | |
140 | | |
141 | | $q_L(\Delta t/2) = q_L(0) + \displaystyle \int_0^{t/2} \left [ \sum R_i \lambda_i \frac{\partial (q(x_L,t'))}{\partial x} + s(q(x_L,t')) \right ] $ |
142 | | |
143 | | |
144 | | |
145 | | and then we can use the solution for the characteristics $q(x_L,t')=q(x_L-\lambda t',0)$ to turn the time integral into a spatial integral |
146 | | |
147 | | $q_L(\Delta t/2) = q_L(0) + A \displaystyle \int_0^{t/2} \left [ \frac{\partial (q(x_L-\lambda t',0))}{\partial x} + s(q(x_L-\lambda t',0)) \right ] $ |
148 | | |
149 | | which we can do a variable substitution $x' = \lambda t'$ |