-4

I would like to produce some 3D distributions of points using Mathematica 7.0, that look like the picture below :

A 3D fractal-like/percolation-like distribution of random points in 3D space

How could I do that ? What are your suggestions ? What Mma 7 codes could do a distribution of points which is approximately like this ? I'm not specifically looking for a diffusion limited aggregation method. Any other method would be interesting, if it's reasonably fast.

Ideally, the code should compile very fast.

Adding some color shades to the distribution would be a nice option.

EDIT : The picture above was taken from this topic : Distribution of random points in 3D space to simulate the Crab Nebula which is about filaments and sub-structures in the Crab nebula.

Community
  • 1
Cham
  • 4,093
  • 21
  • 36
  • Since the image looks MMA-made, do you happen to have a reference or (even better) some code to show for it? – Yves Klett Feb 04 '14 at 17:10
  • The picture was grabbed from another topic, that gives an answer to a similar question. But the code isn't compatible with Mma 7.0 and the results aren't the same at all (it was about filaments and random walks). I'm restarting my project from scratch. – Cham Feb 04 '14 at 17:15
  • 6
    You could have at least linked to where the picture came from, that is here. I did post the code to do the simulation in that very post. (!!) You can't expect people to just guess that this is diffusion limited aggregation. – Szabolcs Feb 04 '14 at 17:15
  • The code doesn't work at all with Mma 7.0. I've added a reference to the picture's origin at the end of my question. – Cham Feb 04 '14 at 17:16
  • @Szabolcs you got my +1 there long ago. Not sure how to treat this Q if is just about back-porting. – Yves Klett Feb 04 '14 at 17:22
  • It's not back porting. I'm restarting my project from scratch. The code given by Szabolcs (that produced the picture above) doesn't work with Mma 7.0. – Cham Feb 04 '14 at 17:24
  • @Szabolcs I'm not looking specificaly for a diffusion limited aggregation method. Any method that could generate random 3D fractal-looking distributions of points, that look-like the picture above, would do (if it works in Mma 7.0). – Cham Feb 04 '14 at 18:01
  • @Cham "Backporting to 7.0" means "make the code work with v 7.0". If, as you say, that is not what you mean, then what are you looking for? – Szabolcs Feb 04 '14 at 18:02
  • Ah, sorry, I interpreted "back-porting" differently. Sorry for the confusion. Yes, I need the code to be compatible with Mma 7. – Cham Feb 04 '14 at 18:02
  • @Cham OK, then if you be useful to spell that out very clearly in the question to avoid people misunderstanding it. Regarding DLA: my suggestion is to implement it in C for performance, and read the data back into Mathematica for plotting. I think it's easier to implement in C than to make it perform well in Mma. This applies to the DLA algorithm only, not to any other possible ideas you might come up with. – Szabolcs Feb 04 '14 at 18:03
  • Hmm, I don't know any C programming. And unfortunately, I only have access to Mma 7.0. The place where I work can't afford an upgrade (too costly). – Cham Feb 04 '14 at 18:04
  • Did you want to translate Szabolcs' answer specifically or any answer on that page? – rm -rf Feb 04 '14 at 18:25
  • What do you mean by "translate" ? I don't want to translate any answer. I need to restart the project from a fresh basis, without reference to "filaments" or "sub-structures". Most of the codes from the "Crab" topic doesn't work with Mma 7.0, or are specifically targeting "filaments". That is not what I'm looking for here. This topic is about a fractal-like distribution or diffusion of points in 3D. – Cham Feb 04 '14 at 18:34
  • 1
    We don't know what you mean by filaments or sub-structures. That answer referred to filaments because you asked for it! Restarting your project is not really relevant to the question here (it's your research problem). To get an answer here, you should be more specific about what exactly you want, specifics about the distribution and references (if possible). Without any of that, it looks like a "Hail Mary" question where you just throw out a complex problem and ask for code specifically for version 7, hoping that someone out there will understand exactly what you want and code it for you... – rm -rf Feb 04 '14 at 19:49
  • @rm -rf : Again, the question here is not about filaments and sub-structures (this was pertinent for the Crab-like model, in the other topic). What I'm asking here is MUCH simpler : just a distribution of random points in space, which has the approximate shape of the picture shown above. I'm asking for any way of doing this, so I can't be more specific about how to do it, since I don't know. It's exactly what I'm asking for, how ? What would be the possible codes to do a distribution like this ? – Cham Feb 04 '14 at 20:42
  • Szabolcs already linked to the process that will lead to such a distribution: diffusion limited aggregation (you can play with the parameters). You haven't explained why you don't want to use a DLA, despite the fact that it produces something very close to what you want. Also, the code is there for v8+, so perhaps it will be useful to others if you show how far you've gotten to writing a v7 implementation. We really really do not do coding requests on this site... I'm just trying to help you improve your question. – rm -rf Feb 04 '14 at 20:57
  • @rm -rf : the DLA solution is NOT working on Mma 7. How many times should I repeat this ? So please, don't refer to that "solution", it isn't a solution to the question ABOVE. – Cham Feb 04 '14 at 21:00
  • 1
    Would a collection of random walks work? Graphics3D[Table[Point@Accumulate[RandomReal[{-1, 1}, {300, 3}]], {10}], BoxRatios -> 1] – Simon Woods Feb 04 '14 at 21:11
  • @Simon Woods : Wow ! This is already an answer to the question above. Why don't you post it as an answer ? It's similar to the answer you already gave to the "Crab" question (filaments), but much simpler. In your answer, could you explain what the Point@Accumulate is doing ? – Cham Feb 04 '14 at 21:15
  • I haven't been following this comment discussion, I'd just like to add one thing regarding the DLA code: if you look carefully, when I wrote it I annotated it, showing what should and what shouldn't work in v7. I looked at it again, and I only found one thing I neglected to mention: you need to remove , CompilationTarget -> "C", RuntimeOptions -> "Speed". The rest should work in v7 correctly. I cannot test this, because I do not have access to v7 any more, I can only say that it works fine in v8. If it really does not work in v7, you should make it clear which parts you're ... – Szabolcs Feb 04 '14 at 21:37
  • ... having trouble with, otherwise I cannot give advice. – Szabolcs Feb 04 '14 at 21:38
  • Cham, try varying the numbers, 300 is the number of points in each "arm" and 10 is the number of arms. More arms with fewer points will give you a less, uh, armlike appearance. – Simon Woods Feb 04 '14 at 21:53
  • @Simon : ok for these two parameters, but what if you want to increase the "density" of the whole thing ? I mean shorter "arms" with more points ? – Cham Feb 04 '14 at 23:32
  • @Cham You didn't read the notes within the code, which I mentioned above. The part that gives the error only works with v9, however it's just an alternative way to display the result and it's non-essential. Regarding the speed, this can certainly be improved somewhat, but I don't believe it's possible to write an efficient 3D DLA in Mathematica. This particular algorithm is very well suited for procedural languages and won't be much more difficult to implement in C/Pascal/FORTRAN/Java/etc. than in Mathematica. It will however be considerably faster. – Szabolcs Feb 04 '14 at 23:51
  • Otherwise you might look at other solutions (i.e. not DLA) that produce a similar look but can run faster, like some of the suggestions above. – Szabolcs Feb 04 '14 at 23:54
  • @Szabolcs : your code gives me some kind of exploded shell, similar to the color picture you have shown for the "Crab" question. I don't get the distribution shown above it. How should I modify it to get the other distribution ? – Cham Feb 05 '14 at 03:01
  • I still don't understand how to get the picture above with Szabolcs's code, which gives a very different beast (exploded shell). Szabolcs, please, how did you got that distribution ? – Cham Feb 05 '14 at 23:14

1 Answers1

3

I found the code I used to generate the graphic you are referencing. (I thought I had deleted it.)

The code from this post is somewhat optimized, and it's made specifically for an inward growth: the particles are always started form the origin for the DLA simulation, not from a random outer position. Also, the "seed" of existing particles is a spherical shell, which made it unnecessary to handle escaping particles. Instead of rewriting this code to work with a central seed, I'm going to give you the original code I used for the image you referenced, with some caveats:

  • It is very slow. It was a quick experiment not intended to run fast. I think that Mathematica is not the right tool for doing a DLA simulation. It's good for prototyping, but the performance will be awful. Also due to the procedural nature of the algorithm, it's not much more difficult to implement this say, in C++ or Java, than in Mathematica. I recommend that if you are serious about making a DLA simulation, use a low level procedural language such as C/C++/Java/FORTRAN/Pascal/etc.

  • I tried to use only v7-compatible functions but I don't have v7. If something doesn't work, let me know which part and I'll see if there's a quick fix.

  • WARNING: If you include the Dynamic part (i.e. the Graphics3D part), it might use up all the memory and crash the front end after a while. I do not know if this is a front end bug or not. It is okay to just not evaluate that cell. In that case you won't see the structure growing in real time but the results will still be recorded.

  • If you prefer the original Point-look to the Sphere-look, just replace Translate[Sphere[{0, 0, 0}, 0.4], Dynamic@points] by Dynamic@Point[points].

The code:

points = N@{{0, 0, 0}};

nf = Nearest[points];

nd[p_] := EuclideanDistance[First@nf[p], p]

(* put this in a separate cell an evaluate on its own *)
Graphics3D[Translate[Sphere[{0, 0, 0}, 0.4], Dynamic@points], 
 Axes -> True, BoxRatios -> Automatic, 
 PlotRange -> (15 {{-1, 1}, {-1, 1}, {-1, 1}}), 
 PlotLabel -> Dynamic@Length[points]]

(* this goes in a separate cell again *)
Do[
 r0 = 2 Max[Max[Norm /@ points], 5];
 r1 = 1.5 r0;
 pt = r0 Normalize@RandomReal[NormalDistribution[0,1], 3];
 While[Norm[pt] < r1 && nd[pt] > 1, pt += RandomReal[NormalDistribution[0,1], 3]];
 If[Norm[pt] < r1,
  AppendTo[points, pt];
  nf = Nearest[points]
 ],

 {100000}
]

Stop the simulation when you like using Alt-.. The results computed so far will be stored in points.

Community
  • 1
Szabolcs
  • 234,956
  • 30
  • 623
  • 1,263
  • Hey ! That's now working like a charm! (using the Dynamic@Point[points]). Thanks a lot for sharing this ! Once the simulation has been stoped, how do I extract the points {x, y, z} coordinates ? – Cham Feb 07 '14 at 00:36
  • @Cham As I said at the very end, the results are stored in points as {{x1,y1,z1}, {x2,y2,z2}, ...}. Should have put that line at the top together with the rest of the text. I can see how one could miss it. – Szabolcs Feb 07 '14 at 00:37
  • Yes, the process is very slow. But the animation is VERY nice, and allows the user to stop the process once the "sculpture" has a nice appearance. I think this is a great piece of code, even if it's very slow. Thanks again ! – Cham Feb 07 '14 at 00:43
  • Oh, just a stupid question : I guess that the process is really random ? I mean that restarting the process will not give the same shape, again and again, isn't ? – Cham Feb 07 '14 at 00:45
  • @Cham Just be careful! Mathematica will crash if you leave it running (with the live animation) for a long time. This is one reason I didn't post this code originally, the other being that it's so slow I couldn't make the filaments grow inwards from a large sphere (as in the other post). Actually I took this code from a bug report I sent about a crash ... – Szabolcs Feb 07 '14 at 00:46
  • @Cham If you use SeedRandom[123] before running the code, it'll give the same result every time. This is often useful for generative art. Sometimes to want to be able to reproduce an "accidentally good looking" result. Play with different seeds (numbers in place of 123) for different results. – Szabolcs Feb 07 '14 at 00:48
  • Ok, thanks again, Szabolcs. Currently, Mma has done 3315 particles without crashing (after about 5 minutes). I'll experiment a lot with this. – Cham Feb 07 '14 at 00:50
  • Szabolcs, in your code above, there are two NormalDistribution[0, 1]. What are they doing ? – Cham Feb 07 '14 at 14:44
  • @Cham Look it up in the docs and check under "Applications". That's what I did to find out how to generate normally distributed random numbers in v7. – Szabolcs Feb 07 '14 at 15:37
  • Szabolcs: well, of course I know what is a Normal distribution ! I just asked what it is doing in your code. – Cham Feb 07 '14 at 15:53
  • @Cham It doesn't do anything. It represents a normal distribution. Then RandomReal draws values from it. It's documented for v7 in the section I referenced. – Szabolcs Feb 07 '14 at 16:05
  • Szabolcs: You didn't understood the question (was it so obscure ?). I mean why are you using a normal distribution in the percolation code ? What is its function there ? – Cham Feb 07 '14 at 16:08
  • Or if you prefer, why the random declaration in the code ? I'm just trying to understand how the code works. – Cham Feb 07 '14 at 16:12
  • @Cham Do you mean, why normal distribution instead of a simple uniform distribution as in RandomReal[{-1,1}]? I used a normal distribution because it's spherically symmetric and I wanted the particle to be able to take a step in any spatial direction with the same probability. RandomReal[{-1,1}] samples from a cube so the particle is slightly more likely to go towards the diagonal of the cube. Actually it makes little difference for the visual appearance of the final result. – Szabolcs Feb 07 '14 at 16:42
  • RandomReal[NormalDistribution[0,1], 3]] simply gives you 3 normally distributed random numbers with mean 0 and spread 1. – Szabolcs Feb 07 '14 at 16:43
  • Ok, thanks. A last question : how can we control (change) the "width" of the arms in the distribution ? Can we have thin percolation arms (a bit like a random walk) ? – Cham Feb 07 '14 at 16:49