(star formation) added star formation setup and bug fix with Bonnor e…

[raomoomoo]

…bert sphere boundary condition

[jameswurster]

What is the necessity of these changes? The new setup you changed is now the same as SETUP=jet, with the exception that JET has MHD=yes. It can be dangerous to change the name of a pre-existing setup block since this may have unintended consequences for other users (likely not in this case since this was a generic block)

The original definition of cs_medium is correct. This was chosen to have continuous pressure (i.e. the pressure at the edge of the sphere is the same as the pressure in the medium).

The value of rhocrit0cgs needs to be between the edge_density and medium density so that the sphere and background are properly controlled by the barotropic EOS.

Your changes to setup_sphereinbox are not bug fixes.

Note: There may be a bug in that the stated central and edge densities are different than the actual densities due to the particle nature of the code. I will look into this shortly.
Update: there is no bug. I was using incorrect splash defaults when plotting.

[danieljprice]

this is mainly a bug fix - we spent about 3 weeks trying to find the bug with the outer boundary, it was definitely not correct as written. The "sphereinbox" setup I suggested to rename as "sphereinbox" is not hugely informative. We're ultimately trying to do dust...

[danieljprice]

but the test failure here was as you predicted @jameswurster from Sedov test going wrong with the recent merge of the du/dt timestep condition to try to prevent the u < 0 failure

[jameswurster]

Re the boundaries: Matthew and I extensively tested the boundary interface for Wurster & Bate (2018). We did find some issues if the cloud was too hot or did not have a consistent BE profile. During later tests with Ian, we found issues if the rotational profile was too large (to be fair, the rotational profile in that case was fast enough that it promoted expansion rather than contraction). As long as the profile is consistent and physically motivated, I found the previous boundaries work well.
Note that creating a BE sphere when setting mass and radius (option 6) works only for very selected pairs of values. I have patched this, but have not yet committed this. If you want to keep your revised boundaries, can we please add this as an option/prompt so that the user can decide if they want pressure equilibrium or not; this will also make my life easier since I am running simulations where the sphere's edge and medium are in pressure equilibrium as originally coded.

Re thermal timestep: Steven and I tested this extensively with no luck. Sedov is a very special case since most of the particles start with u = 0, and as soon as it increases slightly, the timestep drops.
Given the nature of the Verlet leapfrog method, if the energy is changing rapidly enough then it is possible & unavoidable to get negative energy since the update uses an old timestep, an old du/dt & a new du/dt, where the timestep reflects the old but not the new du/dt. This also means that u < 0 can be achieved for both implicit and explicit energy updates since PdV & shock heating are always explicit.
I've commented on this in a paper I'm starting to write with Ian; although the description is not overly detailed, I can forward it to you if it would be helpful.

[danieljprice]

to be clear, we just corrected a bug here, we did not change the boundary condition.

Actually there are two lines here: i) the cs line, we fixed a bug ii) the rho0 line, where we set the threshold density for the barotropic EOS slightly differently. This was because if it is set to exactly the exterior density then any small change in rho causes the temperature to immediately plunge

[jameswurster]

I still don't understand what the bug is that you fixed... i) the sound speed in the medium is currently set to ensure pressure equilibrium (assuming a barotropic EOS); the updated version has a higher background pressure, such that the outer regions of the BE sphere have less gas pressure than the background, meaning that the background will promote collapse rather than hinder expansion. ii) the original rho0 was 2x the background density, so we did account for fluctuations in density without an immediate temperature plunge. In my personal working version, I'm using 4x (assuming a density contrast of 30, although I suspect the exact factor does not matter as long as it's > 1). By default, the barotropic EOS includes a smooth 0.5dex transition across rho0 to prevent a step-function change in temperature.

Clarification on your changes would be much appreciated. Thanks