Meeting Update
- Worked on cleaning up threading algorithm and improved the load balancing algorithm. Early results are quite encouraging. See Ticket #191
- Ran 512^{3} turbulence sims on Kraken using IICool and Isothermal. Velocity spectra seem very similar - but still analyzing other statistics…
- Working on mapping data from 512^{3} simulations onto gravitational collapse simulations. See blog:johannjc03062012b for more info..
- Managed to get code to compile and run with gfortran. Had to build hypre, openmpi, and pgplot with gfortran. They have been added to the modules but some of them are the new default on grass, alfalfa, clover… so check your modules and your paths. See ticket #195
- Updated CollidingFlows page
- Added CameraObjects to allow for projections from any vantage point.
- Added IO routines for layouts to make it easier to transfer data from one simulation to another.
Baowei's Meeting Update 04/24/12
- Add a new global Flag lUseOriginalNewSubGrids to choose the old/new Subgrids-generating algorithm. testing and checking in the code
- Work on running AstroBEAR on Ranger's normal queue:
https://clover.pas.rochester.edu/trac/astrobear/ticket/188
https://clover.pas.rochester.edu/trac/astrobear/ticket/197
Meeting Update 04/24/2012 - Eddie
I successfully added MHD to the 1D radiative shock simulations. The magnetic field is only in the y-direction (perpendicular to the fluid flow). This was initially trickier than I thought it would be. Having magnetic fields changes the shock jump equations and the ODE within the cooling region. I made a lot of changes to my problem.f90 so that the module as a whole could easily switch between hydro and MHD. As expected, the presence of the magnetic field lengthens the cooling region. Here is a movie that compares the temperature profile for hydro and MHD (both simulations use DM cooling). DM_compare_temp.gif
I'm still working out some kinks in using the next type of cooling (NEQ cooling). The tricky part here is now I have to create initialization profiles for hydrogen H and ionized hydrogen HII (I'm not even considering molecular hydrogen, or helium and its ionized species). The cooling routines calculate these ionization and recombination rates, and these can be used to update my initial quantities to create an initial profile. These rates have to be self-consistent with the cooling rate that NEQ cooling gives. So the problem essentially becomes solving two different ODEs simultaneously. I think I have this done correctly, and I am now just running into issues with the way the main cooling routine interfaces with NEQ cooling.
Disk formation in binaries using a co-rotating frame
The simulations in this blog are based on these ideas:
^{+}https://clover.pas.rochester.edu/trac/astrobear/blog/johannjc04032012
^{++}https://clover.pas.rochester.edu/trac/astrobear/blog/johannjc03312012
8May '12
- 30AU, r_{Bondi}=6, vel_{w}=10km/s, t_{w}=1000K, with 96x96x144 + 2amr res, 4hrs, 40 afrank procs, bhive
- 30AU, t_{w}=1000K,
- Movies of logarithmic density (grayscale)and linear velocity [Mach]:
- orbital plane http://www.pas.rochester.edu/~martinhe/2011/binary/30au-1000k.gif
- orbital plane, zoom in http://www.pas.rochester.edu/~martinhe/2011/binary/30au-1000k-z.gif
- longitudinal plane, zoom in http://www.pas.rochester.edu/~martinhe/2011/binary/30au-1000k-z-l.gif
with 100x100x160 + 4amr res. No disk formed after .07 orbits (3days, 1024 ranger procs).
- Movies of logarithmic density (grayscale)and linear velocity [Mach]:
27 Apr '12
- 40AU, t_{w}=3000K, with 168^{3} + 4amr res is running at ranger's long (1024p, 2days) queue.
- 10AU, t_{w}=3000K, with 112x144x160 + 3amr res. No disk formed after .21 orbits.
Orbital plane Longitudinal plane
- Understanding/solving the dynamics of the 40au,
- Here's a temperature map of the whole orbital plane http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1e.gif.
- Hot gas (10^{5} K) seem to discontinually emanate from the particle, specially early in the sim.
- Now running a short test, same parameters, with 10 times more frames.
- Temperature movie 2. 20 times faster frame rate than the temp movie above and zoomed in to the particle, http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1f.gif.
, t_{w}=300K sim (below).
- Here's a temperature map of the whole orbital plane http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1e.gif.
26 Apr '12
- 10AU, t_{w}=3000K, with ~64^{3} + 2amr. No disk formed after 2 orbits (16hrs, 256 ranger procs). Trying 2 times more res.
- 40AU, t_{w}=3000K, with ~128^{3} + 2amr. No disk formed after .44 orbits (1day 256 ranger procs). Images show a zoom in to the particle.
Orbital plane Longitudinal plane
- 40AU, t_{w}=300K, ideal gas solver,
- ~128^{3}+2amr, 0.2 orbits / 1day, 256 ranger procs.
- Movies of logarithmic density (grayscale)and linear velocity [Mach]:
- orbital plane http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1a.gif
- orbital plane, zoom in http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1b.gif
- longitudinal plane http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1c.gif
- longitudinal plane, zoom in http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2-300K-gamma1.1d.gif
- Jet? This corresponds to the last frame of this short test sim.
. Very interesting dynamics.
Though I was expecting differences between the http://adsabs.harvard.edu/abs/2000MNRAS.316..906M. One has to be careful though for their parameters are in general different than ours.
and cases, I'm surprised that they are SO different. This is inconsistent with24 Apr '12. Research meeting. To try:
- 10AU and 40AU low res and t_{wind}=3000K, instead of 300K. Tests ran
- 40AU, t_{wind}=300K, ideal gas solver, . short test ran
- Jonathan's latest thread implementation on the global (2 particle) sims.
23 Apr '12
- separation = 40AU
- m_{prim}= 1.5 M_{o}
- m_{sec} = 1 M_{o}
- temp_{w} = 300 K
- dens_{w} = 1.5x10^{11} gr/cc
- v_{w} = 15 km/s = 15 Mach
- v_{orbit-sec} = 4.5 km/s
- r_{Bondi} = Gm_{sec}/ ( v_{w}^{2} + v_{orbit-sec}^{2} )= 3.6 AU
- r_{soft} = 4 cells
- r_{Bondi}/r_{soft} = 19.3
- r_{orbit}/r_{Bondi} = 6.62
A) low res=56x72x80 + 3 particle Ref. .43 orbits / 26 hrs with ~ 32 afrank nodes in bhive. Running.
Movies:
- dens+vels, secondary zoom in, http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2z.gif
- dens+vels, secondary zoom in, longitudinal plane http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2l.gif
Tilt?. At some point during this simulation I see a rather brief (lasting for ~ 0.005 orbits) tilted disk-like structure similar to the ones I saw in the global (2 particles) sims (https://clover.pas.rochester.edu/trac/astrobear/blog/binaryBondi, top table). |
B) high res= 56x72x80 + 5 particle Ref. .05 orbits / 20 hrs with 256 ranger procs. Running. Note the images below show gas at a point at which the mass of the secondary is increasing to reach its final 1M_{o} value. This was dome for numerical stabilization and, given that we want to run the sim for several orbits, it should not affect the formation/evolution of the disk.
Gas momenta. I'm using this plot trying to understand why gas is not going about the seconday as expected. |
20 Apr '12
We've determined that Jonathan's module template does not have a frame/secondary which follow a circular orbit and has extreme parameters. Presumably, a disk is formed in such setup, with a separation as large as 350AU, because the unstable orbit is subject to a lot shear from the flow which might make rotation about the secondary easier.
The 1st, intermediate-res test with the full wind solution (wind+secondary's gravity) has completed for
- separation = 10AU
- m_{prim}= 1.5 M_{o}
- m_{sec} = 1 M_{o}
- temp_{w} = 300 K
- dens_{w} = 1.5x10^{11} gr/cc
- v_{w} = 15 km/s = 15 Mach
- v_{orbit-sec} = 8.9 km/s
- r_{Bondi} = Gm_{sec}/v_{w}^{2} = 4 AU
- r_{soft} = 4 cells
- r_{Bondi}/r_{soft} = 21
- r_{orbit}/r_{Bondi} = 1.52
- running time ~ 1 orbit/11hrs (16afrank, bhive nodes)
Movies:
- dens+vels, http://www.pas.rochester.edu/~martinhe/2011/binary/corot-10au1.gif
- dens+vels, secondary zoom in, http://www.pas.rochester.edu/~martinhe/2011/binary/corot-10au2.gif
- ram pressure, zoom in, http://www.pas.rochester.edu/~martinhe/2011/binary/corot-10au3.gif
Higher res test running now. Note that we've not used this (new) wind solution for the 40AU case.
19 Apr '12 Trying to understand why Jonathan's setup^{+} forms a disk while mine (below) doesn't. Here's a parameter comparison.
Jonathan | Martin | |
---|---|---|
m_{prim} [M_{o}] | 1.1 | 1.5 |
m_{sec} [M_{o}] | .19 | 1 |
separation [AU] | 350 | 40 |
temp_{w} [K] | 0.01 | 300 |
dens_{w} [part cm^{-3}] | 1^{a} | 1.5x10^{11} |
v_{w}^{b} [km/s] | .9 | 15 |
v_{orbit-sec} [km/s] | .54^{c} | 4.5 |
r_{Bondi} [AU] | 200 | 4 |
r_{soft} [cells] | 4 | 4 |
^{a}
.
^{b}
^{c}
.The combined mass and the separation in the run of the left column seem extreme.
16 Apr '12
Same parameters as the run below, but a slightly larger grid and 3 spatially fixed refinement (particle) levels, http://www.pas.rochester.edu/~martinhe/2011/binary/corot1.gif.
- Discussion and further analysis needed.
- I do not see a disk forming. Will try with 2 more refinements in order to have ~ 40cells/r_{Bondi}.
- The flow shows some strange pulsations not related to protections.
- The accretion tail is not a collimated one, like that in the low res run (below), and it broadens up.
6 Apr '12
I've got Jonathan's new co-rotating frame setup^{+,++} running for a binary system with the parameters below which are representative of my previous lab frame binary runs (see the table below):
- Temp_{wind}=300K
- AGB mass-loss = 10^{-5} M_{o} yr^{-1}
- vel_{wind}=9 km/s (as in Mastrodemos & Morris)
- vel_{wind}/vel_orbit_{sec} = 2.01
- a=40AU
- circular orbit
- q=m_{pri}/m_{sec}=1.5
- r_{soft}=4 cells
- 64x64x128 cells, fixed^{+++}
- grid: 40x40x80 AU.
- t_{run}~.5 orbit/hr in 16 debug bluehive procs
- dx_{co-rotating} ~ 13 dx_{lab frame}
Log(dens) and vel field movie. Pole-on (left) edge-on (right). The green diagonal line in the left panel show the direction of the plane corresponding to the right panel. A disk forms. The tail formed by accretion eventually gets inflow from -Y. The tail is bent then. The edge-on panel shows asymmetries. We are unresolved with respect to the global (lab frame) simulations to fairly compare/conclude anything about disk tilting at this point. |
http://www.pas.rochester.edu/~martinhe/2011/binary/coRot1.gif |
Here's a second try at this sim but with a larger grid. This last sim will run, any time now, at bhive with 512procs and 128^{3}+3amr. |
http://www.pas.rochester.edu/~martinhe/2011/binary/coRot2.gif |
Martin's update, 24 Apr '12
Running in Ranger's normal queue, 1024procs for 2 days.
Binary. Jonathan and I have completed the implementation of the full wind solution for the co-rotating frame. The one we used before was not including the secondary's gravity. See more updates at https://clover.pas.rochester.edu/trac/astrobear/blog/martinhe04232012
Magnetic tower.
- I did a little presentation in a workshop organized by Dan and Manoj last weekend. They are quite interested in our magnetic tower/cooling jet models to explain some of their recent observations.
- Off to Florida on Sunday for HEDLA12.
AGN jet truncation.
- I've registered to the x-ray conference (http://www.brera.inaf.it/xrayastronomy_anniversary/).
- Time to add magnetic fields to the AGN jet and think about production runs.
Meeting Check-in
Hi all, I am running into some resistance from astrobear in finishing up my simulation of the Bannerjee set up. I posted a ticket that describes the sort of trouble I running into. Well, the good news is the code does appear to be running the problem accurately. Up until frame 4/6, everything looks pretty spot on. The spike at the edge of the sphere in the radial velocity plot is likely attributable to poor resolution of the sphere's edge. I bet that spike would average out to 0 in a radially averaged plot. I would like to rerun the sim with ambient rho = sphere's outer edge rho, but am curious to find out what has been happening with my sim before progressing. Any thoughts, greatly appreciated. Thanks E
Some Potential Hedla Movies
Toroidal Only, aligned with the shock sequence:
Click to enlarge.
Movie:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/torxvr.gif
Poloidal Only, aligned with the shock sequence:
Click to enlarge.
During expansion, a shaft shaped "core" formed along the axis due to high field concentration (harder to deform). This core is then gradually deformed (slower than the cloud expansion) due to higher total pressure. At the same time, we see red dots of dense, cold material in the residue, they probably contain tangled magnetic field which make them harder to be further fragmented.
Volume Rendered Movie:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/polxvr.gif
We also study the shocked behavior when the clump contains small scale field (with a length scale much smaller than the clump radius). The following sequence shows a clump with contained magnetic field with flat magnetic power spectrum PB(k) ~ k^{0, initially(left) and after 1 cloud crushing time(right). }
We will also do a more realistic power spectrum with spectral index -5/3 (Kolmogorov type).
Movie:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/tp1vr.gif
Martin's update, 17 Apr '12
More info (images) soon.
-Running in ranger (teragrid)! Limited to 2hrs in 256 procs. though; only the development queue has the updated libraries necessary to run the code.
-Binary.
- Co-rotating frame, new sim with more resolution, https://clover.pas.rochester.edu/trac/astrobear/blog/binaryBondi
- Testing a new wind solution, one that includes the secondary's gravity.
-AGN jet truncation. The RG has crossed the AGN jet's path, http://www.pas.rochester.edu/~martinhe/2011/agn-6apr-dens1.gif
Meeting Update
- Created SpectraObject page as well as PFFT and LayoutObjects
- Worked on creating isotropic turbulence of a given mach number for initial conditions for cloud collappsing cloud sims. See this blog
- Modified IICool routine to use analytic formula that can be adjusted at run time. See this blog
- Started isotropic turbulence runs on Kraken
- Added shapes to various ProcessingObjects. See this blog
Meeting Update 04.10
Toroidal/Poloidal only runs finished.
Toroidal only, perpendicular to the shock:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/dm_z.gif
Toroidal only, aligned with the shock:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/torx.gif
Poloidal only, perpendicular with the shock:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/polp.gif
Poloidal only, aligned with the shock:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/polx.gif
Toroidal + Poloidal running, at frame 42:
Notice it's not symmetric not because noise but the field configuration does not have mirror symmetry when adding toroidal to poloidal components.
Bx, By, Bz amplitudes are plotted below.
I will meet with Ruka tomorrow to help him making the movies.
Poster Suggestions
Unfinished poster http://www.pas.rochester.edu/~blin/poster.pdf
Lack of jet w/o rings data. Perhaps just present poster with rings?
Meeting Update 04/10/2012 - Eddie
DM cooling added to Radiative Shock module
Due to the initial conditions of the problem, the immediate post-shock temperature is near the peak of the DM cooling curve. As a result, the simulations that use DM cooling cool much faster than the analytic simulations where the rate ~ T^{2}. I found that, on average, the DM rate is approximately 250 times greater than the analytic rate within the cooling region. So for the DM simulations, I scaled the length down by a factor of 250. This resulted in a cooling region that was about 10% of the problem domain as it was in the analytic case.
However, even with the same resolution analytic and adiabatic simulations(~40 cells across the cooling region), the DM simulation was not very steady. I ran the DM simulations 10x longer than the analytic and adiabatic simulations that we looked at last week. I also looked at the analytic simulation at this longer runtime, and it was still steady. The snapshots below are all at time 0.12 (the final runtime of the analytic and adiabatic runs). Here is the DM 400 cell simulation:
I therefore doubled the resolution to ~80 cells per cooling region. Here is the 800 cell simulation:
I again doubled the resolution to ~160 cells per cooling region. Here is the 1600 cell simulation:
It looks like I am again running into resolution issues. However, I think that the interpolation of the DM cooling curve might be to blame?
Update
Decreasing the ambient velocity (aka shock velocity) solved the problem. The ambient conditions used for analytic cooling only work for analytic cooling. They lead to an immediate post-shock temperature that is beyond the peak of the DM cooling curve. Therefore, with DM cooling, these radiative shocks would be unstable. Reducing the shock velocity reduces the post-shock temperature. I found that changing the length scale as previously mentioned, and changing the velocity from 100 km/s to 80 km/s keeps the post-shock temperature below the critical value for stability (which is approximately 2e5 K). Here is the simulation for DM cooling (400 cells):
Baowei's Meeting Update 04/10/12
Ticket 185: new grid-generating algorithm
I used test suites to test the performance of the new algorithm. The new algorithm outperforms the old one for most of modules. Results can be found at: https://clover.pas.rochester.edu/trac/astrobear/ticket/185.
Ticket 188: Install AstroBEAR on Ranger
Ranger has a standard environment which is too old for AstroBEAR code. I compiled AstroBEAR on ranger with newer testing environment but got trouble submitting jobs on it. Details can be found at: https://clover.pas.rochester.edu/trac/astrobear/ticket/188
The Afrank Queue on Bluehive
This might not be a big issue for the time being. When I ran a benchmark on bluehive, I found the queue afrank used Ethernet instead of Infiniband which was not the way it's supposed to be. So a lot of time was wasted on waiting for communications when running with multiple nodes. Russell is checking the issue.
Martin's update, 10 Apr '12
See updated in:
Binary https://clover.pas.rochester.edu/trac/astrobear/blog/binaryBondi
AGN jet truncation https://clover.pas.rochester.edu/trac/astrobear/blog/agnJetTruncation
PN wind paper. Waiting for reply from Raghvendra, Adam and Bruce t resubmit to MNRAS.
Histogram funniness
So the smooth colliding flows run showed a peak at .4 as well as 1 in the mass-weighted histograms of the mixing ratio. (Red line in plot). It turns out that the expression for the mixing ratio was to blame. It was set to
to avoid nan's in the region where neither tracer was present instead ofThis means however in regions near the edges of the cylinder, where both flows collide and mix with the ambient, the mixing ratio will be ~ .5 or so… I replaced the expression with
and the second hump went away (Green line)I also remade the histograms using only the data within the colliding region and got the same result regardless of the expression… (dark and light blue lines)
First run on Kraken
Kraken Simulations
res | wall time | advances | restarts | cores | cell updates/s/cpu | notes |
256^{3} | 86410 | 4442 | 124 | 24 | 36934 | Isothermal |
128^{3} | 27475 | 12664 | 194 | 24 | 40893 | Isothermal |
Achieving a target kinetic energy by varying the forcing strength...
Energy dissipates in a crossing time which is
which means that the specific kinetic energy dissipation rate is . If we have a given accleration field , it should inject specific kinetic energy at a rate proportional toSo assume we have a system
where
and are unknown constantsand we want to adjust
so that or so thatwe could start by assuming that
and and set and then run until we reach a steady state energy. At that point we can evaluate and then adjust the accelerationThe kinetic energy could potentially be widely varying on short time scales - which makes it tricky to determine a steady state energy. Additionally any change in the forcing will take of order a dynamical time for the system to relax… So each adjustment to the forcing should be followed by 1 dynamical time of allowing the system to relax - followed by a half dynamical time of averaging. The initial state should probably be allowed to evolve for 2 dynamical times before averaging.
Meeting Summary
Spent most of the past week sorting out sources of asymmetry in the code…
The sources of asymmetry could be broken down into 4 categories…
- Non-symmetric addition [ a+b-c-d /= (a-c)+(b-d) ]
- Diffusion algorithm's calculation of the divergence
- Characteristic limiters that do matrix multiplication from slowest to fastest wave speed
- Calculation of back forces on particles
- Upwind type calculations that expect quantities to be < or > 0…
- Riemann fans that have an even number of states (odd number of waves) such as HLLC, HLLD, ExactRS
- MHD eigenvectors which depend on the sign of Bx or the 'polarization angle' [atan(Bz/By)] which break down when Bx = 0 or By=Bz=0
- Bugs in the code
- Extrapolated MHD physical boundaries not updating cell centered b fields using aux
- Clump object's default behavior to persist in boundaries
- Wind object's failure to zero out other fields - miscellaneous tracers etc…
- Machine optimization
- Apparently optimization can lead to a breakdown in accuracy (perhaps because it reorders additions/subtractions to be faster?). Fortunately this behavior can be corrected by the addition of -fp-model precise to the compiler options.
Developed co-rotating reference frame source terms
or if we assume that
then we havewhere
movie |
movie |
Meeting update
I don't have that much to post currently as coursework is particularly heavy this week.. :/ I did confirm that periodic boundaries reduced the flow of ambient gas of my sphere set up. I am waiting on the output from my high res run as of now — running with the periodic boundaries and at the correct resolution. Will post the results soon.
Meeting Update 04/03/2012 - Baowei
1. Ticket #169 links to papers using Astrobear on the wiki
I created the wiki page at:
https://clover.pas.rochester.edu/trac/astrobear/wiki/AstroBearPublication
The ticket was closed but anyone find a new paper can send to me or edit the page directly.
2. Ticket # 185 Compare new grid-generating algorithm and the old one using AstroBEAR modules
Results can be found at the page: https://clover.pas.rochester.edu/trac/astrobear/ticket/185
The new algo works as expected and is faster than the general case of the old algo.
Martin's update, 2 Apr '12
I plan to show images on some of the projects during the meeting.
Binary.
-Embedded disks in an AGB wind. The 1st test has been running in bgene since Sunday morning. I do not see a strong tilt so far, but gas in the innermost orbits does show some tilt.
No. density iso-surfaces for the case with a disk tilted 30^{o}, facing the primary. r_{disk}=r_{Bondi}/3 to make the test run faster. |
-40au tests, see the table of 29th March at https://clover.pas.rochester.edu/trac/astrobear/blog/binaryBondi. Only run 7 is still running. The others were killed during a weekend bgene maintenance session. I've requested another reservation to continue working on these runs.
-Velocity field 2-plane evolution. I'll soon start doing the movie we discussed last friday.
CRL618
-See https://clover.pas.rochester.edu/trac/astrobear/blog/crl618Figures for new movies.
-Bin will present a poster about this in about a week time.
AGN jet truncation. See https://clover.pas.rochester.edu/trac/astrobear/blog/agnJetTruncation for new movie.
Magnetic tower paper submitted!, ApJ page charges? HEDLA talk?'
Teragrid. Baowei and I are doing the final few installation steps needed to run the code at ranger.
Meeting Update 04/03/2012 - Eddie
The Radiative Shock simulation works!
I've updated my wiki page with a description of the problem, some relevant equations, plots, and movies. Scroll to bottom of zcooling for details.