328 | | [[latex(\frac{\partial e}{\partial t} + \nabla \cdot \left [ \left ( e + P \right ) \mathbf{v} \right ] = \color{red}{-\kappa_{0P}(4 \pi B-cE)} \color{orange}{+\lambda \left ( 2 \frac{\kappa_{0P}}{\kappa_{0R}}-1 \right ) \mathbf{v} \cdot \nabla E} \color{blue}{-\frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E})]] |
329 | | [[latex(\frac{\partial E}{\partial t} \color{red}{ - \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E} = \color{red}{\kappa_{0P} (4 \pi B-cE)} \color{orange}{-\lambda \left(2\frac{\kappa_{0P}}{\kappa_{0R}}-1 \right )\mathbf{v}\cdot \nabla E} \color{green}{-\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right )} \color{blue}{+\frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E} )]] |
330 | | |
331 | | For static diffusion, the terms in blue with v^2^/c can be dropped and the system can be split into the usual hydro update (black), radiative source terms (green), and a coupled implicit solve (red) for the radiation energy density and thermal energy density (ie temperature). And the terms in orange can be solved semi-implicitly (velocity terms are explicit, while radiation energy terms are implicit) |
| 328 | [[latex(\frac{\partial e}{\partial t} + \nabla \cdot \left [ \left ( e + P \right ) \mathbf{v} \right ] = \color{red}{-\kappa_{0P}(4 \pi B-cE)} \color{orange}{+\lambda \left ( 2 \frac{\kappa_{0P}}{\kappa_{0R}}-1 \right ) \mathbf{v} \cdot \nabla E} \color{orange}{-\frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E})]] |
| 329 | [[latex(\frac{\partial E}{\partial t} \color{red}{ - \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E} = \color{red}{\kappa_{0P} (4 \pi B-cE)} \color{orange}{-\lambda \left(2\frac{\kappa_{0P}}{\kappa_{0R}}-1 \right )\mathbf{v}\cdot \nabla E} \color{orange}{-\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right )} \color{orange}{+\frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E} )]] |
| 330 | |
| 331 | It is possible that with included the terms in orange in a semi-implicit method, the dynamic diffusion regime may be stable... In any event, it costs very little to add all of the terms in orange to the implicit solve (using the old velocity). Then the momentum update can be done explicitly - though using a time centered radiation field. |
346 | 346 | * Update Edot,,0,, using pre IR and post IR E,,0,, |
347 | 347 | * Ghost e,,2*mbc,,, E,,mbc+1,,, Edot,,mbc+1,, |
351 | 351 | * Step 2 |
352 | 352 | * Overlap d, p, e, E, and do physical BCs |
353 | 353 | * Do IR which updates e,,0,,, and E,,0,, using d,,1,,, e,,1,,, E,,1,, |
| 354 | * Update E,,mbc,, using Edot,,mbc,, |
| 355 | * Update e,,2*mbc,, using E,,2*mbc,,, Edot,,2*mbc,,, and e,,2*mbc,, |
358 | | * Do ER,,0,, --- Terms with grad E can be done without ghosting since EH did not change E. The del dot vE term needs time centered face centered velocities which can be stored during the hydro update. |
359 | | |
360 | | |
361 | | However with AMR is we need coarse boundary values for E at t=0, and t=dt, and t = 2dt and we would like coarse boundary values at t=2dt to match solution from coarse grid... |
362 | | |
363 | | But the coarse update involves one implicit and one explicit solve. Is there a way to interpolate the fine grid ghost zones in time without storing two time derivatives? If we store dE_i+dE_e then we could do an explicit radiative update internally (E^n^ -> E*), then an implicit radiative update (E* -> E^n+1^ using the fully updated ghost zones (E^n+1^,,g,, = E^n^,,g,,+dE,,g,,) and then store the new time derivative (dE=E^n+1^-E^n^) |
364 | | |
365 | | This seems like it would work for the radiative field, but what about the internal energy terms? We have an equation for e^n+1^ that is a function of E^n+1^, e* and E* but we have no way of getting E* in the ghost zones... |
366 | | |
367 | | So the solution is perhaps to store the contributions from the implicit update and the non-conservative heating terms that appear in the internal energy equation. Then we can update E in the ghost zones for the implicit and the first explicit term using the time derivative - and then use the new E to update e... Then the conservative explicit RadEnergy term can be calculated after the hydro step as well as the momentum explicit term - and then this flux can be coarsened - to keep the value of E in the coarse cells consistent with the value of E in the fine cells... Since both will have been updated by the same eDot and fluxes... |
368 | | |
| 360 | * Do ER,,0,, --- Only momentum terms. |
381 | | [[latex(p^{n+1}_i=p^n_i - \frac{1}{2}\frac{\Delta t}{\Delta x} \lambda_{i} \left ( E^n_{i+1}-E^n_{i-1} \right ) )]] |
382 | | |
383 | | [[latex(E^{n+1}_i=E^n_i + \frac{\Delta t}{\Delta x} \left ( F_{i+1/2}-F_{i-1/2} \right ) )]] |
384 | | |
385 | | where |
386 | | [[latex(F_{i+1/2}=\frac{3-R_{2,i+1/2}}{8} \left(v_{i}^n+v_{i+1}^n \right ) \left ( E^n_i+E^n_{i+1} \right ) )]] |
387 | | |
388 | | where |
389 | | |
390 | | [[latex(R_{2,i+1/2} = \lambda_{i+1/2}+\lambda_{i+1/2}^2 R_{i+1/2}^2)]] |
391 | | |
392 | | and |
393 | | [[latex(R_{i+1/2} = \frac{\left | E^n_{i+1}-E^n_{i} \right | }{2 \kappa_{0R,i+1/2} \left ( E^n_i+E^n_{i+1} \right )})]] |
394 | | |
395 | | and |
396 | | [[latex(\lambda_{i+1/2} = \frac{1}{R_{i+1/2}} \left ( \coth R_{i+1/2} - \frac{1}{R_{i+1/2}} \right ) )]] |
397 | | |
398 | | and |
399 | | [[latex(\kappa_{0R,i+1/2} = \frac{\kappa^n_{0R,i}+\kappa^n_{0R,i+1}}{2} \mbox{ and } \lambda_{i+1/2} = \frac{1}{R_{i+1/2}} \left ( \coth R_{i+1/2} - \frac{1}{R_{i+1/2}} \right ) )]] |
400 | | |
401 | | and |
| 372 | [[latex(p^{n+1}_i=p^n_i - \frac{1}{4}\frac{\Delta t}{\Delta x} \lambda_{i} \left ( \left ( E^{n}_{i+1}+E^{n+1}_{i+1} \right ) - \left ( E^n_{i-1}+E^{n+1}_{i-1} \right ) \right ) )]] |
415 | | [[latex(\frac{\partial e}{\partial t} = -\kappa_{0P}(4 \pi B-cE) + \lambda \left ( 2 \frac{\kappa_{0P}}{\kappa_{0R}}-1 \right ) \mathbf{v} \cdot \nabla E )]] |
416 | | [[latex(\frac{\partial E}{\partial t} = \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E + \kappa_{0P} (4 \pi B(T)-cE) -\lambda \left(2\frac{\kappa_{0P}}{\kappa_{0R}}-1 \right )\mathbf{v}\cdot \nabla E )]] |
| 386 | [[latex(\frac{\partial e}{\partial t} = -\kappa_{0P}(4 \pi B-cE) + \lambda \left ( 2 \frac{\kappa_{0P}}{\kappa_{0R}}-1 \right ) \mathbf{v} \cdot \nabla E -\frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E )]] |
| 387 | [[latex(\frac{\partial E}{\partial t} = \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E + \kappa_{0P} (4 \pi B(T)-cE) -\lambda \left(2\frac{\kappa_{0P}}{\kappa_{0R}}-1 \right )\mathbf{v}\cdot \nabla E + \frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E -\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right ) )]] |
| 388 | |
| 389 | which we can also write as |
| 390 | |
| 391 | [[latex(\frac{\partial e}{\partial t} = f \left ( e,E,\nabla E \right ) |
| 392 | [[latex(\frac{\partial E}{\partial t} = \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E -\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E - f \left (e,E, \nabla E \right ) \right ) )]] |
| 393 | where |
| 394 | [[latex( f \left (e,E, \nabla E \right ) = -\kappa_{0P}(4 \pi B-cE) + \lambda \left ( 2 \frac{\kappa_{0P}}{\kappa_{0R}}-1 \right ) \mathbf{v} \cdot \nabla E -\frac{3-R_2}{2}\kappa_{0P}\frac{v^2}{c}E )]] |
| 395 | |
| 396 | Now we can linearize f about e,,0,, |
| 397 | [[latex( f \left (e, E, \nabla E \right ) = f (\left e_0, E, \nabla E \right ) + \phi e )]] |
| 398 | |
| 399 | so that the first equation can be written as |
| 400 | |
| 401 | [[latex(\frac{\partial e}{\partial t} = f \left ( e_0,E,\nabla E \right ) + \phi e)]] |
| 402 | |
| 403 | and then discretized as |
| 404 | |
| 405 | [[latex(e^{n+1}_i-e^{n}_i = \frac{\Delta t}{\Delta x} f \left ( e_0,E*,\nabla E* \right ) \right ) + \bar{\psi} \phi e^n_i + \psi \phi e^{n+1}_i \right )]] |
| 406 | |
| 407 | which can be solved for [[latex(e^{n+1} = \frac{1}{1-\psi \phi} \left ( \bar{\psi} \phi e^n_i + f \left ( e_0,E*,\nabla E* \right ) \right ) )]] |
| 408 | |
| 409 | Then if we take the semi-discretized equation for E |
| 410 | |
| 411 | [[latex(\frac{\partial E}{\partial t} = \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E -\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right ) - f \left (e_0,E, \nabla E \right ) - \bar{\psi}\phi e^n_i -\psi \phi e^{n+1}_i )]] |
| 412 | |
| 413 | and then plugin the solution for e^n+1^,,i,, we get |
| 414 | |
| 415 | [[latex(\frac{\partial E}{\partial t} = \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E -\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right ) - f \left (e_0,E, \nabla E \right ) - \bar{\psi}\phi e^n_i - \frac{\psi \phi }{1-\psi \phi} \left ( \bar{\psi} \phi e^n_i + f \left ( e_0,E*,\nabla E* \right ) \right ) )]] |
| 416 | |
| 417 | which simplifies to |
| 418 | |
| 419 | [[latex(\frac{\partial E}{\partial t} = \nabla \cdot \frac{c\lambda}{\kappa_{0R}} \nabla E -\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right ) - \frac{1}{1-\psi \phi} \left ( \bar{\psi} \phi e^n_i + f \left ( e_0,E*,\nabla E* \right ) \right ) )]] |
| 420 | |
| 421 | |