pyPRISM¶
pyPRISM is a Python-based, open-source framework for conducting Polymer Reference Interaction Site Model (PRISM) theory calculations. This framework aims to simplify PRISM-based studies by providing a user-friendly scripting interface for setting up and numerically solving the PRISM equations.
PRISM theory describes the equilibrium spatial-correlations of liquid-like polymer systems including melts, blends, solutions, block copolymers, ionomers, liquid crystal forming polymers and nanocomposites. Using PRISM theory, one can calculate thermodynamic (e.g., second virial coefficients, Flory-Huggins \(\chi\) interaction parameters, potentials of mean force) and structural (e.g., pair correlation functions, structure factors) information for these macromolecular materials. See the Frequently Asked Questions section for examples of systems and calculations that are available to PRISM theory.
pyPRISM provides data structures, functions, and classes that streamline PRISM calculations, allowing pyPRISM to be extended for use in other tasks such as the coarse-graining of atomistic simulation force-fields or the modeling of experimental scattering data. The goal of this framework is to reduce the barrier to correctly and appropriately using PRISM theory and to provide a platform for rapid calculations of the structure and thermodynamics of polymeric fluids and nanocomposites.
Citations¶
If you use pyPRISM in your work, we ask that you please cite both of the following articles
- Martin, T.B.; Gartner, T.E. III; Jones, R.L.; Snyder, C.R.; Jayaraman, A.; pyPRISM: A Computational Tool for Liquid State Theory Calculations of Macromolecular Materials, Macromolecules, 2018, 51 (8), p2906-2922 [link]
- Schweizer, K.S.; Curro, J.G.; Integral Equation Theory of the Structure of Polymer Melts, Physical Review Letters, 1987, 58 (3), p246-249 doi:10.1103/PhysRevLett.58.246 [link]
pyPRISM Example¶
Below is an example python script where we use pyPRISM to calculate the pair correlation functions for a nanocomposite (polymer + particle) system with attractive polymer-particle interactions. Below the script is a plot of the pair correlation functions from this calculation. See Quickstart Guide for a more detailed discussion of this example.
import pyPRISM
sys = pyPRISM.System(['particle','polymer'],kT=1.0)
sys.domain = pyPRISM.Domain(dr=0.01,length=4096)
sys.density['polymer'] = 0.75
sys.density['particle'] = 6e-6
sys.diameter['polymer'] = 1.0
sys.diameter['particle'] = 5.0
sys.omega['polymer','polymer'] = pyPRISM.omega.FreelyJointedChain(length=100,l=4.0/3.0)
sys.omega['polymer','particle'] = pyPRISM.omega.InterMolecular()
sys.omega['particle','particle'] = pyPRISM.omega.SingleSite()
sys.potential['polymer','polymer'] = pyPRISM.potential.HardSphere()
sys.potential['polymer','particle'] = pyPRISM.potential.Exponential(alpha=0.5,epsilon=1.0)
sys.potential['particle','particle'] = pyPRISM.potential.HardSphere()
sys.closure['polymer','polymer'] = pyPRISM.closure.PercusYevick()
sys.closure['polymer','particle'] = pyPRISM.closure.PercusYevick()
sys.closure['particle','particle'] = pyPRISM.closure.HyperNettedChain()
PRISM = sys.solve()
pcf = pyPRISM.calculate.prism.pair_correlation(PRISM)
External Resources¶
Source Code Repository | |
Question/Issue Tracker | |
Interactive Binder Tutorial | |
Anaconda Cloud | |
Python Package Index |
Table of Contents¶
- Frequently Asked Questions
- What systems can be studied with PRISM?
- What thermodynamic and structural quantities can PRISM calculate?
- For what systems is PRISM theory not applicable for?
- What are the benefits of using PRISM over other simulation or theory methods?
- Why can’t I import pyPRISM?
- Can pyPRISM be used outside of Jupyter?
- Can pyPRISM handle anisotropic systems?
- How do I set up my specific system?
- What do I do if the solver isn’t converging?
- Why doesn’t pyPRISM doesn’t have the feature I need?
- How to file a bug report, suggest a feature or ask a question about pyPRISM?
- Self-Consistent PRISM Method
- Convergence Tips