Posts for the month of April 2012

Meeting Update

  • Worked on cleaning up threading algorithm and improved the load balancing algorithm. Early results are quite encouraging. See Ticket #191
  • Ran 5123 turbulence sims on Kraken using IICool and Isothermal. Velocity spectra seem very similar - but still analyzing other statistics…
  • Working on mapping data from 5123 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

  1. Add a new global Flag lUseOriginalNewSubGrids to choose the old/new Subgrids-generating algorithm. testing and checking in the code
  1. 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, rBondi=6, velw=10km/s, tw=1000K, with 96x96x144 + 2amr res, 4hrs, 40 afrank procs, bhive

http://www.pas.rochester.edu/~martinhe/2011/binary/8may12.png

27 Apr '12

  • 40AU, tw=3000K, with 1683 + 4amr res is running at ranger's long (1024p, 2days) queue.
  • 10AU, tw=3000K, with 112x144x160 + 3amr res. No disk formed after .21 orbits.
Orbital plane Longitudinal plane
http://www.pas.rochester.edu/~martinhe/2011/binary/corot-10au-sol2-3000K-3amr.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-10au-sol2-3000K-3amr-long.png

26 Apr '12

  • 10AU, tw=3000K, with ~643 + 2amr. No disk formed after 2 orbits (16hrs, 256 ranger procs). Trying 2 times more res.
  • 40AU, tw=3000K, with ~1283 + 2amr. No disk formed after .44 orbits (1day 256 ranger procs). Images show a zoom in to the particle.
Orbital plane Longitudinal plane
http://www.pas.rochester.edu/~martinhe/2011/binary/40au-sol2-3000K-2amr-chombo146.png http://www.pas.rochester.edu/~martinhe/2011/binary/40au-sol2-3000K-2amr-chombo146B.png

Though I was expecting differences between the and cases, I'm surprised that they are SO different. This is inconsistent with http://adsabs.harvard.edu/abs/2000MNRAS.316..906M. One has to be careful though for their parameters are in general different than ours.

24 Apr '12. Research meeting. To try:

  • 10AU and 40AU low res and twind=3000K, instead of 300K. Tests ran
  • 40AU, twind=300K, ideal gas solver, . short test ran
  • Jonathan's latest thread implementation on the global (2 particle) sims.

23 Apr '12

  • separation = 40AU
  • mprim= 1.5 Mo
  • msec = 1 Mo
  • tempw = 300 K
  • densw = 1.5x1011 gr/cc
  • vw = 15 km/s = 15 Mach
  • vorbit-sec = 4.5 km/s
  • rBondi = Gmsec/ ( vw2 + vorbit-sec2 )= 3.6 AU
  • rsoft = 4 cells
  • rBondi/rsoft = 19.3
  • rorbit/rBondi = 6.62

A) low res=56x72x80 + 3 particle Ref. .43 orbits / 26 hrs with ~ 32 afrank nodes in bhive. Running.

Movies:

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). http://www.pas.rochester.edu/~martinhe/2011/binary/question0000.png

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 1Mo 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.

http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b1.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b2.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b3.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b4.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b5.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b6.png
http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b7.png http://www.pas.rochester.edu/~martinhe/2011/binary/corot-40au-sol2b8.png
Gas momenta. I'm using this plot trying to understand why gas is not going about the seconday as expected. http://www.pas.rochester.edu/~martinhe/2011/binary/idea1.png

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
  • mprim= 1.5 Mo
  • msec = 1 Mo
  • tempw = 300 K
  • densw = 1.5x1011 gr/cc
  • vw = 15 km/s = 15 Mach
  • vorbit-sec = 8.9 km/s
  • rBondi = Gmsec/vw2 = 4 AU
  • rsoft = 4 cells
  • rBondi/rsoft = 21
  • rorbit/rBondi = 1.52
  • running time ~ 1 orbit/11hrs (16afrank, bhive nodes)

Movies:

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
mprim [Mo] 1.1 1.5
msec [Mo] .19 1
separation [AU] 350 40
tempw [K] 0.01 300
densw [part cm-3] 1a 1.5x1011
vwb [km/s] .9 15
vorbit-sec [km/s] .54c 4.5
rBondi [AU] 200 4
rsoft [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/rBondi.
  • 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):

  • Tempwind=300K
  • AGB mass-loss = 10-5 Mo yr-1
  • velwind=9 km/s (as in Mastrodemos & Morris)
  • velwind/vel_orbitsec = 2.01
  • a=40AU
  • circular orbit
  • q=mpri/msec=1.5
  • rsoft=4 cells
  • 64x64x128 cells, fixed+++
  • grid: 40x40x80 AU.
  • trun~.5 orbit/hr in 16 debug bluehive procs
  • dxco-rotating ~ 13 dxlab 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 1283+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.

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

http://www.pas.rochester.edu/~erica/Bannerjee_vrad.png

http://www.pas.rochester.edu/~erica/densityBannerjee.png

Some Potential Hedla Movies

Toroidal Only, aligned with the shock sequence:

Click to enlarge.
http://www.pas.rochester.edu/~shuleli/HedlaPics/torxvrshow.png

Movie:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/torxvr.gif



Poloidal Only, aligned with the shock sequence:

Click to enlarge.
http://www.pas.rochester.edu/~shuleli/HedlaPics/polsequence.png
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) ~ k0, initially(left) and after 1 cloud crushing time(right).

http://www.pas.rochester.edu/~shuleli/HedlaPics/shockedsmallscale_0.png

We will also do a more realistic power spectrum with spectral index -5/3 (Kolmogorov type).



http://www.pas.rochester.edu/~shuleli/HedlaPics/poloidalcomponenteffect.png

Movie:
http://www.pas.rochester.edu/~shuleli/HedlaMovie/tp1vr.gif

http://www.pas.rochester.edu/~shuleli/HedlaPics/fieldscale2.png

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.

-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

Density lineout - matches Bannerjee

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:
http://www.pas.rochester.edu/~shuleli/HedlaPics/screencut.png

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.
http://www.pas.rochester.edu/~shuleli/HedlaPics/screencut2.png
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 ~ T2. 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:

movie

I therefore doubled the resolution to ~80 cells per cooling region. Here is the 800 cell simulation:

movie

I again doubled the resolution to ~160 cells per cooling region. Here is the 1600 cell simulation:

movie

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):

movie

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 of

This 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
2563 86410 4442 124 24 36934 Isothermal
1283 27475 12664 194 24 40893 Isothermal

IICool modications

Fabian suggested a modified cooling curve that gives equilibrium temperatures of 10 K at densities > 1000… Here are two plots showing various quantities using those two equilibrium curves


Current Cooling Curve


Modified Cooling Curve


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 to

So assume we have a system

where and are unknown constants

and we want to adjust so that or so that

we 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 acceleration

The 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 Update 04.03

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 have

where


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 30o, facing the primary. rdisk=rBondi/3 to make the test run faster. http://www.pas.rochester.edu/~martinhe/2011/binary/binaryNdisk1.png

-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.