udpate 20160229
Erica: I tried putting all three CALL CreateSpectra function together but the data is not separating the three boxes. So I separated them and ran each one. Once I increased the GmX to (640,64,64) and the final time to 100, the code seems to run extremely slowly and it chocked at the 11th output. I'm not sure if I can revise the code or change a machine to resolve this problem.
Here are some of the histograms: The first one is more like a initial try, I used GmX = 160,16,16 and a very small final time.
The second one is with high resolution and large final time, and I plotted in the log scale.
Update 2/29
Read the Stone-Proga paper - a comparison of 2D simulations to spherically symmetric simulations run by others. They characterize escape from the planet with the hydrodynamic escape parameter, use no magnetism (hydrodynamic rather than MHD), and ignore heating and small-scale effects at the base of the wind.
In models with no stellar wind and an isothermal outflow, they find that the sonic surface of the wind is closer to the planet, with a slower radial velocity on the night side of the planet, and a very evident shock in the |v|/vs plot at various angles. For larger values of gamma, this shock produces a delta T, and the sonic surface is farther from the planet (but still nearer than in the spherically symmetric models). Introducing a stellar wind creates a back-swept profile, but has little effect on the sonic surfaces or mass-loss rate.
Should have more time to read Murray-Clay and a few other papers this week.
New sims for HD209458b 102% done
Here are the sims to frame 50 - though they are completed to frame 100
Neutral hydrogen |
Density |
volume rendering of log(rho) |
Plan is to follow two more orbits with maximally refined around planet out to 7 planet radii
t=10 orbits | t=10.1 orbits | t=10.2 orbits |
Also, new runs fix issues with odd flow patterns near star
link to krumholz paper
Equations 18-22 and surrounding text discuss various opacities…
Update 02/18/16 - Eddie
- Monday, I was busy gathering/generating some old data for Pat for the Mach stem paper. Revisions according to referee report should have been submitted yesterday.
- Been putting my thesis together and writing the cooling paper which will also be a chapter in my thesis.
- Running the 3-D pulsed jets on davinci. Might move to bluehive since the queue on davinci is taking forever. Below is the beta = 1 run. I changed the scaling on the emission to make it brighter. It doesn't look like there is enough resolution to separate the H-alpha and [S II] lines within a clump, but we do see H-alpha (green) leading the clumps (yellow).
Meeting Update
Last week:
- I spent the week so far debugging and testing the rad feedback module, see blog post. I began writing a page on this module as well, here (figured this was useful for new users, as well as my thesis — which will likely incorporate this as a section).
To still do this week:
- Email people back who have gotten in touch with me about paper
- Email people links to my paper (i.e. do some networking)
- Look into Marshak waves and how they pertain to my sink tests
- Address the referee report for the CF paper
- Set up a 2D collapse problem, motivated by different timescales for the radiation equations.
Things I am thinking of putting on the back burner for now:
- Getting Rad Feedback working and tested in 3D, as well as AMR
Things I want to begin working on next week
- Next week I am thinking of switching gears and beginning on outflow feedback. This entails reading the paper again, planning the code revisions, and then of course, coding them.
- Reading papers on radiation feedback.
- Applying for conference talks/posters.
Rad Feedback from Sink Particles Update
This last week I spent testing the feedback code. To do this, I started with a sink particle at the center of the grid and prescribed some energy to be injected into the grid each timestep. Here are some results.
Run 1.
Final time - 1e-4
Energy injection rate - 28052.48 (chosen so that by the end of the simulation a total of 2 Erad will be injected into the grid, where Erad is the total Erad of the simulation at t=0)
Opacities, both set to 1d-6.
This run produced the following graph of Erad vs t, verifying the code is behaving as expected:
That is, we recover the energy injection rate, as given by the 'predicted' curve.
Here is how the radiative, thermal, and total energy behave over time.
A pseudo-color plot of Erad over time seems a little odd though, and I am trying to understand the behavior here. As this movie shows, the min and max remain very close over the simulation (a query shows they only differ in their 7th decimal place), and that by frame 1 there is already a large sphere on the grid. The legend isn't fixed (if it were, would be hard to capture the sphere at all), and so the breathing may be an artifact of the plotting. Why is the min/max so close like this?
I then lowered the final simulation time while keeping the energy injection rate the same. Here is a psuedo plot movie of that, but I see the same behavior. Perhaps the rate is small compared to the diffusion rate so that it all smoothes out within a timestep? Need to look at the different rates in more detail.
In the meantime, I kept the final simulation time small but increased the energy injection rate to play with different opacities. These results I will go through next.
Run 2.
tfinal - 1e-9
energy injection rate - 5,610,800,000 (such that 4*Erad is injected over course of simulation, where Erad is initial total Erad in the grid)
opacities - both set to 1d-6
So as you can see, AstroBEAR is losing some of the energy (supposedly the rad transfer routines aren't conservative??) over the course of the simulation, but performing quite well still.
For this case, the pseudo plot looks quite nice. It is shown in my previous blog post under 'Run 1'. I made a movie that shows how the q variables behave over time for a radial cut of the data. As the pseudo plot suggests, seeing little coupling between the gas and the radiation. For the next run, I will increase the roseland opacity to try and confine the radiation to a small region around the sink particle, but not change the planck opacity so that the gas will remain uneffected by the radiation.
Run 3.
tfinal - 1e-9
energy injection rate - 5,610,800,000
opacities - roseland=1d-3, planck=1d-6
The q-variables are shown in this movie here. For kicks, here are some images of the kernel:
All looks as expected.
Run 4.
Now am trying to couple the gas to the radiation, so boosted the Planck opacity. I found it had to be increased many orders above the Roseland so that it would couple.
tfinal - 1e-9
energy injection rate - 5,610,800,000
opacities - roseland=1d-3, planck=1000
Here is a movie of the q-variables again. And a movie of just Erad here. The pressure, temperature, and total energy grow at the same rate, but density remains flat. I am curious why density is not being affected at all.
mean free paths
cross section for hydrogen hydrogen collisions can be estimated as
where is approximately the bohr radius… This gives
The MFP is then
The material leaving the planet is at
so the MFP is then 1e-7 AU << dxUpdate 2/15
Ran OutflowWind? simulations with varying parameters. Thickness and wind velocity don't appear to have any effect - perhaps an effect of rescaling? Velocity speeds up the simulation as a whole. Radius increases size of wind, but otherwise has little effect, at least at small values. Increasing density appears to decrease relative density at the center of the front of the outflow. Increasing temperature creates a more diffuse cloud at the end. See page for comparison .gifs.
Plan to read some papers this week.
Updates on rad feedback
Run 1.
2D, periodic boundary conditions, 1282 resolution. There is a sink particle in the middle of the box that is injecting a constant 'accretion energy' into the center zones of the grid.
Plot is of Etherm (left) and Erad (right) over time,
Thermal energy is derived from subtracting the specific kinetic energy from the total energy (as there are no magnetic fields or gravity in this simulation). The total radiative energy in the box at the start of the sim is, Erad=90. By the end, Erad=10,038.
Both the planck opacity and the roseland are set to =1d-6. The right plot shows that the radiative energy diffuses through the grid much faster than it is exchanging energy with the gas (i.e. changes in thermal energy happen slower). Hardly any change in thermal energy is happening at all, in fact. Unless I made an error with my expression, not sure how the radiative energy can be increasing, yet the thermal and total energy are equal..
There is a slight 'blip' noticeable between the 1st and 3rd frame. Not immediately sure what this is being caused by either..
I would like to check that the total energy is being conserved. To do this, should I make an expression that sums E and Erad and check that it is constant over the course of the sim?
Overall though, the run seems to proceed nicely for these parameters (need to check my thermal energy expression still). I am trying to set up a simulation where I control the rate of energy injection, such that the total Erad doubles over the course of the simulation. In the meantime, experimented with differing opacities, as will show next.
Run 2.
High roseland opacity (radiation trapped, k=.001), low planck opacity (no coupling of radiation and gas, k=.1d-6).
Note, no blip in this run. I see that the thermal energy is not changing much here either. The max Erad is much greater than compared to above. The radiation being confined to a smaller space explains this.
Run 3.
High roseland opacity (radiation trapped, k=.001), High planck opacity (coupling of radiation and gas, k=10).
Now am starting to see an increase in thermal energy. Need to look at the equations and see how the opacities are coming in and how they are related to each other.
Update 2/8
Not much progress since last week, so plans for this week: Read at least two of my growing library of papers. Experiment with parameter space of OutflowWind code.
20160208
- Velocity Log Spectrum
- Attached are an old plot from last semester with low resolution, and a new one with GmX = (640, 64, 64) resolution.
- My intuition is that the curve should smooth out when resolution is higher but the pattern is not quite obvious.
- Another curiosity is why the curve for Spectra_00100 is smoother than Spectra_00200. Is it because as time progresses, the data has more details thus more bumps?
- Histogram of density
- Still working on this. There is a segmentation fault in the code after I added the histogram function, and Baowei is helping me debugging it. Error message shows as :
forrtl: severe (408): fort: (2): Subscript #1 of the array Q has value 6 which is greater than the upper bound of 5
- I also adjusted "spectra.f90" in astrobear/src/processing, below is the code:
CASE(COSINE_WINDOW)
DO i=layout%mB(MPI_ID,1,1), layout%MB(MPI_ID,1,2)
x=(REAL(i,KIND=qPREC)-origin(1))/radius(1) DO j=layout%MB(MPI_ID,2,1), layout%MB(MPI_ID,2,2)
y=(REAL(j,KIND=qPREC)-origin(2))/radius(2) DO k=layout%MB(MPI_ID,3,1), layout%MB(MPI_ID,3,2)
z=(REAL(k,KIND=qPREC)-origin(3))/radius(3) r=x IF (r > 1d0) THEN
weight=0d0
ELSE
weight=cos(half*Pi*r)
END IF data(i,j,k,=data(i,j,k,*weight
! Spectra%WindowFunction=COSINE_WINDOW
WindowFunction=COSINE_WINDOW
- A general inquiry for Jonathan:
If I want to divide the whole space into three sections to study the spectrum in each one, which part of the code should I adjust, is it line 296 in "spectra.f90" (r=x) again?
Update 02/08/16 - Eddie
- Clumpy paper is done. Time to submit?
- Writing cooling paper, and thesis.
- Running 3-D jet sims. The only thing I don't like so far is that the magnetized runs seem to have pressure (probably magnetic pressure) pushing radially outwards at the injection region. See images and movies below…
Make high-res 3D volume-rendering image using VisIT2.8.2 on Bluehive
Jonathan installed VisIT 2.8.2 on Bluehive with SLIVR (a volume rendering library supported by GPU, http://www.visitusers.org/index.php?title=Volume_Rendering#SLIVR ). This makes doing 3D volume rendering images a lot faster and fancier. You can use 2D transfer functions and you can manipulate the pictures in real time since it's accelerated by GPU. Here's a short introducing movie.
To set the limits (maximum and minimum value) of the variable and to change the opacity, it has to switch to 1D function, as shown in the following two images:
And Here are some AstroBEAR results
| |
Update 2/1
- Read Matsakos paper. They use static rather than adaptive mesh, with flows modelled isothermally. 16 models (4 parameters, 2 values each - high/low flux, high/low escape velocity, large/small orbital radius, and weak/strong magnetic field). Some showed tail, some showed accretion stream - in the end, organized into 4 categories, based on relative strengths of planetary wind, gravity, and magnetic field - 2 categories with tails, 2 with accretion at 90 degrees or directly under planet.
- Ran two more test runs of OutflowWind with longer simulation time and more frames.