Commit b2afa55e by Laurien

### Overwritten

parent 12260e14
 ############################################### IMPORT MODULES ############################################### import numpy as np np.seterr(divide='ignore', invalid='ignore') import matplotlib.pyplot as plt #from anim import make_3d_animation import scipy as sci from scipy import integrate from scipy import optimize from scipy.optimize import curve_fit import math from matplotlib import pyplot as plt from mpl_toolkits.mplot3d import Axes3D from matplotlib.colors import cnames from matplotlib import animation from IPython.display import HTML ################################# FUNCTIONS ############################################# def roulette(i, d, Ntheta, position, theta): """ function: determines the possible new positions of the particle, based on a distribution of theta input: previous position, theta output: new position, distance between the two particles(x,y), absolute distance between the particles (radius) """ alpha = 2*np.pi * np.random.rand(1) possiblestep = np.array([d*np.cos(theta+alpha), d*np.sin(theta+alpha)]).T pos_position = np.tile(position[i-1, :],(Ntheta,1)) + possiblestep return pos_position 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) return E_theta def extrabead(i, k_b, T, Ntheta, E_theta, sumW_l, w_lfinal, position, pos_position): """ 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_l = np.exp(-E_theta/(k_b*T)) sumW_l = np.sum(w_l) probability = w_l/sumW_l choices = np.arange(Ntheta) n = int(np.round(np.random.choice(a=choices, p=probability))) #Chooses a value between 1 and 5 according to the assigned weights w_lfinal = w_l[n] position[i,:] = pos_position[n,:] return w_lfinal, sumW_l, position[i,:] def upperlowerlimit(alpha1, alpha2, w_lfinal, w_ave, lowlimit, uplimit): w_ave=np.average(w_lfinal) weight3=w_lfinal[2] lowlimit=alpha1*(w_ave/weight3) uplimit=alpha2*(w_ave/weight3) return lowlimit, uplimit def crossing(i, position): cross = False dx = np.zeros(shape=(L, L-1)) for j in range(0,i): dx[i,j] = np.linalg.norm(position[i,:] - position[j,:], axis=0) check = dx[i, 0:i] if np.min(check) < treshold: cross = True return cross def endtoend(position, maxbead): e2e = position[maxbead-1,:] - position[0,:] e2e = np.linalg.norm(e2e) return e2e #ef includeweight(obs, weight): #weightedobs = weight * obs #return weightedobs def calc_mean(obs, weight): weightedobs = obs * weight meanvalue = np.sum(weightedobs, axis=1) / np.sum(weight, axis=1) return meanvalue def calc_std(obs, weight): variance = calc_mean(obs**2, weight) - calc_mean(obs, weight)**2 std = np.sqrt(variance) return std def function(x,a): return a* (x)**(0.75)
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