Commit f0bc902a authored by Laurien's avatar Laurien

added docstrings

parent a1851239
......@@ -37,7 +37,7 @@ def lennardjones(eps, sigma, Ntheta, pos_position, position):
"""
positionmatrix = position[:, np.newaxis, :]
radius = np.linalg.norm(np.repeat(positionmatrix, Ntheta, axis=1) - pos_position + 0.00001, axis=2)
E_theta = 4*eps * np.sum((radius/sigma)**-12 - (radius/sigma)**-6, axis=0)
E_theta = 4*eps * np.sum((radius/sigma)**-12 - (radius/sigma)**-6, axis=0)
return E_theta
......@@ -55,16 +55,16 @@ def extrabead(i, k_b, T, Ntheta, E_theta, sumW_l, w_lfinal, position, pos_positi
w_lfinal = w_l[n]
position = pos_position[n,:]
return w_lfinal, sumW_l, position
return w_lfinal, sumW_l, position, n
def upperlowerlimit(alpha1, alpha2, w_lfinal, w_ave, lowlimit, uplimit):
def upperlowerlimit(i, k_b, T, alpha1, alpha2, w_lfinal, w_ave, lowlimit, uplimit):
"""
function: calculate weights corresponding to the functions
input: potential energy, zero matrices
output: weight corresponding to the chosen position, sum of all weights, chosen position of theta: n
"""
w_ave=np.average(w_lfinal)
w_ave=np.average(w_lfinal[:i])
weight3=w_lfinal[2]
lowlimit=alpha1*(w_ave/weight3)
uplimit=alpha2*(w_ave/weight3)
......@@ -93,43 +93,86 @@ def crossing(i, position):
return cross
def endtoend(position, maxbead):
"""
function: calculates the end-to-end length of a polymer
input: maximum bead
output: R
"""
e2e = position[maxbead-1,:] - position[0,:]
e2e = np.linalg.norm(e2e)
return e2e
def calc_mean(obs, weight):
"""
function: calculates the mean of an observable
input: observable, weights
output: mean value
"""
weightedobs = np.multiply(obs, weight)
meanvalue = np.sum(weightedobs, axis=1) / np.sum(weight, axis=1)
return meanvalue
def calc_std(obs, weight):
"""
function: calculates the standard deviation of an observable
input: observable, weights
output: standard deviation
"""
variance = calc_mean(obs**2, weight) - calc_mean(obs, weight)**2
std = np.sqrt(variance)
return std
def function(x,a,b):
return a * x**b
"""
function: theoretical function a * x^b to perform the least-square fit on
input: x data, parameters a and b
output: function y = a * x^b
"""
return a * (x-1)**b
def enrich(PolWeight, NewWeight):
def enrich(PolWeight, NewWeight, position):
"""
function: performs enriching of the polymer ensemble (PERM)
input: positions, polymer weight
output: New weight, positions duplicate
"""
NewWeight1 = 0.5 * PolWeight
newchain = position
NewWeight = 0.5 * PolWeight
#print("Weight is halved!")
NewWeight = NewWeight1
return NewWeight
return NewWeight, newchain
def prune(PolWeight, NewWeight):
def prune(PolWeight, NewWeight, stopcondition):
"""
function: performs pruning of the polymer ensemble (PERM)
input: polymer weight
output: New weight, stop condition
"""
R=np.random.rand(1)
if(R<0.5):
NewWeight = 2 * PolWeight
else:
stopcondition = True
return NewWeight
return NewWeight, stopcondition
def fitting(x, obs, min, max):
"""
function: performs a least-squares fit on the data points, with respect to the theoretical exponential function
Input: R
Output: optimum values and (co)variance
"""
x = np.linspace(2,max-min+2,max-3)
opt, cov = curve_fit(function, x, obs[min:max])
perr = np.sqrt(np.diag(cov))
return opt, cov, perr
......
Markdown is supported
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