dcdfort
Modern Fortran library for analyzing DCD trajectories
Functions/Subroutines
dcdfort_utils Module Reference

Module that contains some useful utilities. More...

Functions/Subroutines

pure real(kind=real64) function, dimension(3) pbc (a, box)
 Corrects for the periodic boundary condition. More...
 
pure real(kind=real64) function, dimension(3) cross (a, b)
 Performs cross product between two vectors. More...
 
pure real(kind=real64) function distance2 (a, b, box)
 Calculates the distance squared between two points. More...
 
pure real(kind=real64) function distance (a, b, box)
 Calculates the distance between two points. More...
 
pure real(kind=real64) function, dimension(3) bond_vector (a, b, box)
 Calculates the bond vector between two points. More...
 
pure real(kind=real64) function magnitude (a)
 Calculates the magnitude of a vector. More...
 
pure real(kind=real64) function bond_angle (a, b, c, box)
 Calculates the bond angle between two vectors. More...
 
pure real(kind=real64) function dihedral_angle (i, j, k, l, box)
 Calculates the dihedral angle between two planes formed by four atoms. More...
 

Detailed Description

Module that contains some useful utilities.

Function/Subroutine Documentation

◆ bond_angle()

pure real(kind=real64) function dcdfort_utils::bond_angle ( real(kind=real64), dimension(3), intent(in)  a,
real(kind=real64), dimension(3), intent(in)  b,
real(kind=real64), dimension(3), intent(in)  c,
real(kind=real64), dimension(6), intent(in)  box 
)

Calculates the bond angle between two vectors.

Calculates the angle between the vector formed by a-b and b-c.

Parameters
[in]afirst point
[in]bmiddle point
[in]cthird point
[in]boxbox, if pbc to be accounted for
Returns
bond angle between a-b-c
227 
228  implicit none
229  real(kind=real64) :: bond_angle
230  real(kind=real64), intent(in), dimension(3) :: a, b, c
231  real(kind=real64), intent(in) :: box(6)
232  real(kind=real64), dimension(3) :: bond1, bond2
233 
234  bond1 = bond_vector(b, a, box)
235  bond2 = bond_vector(b, c, box)
236 
237  bond_angle = dacos(dot_product(bond1, bond2)/(magnitude(bond1)*magnitude(bond2)))
238 

◆ bond_vector()

pure real(kind=real64) function, dimension(3) dcdfort_utils::bond_vector ( real(kind=real64), dimension(3), intent(in)  a,
real(kind=real64), dimension(3), intent(in)  b,
real(kind=real64), dimension(6), intent(in)  box 
)

Calculates the bond vector between two points.

Parameters
[in]afirst point
[in]bsecond point
[in]boxbox, if pbc to be accounted for
Returns
bond vector pointing from a to b
196 
197  implicit none
198  real(kind=real64) :: bond_vector(3)
199  real(kind=real64), intent(in), dimension(3) :: a, b
200  real(kind=real64), intent(in) :: box(6)
201 
202  bond_vector = pbc(a-b, box)
203 

◆ cross()

pure real(kind=real64) function, dimension(3) dcdfort_utils::cross ( real(kind=real64), dimension(3), intent(in)  a,
real(kind=real64), dimension(3), intent(in)  b 
)

Performs cross product between two vectors.

Parameters
[in]afirst vector
[in]bsecond vector
Returns
resulting cross product
137 
138  implicit none
139  real(kind=real64) :: cross(3)
140  real(kind=real64), intent(in), dimension(3) :: a, b
141 
142  cross(1) = a(2) * b(3) - a(3) * b(2)
143  cross(2) = a(3) * b(1) - a(1) * b(3)
144  cross(3) = a(1) * b(2) - a(2) * b(1)
145 

◆ dihedral_angle()

pure real(kind=real64) function dcdfort_utils::dihedral_angle ( real(kind=real64), dimension(3), intent(in)  i,
real(kind=real64), dimension(3), intent(in)  j,
real(kind=real64), dimension(3), intent(in)  k,
real(kind=real64), dimension(3), intent(in)  l,
real(kind=real64), dimension(6), intent(in)  box 
)

Calculates the dihedral angle between two planes formed by four atoms.

Calculates the dihedral angle between the vectors formed by i-j, j-k, k-l

Parameters
[in]ifirst point
[in]jmiddle point
[in]kthird point
[in]lfourth point
[in]boxbox, if pbc to be accounted for
Returns
dihedral angle forms by i-j-k-l
250 
251  implicit none
252  real(kind=real64) :: dihedral_angle
253  real(kind=real64), intent(in), dimension(3) :: i, j, k, l
254  real(kind=real64), intent(in), dimension(6) :: box
255  real(kind=real64) :: A_mag, B_mag, G_mag
256  real(kind=real64), dimension(3) :: H, G, F, A, B, cross_BA
257 
258  h = bond_vector(k, l, box)
259  g = bond_vector(k, j, box)
260  f = bond_vector(j, i, box)
261  a = cross(f, g)
262  b = cross(h, g)
263  cross_ba = cross(b, a)
264  a_mag = magnitude(a)
265  b_mag = magnitude(b)
266  g_mag = magnitude(g)
267 
268  dihedral_angle = atan2(dot_product(cross_ba, g) / (a_mag*b_mag*g_mag), dot_product(a, b) / (a_mag*b_mag))
269 

◆ distance()

pure real(kind=real64) function dcdfort_utils::distance ( real(kind=real64), dimension(3), intent(in)  a,
real(kind=real64), dimension(3), intent(in)  b,
real(kind=real64), dimension(6), intent(in), optional  box 
)

Calculates the distance between two points.

Parameters
[in]afirst point
[in]bsecond point
[in]boxbox, if pbc to be accounted for
Returns
distance
176 
177  implicit none
178  real(kind=real64) :: distance
179  real(kind=real64), intent(in), dimension(3) :: a, b
180  real(kind=real64), intent(in), optional :: box(6)
181 
182  if (present(box)) then
183  distance = dsqrt(distance2(a, b, box))
184  else
185  distance = dsqrt(distance2(a, b))
186  end if
187 

◆ distance2()

pure real(kind=real64) function dcdfort_utils::distance2 ( real(kind=real64), dimension(3), intent(in)  a,
real(kind=real64), dimension(3), intent(in)  b,
real(kind=real64), dimension(6), intent(in), optional  box 
)

Calculates the distance squared between two points.

Parameters
[in]afirst point
[in]bsecond point
[in]boxbox, if pbc to be accounted for
Returns
distance squared
154 
155  implicit none
156  real(kind=real64) :: distance2
157  real(kind=real64), intent(in), dimension(3) :: a, b
158  real(kind=real64) :: c(3)
159  real(kind=real64), intent(in), optional :: box(6)
160 
161  if (present(box)) then
162  c = pbc(a - b, box)
163  else
164  c = a - b
165  end if
166  distance2 = dot_product(c, c)
167 

◆ magnitude()

pure real(kind=real64) function dcdfort_utils::magnitude ( real(kind=real64), dimension(3), intent(in)  a)

Calculates the magnitude of a vector.

Parameters
[in]avector
Returns
magnitude of a
210 
211  implicit none
212  real(kind=real64) :: magnitude
213  real(kind=real64), intent(in) :: a(3)
214 
215  magnitude = dsqrt(dot_product(a, a))
216 

◆ pbc()

pure real(kind=real64) function, dimension(3) dcdfort_utils::pbc ( real(kind=real64), dimension(3), intent(in)  a,
real(kind=real64), dimension(6), intent(in)  box 
)

Corrects for the periodic boundary condition.

Moves particle (or vector) the distance of half the box if it is more than half the distance of the box

Parameters
[in]aoriginal coordinates
[in]boxsimulation box
Returns
the shifted coordinates
55 
56  implicit none
57  real(kind=real64), intent(in) :: a(3), box(6)
58  real(kind=real64) :: pbc(3), tbox(3,3)
59  integer(kind=int32) :: shift
60 
61  tbox = 0.0d0
62 
63  ! A = box(1), length of unit cell vector along x-axis;
64  ! B = box(2), length of unit cell vector in xy-plane;
65  ! C = box(3), length of unit cell vector in yz-plane;
66  ! alpha = box(4), cosine of angle between B and C
67  ! beta = box(5), cosine of angle between A and C
68  ! gamma = box(6), cosine of angle between A and B
69 
70  ! New matrix made up of box vectors:
71 
72  ! ax bx cx
73  ! 0 by cy
74  ! 0 0 cz
75 
76  ! convert angles to box vectors
77 
78  ! A**2 = ax**2 + ay**2 + az**2
79  ! A**2 = ax**2 + (0)**2 + (0)**2
80  ! ax = A
81  tbox(1,1) = box(1)
82 
83  ! dot(A_vec,B_vec) = ax*bx + ay*by + az*bz
84  ! dot(A_vec,B_vec) = ax*bx + (0)*by + (0)*bz
85  ! A*B*gamma = ax*bx
86  ! bx = B*gamma
87  tbox(1,2) = box(2)*box(6)
88 
89  ! dot(A_vec,C_vec) = ax*cx + ay*cy + az*cz
90  ! dot(A_vec,C_vec) = ax*cx + (0)*cy + (0)*cz
91  ! A*C*beta = ax*cx
92  ! cx = C*beta
93  tbox(1,3) = box(3)*box(5)
94 
95  ! B**2 = bx**2 + by**2 + bz**2
96  ! B**2 = bx**2 + by**2 + (0)**2
97  ! by = sqrt(B**2 - bx**2)
98  tbox(2,2) = dsqrt(box(2)**2-tbox(1,2)**2)
99 
100  ! dot(B_vec,C_vec) = bx*cx + by*cy + (0)*cz
101  ! B*C*alpha = bx*cx + by*cy
102  ! cy = (B*C*alpha - bx*cx) / by
103  tbox(2,3) = (box(2)*box(3)*box(4) - tbox(1,2)*tbox(1,3))/tbox(2,2)
104 
105  ! C = cx**2 + cy**2 + cz**2
106  ! cz = sqrt(C**2 - cx**2 - cy**2)
107  tbox(3,3) = dsqrt(box(3)**2-tbox(1,3)**2-tbox(2,3)**2)
108 
109  pbc = a
110 
111  ! Now perform PBC calculation
112  shift = nint(pbc(3) / tbox(3,3))
113  if (shift .ne. 0) then
114  pbc(3) = pbc(3) - tbox(3,3) * shift
115  pbc(2) = pbc(2) - tbox(2,3) * shift
116  pbc(1) = pbc(1) - tbox(1,3) * shift
117  end if
118 
119  shift = nint(pbc(2) / tbox(2,2))
120  if (shift .ne. 0) then
121  pbc(2) = pbc(2) - tbox(2,2) * shift
122  pbc(1) = pbc(1) - tbox(1,2) * shift
123  end if
124 
125  shift = nint(pbc(1) / tbox(1,1))
126  if (shift .ne. 0) then
127  pbc(1) = pbc(1) - tbox(1,1) * shift
128  end if
129 
Generated by   doxygen 1.8.14