Answered step by step
Verified Expert Solution
Link Copied!

Question

1 Approved Answer

Use python to compute step1-4 (I already compute step1-3 and it is correct, the code of step 1-3 are shown at the picture below, only

Use python to compute step1-4 (I already compute step1-3 and it is correct, the code of step 1-3 are shown at the picture below, only need help in step 4)PLEASE HELP ME!!! U need to numpy package i imported?How can i send u the file?where can i find the numpy code?

image text in transcribed

  • Step 1: Build a 2D square box (edge length = 2 nm) containing a random arrangement of particles, each separated by a minimum distance of 0.24 nm.
  • Step 2: Calculate the g(r) of the configuration built in Step 1, without normalising it with respect to the ideal gas.
  • Step 3: Build a number of configurations (such as the one built in Step 1) sufficient to converge the g(r) of the system to an acceptable degree of accuracy. Compute the g(r) of the system and comment on why do you think you have generated enough configurations.
  • Step 4: Find the largest portion of empty space (i.e. void) within each of the configurations you have generated in Step 3. Calculate and visualize the probability density function of the largest void in the system as a function of its size.

image text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribedimage text in transcribed

o jupyter bridson Dact dat Dave jpeg 0 colors css D confpdb d _design.png data_coding.jpg geek.jpg google png Din_logo.png 0 kb_1.jpg Omd_data dat O m d_trajectory gro poisson disc.py python_logo png D README.md temperature_ensemble dat temperature_MD.dat ten_c jpeg Dtmp water_box png water_box tga Type here to search Edit View Insert Help + Cell Kernel N Run O C Widgets Code Validate L. In [34]: # Import a number is useful packages... import numpy as np import matplotlib.pyplot as plt from scipy.interpolate import interpid from scipy.signal import argrelextrema from mpl_toolkits.mplot 3d import Axes 3D from STUFF.poisson_disc import Grid import random In [35]: #step 1 # Set the Length and width (the x and y dimensions) of the box length = 20.0 # can be any units you like. A meaningful choice: [A] width = length # This is a square 2D box # Set the minimum distance between the particles (as you can see, this is not really a hard sphere Liquid! why?) min_r = 2.4 # [A]rad_py # Generate a random arrangement of particles - according to an algorithm called Poisson sampling. # Note that every time you run this cell, the algorithm will generate a DIFFERENT configuration! grid = Grid(min_r, length, width) # Random seed rand = (random.uniform(e, length), random.uniform(e, width)) # Do the sampling data = grid.poisson(rand) In [36]: # We use "def" to tell Python we are defining a function here to search Edit View Insert Cell Kemel Widgets Help Trusted 8 B Run I C Code Validate w OLO POLODUIT anu) In [36]: # We use "def" to tell Python we are defining a function # "unzip" is the name of the function - I chose that. You can be creative, but long names are cumbersome to be called # Later on within the Notebook, and if the name of a function is the same as the name of a variable, things get ugly. def unzip(items): # The list of arguments given as input is passed in the brackets. In this case,"items # is the one and only argument required as input # We use "return" to tell Python this is the quantity we want to get as output anytime we call this function return ([item[i] for item in items] for i in range(Ten(items[0]))) In [37]: # Let's visualize the box and the particles - i.e. our hard sphere-ish Liquid plt. figure(figsize=(6, 6), dpi= 80, facecolor='w', edgecolor='k') plt.scatter(*unzip(data), s=250, facecolors='orange', edgecolors='blue',linewidth=2.5) plt.tick_params (axis='both', which='both'; length=0, labelleft-False, labelbottom=False) plt.xlim(e, length) plt.ylim(e, length) plt.gca().set_aspect('equal', adjustable='box') pit.rcParams ["figure, figsize"] = (8,8) plt.show() here to search 1 X + N Run OC Code P sa _aspeLLCyuda aujusLOVAC VUA plt.rcParams ["figure. figsize"] = (8,8) plt.show() Validate all 3061-45....jpg 880682d8-d667-4....jpg a71c99e5-68eb-4b....jpg 2ddde91d-0192-4....jpg pe here to search Et S29 + 8 Run I C Code - Validate will In [38]: #step 2 def py_rdf(r, s, dr, dim): from numpy import zeros, sqrt, where, pi, mean, arange, histogram, absolute num_particles = len(r) rMax = S/2.0; edges = arange(e., rMax + dr, dr) num_increments = len(edges) - 1 = zeros(num_increments) radii = zeros(num_increments) numberDensity = len(r) / S**dim # Compute pairwise correlation for each particle for index in range(num_particles): d = 0.0 for i in range (dim): dp = absolute(r[index,i] - r[:,i]) mask = dp>S/2.0 dp[mask] = 5 - dp[mask] d += dp dp d = sqrt(d) d[index] = 2 * rMax .(result, bins) = histogram(d, bins-edges, density False) 8 += result 8 = g/(num_particles * numberDensity) Upe here to search + Validate Run IC Code dp[mask] = S - dp[mask] d += dp dp d = sqrt(d) d[index] = 2 * rMax (result, bins) = histogram(d, bins=edges, density=False) 8 += result 8 = g/(num_particles . numberDensity) # Normalize the g(r) dividing by the g(r) of an ideal gas - in 20! if dim == 2: for i in range(num_increments): radii[i] = (edges[i] + edges[i+1]) / 2. router = edges [i + 1] rInner = edges[i] g[i] = g[i] / (2.0 - pi * (router-rInner)* radii[i]) # Needed to compute the 3D g(r) (blue box) # Normalize the g(r) divinding by the g(r) of an ideal gas - in 3D! if dim == 3: for i in range(num_increments): radiii] = (edges[i] + edges [i+1]) / 2. router = edges [i + 1] rInner = edges[i] g[i] = g[i] / (4.0 * pi . (router-rInner)* radii[i] * radii[i] ) return (radii, g) Type here to search ge B Run O C Code - Validate wil return (radii, B) In [39]: # Compute the g(r) res_dr = 0.1 # resolution n_data-np.array(data) # We store the atomic positions into a numpy array # call the py_rdf function. Note that the name of the arguments can be different, but not their order. # For instance, we have called n_data the first argumentwhich is indicated as r in the definition # of py_rdf above. Which is fine, Python will take care of that. But, we are not allowed to swap, say, #n_data and Length!! rad_py, 8_r_py = py_rdf(n_data, length, res_dr, 2) # We interpolate the result to get a smooth Line connecting the dots res = 200 # How many points do we want in the interpolated Line safe = le-1 # Don't you worry about this one... # We interpolate the g(r) using a cubic spline f_cub = interpid(rad_py, g_r_py, kind='cubic') xnew = np.linspace(safe, (length/2.0)-safe, num=res, endpoint=True) In [40]: # Plot plt. figure(figsize=(6, 6), dpi= 80, facecolor='w', edgecolor='k') plt.tick_params (axis='both', which='both', length=10.0, labelleft=True, labelbottom=True, labelsize=2e.e) # Interpolated g(r) plt.plot(xnew, f_cub(xnew), color='purple', linestyle='solid', linewidth=2.0, label='Interpolation') # Actual alc) points Type here to search t 529 * Run C Code Validate will In [40]: # Plot pit.figure(figsize=(6, 6), dpi= 8e, facecolor='w', edgecolor='k') plt.tick_params (axis='both', which='both', length=10.0, labelleft=True, labelbottom=True, labelsize=20.0) # Interpolated g(r) plt.plot(xnew, f_cub(xnew), color='purple', linestyle='solid', linewidth=2.e, label='Interpolation') # Actual g(r) points plt.plot(rad_py, B_r_py, 'o', markerfacecolor='green', markersize=8, markeredgecolor='blue', linewidth=2.5, label='Actual 8(r)) plt.legend() plt.xlabel('r [$\AA$]', fontsize=26) plt.ylabel('s(r)', fontsize=26) plt.show() Interpolation Actual gir) 2.5 2.0 - g(r) Type here to search BE S29 + 8 B A M Run IC Code Validate will PALUEILLIMO) DILE- plt.ylabel('g(r)', fontsize=26) plt.show() Interpolation Actual g(n) g(r) 2.5 5.0 r [8] 7.5 10.0 Type here to search Trusted + X E M Run C Code Validate Lalit In [41]: #step 3 # Needed for one of the blue boxes below... # This may take a while, check the - in the [] on the left of the cell, when that's gone, you are good to go n_confs=25 # Number of configurations we are going to take into account g_r=np.zeros((len(8_r_py), n_confs)) 8_r_ave=np.zeros(len(_r_py)) n_ave- for i in range(e,n_confs): data = grid.poisson(rand) n_data=np.array(data) rad_py, B_r_py = py_rdf (n_data, length, res_dr, 2) g_r[:,il-g_r_py n_ave=n_ave+len(n_data) # Ensemble average of the gir) 8_r_ave=(np. sum(_r, axis-1))_confs # Average number of particles n_ave=n_ave_confs In [42]: # Needed for one of the blue boxes below... them tobe a while, check the in the 1 on the left of the cell, when that's gone, you are good to go Type here to search Bt S2 9 Trusted + * B A Run C Code 5_r_ave-inp. sum gr, axis=1))_conts 2 Validate # Average number of particles n ave=n_ave confs In [42]: # Needed for one of the blue boxes below... # This may take a while, check the - in the [] on the left of the cell, when that's gone, you are good to go n_confs=5 # Number of configurations we are going to take into account g_r=np.zeros((len(8_r_py),n_confs)) 8_r_ave=np.zeros(len(8_r_py)) n_ave=0 for i in range(e,n_confs): data = grid.poisson(rand) n_data=np.array(data) rad_py, _r_py = py_rdf(n_data, length, res_dr, 2) g_r[:,i]=g_r_py n_ave=n_ave+len(n_data) # Ensemble average of the g(r) 8_r_ave=(np.sum(@_r, axis-1)_confs # Average number of particles n_ave=n_ave_confs In [43]: # Comparison of the g(r) for just one configuration and that obtained upon taking an ensemble average Type here to search + 8 N Run O C Code Validate will In [43]: # Comparison of the g(r) for just one configuration and that obtained upon taking an ensemble average f_cub_ave = interpid(rad_py, g_r_ave, kind='cubic') xnew = np.linspace(safe, (length/2.0)-safe, num=res, endpoint=True) pit.figure(figsize=(6, 6), dpi= 8e, facecolor='w', edgecolor='k') # figsize determine the actual size of the figure plt.tick_params (axis='both', which='both', length=10., labelleft=True, labelbottom=True, labelsize=20.8) plt.plot(xnew, f_cub(xnew), color='blue', linestyle='solid', linewidth=1.0, label='One configuration) plt.plot(xnew, f_cub_ave(xnew), color='purple', linestyle='solid', linewidth=3.0, label='Ensemble average') plt.xlabel('r [S\AA$]', fontsize=26) plt.ylabel('g(r)', fontsize=26) plt.legend() plt.show() One configuration Ensemble average g(r) 0.5 L Type here to search Run Ic Code Validate One configuration Ensemble average 2.5 2.0 1.5 g(r) 0.5 l woman 0.0 2.5 7.5 10.0 5.0 r[] Type here to search

Step by Step Solution

There are 3 Steps involved in it

Step: 1

blur-text-image

Get Instant Access to Expert-Tailored Solutions

See step-by-step solutions with expert insights and AI powered tools for academic success

Step: 2

blur-text-image

Step: 3

blur-text-image

Ace Your Homework with AI

Get the answers you need in no time with our AI-driven, step-by-step assistance

Get Started

Recommended Textbook for

Students also viewed these Databases questions