​https://astrobear.pas.rochester.edu/trac/wiki/u/erica/scratch10

Work in progress

Simulation Status

Attempting to run Run4 until transients die down. After 20 hours, timestep is still about 1x10^-7 - doesn't look like it's going to be successful. It's almost certainly related to the high ambient temperature of ~10¹³ K (not sure how much the transient flows are important in the final analysis - sound speed is very high).

Could try lowering density ratio based on recombination timescale (make sure it's large). Could also up resolution so planet gets closer to zero density, but this wouldn't be particularly effective in saving computational time (because of line transfer step - planet itself is only small fraction of grid, so high resolution hydro there is not a problem - which brings an idea - limit line transfer to something higher than max level?).

Movies in progress for corotating Run2.

Radiation pressure

Merged charge exchange and radiation pressure into local line_transfer branch. Need to test to make sure merge was successful, then can queue radiation pressure version of Run2.

WASP-12b Hill vs. bow shock radii

Using Matlab script for planetary parameters, we would need a stellar wind density (at base) of about the same as the planetary wind density:

\rho_\star = 4.95\times 10^{-11} \text{g/cm}^3

with

\rho_P = 5.881\times 10^{-11} \text{g/cm}^3

Flux to ionize Jupiter

An O-type star would need to provide about 1.2x10³⁸ ionizing photons per second to a planet at Jupiter's orbit for the equivalent flux. An O star actually provides about 10⁵⁰ phot/s, resulting in a flux of 1.3x10²¹ phot/s/cm² (or about 4 orders of magnitude more than our high-flux case above).

New Work

1. We solved the problem of mixing variables from different meshes (fluid and particle), which enables plotting e.g. fluid velocities relative to the companion.
2. Extended simulation 143 to 40 days.
3. Figured out how to make face-on movies in rotating frame (e.g. rotating along with orbit), and made movie of density.
4. Made movie of density in frame of orbit with velocity vectors relative to the companion.
5. Plotted circular velocity around the companion, as well as its ratios with the Keplerian circular velocity and with the sound speed.

Results (of new work, point by point from above)

1. Apparently there exists a more robust way to do this involving Cross Mesh Field Evaluation or CMFE functions, but Baowei's "quick fix" works!
2. Will not go beyond 40 days for now as it is a bit too expensive. Plan is to try to get away with less refinement in ambient region, and also to make central refinement zone smaller in volume. If results are the same, can then get away with this to do more runs.
3. This allows us to avoid getting quite as dizzy when looking at the movies!
4. We can now see the velocity field in the frame of the companion, rather than just in the lab frame.
5. These are key quantities that we are trying to investigate. The ratios are less than unity near the companion but as high as ~0.5.

Run 143:

• Similar to run 132 but with subgrid accretion turned off, twice larger box, resolution and softening length that evolve with time, less aggressive refinement, i.e. larger max refinement zones, and no relaxation (damping) run initially
  

Relaxation run: no relaxation run
First frame: 0
Last frame: 173
Total simulation time: 40 sim-days
Machine and partition: Stampde 2 normal
Number of nodes: 128 or 64, each with 68 (standard nodes) or 96 (skx nodes) cores
Total wall time: TBD
Hydro BCs: extrapolated
Poisson BCs: multipole expansion
Box size: L=8e13 cm (1150 Rsun)
Max refinement level: 4 (frames 0 to 72), then 5 (frames 72 to 117) Base resolution: 2.25 Rsun (256³ cells)
Highest resolution: 0.14 Rsun (4096³ cells, 4 levels AMR, frames 0 to 72), and 0.07 Rsun (8192³ cells, 5 levels AMR, frames 72 to 173)
AMR implementation: set by hand to have max level around RG core (frames 0 to 72) or companion (frames 72 to 173)
Max resolution zone: sphere around RG core, radius 5d12 cm (frames 0 to 46), radius 4d12 cm (frames 46 to 72), radius 3d12 (frames 72 to 103), radius 2.5d12 (frames 103 to 161), radius 1.75d12 (frames 161-173)
Buffer zones: 16 cells
Softening length: 2.4 Rsun (frames 0 to 72), 1.2 Rsun (frames 72 to 173)
Ambient density: 6.7e-9 g/cc
Ambient pressure: 10⁵ dyne/cm²
DefaultAccretionRoutine=0 (no accretion)

New Movies of run 143

1. High zoom-in face-on slice, centered on companion, in reference frame of companion that is corotating with the instantaneous orbital angular velocity, with RG core always on the left.
   

​face-on slice, movie rotating with orbit

2. As above, density in the frame rotating with orbit, but now including velocity vectors drawn in the reference frame of the companion, but in the non-rotating frame:

​face-on slice, movie rotating with orbit, with velocity vectors in frame of companion (but not in any rotating frame)

New Snapshots of run 143
Left: v_phi relative to the companion. Face-on, slice through companion, view centered on companion with side of 4e12 cm. RG core is visible at the left. Right: same but zoomed out by factor of 4.

As above but now the ratio of v_phi to the Keplerian circular velocity.

Same as the zoomed-in plot above but different color scheme and range.

As above but now the ratio of v_phi to the local sound speed.

Inter-particle separation (up to 40 days)
Red horizontal lines show radius of maxLevel refinement region. Green horizontal line shows softening length (solid) or five times softening length (dashed) for each particle. Finest resolution is proportional to softening length, i.e. resolution improves by factor of 2 when softening length is reduced by factor of 2. Vertical lines show transition from one refinement radius to the next (dotted red, 5e12 to 4e12 to 3e12 to 2.5e12 cm), and transition from one softening radius to the next (green dashed, ~2.4 Rsun to ~1.2 Rsun)

Examples of streamline (left) and integral curve (right) plots for smaller run 125

Next steps

• Line Integral Convolution using Jonathan's MatLab code
• Plot streamlines in VisIt
• Countours of energy to see if material is bound

The following is from Bruce's email. Just want to put here and see if any comments/ideas:'

Thin knots seem to arise in many outflow models along the y-axis shortly after the launch of a jet. In brief, I’m convinced that the biggest cause of such knots is the shape of the nozzle’s surface (a sphere). A flat or highly conical nozzle will suppress the knots.

The simplest flow is that of a cylindrical jet at the origin moving into an ambient medium of constant density on a Cartesian grid. In principle, such a flow has no way to deviate from a simple cylindrical flow unless shears (at the edge) or kink instabilities develop (they don’t).

Heavy flows: This is obviously the case if the flow density > ambient density. The flow is simply a telephone pole flying through something like a vacuum.

Light flows: If the flow density < ambient density then the flow will interact strongly with the dense medium through which it pushes. Even so, there is no apriori expectation that a dense, thin knot will develop almost immediately along the y axis. But it does: that’s what I find in the sims using the present version of AstroBEAR. See the attached figure where I move the viewing window at the same speed as the head of the flow.

Notes:

1. the spatial units in the graph should be multiplied by two if the basic cell size = 500 AU. I had to mess with the scaling factors in VisIt (0.25 instead of 0.5) to get a good display. That is, the basic cell in the figure wiull have dimensions of 250 AU.
2. I used Nlevel=5 in these sims. Changing it by ± 1 has no effect.

The panels show a light flow of density 10^{2 and speed 200 km/s moving into a uniform ambient medium of density 10}4. The bottom panel shows the geometry at t=0. You are looking at the nozzle (round) and (unit) flow vectors that will emrge through its surface at t=0+. The vectors are perfectly vertical. The nozzle’s surface isn’t a perfect shpere, but that doesn’t matter much.

The vectors along the inside edges of the gas displaced by the round jet (the “swept-up, compressed rim”) almost immediately start to curve towards the y axis. This is exactly what should happen when the flow strikes the inner edge of the rim of displaced gas obliquely. The flow along the rim starts to converge towards the y axis. This convergence forms an incipient knot in 100 y (the nozzle crossing time). The knot rapidly becomes longer and denser as mass continues to continues to flow into it.

It’s what you will get if you put a squishy ball bearing between the jaws of a closing scissors.

My point is that artifical knots are inevitable using spherical nozzles. The formation of this axial knot can be suppressed if the nozzle were a flat surface or a long and thin cone. The flow from a flat nozzle would displace and sweep up a flat plug (a disc) whose speed decreases as ambient gas is incorporated into it. The only wat to completely avoid any axial knot is to introduce a flow with a sharply conical head, like the nose cone of a rocket.

Baowei knows through bitter experience that forming a flat nozzle is difficult in AstroBEAR. It’s even more difficult to make a nozzle shaped like a nose cone. But you might think about it. (Of course, some axial knots might form after a simple jet starts to break up or become unstable and pinch. Such knots are ‘real’, not artificial.)

Of course, no one has any idea what a nozzle looks like on large sice scales. Zhou’s sims might provide some guidance on this. They look highly conical to me.

This email sounds like its just about details of flow geometries. It’s really more about model outcomes. There’s potentially important science at stake.

Referee report, ​WASP-12 paper

Moderate revision - only one real concern, that being the effect of stellar rotation on the stellar wind-planet wind interaction. Eric's comment that "If the result depended on the distribution of mass around the planet, then the point of the referee would lead to very important differences" is accurate, I think, and I do believe that's what we were attempting to get across in that section - so we could either revise that section for clarity, or we could run a new simulation (already set up) where the stellar rotation rate isn't locked to the orbital rotation rate.

The other concern was Fig 2, the 3D rendering. I feel it adds something, but I'm not sure I can make a strong argument for its continued inclusion.

Simulation Status

I plan to prioritize finishing rot_Run2 over the weekend (Jonathan - change frame numbers before restart?).

Now that I've looked at a few frames of Run3, the problem is apparent. The ambient medium is too dense, so no flux is making it to the planet. I believe there's a quick fix for this (extend the planet profile a bit farther, so that we get closer to zero density at the pressure-matching condition), although it may require upping the resolution at the outer boundary as well.

I do have movies for rot_run1 and noRot_run1/run2:

Rot, Run 1 (side view)

Rot, Run 1 (top view)

noRot, Run 1 (side view)

noRot, Run 2 (side view)

Rad transfer module test

Re-running the ionization-only, radiation pressure-only, and combined tests, they look good. The secondary effects in the ionization-only test are pretty strong still, but in general it appears to act like we expect. In the combined test, the effects of ionization dominate radiation pressure at the current fluxes (2x10¹³ ionizing flux, 5x10¹³ Lyman-alpha flux), but you can still see the effects of radiation pressure (particularly after the clump has been blown away).

Fellowships

NESSF - not sure about my chances, but worth applying? Mostly not sure how strong of a proposal I can write.

Simulation Status

Still need to look at 3&4 - I suspect the excessive cooling is related to the higher density in those runs, though.

Rad transfer module test

I've fixed the bug causing the 1D ionization test to fail:

However, the bug was such that it wouldn't have affected the clump test - and in fact the tests are identical pre- and post-fix:

New Work

1. New run (~half completed) with evolving softening length and resolution, larger box, and subgrid accretion turned off
2. Enabled parallel HDF5 output and optimized code to get a bit of speed-up
3. Started to analyze accretion onto the particles.

Results

1. With accretion turned off, the accretion disk (torus) morphology is no longer present.
2. It currently takes about 7000 node-hours on stampede2 skx nodes to complete 43 frames (10 days). But the speed actually varies because the radius of the refinement region as well as highest resolution are made to vary during the run. For comparison, our allocation for 2018 is for 166,000 node-hours. So if each run lasts for ~130 days, we are talking about no more than ~2 runs.
3. The accretion rate increases montonically with radius of the control region, as might be expected. It approaches a state where Mdot oscillates around 0 for both the RG core and companion (the mass aroudn each particle stabilizes) BUT when the softening length is reduced, this seems to (artificially) cause more accretion to happen.

Run 143:

• Similar to run 132 but with subgrid accretion turned off, twice larger box, resolution and softening length that evolve with time, less aggressive refinement, i.e. larger max refinement zones, and no relaxation (damping) run initially
  

Relaxation run: no relaxation run
First frame: 0
Last frame: 117 (so far)
Total simulation time: 27 sim-days
Machine and partition: Stampde 2 normal (completed up to frame 117, or 27 sim-days)
Number of nodes: 128 or 64, each with 68 (standard nodes) or 96 (skx nodes) cores
Total wall time: TBD
Hydro BCs: extrapolated
Poisson BCs: multipole expansion
Box size: L=8e13 cm (1150 Rsun)
Max refinement level: 4 (frames 0 to 72), then 5 (frames 72 to 117) Base resolution: 2.25 Rsun (256³ cells)
Highest resolution: 0.14 Rsun (4096³ cells, 4 levels AMR, frames 0 to 72), and 0.07 Rsun (8192³ cells, 5 levels AMR, frames 72 to 117)
AMR implementation: set by hand to have max level around RG core (frames 0 to 72) or companion (frames 72 to 117)
Max resolution zone: sphere around RG core, radius 5d12 cm (frames 0 to 46), radius 4d12 cm (frames 46 to 72), radius 3d12 (frames 72 to 103), radius 2d12 (frames 103 to 117)
Buffer zones: 16 cells
Softening length: 2.4 Rsun (frames 0 to 72), 1.2 Rsun (frames 72 to 117)
Ambient density: 6.7e-9 g/cc
Ambient pressure: 10⁵ dyne/cm²
DefaultAccretionRoutine=0 (no accretion)

Movies of run 143
1) Edge-on slice, centered at companion. Left: slice through both particles. Right: slice through P2 as viewed from direction of P1.
​edge-on
2) Face-on slice, centered at companion.
​face-on slice
3) Zoom-in face-on slice, inertial frame.
​face-on slice
4) Full box face-on slice, inertial frame.
​face-on slice, full box

Snapshots of run 143
Left: edge-on, slice through both particles (view centered on companion with side of 4e12 cm). Right: same but zoomed in by 4x.

Inter-particle separation
Red horizontal lines show radius of maxLevel refinement region. Green horizontal line shows softening length (solid) or five times softening length (dashed) for each particle. Finest resolution is proportional to softening length, i.e. resolution improves by factor of 2 when softening length is reduced by factor of 2. Vertical lines show transition from one refinement radius to the next (dotted red, 5e12 to 4e12 to 3e12 to 2.5e12 cm), and transition from one softening radius to the next (green dashed, ~2.4 Rsun to ~1.2 Rsun)

Mass accretion
Accreted mass with time (blue) and accretion rate with time (red) inside spheres around RG core (top) and companion (bottom) of different control radii. Plot of orbital separation is also shown for comparison. The bottom two plots are the same as the upper two plots except that data from an additional control radius is added, and individual points are not overplotted on the curve.


AA

Comments

• Accretion onto companion takes place within ~10 days and then the mass becomes fairly steady.
• Which control radius is most appropriate? The total mass accreted seems to be of order 0.01 Msun or 1% of the companion mass.
• There is an anticorrelation between the mass inside the sphere and the particle separation, which makes sense.
• There is a sudden increase in the accreted mass when the softening length is reduced by a factor of 2, especially for the RG core. This makes sense because the particle's gravity becomes stronger inside of the original softening radius. It's not clear how much of a role this is playing, and it will be interesting to see what happens when we halve the softening length again later in the run.
• Accretion rates onto the companion are of order 0.1 to 1 Msun/yr, which is super Edington, though not as high as what we were getting with Krumholz accretion turned on.

Next steps

• Complete the simulation up to at least 300 frames (about 70 days), using stampede2 skx nodes
• Analysis of:
  • mass accretion with time
  • angular momentum accretion with time
• Clean up presentation of movies/snapshots
  • reference frames/views
  • velocity vectors in particle frame
  • units
  • labels

To think about

• Had we been able to use an arbitrarily small softening radius from t=0, would the steady-state accreted mass be higher? If yes, then would the amount of gas left over in the envelope have been significantly affected? Would the orbit of the particles have been significantly affected?
• As a reminder, the strategy is to keep the number of resolution cells per softening length constant (~17), and to keep the number of softening lengths spanning the separation between the two particles at >5 (as done by Ohlmann 2016).