Disc

[giovanni-rosotti]

Implementation of the disc setup. I wanted to describe it in the userguide, but I can't find a section about ICs. We do have something very short in section 5.1 (description of the parameters), but it's far too short. Shall I create a new section (which for the moment will contain only the description of the disc setup)?

[rbooth200]

Looks good so far, but I haven't been through in detail. Just a couple of comments:

1) Can we have the default disc and planet parameters such that someone can just turn on the planet and get a 'sensible' simulation. Currently the planet would be at the disc's inner edge. Maybe use the De Val-Borro setup?

2) Is the switch in the implementation of the locally isothermal eos rubust? Maybe we should discuss what the desired behaviour actually should be?

3) Can we add parameters to set the eccentricity and inclination of the planet

4) Could multiple planets be handled naturally, or should we solve this by adding a multiple planet setup in python, perhaps?

[giovanni-rosotti]

1- oh right i didn't notice i left the planet at 1, will change that. i'll see what the del-borro setup does for the disc surface density and temperature slopes 2-yes let's chat about that 3-will do 4-no, it would need a vector of parameters and there's no way to do it from the parameter file. The solutions are either python as you say, or reading another file (as FARGO does)

[rbooth200]

Ok sounds good. Re 4), should we open an issue and leave it for now? People might want this, but I'd agree that the best solution is not obvious right now.

[dhubber]

I wanted to describe it in the userguide, but I can't find a section about ICs. We do have something very short in section 5.1 (description of the parameters), but it's far too short. Shall I create a new section (which for the moment will contain only the description of the disc setup)?

Hmmmm, yes, this is perhaps important, at the very least a small description. We don't want to go into too much detail or it could bloat the userguide considerably. If we do have a section like this, then perhaps it should be limited to a concise description of the setup, a reference to any paper (if one exists) that perhaps describes it in more detail, and finally a list of the essential parameters that are required. Hopefully each one would be on average about half a page in length so it will only add 10 pages or so. I think that's not too much. Maybe we don't want to concentrate on this right now since the paper is the priority, but at least having a chapter/section in place with stubs for each test so we have it ready to be filled-in in the future?

Regarding point 4 above, if you're talking about general configurations, then sure it's best to leave these to the user to make their own IC files or create them via python. If you're talking about some 'standard' multi-planet simulation that is done often for comparison (and it's only one of a handful of cases), then this can also be hard-wired into the IC code. But since planets is not in my jurisdiction, I'll leave you guys decide this ;-)

[giovanni-rosotti]

All right, so I think what I will do is I will add a section where I explain the disc setup (better to do it now that I remember what I did). If later on (after paper submission!) we want to describe other setups, we can do it.

[giovanni-rosotti]

This should now be ready to be merged

[rbooth200]

Looks fine to me, apart from the one thing I mentioned above.

[giovanni-rosotti]

Should now be really ready to be merged. I've also fixed the rotation curve for the 2d case, for what is worth. Finally, I discovered a horrible bug with integer division, since the only standard function for random number generation in C returns integers...

[dhubber]

Okay, I'll have a quick look sometime this weekend and hopefully will be fine to merge. One point though, I'd recommend using the RandomNumber class for generating ICs like how it's done above for two reasons (i) it helps stop the kind of horrendous/horrible (any other adjectives?? :-) ) bug you found above (I've done the same myself btw) (ii) it allows reproducible ICs using, for example the xorshift random number generator (or any other reproducible generator we care to implement in the future). On this subject, I see this code uses std::rand whereas I've used the C random number function for the default system generator. Maybe I should 'upgrade' this to the C++ std version?

[rbooth200]

I see this code uses std::rand whereas I've used the C random number function for the default system generator. Maybe I should 'upgrade' this to the C++ std version?

I would not bother. Unless we switch to the fully flexible C++ random number generators that offer several different algo's. But what we have works and its not like we're doing MCMC, so who cares?

[giovanni-rosotti]

Done...

I think that xorshift for this application is fine, but I wouldn't want to do serious Monte Carlo stuff with it. If we needed it, I would go for an external library, rather than rolling our own.

One small thing: I saw that when using the standard library generator, there is a modulus operation ((rand()%RAND_MAX)/(FLOAT)RAND_MAX). Why is this needed? Isn't the division enough? The modulus is essentially a no-op since rand() is guaranteed to be smaller than RAND_MAX...

[dhubber]

I think that xorshift for this application is fine, but I wouldn't want to do serious Monte Carlo stuff with it. If we needed it, I would go for an external library, rather than rolling our own.

Yes, I was just mainly talking about reproducibility of ICs. People have done MCRT with worse random number generators though ;-) But we're probably not doing that anyway (at least not anymore) so nothing to worry about.

One small thing: I saw that when using the standard library generator, there is a modulus operation ((rand()%RAND_MAX)/(FLOAT)RAND_MAX). Why is this needed? Isn't the division enough? The modulus is essentially a no-op since rand() is guaranteed to be smaller than RAND_MAX...

Yes, good point. I think I saw several times the old C-style functions written like this (maybe just as a pointer of how to limit the random number range using the modulus?) so just copied it but guess it probably doesn't make any difference like you mention.

One last point regarding the code; it seems like the particle properties are set up with some Monte-Carlo sampling. Wouldn't this be a good candidate for the regularised IC approach? An analytical density function can be provided plus the SetParticleProperties function for setting velocities and other properties after regularisation. I'm not saying to do it now (since time is running short) but maybe later.

Otherwise, I'm happy with merging if Richard is (assuming your previous concerns were resolved?)

[rbooth200]

I'm happy with it so I will check with Giovanni and merge.

RE. the regularization: I agree that it would be a good thing, but it probably needs some tuning. How well would it work with sharp edges to the disc?

[dhubber]

How well would it work with sharp edges to the disc?

I had problems with sharp edges in general before, but this is why I was experimenting with this 'GetSmoothedDensity' function (or whatever it was called) which effectively convolves your given density field with the kernel and uses that smoothed density field. Then it works perfectly well. In any case, I was thinking this is a case where it only needs to be run for a handful (5 - 10) of iterations to remove the worst of the start-up noise. We can try it later once the paper is out of the way.