Changes between Version 118 and Version 119 of FluxLimitedDiffusion


Ignore:
Timestamp:
03/30/13 22:05:59 (12 years ago)
Author:
Jonathan
Comment:

Legend:

Unmodified
Added
Removed
Modified
  • FluxLimitedDiffusion

    v118 v119  
    326326
    327327  [[latex(\frac{\partial }{\partial t} \left ( \rho \mathbf{v} \right ) + \nabla \cdot \left ( \rho \mathbf{vv} \right ) = -\nabla P\color{green}{-\lambda \nabla E})]]   
    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
     331It 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.
    332332
    333333
     
    346346 * Update Edot,,0,, using pre IR and post IR E,,0,,
    347347 * Ghost e,,2*mbc,,, E,,mbc+1,,, Edot,,mbc+1,,
     348 * Store EDot in child arrays to be prolongated
    348349 * Do first EH,,mbc,,
    349  * Do ER,,mbc,, --- 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.
    350  * Store Edot in child arrays to be prolongated
     350 * Do ER,,mbc,, --- Only momentum terms.
    351351* Step 2
    352352 * Overlap d, p, e, E, and do physical BCs
    353353 * 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,,
    354356 * Update Edot,,0,, using pre IR and post IR E,,0,,
    355  * Update E,,1,, using Edot,,1,,
    356357 * Ghost e,,mbc,,, E,,1,,, Edot,,1,,
     358 * Store EDot in child arrays to be prolongated
    357359 * Do second EH,,0,,
    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.
    369361
    370362[[CollapsibleEnd()]]
     
    375367
    376368 [[latex(\frac{\partial }{\partial t} \left ( \rho \mathbf{v} \right ) =-\lambda \nabla E)]]   
    377  [[latex(\frac{\partial E}{\partial t}  = -\nabla \cdot \left ( \frac{3-R_2}{2}\mathbf{v}E \right ) )]] 
    378369
    379370These can be discretized as follows:
    380371
    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 ) )]]
    402373
    403374[[latex(\lambda_{i} = \frac{1}{R_{i}} \left ( \coth R_{i} - \frac{1}{R_{i}} \right ) )]]
     
    413384For now we will assume that [[latex(\kappa_{0P})]] and [[latex(\kappa_{0R})]] are constant over the implicit update and we will treat the energy as the total internal energy ignoring kinetic and magnetic contributions.  In this case we can solve the radiation energy equations:
    414385
    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
     389which 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
     396Now 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
     399so 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
     403and 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
     407which 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
     409Then 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
     413and 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
     417which 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
    417422
    418423=== Expanding about e,,0,, ===