radonpy / radonpy

RadonPy is a Python library to automate physical property calculations for polymer informatics.

Overview

RadonPy is the first open-source Python library for fully automated calculation for a comprehensive set of polymer properties, using all-atom classical MD simulations. For a given polymer repeating unit with its chemical structure, the entire process of the MD simulation can be carried out fully automatically, including molecular modelling, equilibrium and non-equilibrium MD simulations, automatic determination of the completion of equilibration, scheduling of restarts in case of failure to converge, and property calculations in the post-process step. In this release, the library comprises the calculation of 15 properties at the amorphous state.

Requirement

• Python 3.7, 3.8, 3.9, 3.10, 3.11, 3.12
• LAMMPS >= 3Mar20
• rdkit >= 2020.03
• psi4 >= 1.5
• resp
• dftd3
• mdtraj >= 1.9
• scipy
• matplotlib

Installation and usage

User manual and conda packages are currently in preparation.

PyPI package is available, but Psi4 can not be installed by pip install.

PDF file of RadonPy tutorial is available.

Installation for conda (for Psi4 >= 1.8):

1. Create conda environment

conda create -n radonpy python=3.11
conda activate radonpy

2. Installation of requirement packages by conda

conda install -c conda-forge/label/libint_dev -c conda-forge -c psi4 rdkit psi4 resp mdtraj matplotlib

3. Installation of LAMMPS by conda

conda install -c conda-forge lammps

or manually build from source of LAMMPS official site. In this case, the environment variable must be set:

export LAMMPS_EXEC=<Path-to-LAMMPS-binary>

4. Installation of RadonPy

pip install radonpy-pypi

Installation for conda (for Psi4 <= 1.7):

1. Create conda environment

conda create -n radonpy python=3.9
conda activate radonpy

2. Installation of requirement packages by conda

conda install -c psi4 -c conda-forge rdkit psi4 resp mdtraj matplotlib

3. Installation of LAMMPS by conda

conda install -c conda-forge lammps

or manually build from source of LAMMPS official site. In this case, the environment variable must be set:

export LAMMPS_EXEC=<Path-to-LAMMPS-binary>

4. Installation of RadonPy

pip install radonpy-pypi

Installation from PyPI

RadonPy can be also installed by using only pip install. However, this intallation method can not install Psi4.

• Without LAMMPS installation

pip install radonpy-pypi

This is minimal installation of RadonPy. Many functions, such as polymer structure builder, force field assignment, force field descriptor, and tools for polymer informatics, are available, but automated DFT and MD simulations are not available.

• With LAMMPS installation

pip install radonpy-pypi[lammps]

MD simulations are available in this installation, but DFT calculations (conformation search, cherge calculation, and electronic property calculation) are not available.

Features

• Fully automated all-atom classical MD calculation for polymeric materials
  • Conformation search
  • Cherge calculation (RESP, ESP, Mulliken, Lowdin, Gasteiger)
  • Electronic property calculation (HOMO, LUMO, dipole moment, polarizability)
  • Generation of a polymer chain
    • Homopolymer
    • Alternating copolymer
    • Random copolymer
    • Block copolymer
  • Generation of a simulation cell
    • Amorphous
    • Polymer mixture
    • Polymer solution
    • Crystalline polymer
    • Oriented structure
  • Run for equilibration MD
  • Checking archivement of equilibrium
  • Run for non-equilibrium MD (NEMD)
  • Calculation of physical properties from the MD calculation results
    • Thermal conductivity
    • Thermal diffusivity
    • Density
    • Cp
    • Cv
    • Linear expansion coefficient
    • Volumetric expansion coefficient
    • Compressibility
    • Bulk modulus
    • Isentropic compressibility
    • Isentropic bulk modulus
    • Static dielectric constant
    • Refractive index
    • Radius of gyration
    • End-to-end distance
    • Nematic order parameter
  • Using LAMMPS and Psi4 as calculation engines of MD and DFT calculations
• Implementation of add-on like presets to allow for proper and easy execution of polymer MD calculations
  • Equilibration MD
  • Calculation of thermal conductivity with NEMD
• Easy installation
  • Only using open-source software
• Tools for polymer informatics
  • Force field descriptor (How to use)
  • Generator of macrocyclic oligomer for descriptor construction of polymers
  • Full and substruct match function for polymer SMILES
  • Extractor of mainchain in a polymer backbone
  • Monomerization of oligomer SMILES
  • Emulator of polymer classification in PoLyInfo

MD calculated data

• 1070 amorphous polymers

Publications

1. Y. Hayashi, J. Shiomi, J. Morikawa, R. Yoshida, "RadonPy: Automated Physical Property Calculation using All-atom Classical Molecular Dynamics Simulations for Polymer Informatics," npj Comput. Mater., 8:222 (2022) [Link]
2. M. Kusaba, Y. Hayashi, C. Liu, A. Wakiuchi, R. Yoshida, "Representation of materials by kernel mean embedding", Phys. Rev. B, 108:134107 (2023)[Link]

Contributors

• Yoshihiro Hayashi (The Institute of Statistical Mathematics)

Related projects

• XenonPy (Machine learning tools for materials informatics) [Link]
• SMiPoly (Polymerization rule-based virtual polymer generator) [Link]

Acknowledgements

The development of RadonPy was financially supported by the following grants

• Japan Science and Technology Agency (JST) CREST (Grant Number: JPMJCR19I3)
• Ministry of Education, Culture, Sports, Science and Technology (MEXT) as “Program for Promoting Researches on the Supercomputer Fugaku” (Project ID: hp210264)
• The Japan Society for the Promotion of Science (JSPS) as the Grant-in-Aid for Scientific Research (A) (Grant Number: 19H01132)
• JSPS as the Grant-in-Aid for Scientific Research (C) (Grant Number: 22K11949)

The numerical calculations were conducted on the following supercomputer systems

• Fugaku at the RIKEN Center for Computational Science, Kobe, Japan (Project ID: hp210264, hp210213)
• The supercomputer at the Research Center for Computational Science, Okazaki, Japan (Project ID: 21-IMS-C126, 22-IMS-C125, 23-IMS-C113)
• The supercomputer Ohtaka at the Supercomputer Center, the Institute for Solid State Physics, the University of Tokyo, Tokyo, Japan
• The supercomputer TSUBAME3.0 at the Tokyo Institute of Technology, Tokyo, Japan
• The supercomputer ABCI at the National Institute of Advanced Industrial Science and Technology, Tsukuba, Japan

Copyright and licence

©Copyright 2023 The RadonPy developers, all rights reserved. Released under the BSD-3 license.

radonpy's Issues

Calculate FF descriptor for small molecule?

In the present implementation, is it possible to calculate FF descriptor for small molecules (without connection atoms)? If not, will it be expanded to do so?

psi4のマルチプロセスでの動作について

@yhayashi1986
下記フォークにてRISを実装しています
https://github.com/ryuichi-sakashita/RadonPyRIS

Slackにも報告したとおり、
二面角のスキャンはできました
その後の二面角×3の全組み合わせを計算するときにPsi4がプロセスごと落ちてしまいます

別プロセスとすることでエラーをキャッチしようとしましたが、
今度はサブプロセスとしてPsi4を動かすときにPsi4に下記の原因不明のエラーが出るようです

RadonPy info: Psi4 optimization is running...
forrtl: severe (174): SIGSEGV, segmentation fault occurred

変更点：
RISメイン処理：sim/preset/ris.py
sim/psi4_wrapper.pyに　localminima_dihedrals_scan(self, dihedral, step=15, **kwargs)を追加
実行スクリプト：sample_script/ris_test.py　←これを実行したら再現できると思います

お手すきのときに見ていただけたら幸いです

KeyError: 'AtomicCharge' when running equilibrium md functions

Hello,

Thank you for writing this package. It's a very exciting tool!

I have a question about forcefield assignment at the polymer level.

I'm running the following test on polystyrene now following the tutorial pdf:

from radonpy.core import poly, utils
from radonpy.sim import qm
from radonpy.ff.gaff2_mod import GAFF2_mod
from radonpy.sim.preset import eq
import os

cpus = utils.cpu_count()

ps_smiles = 'CC()c1ccccc1'
mol = utils.mol_from_smiles(ps_smiles)

work_dir = 'path/to/work/dir'

ff = GAFF2_mod()
mol, energy = qm.conformation_search(mol,
ff=ff,
work_dir = work_dir
)

qm.assign_charges(mol,
charge='RESP',
work_dir=work_dir
)

repeat_units = 10
atom_nubmer = monomer_atom_num*repeat_units
print(atom_nubmer)

atom_number = 1000
ter = utils.mol_from_smiles('*C')
n = poly.calc_n_from_num_atoms(mol,
atom_number ,
terminal1=ter
)
homopoly = poly.polymerize_rw(mol,
n,
tacticity='atactic'
)
homopoly = poly.terminate_rw(homopoly,
ter)

result = ff.ff_assign(homopoly)

It seems that the forcefield is assigned at the conformation search level. Is there a reason it needs to be reassigned after polymerization? Is there a way to pass the forcefield information from the qm calculations to the polymerize function? More importantly, the partial charges from the charge calculation are not transferred from the monomers to the polymer. This leads to the following key error when you try to run the md equilibration.

n_chains = 10
ac = poly.amorphous_cell(homopoly,
n_chains
)

temp = 300 # K
press = 1.0 # bar
eqmd = eq.EQ21step(ac,
work_dir = work_dir)
ac = eqmd.exec(temp = temp,
press = press,
omp = 1,
gpu = 2,
mpi = utils.cpu_count())

KeyError Traceback (most recent call last)
Cell In[16], line 5
2 press = 1.0 # bar
3 eqmd = eq.EQ21step(ac,
4 work_dir = work_dir)
----> 5 ac = eqmd.exec(temp = temp,
6 press = press,
7 omp = 1,
8 gpu = 2,
9 mpi = utils.cpu_count())

**File anaconda3/envs/radonpy_py311/lib/python3.11/site-packages/radonpy/sim/preset/eq.py:351, in EQ21step.exec(self, confId, f_density, max_temp, temp, press, max_press, step_list, press_ratio, time_step, eq_step, omp, mpi, gpu, intel, opt, **kwargs)
349 utils.MolToPDBFile(self.mol, os.path.join(self.work_dir, self.pdb_file))
350 lmp = lammps.LAMMPS(work_dir=self.work_dir, solver_path=self.solver_path)
--> 351 lmp.make_dat(self.mol, file_name=self.dat_file1, confId=confId)
353 dt1 = datetime.datetime.now()
354 utils.radon_print('Packing simulation (eq1) by LAMMPS is running...', level=1)

File anaconda3/envs/radonpy_py311/lib/python3.11/site-packages/radonpy/sim/lammps.py:333, in LAMMPS.make_dat(self, mol, confId, file_name, dir_name, velocity, temp, drude)
331 dir_name = dir_name if dir_name is not None else self.work_dir
332 dat_file = os.path.join(dir_name, file_name)
--> 333 MolToLAMMPSdata(mol, dat_file, confId=confId, velocity=velocity, temp=temp, drude=drude)
335 return dat_file
...
2755 coord[i][0], coord[i][1], coord[i][2]))
2757 if velocity:
2758 lines.append('')

KeyError: 'AtomicCharge'**

The work around to this is to recalculate the charge after polymerization by running result = ff.ff_assign(homopoly, charge='RESP'), but redoing the qm calculations of partial charges at the polymer level becomes prohibitively expensive as the number of monomers increases. Is there a way to pass the ff and partial charges from the monomer to the polymer mol structure so that these calculations do not have to be redone?

Thanks!
Mike Boyle

Usage of package OPENMP

Hello
Thanks for your repo! it 's a really useful tool for us who don't understand simulation very much
But i came across a situation when trying the tutorial as

''RadonPy warning: OPENMP package is not available. Parallel number of OPENMP is changed to zero.''

I have tired command

'conda install -c conda-forge openmpi '

but it doesn't work, it seemd like i still fail to install some of the necessary packages, so how can i solve this problem?
Looking forward for your reply! Thanks!

About the parameters in sample_script

Hi there!
I am trying to understand the code of sample_script repo.
However, I find it out that the data dictionary has a piece of code os.environ.get('RadonPy_DBID').

Should I have a database preparation before I start the code?
Otherwise, it may report
smi_list = data['smiles_list'].split(',') AttributeError: 'NoneType' object has no attribute 'split'

Hope you may drop me a line.

TypeError in quick_rw when using poly.polymerize_rw with LAMMPS

I find a TypeError when using the poly.polymerize_rw function. The error occurs in the quick_rw function during the invocation of md.add_md():
TypeError: add_md() missing 1 required positional argument: 'step'
I think the issue lies in the quick_rw function, in this line:
md.add_md('nve', step, time_step=time_step, shake=shake, nve_limit=limit, **kwargs)
The step parameter should be passed as a keyword argument to avoid confusion, like so:
md.add_md('nve', step=step, time_step=time_step, shake=shake, nve_limit=limit, **kwargs)
Thank you very much

Adding new atom types

Radonpy is an easy tool for creating MD input. I could generate ff for molecules and polymers. I'm interested in generating a force field for a metal oxide cluster. Is it feasible to load force fields for this system and combine them with polymers/molecules?