Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

v0.1.14 #16

Merged
merged 1 commit into from
May 24, 2024
Merged

v0.1.14 #16

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
59 changes: 29 additions & 30 deletions fast_grid/libs/distance_matrix.pyx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,51 +8,50 @@ cimport numpy as np
from libc.math cimport sqrt


def minimum_image_triclinic(np.ndarray[np.float64_t, ndim=1] dx, np.ndarray[np.float64_t, ndim=2] box):
cdef int ix, iy, iz
cdef float rx, ry0, ry1, rz0, rz1, rz2, dsq
cdef float dsq_min = np.finfo(np.float64).max
cdef np.ndarray[np.float64_t, ndim=1] dx_min = np.zeros(3, dtype=np.float64)

for ix in range(-1, 2):
rx = dx[0] + box[0, 0] * ix
for iy in range(-1, 2):
ry0 = rx + box[1, 0] * iy
ry1 = dx[1] + box[1, 1] * iy
for iz in range(-1, 2):
rz0 = ry0 + box[2, 0] * iz
rz1 = ry1 + box[2, 1] * iz
rz2 = dx[2] + box[2, 2] * iz
dsq = rz0 * rz0 + rz1 * rz1 + rz2 * rz2
if dsq < dsq_min:
dsq_min = dsq
dx_min[0] = rz0
dx_min[1] = rz1
dx_min[2] = rz2

dx[:] = dx_min


def distance_matrix_triclinic_cython(np.ndarray[np.float64_t, ndim=2] pos1,
np.ndarray[np.float64_t, ndim=2] pos2,
np.ndarray[np.float64_t, ndim=2] box):
cdef int i, j
cdef int n = pos1.shape[0]
cdef int m = pos2.shape[0]
cdef double r2

cdef int ix, iy, iz
cdef float rx, ry0, ry1, rz0, rz1, rz2, dsq
cdef float dsq_min = np.finfo(np.float64).max

cdef np.ndarray[np.float64_t, ndim=2] distances = np.zeros((n, m), dtype=np.float64)
cdef np.ndarray[np.float64_t, ndim=1] diff = np.zeros(3, dtype=np.float64)
cdef np.ndarray[np.float64_t, ndim=1] dx_min = np.zeros(3, dtype=np.float64)

for i in prange(n, nogil=True):
for j in range(m):
diff[0] = pos2[j, 0] - pos1[i, 0]
diff[1] = pos2[j, 1] - pos1[i, 1]
diff[2] = pos2[j, 2] - pos1[i, 2]

with gil:
minimum_image_triclinic(diff, box)
r2 = diff[0] * diff[0] + diff[1] * diff[1] + diff[2] * diff[2]

dsq_min = 10000.
dx_min[0] = 0
dx_min[1] = 0
dx_min[2] = 0

for ix in range(-1, 2):
rx = diff[0] + box[0, 0] * ix
for iy in range(-1, 2):
ry0 = rx + box[1, 0] * iy
ry1 = diff[1] + box[1, 1] * iy
for iz in range(-1, 2):
rz0 = ry0 + box[2, 0] * iz
rz1 = ry1 + box[2, 1] * iz
rz2 = diff[2] + box[2, 2] * iz
dsq = rz0 * rz0 + rz1 * rz1 + rz2 * rz2
if dsq < dsq_min:
dsq_min = dsq
dx_min[0] = rz0
dx_min[1] = rz1
dx_min[2] = rz2

r2 = dx_min[0] * dx_min[0] + dx_min[1] * dx_min[1] + dx_min[2] * dx_min[2]
distances[i, j] = sqrt(r2)

return distances
return distances
97 changes: 58 additions & 39 deletions fast_grid/libs/potential.pyx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,34 +8,6 @@ cimport numpy as np
from libc.math cimport exp


def minimum_image_triclinic(np.ndarray[np.float64_t, ndim=1] dx, np.ndarray[np.float64_t, ndim=2] box):
cdef int ix, iy, iz
cdef float rx, ry0, ry1, rz0, rz1, rz2, dsq
cdef float dsq_min = np.finfo(np.float64).max
cdef np.ndarray[np.float64_t, ndim=1] dx_min = np.zeros(3, dtype=np.float64)

for ix in range(-1, 2):
rx = dx[0] + box[0, 0] * ix
for iy in range(-1, 2):
ry0 = rx + box[1, 0] * iy
ry1 = dx[1] + box[1, 1] * iy
for iz in range(-1, 2):
rz0 = ry0 + box[2, 0] * iz
rz1 = ry1 + box[2, 1] * iz
rz2 = dx[2] + box[2, 2] * iz
dsq = rz0 * rz0 + rz1 * rz1 + rz2 * rz2
if dsq < dsq_min:
dsq_min = dsq
dx_min[0] = rz0
dx_min[1] = rz1
dx_min[2] = rz2

dx[:] = dx_min


@cython.boundscheck(False)
@cython.wraparound(False)
@cython.cdivision(True)
def lj_potential_cython(np.ndarray[np.float64_t, ndim=2] pos1,
np.ndarray[np.float64_t, ndim=2] pos2,
np.ndarray[np.float64_t, ndim=2] box,
Expand All @@ -48,6 +20,11 @@ def lj_potential_cython(np.ndarray[np.float64_t, ndim=2] pos1,
cdef double r2, lj6, lj12, inv_r2, inv_r6, inv_r12, e, s, s6, s12, energy
cdef double threshold = 1e-10
cdef double cutoff_squared = cutoff * cutoff

cdef int ix, iy, iz
cdef float rx, ry0, ry1, rz0, rz1, rz2, dsq
cdef float dsq_min = np.finfo(np.float64).max
cdef np.ndarray[np.float64_t, ndim=1] dx_min = np.zeros(3, dtype=np.float64)

cdef np.ndarray[np.float64_t, ndim=1] energy_grid = np.zeros((n), dtype=np.float64)
cdef np.ndarray[np.float64_t, ndim=1] diff = np.zeros(3, dtype=np.float64)
Expand All @@ -58,11 +35,31 @@ def lj_potential_cython(np.ndarray[np.float64_t, ndim=2] pos1,
diff[0] = pos2[j, 0] - pos1[i, 0]
diff[1] = pos2[j, 1] - pos1[i, 1]
diff[2] = pos2[j, 2] - pos1[i, 2]

with gil:
minimum_image_triclinic(diff, box)
r2 = diff[0] * diff[0] + diff[1] * diff[1] + diff[2] * diff[2]

dsq_min = 10000.
dx_min[0] = 0
dx_min[1] = 0
dx_min[2] = 0

# minimum_image_triclinic
for ix in range(-1, 2):
rx = diff[0] + box[0, 0] * ix
for iy in range(-1, 2):
ry0 = rx + box[1, 0] * iy
ry1 = diff[1] + box[1, 1] * iy
for iz in range(-1, 2):
rz0 = ry0 + box[2, 0] * iz
rz1 = ry1 + box[2, 1] * iz
rz2 = diff[2] + box[2, 2] * iz
dsq = rz0 * rz0 + rz1 * rz1 + rz2 * rz2
if dsq < dsq_min:
dsq_min = dsq
dx_min[0] = rz0
dx_min[1] = rz1
dx_min[2] = rz2

r2 = dx_min[0] * dx_min[0] + dx_min[1] * dx_min[1] + dx_min[2] * dx_min[2]

if r2 < cutoff_squared and r2 > threshold:
e = epsilon[j]
s = sigma[j]
Expand All @@ -81,9 +78,6 @@ def lj_potential_cython(np.ndarray[np.float64_t, ndim=2] pos1,
return energy_grid


@cython.boundscheck(False)
@cython.wraparound(False)
@cython.cdivision(True)
def gaussian_cython(np.ndarray[np.float64_t, ndim=2] pos1,
np.ndarray[np.float64_t, ndim=2] pos2,
np.ndarray[np.float64_t, ndim=2] box,
Expand All @@ -97,6 +91,11 @@ def gaussian_cython(np.ndarray[np.float64_t, ndim=2] pos1,
cdef double threshold = 1e-10
cdef double width_squared = width * width
cdef double cutoff_squared = cutoff * cutoff

cdef int ix, iy, iz
cdef float rx, ry0, ry1, rz0, rz1, rz2, dsq
cdef float dsq_min = np.finfo(np.float64).max
cdef np.ndarray[np.float64_t, ndim=1] dx_min = np.zeros(3, dtype=np.float64)

cdef np.ndarray[np.float64_t, ndim=1] energy_grid = np.zeros((n), dtype=np.float64)
cdef np.ndarray[np.float64_t, ndim=1] diff = np.zeros(3, dtype=np.float64)
Expand All @@ -107,10 +106,30 @@ def gaussian_cython(np.ndarray[np.float64_t, ndim=2] pos1,
diff[0] = pos2[j, 0] - pos1[i, 0]
diff[1] = pos2[j, 1] - pos1[i, 1]
diff[2] = pos2[j, 2] - pos1[i, 2]

with gil:
minimum_image_triclinic(diff, box)
r2 = diff[0] * diff[0] + diff[1] * diff[1] + diff[2] * diff[2]

dsq_min = 10000.
dx_min[0] = 0
dx_min[1] = 0
dx_min[2] = 0

# minimum_image_triclinic
for ix in range(-1, 2):
rx = diff[0] + box[0, 0] * ix
for iy in range(-1, 2):
ry0 = rx + box[1, 0] * iy
ry1 = diff[1] + box[1, 1] * iy
for iz in range(-1, 2):
rz0 = ry0 + box[2, 0] * iz
rz1 = ry1 + box[2, 1] * iz
rz2 = diff[2] + box[2, 2] * iz
dsq = rz0 * rz0 + rz1 * rz1 + rz2 * rz2
if dsq < dsq_min:
dsq_min = dsq
dx_min[0] = rz0
dx_min[1] = rz1
dx_min[2] = rz2

r2 = dx_min[0] * dx_min[0] + dx_min[1] * dx_min[1] + dx_min[2] * dx_min[2]

if r2 < cutoff_squared and r2 > threshold:
energy += height * exp(-1 * r2 / width_squared)
Expand Down
4 changes: 2 additions & 2 deletions setup.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,10 @@

setup(
name="fast-grid",
version="0.1.13",
version="0.1.14",
description="Fast grid calculation",
author="Hyunsoo Park",
author_email="hpark@ic.ac.uk",
author_email="phs68660888@gmail.com",
url="https://github.com/hspark1212/fast-grid.git",
install_requires=[
"ase",
Expand Down
47 changes: 9 additions & 38 deletions tutorial.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

Loading