Decomposing gradients into wave strengths
Consider the 2D hydro linearization of the wave equation in each dimension
where
The matrix
has eigenvalues , , , corresponding to left propagating sound wave, contact, shear discontinuity, and right moving sound wave. And we can construct the matrix of left eigenvalues (where each row is a left eigenvector)and the matrix of right eigenvectors
(where each column is a right eigenvector)
and the diagonal matrix of eigenvalues
where each diagonal term is an eigenvalueThe eigenvectors are orthogonal
and we also have
which combining gives
Now since the eigenvectors are orthonormal, we can multiply any left eigenvector by an arbitrary normalization constant, provided we divide the corresponding right eigenvector by the same amount.
Rewriting the equation
Now if we renormalize the left eigenvectors
and the right eigenvectors
then
where
has units of
This gives us the relative strength of each wave at each cell. This looks at gradients in the fluid variables and how they correlate with the various types of waves.
For example, if we look at the Sod Shock Tube test problem whose solution entails a left rarefaction, contact, and right moving shock, we see that the 3 different waves are correlated with the amplitude of the projections of the gradients onto the left eigenvectors associated with each wave.
And if we look at the Brio Wu shock tube which consists of a fast rarefaction, a slow compound wave, a contact, a slow shock and a fast rarefaction, we see the following wave strengths
meeting update
Heavily updated the M2FR and the protocluster formation sections of the paper. Also cleaned up the intro. Am now working the rest of the way through the paper.
HD209458b: PlanetaryWind and Sonic Surface
When trying to produce a high-res picture of HD planetary wind with larger planet, I found the different-looking of sonic surface for L=60 and L40 as shown in blog:bliu10192015 .
I tended to think the L=60 is more correct(?) as
- The SF of L=60 is more like the 2D PW simulations I did before (comparing with 2D sonic surface ).
- In early frames of L=40 has similar sonic surface as shown in this picture
Thought this is resolution related. I did several restarts with higher AMR levels. Here's some of the results
L=40 4-AMR outside the planet | |
L=40 5-AMR outside the planet |
Will try higher-res for L=60 to see if the results could be consistent..
Higher-res results for L=60 shows similar sonic surface as those of L=40. So we conclude the sonic surface of higher resolution is more correct although it looks slightly different from our 2.5D results
L=40 4-AMR outside the planet | movie |
Update 10/19/2015 - Eddie
Below is image/movie showing the hydro 3-D pulsed jet: slice of x-z plane, slice of y-z plane, and column density.
It might be helpful to review how the outflow object is set up: OutflowObjects
My jet uses a radius = Rjet/2, thickness = Rjet/4, base = Rjet, and open_angle = 0. The outflow object steps on some cells, but the rest of the lower boundary is extrapolated; should probably change this to reflecting.
Things To Try
- use different Riemann solver
- adjust CFL
- adjust boundary conditions of outflow object (radius, thickness, fade, etc.)
- turn off pulsations
- have precession fade from 0 to Rjet
UPDATE
Things i have tried thus far:
- For some reason, the exact Riemann solver does not work well. Code stops early due to nan in flux.
- Added velocity fade from 0 to Rjet. This is not a precession fade, just a velocity magnitude fade. Still see waves.
- Turned off characteristic limiters, use primitive limiters instead. Still see waves.
- Turned on ApplyDiffusion.
- Lowered target CFL from 0.3 to 0.1.
- Turned off pulsations.
M2-9: 3D pnStudy with rotating/spinning conical wind
- Test the 3D pnStudy module by adding a rotating velocity to the conical wind:
Non-rotate 4AMR non Rotate movie; 4AMR non Rotate 2D slice movie; Rotate with 150y period 4AMR Rotate movie; 4AMR Rotate 2D slice movie;
- The data file for this run:
tamb = 1d3 ! ambient temp, 1cu = 0.1K (100K=1000cu) namb = 4e4 ! ambient central density cm^-3. Usually 400 for 1/r^2 or torus. stratified = f ! true = add a 1/r^2 background 'AGB stellar wind' torus = f ! true - add torus to the background torusalpha = 0.7 ! alpha and beta specify the geometry torusbeta = 10d0 ! see Frank & Mellema, 1994ApJ...430..800F rings = f ! true - add radial density modulations to AGB wind ! ! FLOW DESCRIPTION SECTION, values apply at origin at t=0 outflowType = 2 ! TYPE OF FLOW 1 cyl jet, 2 conical wind, 3 is clump njet = 4d2 ! flow density at launch zone, 1cu = 1cm^-3 Rjet = 2d0 ! flow radius at launch zone, 1cu = 500AU (outflowType=1 only) vjet = 2e7 ! flow velocity , 1cu = cm/s (100km/s=1e7cu) tjet = 1d3 ! flow temp, 1cu = 0.1K (100K=1000cu) tt = 0.0d0 ! flow accel time, 1cu = 8250y (0.02 = 165y) open_angle = 90d0 ! conical flow open angle (deg) tf = 15d0 ! conical flow Gaussian taper (deg) for njet and vjet; 0= disable sigma = 0d0 ! !toroidal.magnetic.energy / kinetic.energy, example 0.6
HD29:planetary wind and Mass loss
1. Planetary Wind With Temperature Profile
I have been trying to get a high-res picture for planetary wind HD29 with temperature profile. This picture shows there's a shock from the night-side and there's a back-flow.
frame 40 | |
frame 48 |
The shock disappeared as time goes further, although still need to double-check if this comes form the restarting problem
frame 52 (2AMR) | movie for this run | |
frame 76 (3AMR) | movie for this run |
2. Estimate of the Mass Loss Rate
- 1). The mass loss rate is estimated by calculating the flux out of the cubic surfaces as show below:
Where the size of box is 0.4*Orbital separation between the star and the planet, or 0.019 AU and the planet is sitting in the center of the box.
- 2). The mass loss rate for no-temperature-profile case is 3.82E+09 gram per sec, see blog:bliu09222015
- 3). For with-temperature-profile case, I checked the mass loss rate for different frames
frame | Picture | Mass loss rate |
Frame 48 | 1.33E+09 gram per sec | |
Frame 70 3AMR | 1.55E+09 gram per sec | |
Frame 76 3AMR | 1.56E+09 gram per sec |
Mach stem exploration
Running Eddie's setup,
time = 39.7 tau; M = 30.02; d = distance between clumps centers;
Summary: Simulation Extrapolation/Interpolation
d [R_clump]\gamma | 5/3 | 1.4 | 1.2 |
---|---|---|---|
3 | TMS, SB | TMS, SB | TMS, SB |
3.25 | TMS, SB | TMS, SB | SMS |
3.5 | TMS, SB | TMS, SB | RR |
4 | TMS, SB | TMS, SB | RR |
4.5 | TMS, SB | TMS, SB | RR |
5 | TMS, SB | SMS | RR |
5.5 | TMS, SB | RR | RR |
6 | TMS, SB | RR | RR |
6.5 | TMS, SB | RR | RR |
7 | SMS | RR | RR |
7.5 | SMS | RR | RR |
8 | SMS | RR | RR |
8.5 | SMS | RR | RR |
9 | RR | RR | RR |
Snapshots:
gamma = 1.4; d = 6, 5, 4 clump radii:
Update 10/13/2015 - Eddie
2D Mach stem results with low M
low M = 1.55
TMS, SB = Transient Mach Stem, develops into a Single Bow shock
SMS = Stable Mach Stem
RR = Regular Reflection
gamma | d | Result |
---|---|---|
5/3 | 24 | RR |
5/3 | 23 | SMS |
5/3 | 22 | SMS |
5/3 | 21 | SMS |
5/3 | 18 | SMS |
5/3 | 12 | SMS |
5/3 | 11 | SMS |
5/3 | 10 | SMS |
5/3 | 9 | TMS, SB |
1.4 | 21 | RR |
1.4 | 20 | SMS |
1.4 | 19 | SMS |
1.4 | 18 | SMS |
1.4 | 17 | SMS |
1.4 | 15 | SMS |
1.4 | 12 | SMS |
1.4 | 11 | SMS |
1.4 | 10 | SMS |
1.4 | 9 | SMS |
1.4 | 8 | TMS, SB |
1.4 | 6.5 | TMS, SB |
1.4 | 6 | TMS, SB |
1.4 | 5 | TMS, SB |
1.4 | 4.5 | TMS, SB |
1,2 | 18 | RR |
1,2 | 17 | SMS |
1,2 | 16 | SMS |
1,2 | 15 | SMS |
1,2 | 13 | SMS |
1.2 | 12 | SMS |
1.2 | 9 | SMS |
1.2 | 8 | SMS |
1.2 | 7 | TMS, SB |
1.2 | 6 | TMS, SB |
1.2 | 5 | TMS, SB |
1.2 | 4 | TMS, SB |
3D Pulsed Jets
Ran the MHD, beta = 1 model farther at 32 cells/Rjet.
Started a production-quality hydro run at 64 cells/Rjet.
3D Clump/Bow Interaction Paper
- Having some trouble with ds9 and the .fits images. I just want images that look like the ones from Pat's 2011 paper, but I can't seem to get the scaling quite right.
Meeting update
- Continue to work on implementing ThermalConduction solvers.
- Investigated Newton's Method for solving non-linear PDE's
- Currently, Anisotropic & Isotropic, Implicit, & Explicit Thermal diffusion seems to work on Guassian test. johannjc10092015
- Works in parallel and in AMR, with Hydro turned off. Still need to
- Resolve issue with ghost zoning
- Implement flux limiting
- Run various tests - MTI, Ablative RT, ???
- Plugging away at MHD Wire Turbulence sims using ½ of BGQ. Getting a 26GB frame every 8 hours…
- Ran HD… sim out longer - still hasn't reached a steady state…
On the research
This model is still the dense shell + constant stellar wind model.
AGB mass: 2 solar mass with alpha=0.134
Secondary mass: 0.5 solar mass
Mass loss rate in the burst is: 9.3E-6 solar mass / year last for 11.6 years
Mass loss rate of stellar wind is: 1.65E-7 solar mass / year
The speed of burst is between the escape velocity of the AGB star and the binary system. The speed of stellar wind is above the escape velocity of the binary system.
Simulation is in isothermal state at 400 K.
0.65 micron and 0.55 micron respectively. They looks identical…
SED looks more reasonable. The red line is the SED I got in current model. Note the flat spectrum energy signature from 1 micron to 4 micron. This meets the observation in Kervella et al (2015). Also, the 10 micron bump is there.
Dust composition is 0.8 % of olivine and 0.2 % of pyroxene.
Anisotropic Thermal Diffusion
Working on implementing Anisotropic Thermal Diffusion
The 4 columns correspond to Isotropic diffusion, Anisotropic but with
, parallel diffusion, and perpendicular diffusion. And the 3 rows correspond to Backward Euler, Crank Nicholson, and Forward Euler (explicit).The magnetic field is uniform in the +x+y direction and In each case the effective diffusion constant is 1,
, each box is 2 units across, and the total time = 1. Hydro is turned off and these are run using a single core in fixed grid.The pseudocolor is showing Temperature and the initial circle has a temperature ½ that of the surrounding.
movie |