diff --git a/CHANGELOG.md b/CHANGELOG.md index 989000d..5c6b0e0 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -19,12 +19,14 @@ The rules for this file: ## [Unreleased] ### Authors +- ianmkenney ### Added ### Fixed ### Changed +- reorganized Sphinx docs (PR #26) ### Removed diff --git a/docs/source/api.rst b/docs/source/api.rst index a3c9557..76ff448 100644 --- a/docs/source/api.rst +++ b/docs/source/api.rst @@ -1,4 +1,20 @@ API Documentation ================= -.. automodule:: waterdynamics.waterdynamics +This section contains an overview of the `waterdynamics` API. + +.. autoclass:: waterdynamics.AngularDistribution + :members: + :inherited-members: + +.. autoclass:: waterdynamics.WaterOrientationalRelaxation + :members: + :inherited-members: + +.. autoclass:: waterdynamics.MeanSquareDisplacement + :members: + :inherited-members: + +.. autoclass:: waterdynamics.SurvivalProbability + :members: + :inherited-members: diff --git a/docs/source/examples.rst b/docs/source/examples.rst new file mode 100644 index 0000000..f3d3d36 --- /dev/null +++ b/docs/source/examples.rst @@ -0,0 +1,229 @@ +Examples +======== + +``waterdynamics`` operates around class based analyses. +Find examples of their use below. +For an explanation of the classes' outputs, see the `outputs `_ page. + +HydrogenBondLifetimes +~~~~~~~~~~~~~~~~~~~~~ + +To analyze hydrogen bond lifetime, use +:meth:`MDAnalysis.analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalysis.lifetime`. + +.. seealso:: + + :mod:`MDAnalysis.analysis.hydrogenbonds.hbond_analysis` + +WaterOrientationalRelaxation +~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +Analyzing water orientational relaxation (WOR) +:class:`~waterdynamics.WaterOrientationalRelaxation`. In this case we are analyzing "how fast" +water molecules are rotating/changing direction. If WOR is very stable we can +assume that water molecules are rotating/changing direction very slow, on the +other hand, if WOR decay very fast, we can assume that water molecules are +rotating/changing direction very fast:: + + import MDAnalysis + from waterdynamics import WaterOrientationalRelaxation as WOR + + u = MDAnalysis.Universe(pdb, trajectory) + select = "byres name OH2 and sphzone 6.0 protein and resid 42" + WOR_analysis = WOR(universe, select, 0, 1000, 20) + WOR_analysis.run() + time = 0 + #now we print the data ready to plot. The first two columns are WOR_OH vs t plot, + #the second two columns are WOR_HH vs t graph and the third two columns are WOR_dip vs t graph + for WOR_OH, WOR_HH, WOR_dip in WOR_analysis.timeseries: + print("{time} {WOR_OH} {time} {WOR_HH} {time} {WOR_dip}".format(time=time, WOR_OH=WOR_OH, WOR_HH=WOR_HH,WOR_dip=WOR_dip)) + time += 1 + + #now, if we want, we can plot our data + plt.figure(1,figsize=(18, 6)) + + #WOR OH + plt.subplot(131) + plt.xlabel('time') + plt.ylabel('WOR') + plt.title('WOR OH') + plt.plot(range(0,time),[column[0] for column in WOR_analysis.timeseries]) + + #WOR HH + plt.subplot(132) + plt.xlabel('time') + plt.ylabel('WOR') + plt.title('WOR HH') + plt.plot(range(0,time),[column[1] for column in WOR_analysis.timeseries]) + + #WOR dip + plt.subplot(133) + plt.xlabel('time') + plt.ylabel('WOR') + plt.title('WOR dip') + plt.plot(range(0,time),[column[2] for column in WOR_analysis.timeseries]) + + plt.show() + +where t0 = 0, tf = 1000 and dtmax = 20. In this way we create 20 windows +timesteps (20 values in the x axis), the first window is created with 1000 +timestep average (1000/1), the second window is created with 500 timestep +average(1000/2), the third window is created with 333 timestep average (1000/3) +and so on. + +AngularDistribution +~~~~~~~~~~~~~~~~~~~ + +Analyzing angular distribution (AD) :class:`AngularDistribution` for OH vector, +HH vector and dipole vector. It returns a line histogram with vector +orientation preference. A straight line in the output plot means no +preferential orientation in water molecules. In this case we are analyzing if +water molecules have some orientational preference, in this way we can see if +water molecules are under an electric field or if they are interacting with +something (residue, protein, etc):: + + import MDAnalysis + from waterdynamics import AngularDistribution as AD + + u = MDAnalysis.Universe(pdb, trajectory) + selection = "byres name OH2 and sphzone 6.0 (protein and (resid 42 or resid 26) )" + bins = 30 + AD_analysis = AD(universe,selection,bins) + AD_analysis.run() + #now we print data ready to graph. The first two columns are P(cos(theta)) vs cos(theta) for OH vector , + #the seconds two columns are P(cos(theta)) vs cos(theta) for HH vector and thirds two columns + #are P(cos(theta)) vs cos(theta) for dipole vector + for bin in range(bins): + print("{AD_analysisOH} {AD_analysisHH} {AD_analysisDip}".format(AD_analysis.graph0=AD_analysis.graph[0][bin], AD_analysis.graph1=AD_analysis.graph[1][bin],AD_analysis.graph2=AD_analysis.graph[2][bin])) + + #and if we want to graph our results + plt.figure(1,figsize=(18, 6)) + + #AD OH + plt.subplot(131) + plt.xlabel('cos theta') + plt.ylabel('P(cos theta)') + plt.title('PDF cos theta for OH') + plt.plot([float(column.split()[0]) for column in AD_analysis.graph[0][:-1]],[float(column.split()[1]) for column in AD_analysis.graph[0][:-1]]) + + #AD HH + plt.subplot(132) + plt.xlabel('cos theta') + plt.ylabel('P(cos theta)') + plt.title('PDF cos theta for HH') + plt.plot([float(column.split()[0]) for column in AD_analysis.graph[1][:-1]],[float(column.split()[1]) for column in AD_analysis.graph[1][:-1]]) + + #AD dip + plt.subplot(133) + plt.xlabel('cos theta') + plt.ylabel('P(cos theta)') + plt.title('PDF cos theta for dipole') + plt.plot([float(column.split()[0]) for column in AD_analysis.graph[2][:-1]],[float(column.split()[1]) for column in AD_analysis.graph[2][:-1]]) + + plt.show() + + +where `P(cos(theta))` is the angular distribution or angular probabilities. + +MeanSquareDisplacement +~~~~~~~~~~~~~~~~~~~~~~ + +Analyzing mean square displacement (MSD) :class:`~waterdynamics.MeanSquareDisplacement` for +water molecules. In this case we are analyzing the average distance that water +molecules travels inside protein in XYZ direction (cylindric zone of radius +11[nm], Zmax 4.0[nm] and Zmin -8.0[nm]). A strong rise mean a fast movement of +water molecules, a weak rise mean slow movement of particles:: + + import MDAnalysis + from waterdynamics import MeanSquareDisplacement as MSD + + u = MDAnalysis.Universe(pdb, trajectory) + select = "byres name OH2 and cyzone 11.0 4.0 -8.0 protein" + MSD_analysis = MSD(universe, select, 0, 1000, 20) + MSD_analysis.run() + #now we print data ready to graph. The graph + #represents MSD vs t + time = 0 + for msd in MSD_analysis.timeseries: + print("{time} {msd}".format(time=time, msd=msd)) + time += 1 + + #Plot + plt.xlabel('time') + plt.ylabel('MSD') + plt.title('MSD') + plt.plot(range(0,time),MSD_analysis.timeseries) + plt.show() + +.. _SP-examples: + +SurvivalProbability +~~~~~~~~~~~~~~~~~~~ + +Analyzing survival probability (SP) :class:`SurvivalProbability` of molecules. +In this case we are analyzing how long water molecules remain in a +sphere of radius 12.3 centered in the geometrical center of resid 42 and 26. +A slow decay of SP means a long permanence time of water molecules in +the zone, on the other hand, a fast decay means a short permanence time:: + + import MDAnalysis + from waterdynamics import SurvivalProbability as SP + import matplotlib.pyplot as plt + + universe = MDAnalysis.Universe(pdb, trajectory) + select = "byres name OH2 and sphzone 12.3 (resid 42 or resid 26) " + sp = SP(universe, select, verbose=True) + sp.run(start=0, stop=101, tau_max=20) + tau_timeseries = sp.tau_timeseries + sp_timeseries = sp.sp_timeseries + + # print in console + for tau, sp in zip(tau_timeseries, sp_timeseries): + print("{time} {sp}".format(time=tau, sp=sp)) + + # plot + plt.xlabel('Time') + plt.ylabel('SP') + plt.title('Survival Probability') + plt.plot(tau_timeseries, sp_timeseries) + plt.show() + +One should note that the `stop` keyword as used in the above example has an +`exclusive` behaviour, i.e. here the final frame used will be 100 not 101. +This behaviour is aligned with :class:`AnalysisBase` but currently differs from +other :mod:`waterdynamics` classes, which all exhibit +`inclusive` behaviour for their final frame selections. + +Another example applies to the situation where you work with many different "residues". +Here we calculate the SP of a potassium ion around each lipid in a membrane and +average the results. In this example, if the SP analysis were run without treating each lipid +separately, potassium ions may hop from one lipid to another and still be counted as remaining +in the specified region. That is, the survival probability of the potassium ion around the +entire membrane will be calculated. + +Note, for this example, it is advisable to use `Universe(in_memory=True)` to ensure that the +simulation is not being reloaded into memory for each lipid:: + + import MDAnalysis as mda + from waterdynamics import SurvivalProbability as SP + import numpy as np + + u = mda.Universe("md.gro", "md100ns.xtc", in_memory=True) + lipids = u.select_atoms('resname LIPIDS') + joined_sp_timeseries = [[] for _ in range(20)] + for lipid in lipids.residues: + print("Lipid ID: %d" % lipid.resid) + + select = "resname POTASSIUM and around 3.5 (resid %d and name O13 O14) " % lipid.resid + sp = SP(u, select, verbose=True) + sp.run(tau_max=20) + + # Raw SP points for each tau: + for sps, new_sps in zip(joined_sp_timeseries, sp.sp_timeseries_data): + sps.extend(new_sps) + + # average all SP datapoints + sp_data = [np.mean(sp) for sp in joined_sp_timeseries] + + for tau, sp in zip(range(1, tau_max + 1), sp_data): + print("{time} {sp}".format(time=tau, sp=sp)) diff --git a/docs/source/getting_started.rst b/docs/source/getting_started.rst index 299c651..ed49134 100644 --- a/docs/source/getting_started.rst +++ b/docs/source/getting_started.rst @@ -1,23 +1,46 @@ Getting Started =============== -This page details how to get started with waterdynamics. +Install the ``waterdynamics`` package by any of the methods described below: Installation ------------ -PyPi -~~~~ +From conda-forge +~~~~~~~~~~~~~~~~ -The waterdyanmics package can be installed from PyPI with: +Install `waterdynamics from the conda-forge channel `_ +by first adding ``conda-forge`` to you channels with: -.. code-block:: bash +.. code-block:: sh + + conda config --add channels conda-forge + conda config --set channel_priority strict + +Once the ``conda-forge`` channel has been enabled, ``waterdynamics`` can be installed with :program:`conda`: + +.. code-block:: sh + + conda install waterdynamics + +or with :program:`mamba` (if you have it installed): + +.. code-block:: sh + + mamba install waterdynamics + +From the Python Package Index +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +The waterdynamics package can be installed from PyPI with: + +.. code-block:: sh pip install waterdynamics If you would like to test your installation, use the `test` optional dependencies and run the tests: -.. code-block:: bash +.. code-block:: sh pip install "waterdynamics[test]" pytest --pyargs waterdynamics.tests @@ -27,7 +50,7 @@ From source If you want the latest (non-release) version of waterdynamics, you should install the package from source. -.. code-block:: bash +.. code-block:: sh git clone https://github.com/MDAnalysis/waterdynamics cd waterdynamics @@ -35,7 +58,7 @@ If you want the latest (non-release) version of waterdynamics, you should instal If you want to run tests on your installation, install the `test` dependencies and run the tests: -.. code-block:: bash +.. code-block:: sh pip install ".[test]" pytest --pyargs waterdynamics.tests diff --git a/docs/source/index.rst b/docs/source/index.rst index 84cc4f2..b95fc16 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -6,13 +6,19 @@ Welcome to waterdynamics's documentation! ========================================================= -This package provides classes and functions for the analysis of water dynamics and water interactions with other molecules in molecular simulations. The functions in this package are: water orientational relaxation, hydrogen bond lifetimes, angular distribution, mean square displacement, and survival probability. +This package provides classes and functions for the analysis of water dynamics and water interactions with other molecules in molecular simulations. The functions in this package are: water orientational relaxation (WOR) :footcite:p:`Yeh1999`, hydrogen bond lifetimes (HBL) :footcite:p:`Rapaport1983`, angular distribution (AD) :footcite:p:`Grigera1995`, mean square displacement (MSD) :footcite:p:`Brodka1994`, and survival probability (SP) :footcite:p:`Liu2004`. + +For more information about this type of analysis please refer to +:footcite:p:`ArayaSecchi2014` (water in a protein cavity) and +:footcite:p:`Milischuk2011` (water in a nanopore). .. toctree:: :maxdepth: 2 :caption: Contents: getting_started + examples + outputs api @@ -23,3 +29,7 @@ Indices and tables * :ref:`genindex` * :ref:`modindex` * :ref:`search` + +.. rubric:: References + +.. footbibliography:: diff --git a/docs/source/outputs.rst b/docs/source/outputs.rst new file mode 100644 index 0000000..ef8df7c --- /dev/null +++ b/docs/source/outputs.rst @@ -0,0 +1,63 @@ +Outputs +======= + +Here we outline the output data structures for each of the packages classes. + +WaterOrientationalRelaxation +~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +:class:`~waterdynamics.WaterOrientationalRelaxation` (WOR) data is returned per window timestep, +which is stored in :attr:`WaterOrientationalRelaxation.timeseries`:: + + results = [ + [ # time t0 + , , + ], + [ # time t1 + , , + ], + ... + ] + +AngularDistribution +~~~~~~~~~~~~~~~~~~~ + +:class:`~waterdynamics.AngularDistribution` (AD) data is returned per vector, which is stored in +:attr:`AngularDistribution.graph`. In fact, AngularDistribution returns a +histogram:: + + results = [ + [ # OH vector values + # the values are order in this way: + , , ... + ], + [ # HH vector values + , , ... + ], + [ # dip vector values + , , ... + ], + ] + +MeanSquareDisplacement +~~~~~~~~~~~~~~~~~~~~~~ + +:class:`~waterdynamics.MeanSquareDisplacement` (MSD) data is returned in a list, which each element +represents a MSD value in its respective window timestep. Data is stored in +:attr:`MeanSquareDisplacement.timeseries`:: + + results = [ + #MSD values orders by window timestep + , , ... + ] + +SurvivalProbability +~~~~~~~~~~~~~~~~~~~ + +:class:`~waterdynamics.SurvivalProbability` (SP) computes two lists: a list of taus (:attr:`SurvivalProbability.tau_timeseries`) and a list of +the corresponding survival probabilities (:attr:`SurvivalProbability.sp_timeseries`). + + results = [ tau1, tau2, ..., tau_n ], [ sp_tau1, sp_tau2, ..., sp_tau_n] + +Additionally, a list :attr:`SurvivalProbability.sp_timeseries_data`, is provided which contains +a list of all SPs calculated for each tau. This can be used to compute the distribution or time dependence of SP, etc. diff --git a/waterdynamics/waterdynamics.py b/waterdynamics/waterdynamics.py index 033fd6f..4aa0e39 100644 --- a/waterdynamics/waterdynamics.py +++ b/waterdynamics/waterdynamics.py @@ -35,320 +35,6 @@ :footcite:p:`ArayaSecchi2014` (water in a protein cavity) and :footcite:p:`Milischuk2011` (water in a nanopore). - - -Example use of the analysis classes ------------------------------------ - -HydrogenBondLifetimes -~~~~~~~~~~~~~~~~~~~~~ - -To analyse hydrogen bond lifetime, use -:meth:`MDAnalysis.analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalysis.lifetime`. - -See Also --------- -:mod:`MDAnalysis.analysis.hydrogenbonds.hbond_analysis` - - -WaterOrientationalRelaxation -~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -Analyzing water orientational relaxation (WOR) -:class:`WaterOrientationalRelaxation`. In this case we are analyzing "how fast" -water molecules are rotating/changing direction. If WOR is very stable we can -assume that water molecules are rotating/changing direction very slow, on the -other hand, if WOR decay very fast, we can assume that water molecules are -rotating/changing direction very fast:: - - import MDAnalysis - from waterdynamics import WaterOrientationalRelaxation as WOR - - u = MDAnalysis.Universe(pdb, trajectory) - select = "byres name OH2 and sphzone 6.0 protein and resid 42" - WOR_analysis = WOR(universe, select, 0, 1000, 20) - WOR_analysis.run() - time = 0 - #now we print the data ready to plot. The first two columns are WOR_OH vs t plot, - #the second two columns are WOR_HH vs t graph and the third two columns are WOR_dip vs t graph - for WOR_OH, WOR_HH, WOR_dip in WOR_analysis.timeseries: - print("{time} {WOR_OH} {time} {WOR_HH} {time} {WOR_dip}".format(time=time, WOR_OH=WOR_OH, WOR_HH=WOR_HH,WOR_dip=WOR_dip)) - time += 1 - - #now, if we want, we can plot our data - plt.figure(1,figsize=(18, 6)) - - #WOR OH - plt.subplot(131) - plt.xlabel('time') - plt.ylabel('WOR') - plt.title('WOR OH') - plt.plot(range(0,time),[column[0] for column in WOR_analysis.timeseries]) - - #WOR HH - plt.subplot(132) - plt.xlabel('time') - plt.ylabel('WOR') - plt.title('WOR HH') - plt.plot(range(0,time),[column[1] for column in WOR_analysis.timeseries]) - - #WOR dip - plt.subplot(133) - plt.xlabel('time') - plt.ylabel('WOR') - plt.title('WOR dip') - plt.plot(range(0,time),[column[2] for column in WOR_analysis.timeseries]) - - plt.show() - -where t0 = 0, tf = 1000 and dtmax = 20. In this way we create 20 windows -timesteps (20 values in the x axis), the first window is created with 1000 -timestep average (1000/1), the second window is created with 500 timestep -average(1000/2), the third window is created with 333 timestep average (1000/3) -and so on. - -AngularDistribution -~~~~~~~~~~~~~~~~~~~ - -Analyzing angular distribution (AD) :class:`AngularDistribution` for OH vector, -HH vector and dipole vector. It returns a line histogram with vector -orientation preference. A straight line in the output plot means no -preferential orientation in water molecules. In this case we are analyzing if -water molecules have some orientational preference, in this way we can see if -water molecules are under an electric field or if they are interacting with -something (residue, protein, etc):: - - import MDAnalysis - from waterdynamics import AngularDistribution as AD - - u = MDAnalysis.Universe(pdb, trajectory) - selection = "byres name OH2 and sphzone 6.0 (protein and (resid 42 or resid 26) )" - bins = 30 - AD_analysis = AD(universe,selection,bins) - AD_analysis.run() - #now we print data ready to graph. The first two columns are P(cos(theta)) vs cos(theta) for OH vector , - #the seconds two columns are P(cos(theta)) vs cos(theta) for HH vector and thirds two columns - #are P(cos(theta)) vs cos(theta) for dipole vector - for bin in range(bins): - print("{AD_analysisOH} {AD_analysisHH} {AD_analysisDip}".format(AD_analysis.graph0=AD_analysis.graph[0][bin], AD_analysis.graph1=AD_analysis.graph[1][bin],AD_analysis.graph2=AD_analysis.graph[2][bin])) - - #and if we want to graph our results - plt.figure(1,figsize=(18, 6)) - - #AD OH - plt.subplot(131) - plt.xlabel('cos theta') - plt.ylabel('P(cos theta)') - plt.title('PDF cos theta for OH') - plt.plot([float(column.split()[0]) for column in AD_analysis.graph[0][:-1]],[float(column.split()[1]) for column in AD_analysis.graph[0][:-1]]) - - #AD HH - plt.subplot(132) - plt.xlabel('cos theta') - plt.ylabel('P(cos theta)') - plt.title('PDF cos theta for HH') - plt.plot([float(column.split()[0]) for column in AD_analysis.graph[1][:-1]],[float(column.split()[1]) for column in AD_analysis.graph[1][:-1]]) - - #AD dip - plt.subplot(133) - plt.xlabel('cos theta') - plt.ylabel('P(cos theta)') - plt.title('PDF cos theta for dipole') - plt.plot([float(column.split()[0]) for column in AD_analysis.graph[2][:-1]],[float(column.split()[1]) for column in AD_analysis.graph[2][:-1]]) - - plt.show() - - -where `P(cos(theta))` is the angular distribution or angular probabilities. - - -MeanSquareDisplacement -~~~~~~~~~~~~~~~~~~~~~~ - -Analyzing mean square displacement (MSD) :class:`MeanSquareDisplacement` for -water molecules. In this case we are analyzing the average distance that water -molecules travels inside protein in XYZ direction (cylindric zone of radius -11[nm], Zmax 4.0[nm] and Zmin -8.0[nm]). A strong rise mean a fast movement of -water molecules, a weak rise mean slow movement of particles:: - - import MDAnalysis - from waterdynamics import MeanSquareDisplacement as MSD - - u = MDAnalysis.Universe(pdb, trajectory) - select = "byres name OH2 and cyzone 11.0 4.0 -8.0 protein" - MSD_analysis = MSD(universe, select, 0, 1000, 20) - MSD_analysis.run() - #now we print data ready to graph. The graph - #represents MSD vs t - time = 0 - for msd in MSD_analysis.timeseries: - print("{time} {msd}".format(time=time, msd=msd)) - time += 1 - - #Plot - plt.xlabel('time') - plt.ylabel('MSD') - plt.title('MSD') - plt.plot(range(0,time),MSD_analysis.timeseries) - plt.show() - - -.. _SP-examples: - -SurvivalProbability -~~~~~~~~~~~~~~~~~~~ - -Analyzing survival probability (SP) :class:`SurvivalProbability` of molecules. -In this case we are analyzing how long water molecules remain in a -sphere of radius 12.3 centered in the geometrical center of resid 42 and 26. -A slow decay of SP means a long permanence time of water molecules in -the zone, on the other hand, a fast decay means a short permanence time:: - - import MDAnalysis - from waterdynamics import SurvivalProbability as SP - import matplotlib.pyplot as plt - - universe = MDAnalysis.Universe(pdb, trajectory) - select = "byres name OH2 and sphzone 12.3 (resid 42 or resid 26) " - sp = SP(universe, select, verbose=True) - sp.run(start=0, stop=101, tau_max=20) - tau_timeseries = sp.tau_timeseries - sp_timeseries = sp.sp_timeseries - - # print in console - for tau, sp in zip(tau_timeseries, sp_timeseries): - print("{time} {sp}".format(time=tau, sp=sp)) - - # plot - plt.xlabel('Time') - plt.ylabel('SP') - plt.title('Survival Probability') - plt.plot(tau_timeseries, sp_timeseries) - plt.show() - -One should note that the `stop` keyword as used in the above example has an -`exclusive` behaviour, i.e. here the final frame used will be 100 not 101. -This behaviour is aligned with :class:`AnalysisBase` but currently differs from -other :mod:`waterdynamics` classes, which all exhibit -`inclusive` behaviour for their final frame selections. - -Another example applies to the situation where you work with many different "residues". -Here we calculate the SP of a potassium ion around each lipid in a membrane and -average the results. In this example, if the SP analysis were run without treating each lipid -separately, potassium ions may hop from one lipid to another and still be counted as remaining -in the specified region. That is, the survival probability of the potassium ion around the -entire membrane will be calculated. - -Note, for this example, it is advisable to use `Universe(in_memory=True)` to ensure that the -simulation is not being reloaded into memory for each lipid:: - - import MDAnalysis as mda - from waterdynamics import SurvivalProbability as SP - import numpy as np - - u = mda.Universe("md.gro", "md100ns.xtc", in_memory=True) - lipids = u.select_atoms('resname LIPIDS') - joined_sp_timeseries = [[] for _ in range(20)] - for lipid in lipids.residues: - print("Lipid ID: %d" % lipid.resid) - - select = "resname POTASSIUM and around 3.5 (resid %d and name O13 O14) " % lipid.resid - sp = SP(u, select, verbose=True) - sp.run(tau_max=20) - - # Raw SP points for each tau: - for sps, new_sps in zip(joined_sp_timeseries, sp.sp_timeseries_data): - sps.extend(new_sps) - - # average all SP datapoints - sp_data = [np.mean(sp) for sp in joined_sp_timeseries] - - for tau, sp in zip(range(1, tau_max + 1), sp_data): - print("{time} {sp}".format(time=tau, sp=sp)) - -.. _Output: - -Output ------- - -WaterOrientationalRelaxation -~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -Water orientational relaxation (WOR) data is returned per window timestep, -which is stored in :attr:`WaterOrientationalRelaxation.timeseries`:: - - results = [ - [ # time t0 - , , - ], - [ # time t1 - , , - ], - ... - ] - -AngularDistribution -~~~~~~~~~~~~~~~~~~~ - -Angular distribution (AD) data is returned per vector, which is stored in -:attr:`AngularDistribution.graph`. In fact, AngularDistribution returns a -histogram:: - - results = [ - [ # OH vector values - # the values are order in this way: - , , ... - ], - [ # HH vector values - , , ... - ], - [ # dip vector values - , , ... - ], - ] - -MeanSquareDisplacement -~~~~~~~~~~~~~~~~~~~~~~ - -Mean Square Displacement (MSD) data is returned in a list, which each element -represents a MSD value in its respective window timestep. Data is stored in -:attr:`MeanSquareDisplacement.timeseries`:: - - results = [ - #MSD values orders by window timestep - , , ... - ] - -SurvivalProbability -~~~~~~~~~~~~~~~~~~~ - -Survival Probability (SP) computes two lists: a list of taus (:attr:`SurvivalProbability.tau_timeseries`) and a list of -the corresponding survival probabilities (:attr:`SurvivalProbability.sp_timeseries`). - - results = [ tau1, tau2, ..., tau_n ], [ sp_tau1, sp_tau2, ..., sp_tau_n] - -Additionally, a list :attr:`SurvivalProbability.sp_timeseries_data`, is provided which contains -a list of all SPs calculated for each tau. This can be used to compute the distribution or time dependence of SP, etc. - -Classes --------- - -.. autoclass:: WaterOrientationalRelaxation - :members: - :inherited-members: - -.. autoclass:: AngularDistribution - :members: - :inherited-members: - -.. autoclass:: MeanSquareDisplacement - :members: - :inherited-members: - -.. autoclass:: SurvivalProbability - :members: - :inherited-members: - .. rubric:: References .. footbibliography:: @@ -562,8 +248,8 @@ class AngularDistribution(object): vector, HH vector or dipolar vector of water molecules and a vector :math:`\hat n` parallel to chosen axis (z is the default value). The cosine is define as :math:`\cos \theta = \hat u \cdot \hat n`, where :math:`\hat - u` is OH, HH or dipole vector. It creates a histogram and returns a list - of lists, see Output_. The AD is also know as Angular Probability (AP). + u` is OH, HH or dipole vector. The AD is also know as Angular + Probability (AP). Parameters @@ -870,7 +556,7 @@ class SurvivalProbability(object): select : str Selection string; any selection is allowed. With this selection you define the region/zone where to analyze, e.g.: "resname SOL and around 5 - (resid 10)". See `SP-examples`_. + (resid 10)". verbose : Boolean, optional When True, prints progress and comments to the console.