# Create a random polymer chain

I want to create a polymer chain (in 2D) of a given length such that:
1. First monomer is at {0,0}
2. All other monomers are in the positive x half-plane
3. The distance between two bonded monomers is r0
4. No two non-bonded monomers come closer to each other than rc (rc>r0)

Edit: I found a way to implement these rules:
r0 FoldList[AngleVector,{0, 0},RandomReal[{-1, 1} ArcCos[rc/(2 r0)], n - 1]]
However, this solution does not sample all possible configurations, but only a small subset.

• So... a self-avoiding random walk? Did you try searching the site for random walk implementations? – J. M. will be back soon Apr 11 '18 at 13:23
• @J.M.needshelp. It is not exactly a self-avoiding random walk. But I did search for random walk. This is the closest I got to what I want. – pukkandan Apr 11 '18 at 13:38
• Your fourth condition sounds like "self-avoiding" to me... – J. M. will be back soon Apr 11 '18 at 13:40
• Yes, it has to be self avoiding. But rather than just "not crossing" each other, it has to stay within a certain distance from each other. For example, this is self-avoiding, but does not satisfy rule 4 – pukkandan Apr 11 '18 at 13:59

This could get you started. As always, we use Nearest to determine collision, but put the cheaper collision detection with the left half plane in front. This features also a buffer reservoir for future random steps since creating many random numbers at once is usually much more performant.

r0 = 0.1;
rc = 0.15;
dim = 2;
maxtrials = 200;
reservoircounter = 1;
chaincounter = 2;
maxchainlength = 10000;
reservoirLength = 1000;
getReservoir[n_, r0_] := RandomPoint[Sphere[ConstantArray[0., dim], r0], n];
chain = ConstantArray[0., {maxchainlength, dim}];
reservoir = getReservoir[reservoirLength, r0];
chain[[2]] = x = chain[[1]] + (# Sign[#[[1]]]) &[RandomPoint[Sphere[ConstantArray[0., dim], r0]]];

chaincounter = 2;
While[chaincounter < maxchainlength,
nf = Nearest[chain[[1 ;; chaincounter - 1]] -> Automatic];
ncollisions = 1;
iter = 0;
While[ncollisions > 0 && iter < maxtrials,
reservoircounter++;
iter++;
If[reservoircounter > reservoirLength,
reservoircounter = 1;
reservoir = getReservoir[reservoirLength, r0];
];
xnew = x + reservoir[[reservoircounter]];
ncollisions = If[xnew[[1]] >= 0., Length[nf[xnew, {∞, rc}]], 1];
];
If[iter >= maxtrials,
Break[];
,
chain[[chaincounter + 1]] = x = xnew;
chaincounter++;
]
];


And some visualization:

Graphics[{
Line[chain[[1 ;; chaincounter]]], Blue,
Point[chain[[1 ;; chaincounter]]],
Red, Opacity[0.15], Disk[#, rc/2] & /@ chain[[1 ;; chaincounter]],
Darker@Green, PointSize[0.02], Opacity[1], Point[chain[[1]]],
Darker@Red, Point[chain[[chaincounter]]]}
]


The algorithm can easily be made working for dimension 3 by setting dim = 3. This way, one can obtain something like this:

Graphics3D[{
Orange, Specularity[White, 30],
Sphere[chain[[1 ;; chaincounter]], rc/2],
Red, Darker@Green, Sphere[chain[[1]], rc], Darker@Red,
Sphere[chain[[chaincounter]], rc],
}, Lighting -> "Neutral"]


• Awesome! Thanks I will keep the question unanswered for now in case anyone else has a more elegant solution – pukkandan Apr 12 '18 at 7:37
• You're welcome! – Henrik Schumacher Apr 12 '18 at 7:58