Assignment Questions 2

pdf

School

University of Waterloo *

*We aren’t endorsed by this school

Course

334

Subject

Mechanical Engineering

Date

Apr 3, 2024

Type

pdf

Pages

8

Uploaded by ElderWhale4083

Report
2. Water NVT 2 - Water: NVT Use the function simulate_solidwater() to study the dynamics of two rigid water molecules under NVT (canonical) conditions (change ensemble key to "NVT™ ). 1. Report the parameters of your simulation # change your parameters here parms_1 = { # place values here ‘steps’' :327, ‘skip_steps':1, "temperature':71., ‘dt’': 1 * unit.femtoseconds, "ensemble’ :'NVT', } parms_2 = { # place values here 'steps':61, 'skip_steps':3, 'temperature’:160., ‘dt': 2 * unit.femtoseconds, "ensemble’ :'NVT', } # run the sim here (check usage instructions for help): simulation1 = simulate_solidwater(parms_1) simulation2 = simulate_solidwater(parms_2) Mmoo & #'Progress (%)" "Step" "Potential Energy (kJ/mole)" "Temperature (K)" "Speed (ns/day)" "Time Remaining" 0.3% 1 -27.47715442235899 60.93034157024491 0 - 0.6% 2 -27.465133372268145 60.121087525363414 27.2 0:01 0.9% 3 -27.446789493523514 61.143776906171496 28.2 0:00 1.2% 4 -27.42235559206051 62.86486710856393 29 0:00 1.5% 5 -27.393578497498183 65.28576329863644 27 0:01 1.8% 6 -27.362855633828175 67.66828603991559 26.9 0:01 2.1% 7 -27.330173111407312 64.57612770484852 26.4 0:01 2.45% 8 -27.295204255001767 63.574484715423 26.5 9:01 2.8% 9 -27.261444642635347 64.704842632657 27 0:01 3.1% 10 -27.23459837211833 60.70498632251083 27.4 0:00 3.4% 11 -27.210272343022375 56.751211553395784 19 0:01 3.7% 12 -27.189158425631874 55.10331020189065 18.9 0:01 4.0% 13 -27.1691409901416 54.54243423114256 19.1 0:01 4.3% 14 -27.152099150380117 53.6320611318315 19.2 0:01 4.6% 15 -27.137234897731872 54.02581243346301 19.7 0:01 4.9% 16 -27.121603360205707 54.,37869978785319 19.8 0:01 5.2% 17 -27.107305835266224 52.16867082714445 20.3 0:01 5.5% 18 -27.095537235601892 52.21136293661746 20.3 0:01 5.8% 19 -27.08697093034741 51.17564750885838 20.4 0:01 6.1% 20 -27.078160126654023 52.62916445597792 20.6 0:01 6.4% 21 -27.070148400073048 51.55726699974624 20.7 0:01 6.7% 22 -27.063030765800644 48.029187453920386 20.7 0:01 7.0% 23 -27.053092408944526 48.45244735709546 21 0:01 7.3% 24 -27.044842123131797 48.3059011674182 21.3 0:01 7.6% 25 -27.040053982908148 49.9206198666751 21.5 0:01 8.0% 26 -27.03317295589239 50.47359986670718 21.6 0:01 8.3% 27 -27.02667536810909 48.961141392131566 21.8 0:01 8.6% 28 -27.017315037914106 47.77791399868938 22 0:01 8.9% 29 ~27.009074738089357 46.58479166753099 22.1 0:01
2. Plot the kinetic, potential, and total energies (stored in newly created results dictionary). Note your observations. import matplotlib.pyplot as plt # plots # (feel free to separate into more codeblocks) # (just make sure to enable them when publishing, i.e. double check your output pdf) potentiall = simulationi['potE’] potential2 = simulation2['potE’] kinetic1 = simulationi[ 'kinE'] kinetic2 = simulation2[ kinE'] totall = simulationl[ ' totalE'] total2 = simulation2[ ' totalE'] timel = simulationi['time’] time2 = simulation2['time'] plt.plot(timel, potentiall, label = "Potential Energy for First Set of Parameters") plt.plot(time2, potential2, 1label = "Potential Energy for Second Set of Parameters") plt.legend() plt.show() plt.plot(timel, kinetic1, label = "Kinetic Energy for First Set of Parameters") plt.plot(time2, kinetic2, label = "Kinetic Energy for Second Set of Parameters") plt.legend() plt.show() plt.plot(timel1, totall, label = "Total Energy for First Set of Parameters”) plt.plot(time2, total2, label = "TotalL Energy for Second Set of Parameters") plt.legend() plt.show() —— Potential Energy for First Set of Parameters —25.75 1 —— Potential Energy for Second Set of Parameters —26.00 A —26.25 A —26.50 A —26.75 A —27.00 A —27.25 A —27.50 A1 0.00 0.05 0.10 0.15 0.20 0.25 0.30 4.5 4.0 1 3.5 3.0 A 2.5 4 2.0 A 1.5 4 1.0 1 —— Kinetic Energy for First Set of Parameters —— Kinetic Energy for Second Set of Parameters 0.5 1 0.00 0.05 0.10 0.15 0.20 0.25 0.30
-22.5 4 -23.0 4 -23.5 1 -24.0 4 —24.5 4 —25.0 1 —— Total Energy for First Set of Parameters -25.5 —— Total Energy for Second Set of Parameters 0.00 0.05 0.10 0.15 0.20 0.25 0.30 Ans > In this case it looks like the total energy does not stay constant (makes sense because NVT system), and as expected because of that, you do not see the same trend where the kinetic energy is inversely proportional to the potential energy. In fact, it kind of looks like they follow the same general trend (for both systems) where they tend to increase and decrease togetehr (although not exactly at every moment and not perfectly in sync), just generally. 3. Compute the average and variance of the kinetic, potential, and total energies. Vary the length of the simulation for a constant time step ( dt ), what do you observe? For a fixed number of time steps, what happens when you vary the distance between them ( dt )? Which quantity do you expect will be conserved? # space for multiple simulation runs parmsa = { steps':100, skip_steps’':1, temperature' :90., dt': 1 * unit.femtoseconds, ensemble’ :'NVE', } parmsb = { steps’ :200, skip_steps’':1, temperature' :90., dt': 1 * unit.femtoseconds, ensemble’ : 'NVT } parmsc = { steps' :200, skip_steps':1, temperature':90., dt': 3 * unit.femtoseconds, ensemble’ :'NVT', } asimulation = simulate_solidwater(parmsa) bsimulation = simulate_solidwater(parmsb) csimulation = simulate_solidwater(parmsc) 19.0% 10 ~£0.DY40Y/01£0031/ 30.4/91914£3/114 49.4 vive 79.5% 159 -26.80856085961033 44.71658707557548 49.5 0:00 80.0% 160 ~26.949454125762156 47.317529084450236 49.6 0:00 80.5% 161 ~26.942190248680717 48.18427373000055 49.7 0:00 81.0% 162 ~26.841308030642054 48.02594134534121 49.8 0:00 81.5% 163 ~26.696632097897666 48.12853892447746 49.9 0:00 82.0% 164 ~26.583042422618476 51.012656464962824 49.9 0:00 82.5% 165 ~26.55421441734677 52.28544068898827 50.1 0:00 83.0% 166 -26.635237383723826 51.78444287543882 50.1 0:00 83.5% 167 -26.792836149510425 52.94554962078794 50.2 0:00 84.0% 168 -26.953177261812428 60.317040687855396 50.3 0:00 84.5% 169 -27.023279714871997 67.48051011745048 48.3 0:00 85.0% 170 -26.983689462971135 60.42995599732824 48.3 0:00 85.5% 171 -26.83876803769347 55.64320449862518 48.3 0:00 86.0% 172 -26.640051443656674 50.10106784917126 48.3 0:00 86.5% 173 -26.479206153180456 51.41216066000437 48.4 0:00 87.0% 174 -26.407238667656905 48.4028239073193 48.4 0:00 87.5% 175 -26.409076106610527 50.60176947760392 48.4 0:00 88.0% 176 -26.428505688003966 54.75716585974029 48.4 0:00 88.5% 177 ~26.37790930655158 54.980105723780355 48.4 0:00 89.0% 178 ~26.243618744069487 48.780336516612365 48.5 0:00 89.5% 179 ~26.03062005228424 39.90207617696733 48.6 0:00 90.0% 180 ~25.77532840623163 36.29123259412879 48.7 0:00 90.5% 181 ~25.618700443117937 29.76633262821862 48.8 0:00
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
  • Access to all documents
  • Unlimited textbook solutions
  • 24/7 expert homework help
from statistics import mean from statistics import variance potentiala = asimulation['potE’] potentialb = bsimulation['potE'] potentialc = csimulation['potE'] kinetica = asimulation['kinE’ kineticb = bsimulation['kinE’ kineticc = csimulation[ 'kinE’ —— totala = asimulation['totalE’ totalb = bsimulation['totalE’ totalc = csimulation['totalE’ —_—— timea = asimulation['time'] timeb = bsimulation['time'] timec = csimulation['time'] meanpa = mean(potentiala) meanpb = mean(potentialb) meanpc = mean(potentialc) variancepa = variance(potentiala, meanpa) variancepb = variance(potentialb, meanpb) variancepc = variance(potentialc, meanpc) print('The mean potential energy of System A, System B and System C is:\n', meanpa, ',' , meanpb, ', &', meanpc, ', respectively') print('The variance of the potential energy of System A, System B and System C is:\n', variancepa, ',' , variancepb, ', &', variancepc, meanka = mean(kinetica) meankb = mean(kinetichb) meankc = mean(kineticc) varianceka = variance(kinetica, meanka) variancekb = variance(kineticb, meankb) variancekc = variance(kineticc, meankc) print('The mean kinetic energy of System A, System B and System C is:\n', mean_ka, '," , mean_kb, ', &', mean_kc, ', respectively') print('The variance of the potential energy of System A, System B and System C is:\n', varianceka, ',' , variancekb, ', &', variancekc, meanta = mean(totala) meantb = mean(totalb) meantc = mean(totalc) varianceta = variance(totala, meanta) variancetb = variance(totalb, meantb) variancetc = variance(totalc, meantc) print('The mean total energy of System A, System B and System C is:\n', meanta, ',' , meantb, ', &', meantc, ', respectively') print('The variance of the potential energy of System A, System B and System C is:\n', varianceta, ',' , variancetb, ', &' The mean potential energy of System A, System B and System C is: -25.174706 , -25.449154 , & -26.28153 , respectively The variance of the potential energy of System A, System B and System C is: 0.9260010139365807 , 0.8931646626936449 , & 0.24164393198324016 , respectively The mean kinetic energy of System A, System B and System C is: 2.0674589 , 2.6799562 , & 2.337038 , respectively The variance of the potential energy of System A, System B and System C is: 0.9231701041780523 , 0.8838388057403463 , & 0.5577864537934012 , respectively The mean total energy of System A, System B and System C is: -21.868355 , -22.499918 , & -23.820663 , respectively The variance of the potential energy of System A, System B and System C is: 3.755164571986016e-06 , ©.36222304130216854 , & 0.3520218639182965 , respectively #graphs #graphs plt.plot(timea, potentiala, 'mo', label = "Potential Energy with length = 188") plt.plot(timeb, potentialb, 'go', 1label = "Potential Energy with length = 268") plt.axhline(y = meanpa, color = 'm', linestyle = '--', label = 'Average Potential Energy wtih length = 108") plt.axhline(y = meanpb, color = 'g', linestyle = '--', label = 'Average Potential Energy wtih length = 268') plt.ylabel( 'Energy’) plt.xlabel( ' Time') plt.title("Average Potential Energy of Two Systems With Different Lengths") plt.legend() plt.show() plt.plot(timeb, potentialb, 'go', 1label = "Potential Energy with dt = 1 femtosecond") plt.plot(timec, potentialc, ‘co', 1label = "Potential Energy with dt = 3 femtoseconds"”) plt.axhline(y = meanpb, color = 'g', linestyle = '--', label = 'Average Potential Energy dt = 1 femtosecond') plt.axhline(y = meanpc, color = 'c', linestyle = '--', label = 'Average Potential Energy dt = 3 femtoseconds') plt.ylabel('Energy’) plt.xlabel('Time') plt.title("Average Potential Energy of Two Systems With Different dt") plt.legend() plt.show()
Energy Energy plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. plt. Average Potential Energy of Two Systems With Different Lengths —23.5 1 -24.0 A —24.5 1 =25.0 1 -25.5 1 —26.0 -26.5 A - @ Potential Energy with length = 100 ~27.0 - @ Potential Energy with length = 200 -== Average Potential Energy wtih length = 100 275 - -=-= Average Potential Energy wtih length = 200 0.000 0.025 0.050 0.075 0.100 0.125 0.150 0.175 0.200 Time Average Potential Energy of Two Systems With Different dt ® Potential Energy with dt = 1 femtosecond -24.0 @ Potential Energy with dt = 3 femtoseconds ; o =" Average Potential Energy dt = 1 femtosecond —24.5 1 o - Average Potential Energy dt = 3 femtoseconds " : o co L -25.0 - . 9 - . e -9-- - —25.5 1 o e ® o ° o.‘ \ > o0 [ ] -2601 o 4B L __® o rs L o —26.5 o ¢ g i « ' o -27.01 o © -27.5 1 ' 0.0 0.1 0.2 0.3 0.4 0.5 0.6 Time plot(timea, kinetica, 'mo’', label = "Kinetic Energy with length = 1088") plot(timeb, kineticb, 'go’', 1label = "Kinetic Energy with length = 280") axhline(y = meanka, color = 'm’, linestyle = ', label = 'Average Kinetic Energy wtih length = 180") axhline(y = meankb, color = 'g’, linestyle = ', label = 'Average Kinetic Energy wtih length = 2080') ylabel( 'Energy’) xlabel( ' Time') title("Average Kinetic Energy of Two Systems With Different Lengths") legend() show() plot(timeb, kineticb, 'go’', label = "Kinetic Energy with dt = 1 femtosecond") plot(timec, kineticc, 'co’, label = "Kinetic Energy with dt = 3 femtoseconds”) axhline(y = meankb, color = 'g’, linestyle = ', label = 'Average Kinetic Energy dt = 1 femtosecond') axhline(y = meankc, color = 'c’, linestyle = ', label = 'Average Kinetic Energy dt = 3 femtoseconds') ylabel( 'Energy’) xlabel('Time') title("Average Kinetic Energy of Two Systems With Different dt") legend() show() Average Kinetic Energy of Two Systems With Different Lengths 6 4 ® Kinetic Energy with length = 100 @ Kinetic Energy with length = 200 —-=-= Average Kinetic Energy wtih length = 100 51 —-—=- Average Kinetic Energy wtih length = 200 ¢ 4o NN ® 0 9 4 oo © K |} 4 - e = = = e e ry ___‘__. ________ 34 --:-.- --------- I T - ® 9 1‘ :'. ® ° ° o e [ ] .| R 4W: %% i/ 0.000 0.025 0.050 0.075 0.100 0.125 0.150 0.175 0.200 Time
Average Kinetic Energy of Two Systems With Different dt 6 K 0 ® Kinetic Energy with dt = 1 femtosecond o ® Kinetic Energy with dt = 3 femtoseconds 3 —-==- Average Kinetic Energy dt = 1 femtosecond 54 °® ~ =~ Average Kinetic Energy dt = 3 femtoseconds =] 1 e > 2 2‘ [ oL b S 4 "i o ° 34 [ o o0 -~ 4 ';' & """'T*" 2T TTTTTe T @ | Yee oA Ty WY S . LJ ° @ ° o ? e ® o . ° L ° o 2 1 % A o ° o 2 ¢ ¢ "5 ] P ,' 1;‘ 1 - 0.0 0.1 0.2 0.3 0.4 0.5 0.6 Time plt.plot(timea, totala, 'mo’', label = "Total Energy with length = 108") plt.plot(timeb, totalb, ‘go’', 1label = "Total Energy with length = 280") plt.axhline(y = meanta, color = 'm’, linestyle = , label = 'Average Total Energy wtih length = 180") plt.axhline(y = meantb, color = 'g’, linestyle = '--', label = 'Average Total Energy wtih length = 280°) plt.ylabel( ' Energy’) plt.xlabel( ' Time') plt.title("Average Total Energy of Two Systems With Different Lengths") plt.legend() plt.show() plt.plot(timeb, totalb, 'go’', 1label = “Total Energy with dt = 1 femtosecond”) plt.plot(timec, totalc, ‘co’', 1label = "Total Energy with dt = 3 femtoseconds") plt.axhline(y = meantb, color = 'g’, linestyle = , label = 'Average Total Energy dt = 1 femtosecond') plt.axhline(y = meantc, color = '‘c’, linestyle = , label = 'Average Total Energy dt = 3 femtoseconds') plt.ylabel( 'Energy’) plt.xlabel( ' Time') plt.title(“Average Total Energy of Two Systems With Different dt") plt.legend() plt.show() Average Total Energy of Two Systems With Different Lengths -21.5 L ] LS ° ® = | M ---------------------- -22.0 § ° o TR ° f o) L J % 225 +----- - - - - ——————- e e e g e o 5 “d —23.0 b 3 @ Total Energy with length = 100 * _2354{ ® Total Energy with length = 200 # " | === Average Total Energy wtih length = 100 v: o —-==- Average Total Energy wtih length = 200 LY 0.000 0.025 0.050 0.075 0.100 0.125 0.150 0.175 0.200 Time Average Total Energy of Two Systems With Different dt -21.5 ® Total Energy with dt = 1 femtosecond ® Total Energy with dt = 3 femtoseconds -22.0 - -== Average Total Energy dt = 1 femtosecond -~~~ Average Total Energy dt = 3 femtoseconds -22.5 T-gillf---PE-======-mmemmmccccccccccccccccneee F B L —23.01 ® ? ® e o Y 7] [ = w o ° -23.5 c 2 ? ' °%w & . ________ S LY B N, S 200 | I A of —24.5 - .‘. ‘3 -25.0 v. 0.0 0.1 0.2 0.3 0.4 0.5 0.6
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
  • Access to all documents
  • Unlimited textbook solutions
  • 24/7 expert homework help
4. Plot the distance between the two oxygen atoms of each water as a function of time (results dictionary, key r00 ). Do you observe oscillations? Report the mean OO distance along with its variance. # plot distance here parms5 = { steps':700, skip_steps':1, temperature':75., dt': 1 * unit.femtoseconds, ensemble’ :'NVT', simm = simulate_solidwater(parms5) import matplotlib.pyplot as plt from statistics import variance from statistics import mean plt.plot(simm| time'], simm['r00'], plt.xlabel( Time") plt.ylabel("Distance Between Oxygen Atoms") plt.title("00 Distance Vs. Time") ‘co') meann = mean(simm[ r00']) variancee = variance(simm['r00'], meann) 95.9% 671 -27.011353228961152 60.25108854487167 14.3 0:00 96.0% 672 -26.963704297132054 59.3543439387123 14.3 09:00 96.1% 673 -26.903180310698648 59.09031440037103 14.3 09:00 96.3% 674 -26.82788124365246 58.632210030902314 14.3 0:00 96.4% 675 ~26.7444643602465 53.99077464026635 14.3 09:00 96.6% 676 -26.654174921853784 48.70976247347535 14.3 9:00 96.7% 677 -26.557167063625432 47.06353151180957 14.3 09:00 96.9% 678 ~26.45608559883964 43.33141428293991 14.3 0:00 97.0% 679 ~26.356452176136678 40.13511221044484 14.3 0:00 97.1% 680 -26.251678763530915 39.4841266995876 14.3 0:00 97.3% 681 -26.158681853163518 35.02387335961625 14.3 0:00 97.4% 682 -26.06544619123832 34,22479602259779 14.3 0:00 97.6% 683 ~25.989035619336832 30.902450550344 14.3 0:00 97.7% 684 -25.931935225742425 31.47440577095352 14.3 09:00 97.9% 685 -25.894315456638903 30.40155937347819 14.3 9:00 98.0% 686 -25.87972231939517 29.515236800625956 14.3 0:00 98.1% 687 ~25.88098285128318 29.64250812993263 14.3 0:00 98.3% 688 -25.893204998326723 29.21621276500466 14.2 09:00 98.4% 689 -25.91575851740292 29.170626846031837 14.2 0:00 98.6% 690 -25.95080593910337 30.619246561257143 14.3 0:00 98.7% 691 ~25.987118756869144 30.209260641022638 14.3 0:00 98.9% 692 -26.029141073312267 33.58173040782551 14.3 0:00 99.0% 693 -26.0697067761699 35.96702505406139 14.3 09:00 99.1% 694 -26.10512160053018 36.89279196822475 14.3 0:00 99.3% 695 -26.13263802953812 35.54299380954451 14.3 0:00 99.4% 696 ~-26.150026855370793 36.398848890391086 14.3 0:00 99.6% 697 -26.153766135148583 39.060961705438785 14.3 09:00 99.7% 698 -26.14376160700697 38.42060879419322 14.3 9:00 99.9% 699 ~26.113717989771907 40.68048875572083 14.3 0:00 100.0% 700 -26.0670145895953 38.2015248719017 14.3 09:00 00 Distance Vs. Time 0.2825 A w £ g 0.2800 + c = > 0.2775 A x o & o 0.2750 A = - ] [ Y 0.2725 c s 2 0O 0.2700 A 0.2675 T T T T T T T T 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
print("“Mean:", meann ) print("Variance:", variancee) Mean: 0.27603397 Variance: 1.905367485407796e-05 Ans > Yep there are oscillations!