Commit fad6dc40 authored by Gilliam van Oudenaren's avatar Gilliam van Oudenaren
Browse files

Added more docstrings, fixed errorbars

parent d4da17ee
......@@ -78,6 +78,7 @@ def direction(a, b, c):
"""
Inputs:
- points a, b, c
Checks what the orienation of 3 given points is.
Output correspond to colinear, anti-clockwise or clockwise.
"""
......@@ -93,6 +94,7 @@ def checker(a, b, c, d):
"""
Inputs:
- points a, b, c d
Checks if the linesegments between a and b, and c and d intersect. Calls
the direction function to check for orientation. Operates on the assumption
that starting/end points of the line segments are never exactly on the
......@@ -112,6 +114,11 @@ def checker(a, b, c, d):
def chainmaker(chains = None, k = None, flipping = True):
"""
Inputs:
- numpy array for storing chain coordinates
- chain iteration k index
- boolean for flipping
Function used to produce k chains with a maximum length N. Calculates the
position to place the next bead based on weights determined by Lennard-Jones
potential and Monte Carlo methods. Flips the chain and starts building at
......@@ -160,6 +167,14 @@ def chainmaker(chains = None, k = None, flipping = True):
return beads, length
def data(storechains = True, flipping = True):
"""
Inputs:
- Boolean for storing chain coordinate data or not
- Boolean for flipping the chain if it gets stuck on one end or not
Stores chain data generated by chainmaker function in lists and numpy
array of chain coordinates.
"""
if storechains == True:
global chains
chains = np.zeros((k, N, 2))
......@@ -173,24 +188,29 @@ def e2edis(chains):
"""
Inputs:
- A set of polymer chains with size N
- The total number of beads in each chain
For a given set of chains, this function calculates the end to end length
of the subchains within a chain from point 0 to j. These lengths are added
to a single list and
of the subchains within a chain from point 0 to j. From this, the mean
length of chains with any number of beads up to the maximum amount is
calculated.
"""
chains = np.where(chains!=0, chains, np.nan)
chains0 = np.repeat(chains[:,0,:],chains.shape[1], axis=0).reshape(chains.shape)
chains = np.where(chains[:,2:]!=0, chains[:,2:], np.nan)
chains0 = np.repeat(chains[:,0,:], chains.shape[1],
axis=0).reshape(chains.shape)
R2 = np.sum((chains - chains0)**2, axis=2)
R2std = np.nanstd(R2, axis=0)
R2std = np.nanstd(R2, axis=0)/np.sqrt(np.sum(~np.isnan(R2), axis=0))
R2std = R2std[~np.isnan(R2std)]
R2mean = np.nanmean(R2, axis=0)
R2mean = R2mean[~np.isnan(R2mean)]
return R2mean[2:], R2std[2:]
return R2mean, R2std, R2
def linefit(beads, R2mean):
"""
Inputs:
- Array with bead numbers
- Array with mean lengths squared
Calculates parameter a for the function y = a*(N-1)^1.5. Takes inputs of
number of beads N and corresponding end to end lengths y.
"""
......
......@@ -5,6 +5,5 @@ import plots
plt.close("all")
beads, lengths, chains = functions.data(storechains = True, flipping = True)
plots.lengthplot(chains)
#plots.chainplot(beads, chains, np.arange(10))
\ No newline at end of file
plots.chainplot(beads, chains, np.arange(10))
\ No newline at end of file
import numpy as np
#Physical constants
kbT = 1e+2
kbT = 1e+8
epsilon = 0.25
sigma = 0.8
......
......@@ -4,33 +4,38 @@ import functions
from parameters import kbT
def lengthplot(chains):
# beadlengths = {}
# for i in range(len(beads)):
# if beads[i] in beadlengths.keys():
# beadlengths[beads[i]].append(lengths[i])
# else:
# beadlengths[beads[i]] = [lengths[i]]
# beads = list(beadlengths.keys())
# lengths = [sum(i)/len(i) for i in list(beadlengths.values())]
R2mean, R2std = functions.e2edis(chains)
N = np.arange(len(R2mean))+2
"""
Inputs:
- Array containing coordinates of generated chains
Calls function to calculate the end to end distance and standard deviation
for each amount of beads. Plots mean distance vs. bead number with
corresponding standard deviation.
"""
R2mean, R2std, R2 = functions.e2edis(chains)
N = np.arange(R2mean.shape[0]) + 1
a = functions.linefit(N, R2mean)
y = a*(N-1)**1.5
print("a=%.2f"%a)
exp = int(np.log10(kbT))
plt.figure()
plt.title("kbT = 10^%d" % exp)
plt.loglog(N, y, c = 'r', linestyle = '--')
plt.xlim((1, 2*max(N)))
plt.ylim((0.5, 2*max(R2mean)))
plt.errorbar(N, R2mean, R2std)
plt.xlim((2, 2*max(N)))
plt.ylim((1, 2*max(R2mean)))
plt.errorbar(N, R2mean, R2std,
linestyle = '', marker = 'o', markersize = 1, capsize = 2)
plt.ylabel("<R²>")
plt.xlabel("Polymer beads")
plt.show()
def chainplot(beads, chains, interval):
# if interval == None:
# interval = range(chains.shape[0])
"""
Inputs:
- Amount of beads in each chain
- Array with chain coordinates
- Array with chain numbers to plot
"""
for i in interval:
plt.figure()
plt.title("Chain %d" % i)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment