
From:  Torsten Stuehn 
Subject:  Re: [ESPResSo] bug in random.tcl in mbtools 
Date:  Thu, 19 Mar 2009 12:44:42 +0100 
Useragent:  Thunderbird 2.0.0.12 (X11/20080213) 
Hi all, to generate Random Walk Polymers (polymercommand), ESPResSo (internally) uses the following algorithm to generate homogeneously distributed points on the surface of a sphere with radius=bond_length: zz = (2.0*d_random()1.0)*bond_length; rr = sqrt(SQR(bond_length)SQR(zz)); phi = 2.0*PI*d_random(); x = rr*cos(phi); y = rr*sin(phi); z = zz; The advantage of this algorithm is, that it needs only two random numbers and you never have to throw away any numbers. Experience has shown, that on modern CPUs this algorithm runs faster than the one Burkhard proposed. It is based on the (not very intuitive) fact, that the projection of homogeneously distributed points on a spheresurface onto any axis gives a normal distribution of points on that axis (a prove of that can be found in math books). Greetings, Torsten Tristan Bereau wrote:
Hello,I'm not sure whether a gaussian random number generator has already been implemented in the TCL shell of ESPResSo, but if not, there is a couple of methods described at the end of the "Old Testament" of molecular simulations (i.e. Allen and Tildesley), and Burkhard's is one of them.Best, Tristan On Mar 19, 2009, at 4:52 AM, Jacob Kirkensgaard wrote:Hi againThanks to Tristan and Burkhard  you are of course right. For my purposes my suggestion is sufficient but I will see if I have time to implement a small routine to do it properly as outlined by Burkhard. Is there a routine for generating Gaussian random numbers? As Tristan mentioned this would be a quick way to solve the issue.Best, Jacob On Mar 19, 2009, at 9:40 AM, address@hidden wrote:Hello, from Tristan's remarks I understand that there is a question of how to get random numbers uniformly in or on a sphere? If so, let me tell you one possible answer: 1. You first draw 3 RNs u_1, u_2, u_3 from [0:1] 2. Then you transform u_1 = 2 * u_1  1 etc., resulting in RNs in [1,1] 3. Calculate r2 = u_1**2 + u_2**2 + u_3**2 4. If r2 > 1 > throw everything away, go to step 1, try again 5. If r2 < 1: Either you are done (supposing that you want RN IN the sphere), or you set norm = 1 / sqrt(r2) u_1 = u_1 * norm etc. By this you project the output onto the surface of the unit sphere, so you get something ON the sphere. If that was trivial: Please accept my apologies for bothering you. Regards Burkhard._______________________________________________ ESPResSo mailing list address@hidden https://fias.unifrankfurt.de/mailman/listinfo/espresso
 Dr. Torsten Stühn MPI für Polymerforschung Ackermannweg 10 55128 Mainz / Germany Tel. +49(0)6131379268Fax +49(0)6131379100 EMail address@hidden
[Prev in Thread]  Current Thread  [Next in Thread] 