SDA  7.1
Simulation of Diffusional Association
 All Classes Files Functions Variables Groups Pages
Functions/Subroutines
compute_distance.f90 File Reference

General functions to compute the distance between 2 proteins. More...

Functions/Subroutines

subroutine compute_dist_pbc (geom, prot1, prot2, dist, r_corrected)
 Compute distance and relative position between 2 solutes.
Always r2-r1. More...
 
subroutine compute_distance (prot1, prot2, dist)
 Compute the distance between two proteins
Take into if a surface is present, but not for pbc ( call compute_dist_pbc instead ) More...
 
subroutine compute_distance_test (geom, rtry, prot1_pos, dist)
 Version used by make_bd_2prot to update dist.
. More...
 
subroutine compute_distance_criteria_v2 (prot1, prot2, pos_relat, n1, n2, dist)
 Compute the distance need by reaction criteria
. More...
 
integer function, dimension(:,:),
allocatable 
get_atoms_min_dist (prot1, prot2, max_contact_dist, min_contact_dist)
 Compute the minimum distance between the surface atoms (i.e., solvent accessible atoms) of the two given proteins (prot1 and prot2). More...
 

Detailed Description

General functions to compute the distance between 2 proteins.

Version
{version 7.1 (2015)}

Copyright (c) 2009, 2010, 2015 Heidelberg Institute of Theoretical Studies (HITS, www.h-its.org) formerly EML Research gGmbH (EML-R ) Schloss-Wolfsbrunnenweg 35 69118 Heidelberg, Germany

Copyright (c) 2000, 2003 European Molecular Biology Laboratory Meyerhofstr. 1, Postfach 10.2209 D-69012, Heidelberg, Germany

Please send your contact address to get information on updates and new features to "mcmsoft@h-its.org". Questions will be answered as soon as possible.

References (see also http://mcm.h-its.org/sda7/doc/doc_sda7/references.html):

Brownian dynamics simulation of protein-protein diffusional encounter.

(1998) Methods, 14, 329-341.

Authors: M.Martinez, N.J.Bruce, J.Romanowska, D.B.Kokh, P.Mereghetti, R.R.Gabdoulline, M. Ozboyaci, S.Richter and R.C.Wade


Function/Subroutine Documentation

subroutine compute_dist_pbc ( type ( geometry )  geom,
type ( protein ), intent(in)  prot1,
type ( protein ), intent(in)  prot2,
real ( kind=4 ), intent(out)  dist,
real ( kind=8 ), dimension ( 3 ), intent(out)  r_corrected 
)

Compute distance and relative position between 2 solutes.
Always r2-r1.

Applied with sdamm (should work with sda_2prot). Take into account sphere/box, PBC, surface

It returns the correct relative position (after PBC)
but the distance (dist) is only along the z-axis if one of the solute is a surface.
So sqrt( dot(r_corrected,r_corrected) ) is not always equal to dist.

This function should always be called in force or energy calculation in the case of sdamm. With sda_2proteins, pos_relat is always = prot2 % position, and dist is setup correctly by compute_distance_test

Parameters
geom: instance of mod_geometry
prot1: solute 1, instance of protein
prot2: solute 2, instance of protein
dist: return distance between the 2 solutes
r_corrected: return corrected relative position

Here is the caller graph for this function:

subroutine compute_distance ( type ( protein ), intent(in)  prot1,
type ( protein ), intent(in)  prot2,
real ( kind=4 ), intent(out)  dist 
)

Compute the distance between two proteins
Take into if a surface is present, but not for pbc ( call compute_dist_pbc instead )

Parameters
prot1: solute 1, instance of protein
prot2: solute 2, instance of protein
dist: return distance between the 2 solutes

Here is the caller graph for this function:

subroutine compute_distance_criteria_v2 ( type ( protein )  prot1,
type ( protein )  prot2,
real ( kind = 4 ), dimension( 3 ), intent(in)  pos_relat,
integer, intent(in)  n1,
integer, intent(in)  n2,
real ( kind=4 ), dimension ( n2-n1+1 )  dist 
)

Compute the distance need by reaction criteria
.

version changed because the movable atoms belong to protein now. checked if prot1 is a surface - then, taking only the z-distance

Parameters
prot1: solute 1
prot2: solute 2
pos_relat: relative position
n1: which reaction criterium to start from
n2: which reaction criterium to end on
dist: array with values of squared distance for each reaction criterium

Here is the caller graph for this function:

subroutine compute_distance_test ( type ( geometry )  geom,
real ( kind = 8 ), dimension ( 3 )  rtry,
real ( kind = 8 ), dimension ( 3 )  prot1_pos,
real ( kind=4 )  dist 
)

Version used by make_bd_2prot to update dist.
.

Used only by sda_2proteins, check if one of the solute is a surface

Here tricky, because move not yet accepted before test escape
Particular to sda_2proteins, rtry is the potential new position

Parameters
geom: instance of mod_geometry
rtry: potential new position of the solute
prot1_pos: position of the solute 1
dist: return the new distance if the move rtry is accepted

Here is the caller graph for this function:

integer function, dimension ( :,: ), allocatable get_atoms_min_dist ( type ( protein ), intent(in)  prot1,
type ( protein ), intent(in)  prot2,
real ( kind = 8 ), intent(in)  max_contact_dist,
real ( kind = 8 ), intent(in)  min_contact_dist 
)

Compute the minimum distance between the surface atoms (i.e., solvent accessible atoms) of the two given proteins (prot1 and prot2).

Parameters
prot1: input, solute 1
prot2: input, solute 2
max_contact_dist: input, the maximum distance between two atoms from the two proteins that is required to define a contact
min_contact_dist: input, the minimum ditance between atoms on the same proteins to define two different contacts
Returns
contact_atoms : output, array of pairs of numbers, pointing to the accessible atoms lists of prot1 and prot2, respectively