CEJet

Working on branch ticket_442.

Questions to answer for jet project:

• How does the jet affect the morphology of the envelope?
• How does the jet affect the ejection of the envelope?
• How does the jet evolve, does it get quenched?
• What is the dependence of these questions on when the jet gets turned on?
• What is the dependence on accretion rate?
• Opening angle?

movies from test run 44

Starting with secondary right against the primary envelope. No secondary mass. Jet mass loss rate at 2 M_sun/year, and radial velocity 430 km/s (Keplerian).

​density-edge-on ​temperature-edge-on ​density-face-on

3 test runs planned for this week:

4 levels of AMR with 64 base grid resolution. Each run is planed for 10 days in simulation time, and takes about 2 days to complete on bluehive.

1. [Run 45] Secondary starts right against the primary envelope: a0 = 49 R_sun.

Use lower ambient density and pressure profile: rho_amb = 1e-10 g/cc, p_amb = 1961 dyn/cm^{2 .}

2. [Run 46] Secondary starts right against the primary envelope: a0 = 49 R_sun.

Use ambient profile as the CE fiducial run (run 143): rho amb = 6.67e-9 g/cc, p_amb = 1e5 dyn/cm²

This and the first run will compare the effect of ambient material over the jet. Also a comparison with the CE fiducial run.

3. [Run 47] Start with the secondary farther out: a0 = 72.5 R_sun, and use lower ambient profile.

There is also a CE run without jet to be compared.

Some parameters to explore

• Initial separation (a0) (ref. ​Shiber+2019)

49 R_sun, right against the primary envelope. Jet might be quenched after the secondary enters the envelope.

72.5 R_sun, start farther out, longer time for the jet to evolve (although farther means lower mass loss rate if we use accretion rate to determine the mass loss of the jet).

109 R_sun, Roche distance of our system.

• Jet mass loss rate

From CE paper 1, the secondary can have extremely high accretion rates of ~ 0.2 - 2 M_sun/year. This rate corresponds to 100 - 1000 times the Eddington rate for a MS star, or 1e4 - 1e5 times for a WD.

Currently we are using the upper limit (2 M_sun/year). For a more physical value, we plan to compare between 1 vs 100 times Eddington rate of a WD (100 times WD rate would be ~ Eddington rate for a MS).

Check out the Xsede proposal for more suggested tests.

• Jet radial velocity

​Federrath+2014 suggests using Keplerian velocity at the surface of the star.

In the current test runs, this is computed using initial mass of the secondary (0.978 M_sun), and 1 solar radius. So,

jet_vrad = sprt (G*M2/R_sun) = 430.75 km/s

• Jet radius, distance from the secondary particle within which the jet is initialized. 64 cells (70.4 R_sun) for now. Federrath+2014 also tested this effect, and suggests using 16 or 32 cells.

• Refinement radius (compare to the orbit). See CE paper 1 fig. 3, and paper 2 appendix fig. C1.

• Refinement shapes and qTolerance.

Radiation Pressure Paper

Updated version is on the arXiv today. I've emailed it out to a few people.

Charge Exchange

Waiting for job to start on BlueHive to get diagnostics.

AMR line transfer

Pushed two commits with separated local and MPI rays, and keeping an integer pointer to the next free ray and ray to be processed (instead of looping through entire array each time). Not tested yet, but weren't complicated changes, so (fingers crossed) they should work immediately. Now just need to compare timings to see if it was worth it for plane-parallel.

Also pushed commit with output for following ray propagation across the grid in time order (synchronized time, ray initial position, and processor ID). This may require some testing, particularly to determine the best way to process it for visualization.

Misc.

Set up ​GitHub repository of my Euler code. Working on personal website, as well.

Jet in x-direction

This is mostly to address a comment from Eric, if we put in an equatorial outflow, will it still be collimated?

Instead a "disk" outflow, I put a jet in the +- x-directions, with 15-degree open angle, otherwise the same as model A.

​http://www.pas.rochester.edu/~yzou5/PN_xjet_figures.pdf

On momentum (L/c)

In our simulation, the "momentum rate" (L/c) is set by the outflow's density and velocity

L_{mech}/c = 4\pi r_w^2 \rho v_w^3 / c

for high-momentum flux (models A and C), (L/c)_high = 3.835e+34 (g cm s^{-1 year}-1), and

for low-momentum flux (models B and D), (L/c)_high = 1.917e+33 (g cm s^{-1 year}-1)

At the last frame presented in PN paper, the total momentum ejected is then:

P_{sim} = t_{dur}\cdot (L/c)

Observationally, the upper limit on t comes from the size and jet speed: t<R/v, and maximum momentum is estimated as:

P_{measure}< (R/v_w)\cdot(L/c)

Here I use the extend of the top lobe for R, and wind speed v=300 km/s, the momentum calculated in those two ways and their ratio are:

Initial Mass Flow

This seems problematic. I will check my Visit code and come back with another update. Ignore the following for now.

Outer boundary flows inward, floor density raise remove some mass, but both are small.

total change: -0.2 Msun

outer BC: 1.1e-7 Msun

floor raise: -1.4e-3 Msun

Left: including frame 0; Right: excluding frame 0.

Radiation Pressure Paper

Submitted to MNRAS and the arXiv. Also looked at renormalized absorption:

Will update paper on Overleaf soon.

Charge exchange

Ionized and neutral planetary material appear to be interacting and becoming stellar material. Think the problem must be here:

if (q(iHhot) /= 0) then
  prop(i,j,k,lt_iHhot) = q(iHhot)/(q(iHcold) + q(iHhot))
else
  prop(i,j,k,lt_iHhot) = 0d0
end if

if (q(iHcold) /= 0) then
  prop(i,j,k,lt_iHcold) = q(iHcold)/(q(iHcold) + q(iHhot))
else
  prop(i,j,k,lt_iHcold) = 0d0
end if

if (q(iHIIhot) /= 0) then
  prop(i,j,k,lt_iHIIhot) = q(iHIIhot)/(q(iHIIcold) + q(iHIIhot))
else
  prop(i,j,k,lt_iHIIhot) = 0d0
end if

if (q(iHIIcold) /= 0) then
  prop(i,j,k,lt_iHIIcold) = q(iHIIcold)/(q(iHIIcold) + q(iHIIhot))
else
  prop(i,j,k,lt_iHIIcold) = 0d0
end if

data(:,:,:,lt_iHhot)=data(:,:,:,lt_iHhot)-rates(:,:,:,iIoniz)*dt*prop(:,:,:,lt_iHhot) + rates(:,:,:,iChargeExchange)*dt
data(:,:,:,lt_iHcold)=data(:,:,:,lt_iHcold)-rates(:,:,:,iIoniz)*dt*prop(:,:,:,lt_iHcold) - rates(:,:,:,iChargeExchange)*dt
data(:,:,:,lt_iHIIhot)=data(:,:,:,lt_iHIIhot)+rates(:,:,:,iIoniz)*dt*prop(:,:,:,lt_iHIIhot) - rates(:,:,:,iChargeExchange)*dt
data(:,:,:,lt_iHIIcold)=data(:,:,:,lt_iHIIcold)+rates(:,:,:,iIoniz)*dt*prop(:,:,:,lt_iHIIcold) + rates(:,:,:,iChargeExchange)*dt

rates(i,j,k,iChargeExchange) = chargeExchangeRate*(data(i,j,k,lt_iHIIhot)*data(i,j,k,lt_iHcold) - data(i,j,k,lt_iHhot)*data(i,j,k,lt_iHIIcold))

Charge exchange rate should be zero whenever there is no planetary or stellar material. Proportion of stellar material should be zero when there is no stellar material. So it should be impossible for q(iHIIhot/cold) to become nonzero when they're both zero at the start of the line transfer step. Clearly they are, however - will put some diagnostic variables in to see what exactly is nonzero when it shouldn't be.

Also need to check timing for charge exchange runs on non-rotating case of parameter space runs on BlueHive (will need to get the last frames).

AMR line transfer

Will try to run timing this week (after charge exchange diagnostics are complete), in addition to working on a new commit with optimizations.

Upcoming/Possible Projects

In approximate order of importance (though not necessarily priority).

• Literature review
• Spherical line transfer
• Metastable helium line
• Hot Neptune
• Flares

Fall 2019 registration

Last semester to take courses for credit. Would like to take CSC 482 (designing efficient algorithms), I think - will talk with Laura as well.

Update on plan

The OPAL EOS seems to go very deep into thermodynamics, which will take a long time before I figure every detail out. The (tentative) plan for precede this plan is list as follows

1. Run MESA and understand the output: To use OPAL EOS, we need mentality information, which is now lacking.

   Success so far but still working on the understanding part.
   Remark: MESA took (almost) all my computer's remaining memory away

2. Cross compare the OPAL EOS with MESA output by using \rho and T as the independent variables.
3. Get a new initial condition for our simulation
4. Put into ASTROBear

Meanwhile, I will try to understand the physics behind OPAL EOS.

Note update

Note > ​https://www.pas.rochester.edu/~yishengtu/research_files/CEE_gp_meeting_06032019/convection_in_CEE_06032019.pdf
Added section 3 on Statistical mechanics. In Roger's paper he mentioned that Grand Canonical ensemble is where the derivation start, so I'm trying to understand this part first