Skip to content
Snippets Groups Projects
Commit df9e4aa3 authored by Bienchen's avatar Bienchen
Browse files

Merge branch 'mm' into develop

parents 7c55ff7d 851d6e1a
No related branches found
No related tags found
No related merge requests found
Showing
with 1680 additions and 29 deletions
...@@ -72,3 +72,4 @@ build.ninja ...@@ -72,3 +72,4 @@ build.ninja
modules/gui/src/dngr.qrc.depends modules/gui/src/dngr.qrc.depends
/compile_commands.json /compile_commands.json
/modules/bindings/tests/formatdb.log /modules/bindings/tests/formatdb.log
/modules/mol/mm/src/settings.hh
...@@ -51,6 +51,7 @@ option(COMPILE_TESTS "whether unit tests should be compiled by default" OFF) ...@@ -51,6 +51,7 @@ option(COMPILE_TESTS "whether unit tests should be compiled by default" OFF)
option(ENABLE_STATIC "whether static libraries should be compiled" OFF) option(ENABLE_STATIC "whether static libraries should be compiled" OFF)
option(UBUNTU_LAYOUT "whether Debian/Ubuntu's lib and libexec directory layout should be used" OFF) option(UBUNTU_LAYOUT "whether Debian/Ubuntu's lib and libexec directory layout should be used" OFF)
option(HIDDEN_VISIBILITY "on gcc, use -fvisibility=hidden" OFF) option(HIDDEN_VISIBILITY "on gcc, use -fvisibility=hidden" OFF)
option(ENABLE_MM "whether openmm support is added" OFF)
if (CXX) if (CXX)
set(CMAKE_CXX_COMPILER ${CXX}) set(CMAKE_CXX_COMPILER ${CXX})
...@@ -140,6 +141,12 @@ if (PROFILE) ...@@ -140,6 +141,12 @@ if (PROFILE)
else() else()
set(_PROFILE OFF) set(_PROFILE OFF)
endif() endif()
if(ENABLE_MM)
set(_OPENMM ON)
else()
set(_OPENMM OFF)
endif()
if (UBUNTU_LAYOUT) if (UBUNTU_LAYOUT)
set(_UBUNTU_LAYOUT ON) set(_UBUNTU_LAYOUT ON)
...@@ -147,8 +154,6 @@ else() ...@@ -147,8 +154,6 @@ else()
set(_UBUNTU_LAYOUT OFF) set(_UBUNTU_LAYOUT OFF)
endif() endif()
add_definitions(-DEIGEN2_SUPPORT)
if (COMPOUND_LIB) if (COMPOUND_LIB)
set(_COMP_LIB "${COMPOUND_LIB}") set(_COMP_LIB "${COMPOUND_LIB}")
if (NOT IS_ABSOLUTE "${COMPOUND_LIB}") if (NOT IS_ABSOLUTE "${COMPOUND_LIB}")
...@@ -227,13 +232,20 @@ else() ...@@ -227,13 +232,20 @@ else()
set (PNG_LIBRARY "") set (PNG_LIBRARY "")
endif() endif()
find_package(Eigen 2.0.0 REQUIRED) find_package(Eigen 3.2.0 REQUIRED)
find_package(Python 2.4 REQUIRED) find_package(Python 2.4 REQUIRED)
if(USE_NUMPY) if(USE_NUMPY)
find_package(Numpy REQUIRED) find_package(Numpy REQUIRED)
endif() endif()
if(ENABLE_MM)
find_package(OpenMM REQUIRED)
set(_OPENMM_PLUGINS "${OPEN_MM_PLUGIN_DIR}")
else(ENABLE_MM)
set(_OPENMM_PLUGINS "NONE")
endif(ENABLE_MM)
if (ENABLE_IMG) if (ENABLE_IMG)
find_package(FFTW REQUIRED) find_package(FFTW REQUIRED)
find_package(TIFF REQUIRED) find_package(TIFF REQUIRED)
...@@ -287,9 +299,10 @@ endif() ...@@ -287,9 +299,10 @@ endif()
# basic environment # basic environment
include_directories(${Boost_INCLUDE_DIRS} include_directories(${Boost_INCLUDE_DIRS}
${FFTW_INCLUDE_DIRS} ${FFTW_INCLUDE_DIRS}
${EIGEN2_INCLUDE_DIRS} ${EIGEN3_INCLUDE_DIRS}
${TIFF_INCLUDE_DIR} ${TIFF_INCLUDE_DIR}
${PNG_INCLUDE_DIRS} ${PNG_INCLUDE_DIRS}
${OPEN_MM_INCLUDE_DIRS}
) )
if(USE_NUMPY) if(USE_NUMPY)
include_directories(${PYTHON_NUMPY_INCLUDE_DIR}) include_directories(${PYTHON_NUMPY_INCLUDE_DIR})
...@@ -344,6 +357,8 @@ message(STATUS ...@@ -344,6 +357,8 @@ message(STATUS
" Shader support (-DUSE_SHADER) : ${_SHADER}\n" " Shader support (-DUSE_SHADER) : ${_SHADER}\n"
" Numpy support (-DUSE_NUMPY) : ${_NUMPY}\n" " Numpy support (-DUSE_NUMPY) : ${_NUMPY}\n"
" SpaceNav Device support (-DENABLE_SPNAV) : ${_SPNAV}\n" " SpaceNav Device support (-DENABLE_SPNAV) : ${_SPNAV}\n"
" OpenMM support (-DENABLE_MM) : ${_OPENMM}\n"
" OpenMM plugins (-DOPEN_MM_PLUGIN_DIR) : ${_OPENMM_PLUGINS}\n"
" Optimize (-DOPTIMIZE) : ${_OPT}\n" " Optimize (-DOPTIMIZE) : ${_OPT}\n"
" Profiling support (-DPROFILE) : ${_PROFILE}\n" " Profiling support (-DPROFILE) : ${_PROFILE}\n"
" Double Precision (-DUSE_DOUBLE_PRECISION) : ${_DOUBLE_PREC}\n" " Double Precision (-DUSE_DOUBLE_PRECISION) : ${_DOUBLE_PREC}\n"
......
set(FFTW_INCLUDE_PATH "/import/bc2/soft/app/FFTW/current/Linux/include/" CACHE PATH "") set(FFTW_INCLUDE_PATH "/scicore/soft/apps/FFTW/3.3.3-gompi-1.4.10/include/" CACHE PATH "")
set(FFTW_LIBRARIES "/import/bc2/soft/app/FFTW/current/Linux/lib/libfftw3f.so" CACHE PATH "") set(FFTW_LIBRARY "/scicore/soft/apps/FFTW/3.3.3-gompi-1.4.10/lib/libfftw3f.a" CACHE PATH "")
set(Boost_COMPILER "-gcc41" CACHE PATH "") set(PYTHON_ROOT "/scicore/soft/apps/Python/2.7.5-goolf-1.4.10" CACHE PATH "")
set(PYTHON_ROOT "/import/bc2/soft/app/Python/current/Linux/" CACHE PATH "") set(EIGEN3_INCLUDE_DIR "/scicore/soft/apps/Eigen/3.2.1-goolf-1.4.10/include/Eigen" CACHE PATH "")
set(EIGEN2_INCLUDE_DIR "/import/bc2/soft/app/eigen/current/Linux/include/eigen2" CACHE PATH "") set(BOOST_ROOT "/import/bc2/apps/Boost/1.53.0-goolf-1.4.10-Python-2.7.5" CACHE PATH "")
set(BOOST_ROOT "/import/bc2/soft/app/boost/current/Linux" CACHE PATH "") set(QT_QMAKE_EXECUTABLE "/usr/bin/qmake-qt4" CACHE PATH "")
set(QT_QMAKE_EXECUTABLE "/import/bc2/soft/app/Qt/current/Linux/bin/qmake" CACHE PATH "")
set(COMPOUND_LIB "/import/bc2/home/schwede/GROUP/OpenStructure/ChemLib/compounds.chemlib" CACHE PATH "") set(COMPOUND_LIB "/import/bc2/home/schwede/GROUP/OpenStructure/ChemLib/compounds.chemlib" CACHE PATH "")
set (CMAKE_BUILD_TYPE "Release" CACHE PATH "") set (CMAKE_BUILD_TYPE "Release" CACHE PATH "")
set(FFTW_INCLUDE_PATH "/import/bc2/soft/app/FFTW/current/Linux/include/" CACHE PATH "") set(FFTW_INCLUDE_PATH "/scicore/soft/apps/FFTW/3.3.3-gompi-1.4.10/include/" CACHE PATH "")
set(FFTW_LIBRARIES "/import/bc2/soft/app/FFTW/current/Linux/lib/libfftw3f.a" CACHE PATH "") set(FFTW_LIBRARY "/scicore/soft/apps/FFTW/3.3.3-gompi-1.4.10/lib/libfftw3f.a" CACHE PATH "")
set(DL_LIBRARIES "/usr/lib64/libdl.a" CACHE PATH "") set(DL_LIBRARIES "/usr/lib64/libdl.a" CACHE PATH "")
set(PTHREAD_LIBRARIES "/usr/lib64/libpthread.a" CACHE PATH "") set(PTHREAD_LIBRARIES "/usr/lib64/libpthread.a" CACHE PATH "")
set(ZLIB_LIBRARY "/usr/lib64/libz.a" CACHE PATH "") set(ZLIB_LIBRARY "/scicore/soft/apps/zlib/1.2.8-goolf-1.4.10/lib/libz.a" CACHE PATH "")
set(Boost_COMPILER "-gcc41" CACHE PATH "") set(PYTHON_ROOT "/scicore/soft/apps/Python/2.7.5-goolf-1.4.10" CACHE PATH "")
set(PYTHON_ROOT "/import/bc2/soft/app/Python/current/Linux/" CACHE PATH "") set(EIGEN3_INCLUDE_DIR "/scicore/soft/apps/Eigen/3.2.1-goolf-1.4.10/include/Eigen" CACHE PATH "")
set(EIGEN2_INCLUDE_DIR "/import/bc2/soft/app/eigen/current/Linux/include/eigen2" CACHE PATH "") set(QT_QMAKE_EXECUTABLE "/usr/bin/qmake-qt4" CACHE PATH "")
set(QT_QMAKE_EXECUTABLE "/import/bc2/soft/app/Qt/current/Linux/bin/qmake" CACHE PATH "")
set(COMPOUND_LIB "/import/bc2/home/schwede/GROUP/OpenStructure/ChemLib/compounds.chemlib" CACHE PATH "") set(COMPOUND_LIB "/import/bc2/home/schwede/GROUP/OpenStructure/ChemLib/compounds.chemlib" CACHE PATH "")
set (CMAKE_BUILD_TYPE "Release" CACHE PATH "") set (CMAKE_BUILD_TYPE "Release" CACHE PATH "")
set(PNG_INCLUDE_DIR "/usr/X11/include" CACHE PATH "include path for system libpng") set(PNG_INCLUDE_DIR "/usr/X11/include" CACHE PATH "include path for system libpng")
set(EIGEN2_INCLUDE_DIR "/usr/local/include/eigen3" CACHE PATH "Eigen3-compatible include path") set(EIGEN3_INCLUDE_DIR "/usr/local/include/eigen3" CACHE PATH "Eigen3-compatible include path")
...@@ -10,7 +10,7 @@ set(PNG_LIBRARY "${USER_DIR}/lib/libpng13.lib" CACHE PATH "includes for libpng") ...@@ -10,7 +10,7 @@ set(PNG_LIBRARY "${USER_DIR}/lib/libpng13.lib" CACHE PATH "includes for libpng")
set(ZLIB_INCLUDE_DIR "${USER_DIR}/include/" CACHE PATH "zlib includes") set(ZLIB_INCLUDE_DIR "${USER_DIR}/include/" CACHE PATH "zlib includes")
set(TIFF_INCLUDE_DIR "${USER_DIR}/include/" CACHE PATH "tiff includes") set(TIFF_INCLUDE_DIR "${USER_DIR}/include/" CACHE PATH "tiff includes")
set(TIFF_LIBRARY "${USER_DIR}/lib/libtiff_i.lib/" CACHE PATH "tiff includes") set(TIFF_LIBRARY "${USER_DIR}/lib/libtiff_i.lib/" CACHE PATH "tiff includes")
set(EIGEN2_INCLUDE_DIR "${USER_DIR}/include/eigen2" CACHE PATH "gsl includes") set(EIGEN3_INCLUDE_DIR "${USER_DIR}/include/Eigen" CACHE PATH "gsl includes")
set(USE_DOUBLE_PRECISION "OFF" CACHE PATH "double precision") set(USE_DOUBLE_PRECISION "OFF" CACHE PATH "double precision")
set(FFTW_LIBRARIES "${USER_DIR}/lib/libfftw3f-3.lib" CACHE PATH "fftw libs") set(FFTW_LIBRARIES "${USER_DIR}/lib/libfftw3f-3.lib" CACHE PATH "fftw libs")
#set(FFTW_LIBRARIES "${USER_DIR}/lib/libfftw3-3.lib" CACHE PATH "fftw libs") #set(FFTW_LIBRARIES "${USER_DIR}/lib/libfftw3-3.lib" CACHE PATH "fftw libs")
......
# - Try to find Eigen2 lib # - Try to find Eigen3 lib
# Once done this will define # Once done this will define
# #
# EIGEN2_FOUND - system has eigen lib # EIGEN3_FOUND - system has eigen lib
# EIGEN2_INCLUDE_DIRS - the eigen include directories # EIGEN3_INCLUDE_DIRS - the eigen include directories
find_path(EIGEN2_INCLUDE_DIR NAMES Eigen/Core find_path(EIGEN3_INCLUDE_DIR NAMES Eigen/Core
PATH_SUFFIXES eigen2
PATH_SUFFIXES eigen3 PATH_SUFFIXES eigen3
) )
set(EIGEN2_INCLUDE_DIRS ${EIGEN2_INCLUDE_DIR}) set(EIGEN3_INCLUDE_DIRS ${EIGEN3_INCLUDE_DIR})
include (FindPackageHandleStandardArgs) include (FindPackageHandleStandardArgs)
find_package_handle_standard_args (EIGEN DEFAULT_MSG EIGEN2_INCLUDE_DIR) find_package_handle_standard_args (EIGEN DEFAULT_MSG EIGEN3_INCLUDE_DIR)
mark_as_advanced (EIGEN2_INCLUDE_DIR) mark_as_advanced (EIGEN3_INCLUDE_DIR)
......
# Try to find OpenMM
#
# Once done this will define
# - Find OpenMM
# Find the native OpenMM includes and library
#
# OPEN_MM_INCLUDE_DIR - OpenMM include dirs.
# OPEN_MM_LIBRARIES - List of libraries when using OpenMM.
# OPEN_MM_FOUND - True if OpenMM found.
# OPEN_MM_PLUGIN_DIR - Default for ressource plugins
if (OPEN_MM_INCLUDE_DIR)
set(OPEN_MM_FOUND TRUE)
else (OPEN_MM_INCLUDE_DIR)
find_path (OPEN_MM_INCLUDE_DIR OpenMM.h)
endif (OPEN_MM_INCLUDE_DIR)
if (OPEN_MM_LIBRARY)
set(OPEN_MM_FOUND TRUE)
else (OPEN_MM_LIBRARY)
find_library (OPEN_MM_LIBRARY NAMES OpenMM)
endif(OPEN_MM_LIBRARY)
set(OPEN_MM_LIBRARIES ${OPEN_MM_LIBRARY})
set(OPEN_MM_INCLUDE_DIRS ${OPEN_MM_INCLUDE_DIR})
# handle the QUIETLY and REQUIRED arguments and set OpenMM_FOUND to TRUE if
# all listed variables are TRUE
include (FindPackageHandleStandardArgs)
find_package_handle_standard_args (OpenMM DEFAULT_MSG OPEN_MM_LIBRARY OPEN_MM_INCLUDE_DIR)
if (OPEN_MM_LIBRARY)
set(OPEN_MM_FOUND TRUE)
if (NOT OPEN_MM_PLUGIN_DIR)
set(OPEN_MM_PLUGIN_DIR "/usr/local/openmm/lib/plugins")
endif (NOT OPEN_MM_PLUGIN_DIR)
if (NOT EXISTS "${OPEN_MM_PLUGIN_DIR}")
message(FATAL_ERROR "OpenMM plugin directory does not exist: '${OPEN_MM_PLUGIN_DIR}'")
endif (NOT EXISTS "${OPEN_MM_PLUGIN_DIR}")
endif (OPEN_MM_LIBRARY)
mark_as_advanced (OPEN_MM_LIBRARY OPEN_MM_INCLUDE_DIR OPEN_MM_FOUND
OPEN_MM_PLUGIN_DIR)
This diff is collapsed.
ATOM 1 C1 ETH A 1 -0.246 0.000 1.285 1.00 0.00 C
ATOM 2 C2 ETH A 1 0.566 0.000 -0.011 1.00 0.00 C
ATOM 3 O ETH A 1 -0.322 0.000 -1.130 1.00 0.00 O
ATOM 4 H11 ETH A 1 0.376 0.000 2.068 1.00 0.00 H
ATOM 5 H12 ETH A 1 -0.847 0.866 1.302 1.00 0.00 H
ATOM 6 H13 ETH A 1 -0.847 -0.866 1.302 1.00 0.00 H
ATOM 7 H21 ETH A 1 1.142 0.816 -0.043 1.00 0.00 H
ATOM 8 H22 ETH A 1 1.142 -0.816 -0.043 1.00 0.00 H
ATOM 9 OH ETH A 1 0.209 0.000 -1.977 1.00 0.00 H
END
from ost.mol.mm import *
from PyQt4 import QtCore
"""
MM - Demo
!!!Requires OpenStructure to be compiled with MM-Support!!!
This demo should show the procedure of building a custom forcefield.
The script can directly be started from within DNG.
"""
class Anim(QtCore.QTimer):
def __init__(self,sim,go):
QtCore.QTimer.__init__(self)
self.ent=sim.GetEntity()
self.sim=sim
self.go = go
self.ed = ent.EditXCS()
QtCore.QObject.connect(self, QtCore.SIGNAL("timeout()"), self.OnTimer)
def OnTimer(self):
self.sim.Steps(1)
positions = sim.GetPositions()
for a, pos in zip(self.ent.atoms,positions):
self.ed.SetAtomPos(a,pos)
self.go.UpdatePositions()
#Let's create an ethanol without its hydrogens
ent = ost.mol.CreateEntity()
ed = ent.EditXCS()
chain = ed.InsertChain('A')
res = ed.AppendResidue(chain,"ETH")
c1 = ed.InsertAtom(res,"C1",geom.Vec3(-0.246,0.0,1.285),"C")
c2 = ed.InsertAtom(res,"C2",geom.Vec3(0.566,0.0,-0.011),"C")
o = ed.InsertAtom(res,"O",geom.Vec3(-0.322,0.0,-1.130),"O")
ed.Connect(c1,c2)
ed.Connect(c2,o)
ff = Forcefield()
#generate hydrogen constructor and add it to the forcefield
h_constructor = GromacsHydrogenConstructor()
h_constructor.AddHydrogenRule(3,4,["H11","H12","H13"],["C1","C2","O"])
h_constructor.AddHydrogenRule(2,6,["H21","H22"],["C2","O","C1"])
h_constructor.AddHydrogenRule(1,2,["OH"],["O","C2","C1"])
ff.AddHydrogenConstructor("ETH",h_constructor)
#generate buildingblock and add it to the forcefield
building_block = BuildingBlock()
building_block.AddAtom("C1","C",-0.18)
building_block.AddAtom("C2","C",0.145)
building_block.AddAtom("O","O",-0.683)
building_block.AddAtom("H11","H",0.06)
building_block.AddAtom("H12","H",0.06)
building_block.AddAtom("H13","H",0.06)
building_block.AddAtom("H21","H",0.06)
building_block.AddAtom("H22","H",0.06)
building_block.AddAtom("OH","OH",0.418)
bonds = list()
for i in range(8):
bonds.append(Interaction(FuncType.HarmonicBond))
bonds[0].SetNames(["H11","C1"])
bonds[1].SetNames(["H12","C1"])
bonds[2].SetNames(["H13","C1"])
bonds[3].SetNames(["C1","C2"])
bonds[4].SetNames(["C2","H21"])
bonds[5].SetNames(["C2","H22"])
bonds[6].SetNames(["O","C2"])
bonds[7].SetNames(["O","OH"])
for b in bonds:
building_block.AddBond(b)
ff.AddBuildingBlock("ETH",building_block)
#add masses to the forcefield
ff.AddMass("C",12.011)
ff.AddMass("O",15.999)
ff.AddMass("H",1.008)
ff.AddMass("OH",1.008)
#generate bonds and add them to the forcefield
h_c_bond = Interaction(FuncType.HarmonicBond)
c_c_bond = Interaction(FuncType.HarmonicBond)
c_o_bond = Interaction(FuncType.HarmonicBond)
o_h_bond = Interaction(FuncType.HarmonicBond)
h_c_bond.SetTypes(["H","C"])
c_c_bond.SetTypes(["C","C"])
c_o_bond.SetTypes(["C","O"])
o_h_bond.SetTypes(["O","OH"])
h_c_bond.SetParam([0.1111,269449.6])
c_c_bond.SetParam([0.153,186188.0])
c_o_bond.SetParam([0.142,358150.4])
o_h_bond.SetParam([0.096,456056.0])
ff.AddBond(h_c_bond)
ff.AddBond(c_c_bond)
ff.AddBond(c_o_bond)
ff.AddBond(o_h_bond)
#generate angles and add them to the forcefield
h_c_h_angle = Interaction(FuncType.UreyBradleyAngle)
h_c_c_angle = Interaction(FuncType.UreyBradleyAngle)
h_c_o_angle = Interaction(FuncType.UreyBradleyAngle)
c_c_o_angle = Interaction(FuncType.UreyBradleyAngle)
c_o_h_angle = Interaction(FuncType.UreyBradleyAngle)
h_c_h_angle.SetTypes(["H","C","H"])
h_c_c_angle.SetTypes(["H","C","C"])
h_c_o_angle.SetTypes(["H","C","O"])
c_c_o_angle.SetTypes(["C","C","O"])
c_o_h_angle.SetTypes(["C","O","OH"])
h_c_h_angle.SetParam([1.892,297.064,0.1802,4518.72])
h_c_c_angle.SetParam([1.922,289.5328,0.2179,18853.104])
h_c_o_angle.SetParam([1.900,384.0912,0.0,0.0])
c_c_o_angle.SetParam([1.920,633.4576,0.0,0.0])
c_o_h_angle.SetParam([1.850,481.16,0.0,0.0])
ff.AddAngle(h_c_h_angle)
ff.AddAngle(h_c_c_angle)
ff.AddAngle(h_c_o_angle)
ff.AddAngle(c_c_o_angle)
ff.AddAngle(c_o_h_angle)
#generate dihedrals and add them to the forcefield
x_c_c_x_dihedral = Interaction(FuncType.PeriodicDihedral)
c_c_o_h_dihedral_one = Interaction(FuncType.PeriodicDihedral)
c_c_o_h_dihedral_two = Interaction(FuncType.PeriodicDihedral)
c_c_o_h_dihedral_three = Interaction(FuncType.PeriodicDihedral)
x_c_o_x_dihedral = Interaction(FuncType.PeriodicDihedral)
x_c_c_x_dihedral.SetTypes(["X","C","C","X"])
c_c_o_h_dihedral_one.SetTypes(["C","C","O","OH"])
c_c_o_h_dihedral_two.SetTypes(["C","C","O","OH"])
c_c_o_h_dihedral_three.SetTypes(["C","C","O","OH"])
x_c_o_x_dihedral.SetTypes(["X","C","O","X"])
x_c_c_x_dihedral.SetParam([3,0.0,0.66944])
c_c_o_h_dihedral_one.SetParam([1,0.0,5.4392])
c_c_o_h_dihedral_two.SetParam([2,0.0,1.2552])
c_c_o_h_dihedral_three.SetParam([3,0.0,1.75728])
x_c_o_x_dihedral.SetParam([3,0.0,0.58576])
ff.AddDihedral(x_c_c_x_dihedral)
ff.AddDihedral(c_c_o_h_dihedral_one)
ff.AddDihedral(c_c_o_h_dihedral_two)
ff.AddDihedral(c_c_o_h_dihedral_three)
ff.AddDihedral(x_c_o_x_dihedral)
#add lj parameters
c_lj = Interaction(FuncType.LJ)
o_lj = Interaction(FuncType.LJ)
h_lj = Interaction(FuncType.LJ)
oh_lj = Interaction(FuncType.LJ)
c_lj.SetTypes(["C"])
o_lj.SetTypes(["O"])
h_lj.SetTypes(["H"])
oh_lj.SetTypes(["OH"])
c_lj.SetParam([0.37,0.28])
o_lj.SetParam([0.32,0.64])
h_lj.SetParam([0.24,0.09])
oh_lj.SetParam([0.04,0.2])
ff.AddLJ(c_lj)
ff.AddLJ(o_lj)
ff.AddLJ(h_lj)
ff.AddLJ(oh_lj)
#we do not explicitely set the 1,4 interactions,
#we rather decrease it by a certain factor
ff.SetFudgeLJ(0.5)
ff.SetFudgeQQ(0.5)
#construct settings
settings = Settings()
settings.init_temperature = 310
settings.forcefield = ff
#construct integrator with tiny timesteps
settings.integrator = VerletIntegrator(0.0001)
#create the simulation
sim = Simulation(ent,settings)
io.SavePDB(sim.GetEntity(),'ethanol.pdb')
go = gfx.Entity("ethanol",sim.GetEntity())
go.SetRenderMode(gfx.CUSTOM)
scene.Add(go)
anim = Anim(sim,go)
anim.start(2)
from ost.mol.mm import *
from PyQt4 import QtCore
from math import sqrt
from math import sin
"""
MM - Demo
!!!Requires OpenStructure to be compiled with MM-Support!!!
This demo should show the procedure building a topology from scratch and
setting up a simulation. As an additional piece of sugar, it demonstrates
the capabilities of altering force parameters as the simulation runs.
The script can directly be started from within DNG.
"""
class Anim(QtCore.QTimer):
def __init__(self,sim,go,cc_bond_index):
QtCore.QTimer.__init__(self)
self.ent=sim.GetEntity()
self.sim=sim
self.go = go
self.ed = ent.EditXCS()
self.bond_index = cc_bond_index
top = sim.GetTopology()
param = top.GetHarmonicBondParameters(cc_bond_index)
self.bond_length = param[2]
self.force_constant = param[3]
self.steps = 0
self.counter = 0
QtCore.QObject.connect(self, QtCore.SIGNAL("timeout()"), self.OnTimer)
def OnTimer(self):
if self.steps % 20 == 0:
sim.ResetHarmonicBond(self.bond_index,sin(self.counter*3.1416/100)*self.bond_length*0.75 + self.bond_length,self.force_constant)
self.counter += 1
self.sim.Steps(1)
positions = sim.GetPositions()
for a, pos in zip(self.ent.atoms,positions):
self.ed.SetAtomPos(a,pos)
self.go.UpdatePositions()
self.steps+=1
#create topology by only defining masses
prof = io.IOProfile(dialect='PDB', fault_tolerant=False,
quack_mode=False, processor=conop.HeuristicProcessor())
ent = io.LoadPDB('ethanol.pdb', profile=prof)
masses = [12.011,12.011,15.999,1.008,1.008,1.008,1.008,1.008,1.008]
top = Topology(masses)
#all per atom parameters have to be set at once...
charges = [-0.18,0.145,-0.683,0.06,0.06,0.06,0.06,0.06,0.418]
sigmas = [0.37,0.37,0.32,0.24,0.24,0.24,0.24,0.24,0.04]
epsilons = [0.28,0.28,0.64,0.09,0.09,0.09,0.09,0.09,0.2]
top.SetCharges(charges)
top.SetSigmas(sigmas)
top.SetEpsilons(epsilons)
#we can set the fudge parameters for 1,4 interactions also in the topology
top.SetFudgeQQ(0.5)
top.SetFudgeLJ(0.5)
#add the bonds
top.AddHarmonicBond(0,1,0.153,186188.0) #C-C
top.AddHarmonicBond(0,3,0.1111,269449.6) #C-H
top.AddHarmonicBond(0,4,0.1111,269449.6) #C-H
top.AddHarmonicBond(0,5,0.1111,269449.6) #C-H
top.AddHarmonicBond(1,6,0.1111,269449.6) #C-H
top.AddHarmonicBond(1,7,0.1111,269449.6) #C-H
top.AddHarmonicBond(1,2,0.142,358150.4) #C-O
top.AddHarmonicBond(2,8,0.096,456056.0) #O-H
#add the angles
top.AddUreyBradleyAngle(3,0,4,1.892,297.064,0.1802,4518.72) #H-C-H
top.AddUreyBradleyAngle(3,0,5,1.892,297.064,0.1802,4518.72) #H-C-H
top.AddUreyBradleyAngle(4,0,5,1.892,297.064,0.1802,4518.72) #H-C-H
top.AddUreyBradleyAngle(3,0,1,1.922,289.5328,0.2179,18853.104) #H-C-C
top.AddUreyBradleyAngle(4,0,1,1.922,289.5328,0.2179,18853.104) #H-C-C
top.AddUreyBradleyAngle(5,0,1,1.922,289.5328,0.2179,18853.104) #H-C-C
top.AddUreyBradleyAngle(0,1,6,1.922,289.5328,0.2179,18853.104) #H-C-C
top.AddUreyBradleyAngle(0,1,7,1.922,289.5328,0.2179,18853.104) #H-C-C
top.AddUreyBradleyAngle(0,1,2,1.920,633.4576,0.0,0.0) #C-C-O
top.AddUreyBradleyAngle(6,1,2,1.900,384.0912,0.0,0.0) #H-C-O
top.AddUreyBradleyAngle(7,1,2,1.900,384.0912,0.0,0.0) #H-C-O
top.AddUreyBradleyAngle(6,1,7,1.892,297.064,0.1802,4518.72) #H-C-H
top.AddUreyBradleyAngle(1,2,8,1.850,481.16,0.0,0.0) #C-O-H
#add the dihedrals
top.AddPeriodicDihedral(3,0,1,6,3,0.0,0.66944) #H-C-C-H
top.AddPeriodicDihedral(3,0,1,7,3,0.0,0.66944) #H-C-C-H
top.AddPeriodicDihedral(3,0,1,2,3,0.0,0.66944) #H-C-C-O
top.AddPeriodicDihedral(4,0,1,6,3,0.0,0.66944) #H-C-C-H
top.AddPeriodicDihedral(4,0,1,7,3,0.0,0.66944) #H-C-C-H
top.AddPeriodicDihedral(4,0,1,2,3,0.0,0.66944) #H-C-C-O
top.AddPeriodicDihedral(5,0,1,6,3,0.0,0.66944) #H-C-C-H
top.AddPeriodicDihedral(5,0,1,7,3,0.0,0.66944) #H-C-C-H
top.AddPeriodicDihedral(5,0,1,2,3,0.0,0.66944) #H-C-C-O
top.AddPeriodicDihedral(6,1,2,8,3,0.0,0.58576) #H-C-O-H
top.AddPeriodicDihedral(7,1,2,8,3,0.0,0.58576) #H-C-O-H
top.AddPeriodicDihedral(0,1,2,8,1,0.0,5.4392) #C-C-O-H
top.AddPeriodicDihedral(0,1,2,8,2,0.0,1.2552) #C-C-O-H
top.AddPeriodicDihedral(0,1,2,8,3,0.0,1.75728) #C-C-O-H
#we have to manually set the exclusions for the nonbonded interaction...
#1,2 interactions
top.AddExclusion(3,0)
top.AddExclusion(4,0)
top.AddExclusion(5,0)
top.AddExclusion(0,1)
top.AddExclusion(1,6)
top.AddExclusion(1,7)
top.AddExclusion(1,2)
top.AddExclusion(2,8)
#1,3 interactions
top.AddExclusion(3,4)
top.AddExclusion(3,5)
top.AddExclusion(4,5)
top.AddExclusion(3,1)
top.AddExclusion(4,1)
top.AddExclusion(5,1)
top.AddExclusion(0,6)
top.AddExclusion(0,7)
top.AddExclusion(0,2)
top.AddExclusion(6,2)
top.AddExclusion(7,2)
top.AddExclusion(1,8)
#add lj-pairs, note, that the fudge parameters still get applied during system setup
top.AddLJPair(3,6,0.5*(sigmas[3]+sigmas[6]),sqrt(epsilons[3]*epsilons[6]))
top.AddLJPair(3,7,0.5*(sigmas[3]+sigmas[7]),sqrt(epsilons[3]*epsilons[7]))
top.AddLJPair(3,2,0.5*(sigmas[3]+sigmas[2]),sqrt(epsilons[3]*epsilons[2]))
top.AddLJPair(4,6,0.5*(sigmas[4]+sigmas[6]),sqrt(epsilons[4]*epsilons[6]))
top.AddLJPair(4,7,0.5*(sigmas[4]+sigmas[7]),sqrt(epsilons[4]*epsilons[7]))
top.AddLJPair(4,2,0.5*(sigmas[4]+sigmas[2]),sqrt(epsilons[4]*epsilons[2]))
top.AddLJPair(5,6,0.5*(sigmas[5]+sigmas[6]),sqrt(epsilons[5]*epsilons[6]))
top.AddLJPair(5,7,0.5*(sigmas[5]+sigmas[7]),sqrt(epsilons[5]*epsilons[7]))
top.AddLJPair(5,2,0.5*(sigmas[5]+sigmas[2]),sqrt(epsilons[5]*epsilons[2]))
top.AddLJPair(6,8,0.5*(sigmas[6]+sigmas[8]),sqrt(epsilons[6]*epsilons[8]))
top.AddLJPair(7,8,0.5*(sigmas[7]+sigmas[8]),sqrt(epsilons[7]*epsilons[8]))
top.AddLJPair(0,8,0.5*(sigmas[0]+sigmas[8]),sqrt(epsilons[0]*epsilons[8]))
settings = Settings()
settings.init_temperature = 310
settings.integrator = VerletIntegrator(0.0001)
#let's extract the cc_bond_index to have some fun with it
cc_bond_index = top.GetHarmonicBondIndices(0,1)[0]
#create the simulation
sim = Simulation(top,ent,settings)
go = gfx.Entity("ethanol",sim.GetEntity())
go.SetRenderMode(gfx.CUSTOM)
scene.Add(go)
anim = Anim(sim,go,cc_bond_index)
anim.start(2)
from ost.mol import mm
from PyQt4 import QtCore
"""
MM - Demo
!!!Requires OpenStructure to be compiled with MM-Support!!!
This demo should show the setup of a simple simulation using implicit solvent.
The script can directly be started from within DNG.
"""
class Anim(QtCore.QTimer):
def __init__(self,sim,ent,go):
QtCore.QTimer.__init__(self)
self.ent=ent
self.sim=sim
self.go = go
self.ed = ent.EditXCS(ost.mol.BUFFERED_EDIT)
QtCore.QObject.connect(self, QtCore.SIGNAL("timeout()"), self.OnTimer)
def OnTimer(self):
self.sim.Steps(5)
positions = sim.GetPositions()
for a, pos in zip(self.ent.atoms,positions):
self.ed.SetAtomPos(a,pos)
self.ed.UpdateICS()
self.go.UpdatePositions()
prot=ost.mol.CreateEntityFromView(io.LoadPDB('1CRN.pdb').Select('peptide=true'),True)
settings = mm.Settings()
settings.integrator = mm.LangevinIntegrator(310,1,0.002)
settings.add_gbsa = True
settings.forcefield = mm.LoadCHARMMForcefield()
settings.platform = mm.Platform.CPU
settings.nonbonded_method = mm.NonbondedMethod.CutoffNonPeriodic
sim = mm.Simulation(prot,settings)
sim.ApplyLBFGS(tolerance = 1.0, max_iterations = 200)
sim.UpdatePositions()
ent = sim.GetEntity()
go = gfx.Entity("crambin",ent)
scene.Add(go)
anim = Anim(sim,ent,go)
anim.start(1)
File added
This diff is collapsed.
from ost.gui import trajectory_viewer
"""
MM - Demo
!!!Requires OpenStructure to be compiled with MM-Support!!!
This demo demonstrates how to actually look at the trajectory created
in the gb_example_writing_trajectory.py script.
"""
t = io.LoadCHARMMTraj("gb_example_traj.pdb","gb_example_traj.dcd")
go = gfx.Entity("gb",t.GetEntity())
scene.Add(go)
w = trajectory_viewer.TrajWidget([t],[go])
w.show()
from ost.mol.mm import *
"""
MM - Demo
!!!Requires OpenStructure to be compiled with MM-Support!!!
This demo should show the procedure of setting up a simulation and write
the results directly to disk.
"""
prot=io.LoadPDB('1CRN.pdb')
#set up the simulation
settings = Settings()
settings.integrator = LangevinIntegrator(310,1,0.002)
settings.add_gbsa = True
settings.forcefield = LoadCHARMMForcefield()
settings.nonbonded_cutoff = 8.0
settings.nonbonded_method = NonbondedMethod.CutoffNonPeriodic
settings.platform = Platform.CPU
settings.cpu_properties["CpuThreads"] = "2"
sim = Simulation(prot,settings)
#minimize it
sim.ApplySD(tolerance = 1.0, max_iterations = 200)
#create a trajectory observer and register it to the simulation
observer = TrajWriter(10,"gb_example_traj.pdb","gb_example_traj.dcd")
sim.Register(observer)
#run the simulation
sim.Steps(10000)
#Trajectory Observer needs to finalize, otherwise you might get a corrupt dcd file
observer.Finalize()
...@@ -15,6 +15,7 @@ So far, the binding module includes: ...@@ -15,6 +15,7 @@ So far, the binding module includes:
dssp dssp
blast blast
hhblits
msms msms
tmtools tmtools
clustalw clustalw
......
:mod:`~ost.bindings.hhblits` - Search related sequences in databases
================================================================================
.. module:: ost.bindings.hhblits
:synopsis: Search related sequences in databases
Introduction
--------------------------------------------------------------------------------
HHblits is a sequence search tool like BLAST, able to find more distant
homologs This is achieved by performing `profile-profile searches`. In BLAST, a
\query sequence is searched against a sequence database. That is a
`sequence-sequence search`. HHblits works on a profile database, usually that
one is provided, queried with a sequence profile. The latter one needs to be
calculated before the actual search. In very simple words, HHblits is using
per-sequence scoring functions to be more sensitive, in this particular case
Hidden Markov models. The software suite needed for HHblits can be found
`here <http://toolkit.tuebingen.mpg.de/hhblits>`_.
Examples
--------------------------------------------------------------------------------
A typical search: Get an instance of the bindings, build the search profile out
of the query sequence, run the search and iterate results. Since
:class:`~HHblits` works with a :class:`~ost.seq.SequenceHandle` or a FastA
file, both ways are shown.
First query by sequence:
.. code-block:: python
from ost.bindings import hhblits
# get a SequenceHandle
query_seq = seq.CreateSequence('Query',
'MRIILLGAPGAGKGTQAQFIMEKYGIPQISTGDMLRAAVKSGS'+
'ELGKQAKDIMDAGKLVTDELVIALVKERIAQEDCRNGFLLDGF'+
'PRTIPQADAMKEAGINVDYVLEFDVPDELIVDRIVGRRVHAPS'+
'GRVYHVKFNPPKVEGKDDVTGEELTTRKDDQEETVRKRLVEYH'+
'QMTAPLIGYYYYSKEAEAGNTKYAKVDGTKPVAEVRADLEKILG')
# set up the search environment
# $EBROOTHHMINSUITE points to the root of a HHsuite installation
hh = hhblits.HHblits(query_seq, os.getenv('EBROOTHHMINSUITE'))
# now create a search profile for the query sequence against the NR20 db
# provided on the HHblits web page, nr20_12Aug11 is just the prefix common to
# all db files, so `ls hhtools/nr20_12Aug11*` would list all of them
a3m_file = hh.BuildQueryMSA(nrdb='hhtools/nr20_12Aug11')
# search time! we just search against NR20 again, but every HHblits db is
# working here, e.g. one build from all the sequences in PDB
hit_file = hh.Search(a3m_file, 'hhtools/nr20_12Aug11')
# lets have a look at the resuls
with open(hit_file) as hit_fh:
header, hits = hhblits.ParseHHblitsOutput(hit_fh)
for hit in hits:
print hit.aln
Very similar going by file:
.. code-block:: python
from ost.bindings import hhblits
# set up the search environment
# $EBROOTHHMINSUITE points to the root of a HHsuite installation
hh = hhblits.HHblits('query.fas', os.getenv('EBROOTHHMINSUITE'))
# now create a search profile for the query sequence against the NR20 db
# provided on the HHblits web page, nr20_12Aug11 is just the prefix common to
# all db files, so `ls hhtools/nr20_12Aug11*` would list all of them
a3m_file = hh.BuildQueryMSA(nrdb='hhtools/nr20_12Aug11')
# search time! we just search against NR20 again, but every HHblits db is
# working here, e.g. one build from all the sequences in PDB
hit_file = hh.Search(a3m_file, 'hhtools/nr20_12Aug11')
# lets have a look at the resuls
with open(hit_file) as hit_fh:
header, hits = hhblits.ParseHHblitsOutput(hit_fh)
for hit in hits:
print hit.aln
The alignments produced by HHblits are sometimes slightly better than by BLAST,
so one may want to extract them:
.. code-block:: python
from ost.bindings import hhblits
# set up the search environment
# $EBROOTHHMINSUITE points to the root of a HHsuite installation
hh = hhblits.HHblits('query.fas', os.getenv('EBROOTHHMINSUITE'))
# now create a search profile for the query sequence against the NR20 db
# provided on the HHblits web page, nr20_12Aug11 is just the prefix common to
# all db files, so `ls hhtools/nr20_12Aug11*` would list all of them
a3m_file = hh.BuildQueryMSA(nrdb='hhtools/nr20_12Aug11')
# note that ParseA3M is not a class method but a module function
output = hhblits.ParseA3M(open(a3m_file))
print output['msa']
Binding API
--------------------------------------------------------------------------------
.. autoclass:: ost.bindings.hhblits.HHblits
:members:
.. autoclass:: ost.bindings.hhblits.HHblitsHit
.. autoclass:: ost.bindings.hhblits.HHblitsHeader
.. autofunction:: ost.bindings.hhblits.ParseHHblitsOutput
.. autofunction:: ost.bindings.hhblits.ParseA3M
.. autofunction:: ost.bindings.hhblits.ParseHeaderLine
.. autofunction:: ost.bindings.hhblits.ParseHHM
.. autofunction:: ost.bindings.hhblits.EstimateMemConsumption
.. LocalWords: HHblits homologs
...@@ -8,6 +8,7 @@ dssp.py ...@@ -8,6 +8,7 @@ dssp.py
clustalw.py clustalw.py
utils.py utils.py
naccess.py naccess.py
hhblits.py
blast.py blast.py
cadscore.py cadscore.py
kclust.py kclust.py
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please to comment