Posts for the month of June 2017

256^3

The color shows the neutral density in particles/cc. Also shown are the inner, outer, and tau=1 surface (dashed lines) given by John's python script along with the calculated tau=1 from the code (solid black line). Mach surface is shown in red

movie

Comparison with AMR

Both with recombination/cooling, 2e13 flux, no Lyman alpha cooling

323 + 2 levels AMR

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+20000.png

frames

1283 fixed grid

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux0000.png

frames

Lineouts

Circles are 323 + 2 levels AMR, triangles are 1283 fixed grid. The 1283 fixed grid is 200 frames shorter - the last frame is repeated. Day side looks nearly identical - night side differs, probably due to lower resolution.

Density

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+2_comparison_rho0000.png

Ionization fraction

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+2_comparison_X0000.png

Optical depth

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+2_comparison_tau0000.png

Temperature

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+2_comparison_temp0000.png

Heating/cooling time scales

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+2_comparison_timescales0000.png

Heating/cooling rates

http://www.pas.rochester.edu/~adebrech/PlanetIonization/32+2/32+2_comparison_rates0000.png

Optical Depth

Planet To-do

  • Add diagnostics for cooling and recombination. ✓
  • Subcycling time for planet. ✓
  • Add flag for Lyman cooling. ✓
  • Check cooling terms for correctness. ✓
  • Ramp up flux over time. ✓
  • Check recombination on night side. ✓
  • With stored optical depth, calculate ionization and ionization heating rate (can't be done straightforwardly with a diagnostic variable). ✓
  • Jonathan will be adding optical depth variable. ✓

Runs

  • Lyman cooling, recombination/cooling, 2e13 flux ✓

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling_rho.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling_X.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling_tau.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling_temp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling_timescales.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/all_cooling/all_cooling_rates.gif

frames

  • Recombination/cooling, 2e13 flux ✓

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux_rho.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux_X.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux_tau.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux_temp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux_timescales.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_flux/recom_flux_rates.gif

frames

  • Recombination/cooling, no flux ✓

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only_rho.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only_X.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only_tau.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only_temp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only_timescales.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/recom_only/recom_only_rates.gif

frames

  • Lyman cooling, recombination/cooling, 2e10 flux ✓

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux_rho.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux_X.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux_tau.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux_temp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux_timescales.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/low_flux/low_flux_rates.gif

frames

  • Calculate new profile for planet unaffected by Lyman cooling ✓

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet_rho.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet_X.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet_tau.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet_temp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet_timescales.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/new_planet/new_planet_rates.gif

frames

  • Lyman cooling, recombination/cooling, ramping flux 1.2e13 → 2.02e13 ✓

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp_rho.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp_X.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp_tau.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp_temp.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp_timescales.gif

http://www.pas.rochester.edu/~adebrech/PlanetIonization/flux_ramp/flux_ramp_rates.gif

frames

  • Recombination/cooling, ramping flux 3.6e12 → 2e13 - running

Comparison with O+16a, attempt 2

Recap of Last Post

  • In the last post, I performed a run with the same parameters as the Ohlmann+16a run.
  • However the simulation went from 2-body to 3-body when the secondary was mistakenly reintroduced after about 7 days.
  • This was not a bug in the code, just me forgetting to flip a 'switch'.

New Work

  • I flipped the switch and reran the code up to 20 days.

Summary of New Results

  • The results match qualitatively those of O+16a. Differences are likely caused by remaining differences in the numerical setups, including softening radius, ambient density, resolution, and degree of corotation.

Detailed Results

I) Circular orbit as in Ohlmann16a
Damp088) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient dyne/cm
(comet compute 1728 cores, 2 cores per task to increase memory per task)
( cm, , 5 levels AMR)
Restarted from run Damp062, at s, just after Damping stopped, to s
2d density
2d density (zoomed)
2d density and velocity (zoomed)
2d density edge-on slice (zoomed)
2d particle velocity (zoomed)
2d particle position (zoomed)
2d Particle mass (zoomed)

Here are PROJECTIONS, instead of slices, zoomed in with a better color scheme, and extended to go up to 20 days:
2d density projection
2d density projection edge-on

Here is a version of the SLICE in the orbital plane, going up to 20 days, rotated 180 degrees to match the orbit of O+16a, color range shifted to match O+16a, and with a better color scheme:
2d density
…or with blue labels for better visibility: 2d density
2d density (zoomed)
Here is temperature:
2d temperature
…and Mach number:
Mach number

Below I've plotted Fig 1 of Ohlmann+16a and the equivalent figure with this simulation for comparison. In my plot of separation, I've shown the softening radius as a dotted horizontal line. In the inset showing particle trajectories, I've inverted the axes and shifted the origin to allow direct comparison (our simulation is rotated by 180 degrees with respect to theirs). Circles represent the softening radius, while the green dot shows the initial center of mass of the two stars. The very small green square in the upper right shows the smallest resolution element.
http://www.pas.rochester.edu/~lchamandy/Graphics/Ohlmann/fig1.pdf http://www.pas.rochester.edu/~lchamandy/Graphics/RGB/Post-sink_particle/Post-modified_Lane_Emden/Damp088/p20d.gif

Below are snapshots of density in the orbital plane from t = 10 days and t = 20 days. Note that to compare directly, I have rotated each figure by 180 deg. In the O+16a figure, the '+' denotes the primary core particle, while the 'x' denotes the secondary point particle. The core and secondary are denoted as '0' and '1' respectively in my plot.
http://www.pas.rochester.edu/~lchamandy/Graphics/RGB/Post-sink_particle/Post-modified_Lane_Emden/Damp088/rho2d_zoom_blue_0043.png http://www.pas.rochester.edu/~lchamandy/Graphics/RGB/Post-sink_particle/Post-modified_Lane_Emden/Damp088/rho2d_zoom_blue_0086.png
http://www.pas.rochester.edu/~lchamandy/Graphics/Ohlmann/fig3.pdf

Discussion

  • The results match qualitatively those of Ohlmann+16a, but differ in detail.
  • The orbit of the secondary is larger, while that of the primary is slightly smaller, compared to O+16a. The separation does not become as small at the first minimum or as large at the first maximum, compared to O+16a. This may be due to the ~5 times larger gravitational softening length used.
  • At t = 10 days, the density in the orbital plane appears similar to that of O+16a, but the spiral arm is less pronounced and features are less detailed. This could be due to the resolution being lower by a factor of somewhere between 2 and 12. It may be partly caused by the weaker gravitational interaction due to larger softening radius mentioned above.
  • Also, the common envelope has not expanded as much into the ambient medium. This could be due to the >107 times larger ambient density used compared to that of O+16a.
  • Recall the differences between our setup and O+16a, outlined in the Table 2 of the following pdf file. The main differences between the two setups are:
    • 5x larger softening radius
    • 107 times larger ambient density
    • up to 12 times smaller resolution
    • 9x smaller box size
    • no initial rotation vs 95% solid-body corotation with orbit
    • some differences in the relaxation run to prepare the red giant star (detailed in the table and in previous posts)
    • extrapolated instead of periodic BCs

Next Steps

  • The easiest and most obvious thing to do is to reduce the softening radius (= cutoff radius defening RG core. Recall that for r smaller than this radius, we replace the core with a modified Lane-Emden solution). However, reducing this softening radius to 1 Rsun and keeping the same max resolution would mean that we resolve the softening radius by only 3.5 resolution elements, whereas O+16a claim that 10 resolution elements is needed.
  • That suggests that we should increase the max resolution by a factor of 2 or 4, but this would slow down the code considerably. It would be best to first lower the softening length without changing the resolution.
  • It would also be worth testing lower ambient densities. In past tests the code has crashed, but it should be possible to lower the ambient density by at least a factor of 10 or 30.
  • So this suggests the following 3 runs (starting from the beginning of the relaxation run):
    • As 088 but reduce the ambient density by 10 or 30
    • As 088 but reduce the softening radius by 4.8 to equal that of O+16a
    • As 088 but with both of the above changes

Planet with Line Transfer

  • Fixed minor typo in scaling of hydrogen cross section and added floor to line transfer dt. Pushed to line_transfer branch.
  • Now appears to be working as expected (as far as absorbing the flux, anyway). The dynamics become more interesting around frame 50 (~4 hours).

Reproducing Ohlmann+16a results: 1st attempt

Recap of Last Post

  • In the last post, I performed several runs with varying separation and secondary mass.
  • However it was determined that point particles were not feeling the gravity of the gas, so the results were wrong.

New Work

  • I wrote a script in IDL which calculates the Kepler orbit given the initial velocities, and used it to show that the orbit obtained was just what would be expected if the point particles were not feeling the gas gravity, further confirming the above conclusion.
  • This problem with the code was corrected by Baowei and Jonathan.
  • The base case (run 084 from last blog) was rerun now with the code corrected, with the goal of reproducing the basic results from Ohlmann+16a.

Summary of New Results

  • The simulation now looks "reasonable" at the beginning (first ~dynamical time) but then deviates from the results of Ohlmann+16a.
  • I realized while writing this blog that the problem is caused by the mistaken re-introduction of the secondary (hence the introduction of a third particle with mass equal to that of the secondary) on restarts. In fact, there was a switch in the problem.data file called "secondary_already_present" to deal with restarts, but I had forgotten to change it to .true. before restarting.

Results

I) Circular orbit as in Ohlmann16a
Damp088) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient dyne/cm
(comet compute 1728 cores, 2 cores per task to increase memory per task)
( cm, , 5 levels AMR)
Restarted from run Damp062, at s, just after Damping stopped, to s
2d density
2d density and velocity
2d density and velocity (zoomed)
2d temperature (zoomed)
2d particle velocity (zoomed)
2d particle position (zoomed)
2d Particle mass (zoomed)

Particle trajectories, (with circles representing softening radius)
Particle separation with time
Particle trajectories up to mistaken introduction of third particle
Particle separation with time up to mistaken introduction of third particle

Discussion

  • What has happened is that another secondary got introduced on each restart. This happens at frame 107. It also happens at frame 92 but this won't matter much because I had to redo that frame only on bluehive because there had been an I/O problem at that frame on comet. So the simulation is correct up to frame 106, or elapsed time t = 2.12e6 - 1.5e6 = 6.2e5 s = 7.2 days, about 1 dynamical time, but is wrong after that.
  • Results up to the kink seen in the plot of separation vs time are very close to those of Ohlmann+16a (see below). The introduction of the third point particle occurs before the kink. So the results so far are consistent with those of Ohlmann+16a. (If one looks closely there is another point in the plot of separation vs time that is a little off, and this is due to the one frame that was done on bluehive where the 3rd particle also appeared.)

Next Steps
Once I've edited the problem.data file to make "secondary_already_present=.true.", I will redo the simulation from frame 106. I will then again compare the result with Ohlmann+16a Figure 1.

Here are the other figures from that paper for reference:
Figure 2
Figure 3
Figure 4

Line transfer

3D run, 3 levels AMR, 64x64x64, initial ionization fraction 10-5 (0 was crashing, but I believe that's been fixed), 1 computational time = 1 day.

Lineouts at 6.144 hours:

Comparison to Input Profile

Movie goes to 13 hours.

Update on common envelope simulations

Recap of Last Post

  • In the last post I presented a simulation of a RG translating across the grid.
  • I also presented a common envelope simulation but later realized that the initial velocities of the primary and secondary were wrong (about a factor of 2 too small). I corrected the script to calculate the initial velocities and tested it against the Sun-Earth system.

New Work

I corrected the initial velocities and simulated the following cases:

  • RG star of mass 1.956 Msun including 0.369 Msun core (point mass), with secondary point mass of 0.978 Msun (half the RG mass), at initial separation 49 Rsun, just outside of the RG outer radius of about 48 Rsun. Very close to parameters of Ohlmann+16a.
  • As above but larger separation of 60 Rsun instead of 49 Rsun.
  • As above but even larger separation of 98 Rsun instead of 49 Rsun.
  • As the first case but now make secondary mass equal to only 0.001 Msun so that problem reduces to "1-body".
  • Replace RG with a point particle of the same mass and run in low resolution with uniform background to test orbital dynamics.
  • Same but with original point masses of the first case, i.e. the RG core and the secondary, without RG envelope.

Summary of New Results

  • We do not achieve a circular orbit in the CE sims. In each case the orbit is looser than it should be. The binary separation grows with time and the speed of the secondary reduces with time.
  • When the RG is replaced by a point mass then the correct circular orbit is obtained (in low res with a uniform ambient medium).
  • If the RG envelope is excluded and the point masses representing the core and secondary are made to orbit in a uniform ambient medium, we obtain almost the identical (non-circular) orbit as for the full CE simulation.
  • The natural explanation for these results is that the point masses can feel each other's gravity, the gas can feel the gravity of the point masses, but the point masses do not feel gravity from the gas! This feature must obviously be included in the next set of runs.

Results

I) Circular orbit as in Ohlmann16a
Damp084) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient dyne/cm
(comet compute 1728 cores, 2 cores per task to increase memory per task)
( cm, , 5 levels AMR)
Restarted from run Damp062, at s, just after Damping stopped, to s
2d density
2d density and velocity
2d density and particle velocity
2d density edge-on
2d Temperature
2d Particle position
2d Particle mass

II) As (I) above but wider orbit (binary separation of 60Rsun instead of 49Rsun)
Damp086) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient dyne/cm
(bluehive standard 120 cores)
( cm, , 5 levels AMR)
Restarted from run Damp062, at s, just after Damping stopped
2d density
2d density and velocity

III) As (I) above but even wider orbit (binary separation of 98Rsun instead of 49Rsun)
Damp087) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient dyne/cm
(comet compute 1728 cores, 2 cores per task to increase memory per task)
( cm, , 5 levels AMR)
Restarted from run Damp062, at s, just after Damping stopped
2d density and velocity

IV) As (I) above but very low mass secondary (0.001Msun instead of 0.978Msun)
Damp085) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient dyne/cm
(stampede normal 1024 cores)
( cm, , 5 levels AMR)
Restarted from run Damp062, at s, just after Damping stopped
2d density andparticle velocity

V) As (I) above but RG represented by a point mass, so two point masses orbiting on uniform background
Kepler083) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient density and pressure set to 1d-10
( cm, , 0 levels AMR)
(Started from t=0)
2D projected particle positions
particle velocities
2D projected particle positions, edge-on

VI) As (V) above but only mass of RG core is included in primary point mass, so same point masses as in (I)
Kepler083core) Extrapolated hydro BCs, Multipole expansion Poisson BCs, ambient density and pressure set to 1d-10
(stampede normal 1024 cores)
( cm, , 0 levels AMR)
(Started from t=0)
2D projected particle positions
particle velocities

Comparison between (I) and (VI)
Particle velocities

COMMENTS: Clearly, the secondary does not feel the gravity of the RG envelope. Thus, point masses do not feel gravity of the gas. This obviously needs to be corrected for the next set of runs.

Next Steps
Once this problems is resolved, we need to redo run (I) in particular. If everything is in order, we can then plot the orbital separation as a function of time, and compare the result with Ohlmann+16a Figure 1.

Here are the other figures from that paper for reference:
Figure 2
Figure 3
Figure 4

Aside from resolution and other details, the biggest difference between our setup and that of Ohlmann+ is that their RG is initiated with 95% corotation with the orbit, whereas ours is not rotating initially.

Other points worth mentioning

  • My quota on bluehive has been reduced from 24 TB to 16 TB. I will look to reduce it more in the next two weeks. My quota on blue streak has been reduced from 5 TB to 2 TB.
  • There is a nice paper by MacLeod et al. about accretion onto the secondary during the common envelope that I'm reading at the moment. Aside from being interesting, it is quite pedgogically written. I'm doing a literature review to prepare for the upcoming conference.

Update 6/14

  • Implemented multispecies gas for Roe solver in AstroBEAR, but need a solution for density or sound speed for van der Waals gas to implement it. Relevant code snippet:
     ! Entropy fix      !!! Implement a different entropy fix to use van der Waals gas - density equations have complicated/no solution
     lambdaLL = WL(ivx) - PrimSoundSpeed(WL)
     lambdaRR = WR(ivx) + PrimSoundSpeed(WR)
     if (lambdaLL < 0 .or. lambdaRR > 0) then
      call StarPU(WL((/ irho, iE, ivx /)), WR((/ irho, iE, ivx /)), PrimSoundSpeed(WL), PrimSoundSpeed(WR), pstar, ustar)
      ! Passively advected (tracers; y,z momentum)
      ExactL = WL
      ExactR = WR
      ! Solution in interaction region
      ExactL(iE) = pstar
      ExactL(ivx) = ustar
      ExactR(iE) = pstar
      ExactR(ivx) = ustar
      ! Densities differ
      if (mpi_id == 0 .and. iEOS == EOS_vanDerWaals) write (*,*) 'This solution for density will not be accurate. The attraction parameter is not yet included in these calculations.'
      ExactL(irho) = PrimExactDensity(WL,pstar)
      ExactR(irho) = PrimExactDensity(WR,pstar)
      
      lambdaLR = ustar - PrimSoundSpeed(ExactL)
      lambdaRL = ustar + PrimSoundSpeed(ExactR)
     end if
     if (lambdaLL < 0 .and. lambdaLR > 0) then    ! Check for transonic rarefaction on left side
      print *, 'Using entropy fix.'
      lambda(1) = lambdaLL*((lambdaLR - lambda(1))/(lambdaLR - lambdaLL))
     end if
     if (lambdaRL < 0 .and. lambdaRR > 0) then    ! check for transonic rarefaction on right side
      print *, 'Using entropy fix.'
      lambda(nDim+2) = lambdaRR*((lambda(nDim+2) - lambdaRL)/(lambdaRR - lambdaRL))
     end if
  • Wrote WAF scheme for own solver. Working, but results are not as good as I expect based on Toro's. Not yet sure where the problem lies.

http://www.pas.rochester.edu/~adebrech/code/debugging/test1_HLLC.gif

  • Related: Transmissive boundary conditions for second set of ghost cells don't quite make sense yet. Why move back inward to get wave to flow out?
  • Still getting explosion (but better than implosion) with adiabatic atmosphere.

Should be the same as Jonathan's result - not sure why it appears to not be working.

CE

In this simulation, I placed an outflow sphere at the center. The boundary condition is 30000K, 1e-13 g/cc and 300km/s from 0 to 100 day (mass loss rate is 7.2e-6 solar mass / yr) at the radii of 50 R_sun then the outflow speed changes to 10km/s. This outflow is not doing much at this timescale.

I launched another outburst from 1000 to 1100 day with 30000K, 1e-10 g/cc and 1000km/s. Thus the mass loss rate is 0.024 solar mass / yr.

topdensity

sidedensity

Stable Atmosphere

I've added some data and a Matlab script for creating images such as those below. Density is on a linear plot- looks pretty good overall, although there's a strange wiggle near the inner boundary. Pressure is on a semilog plot, and is orders of magnitude off (I believe these should be in the same (computational) units).

EDIT: After running for longer time, not as stable as it first appeared (although I'm not seeing the same values in the lineouts as before, either). Perhaps the outer boundary is not as balanced as it should be?

Density plot (same scale as below):

Lineouts:

Pressure lineout scale changes dramatically after first frame, but can't get VisIt to save window correctly.

Ionization

Time stepping too large currently for parameters.


To mimic Anjali's low flux case but instead with a solid sphere, let's use an ideal gas with the following setup:

, , , , , , .

Meeting Update --6/2/17

  • Disk Space
    1. archived MachStems data. Currently on bluehive, /scratch/afrank_lab has 4.4TB available.

Update 6/2

  • Roe solver implemented in AstroBEAR:

http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test1.gif http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test3.gif http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test3_noP.gif http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test4.gif http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test4_noP.gif http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test5.gif http://www.pas.rochester.edu/~adebrech/code/AstroBEAR_test5_noP.gif

Still need to implement van der Waals gas and figure out how to implement multispecies (and find tests).

  • Testing charge exchange:

http://www.pas.rochester.edu/~adebrech/code/ChargeExchangeTest/ChargeExchange_part.gif

frames

Seems to be working pretty much as expected, although the proportions of cold HII to hot H seem a little wonky - likely because of the other physics present, though.

  • WASP-12 paper nearly there. Got some comments back from Luca, so I'll look at those in the next few days.