From 955ff888b5f80032d70099eaf3e0d349f811f087 Mon Sep 17 00:00:00 2001 From: "Igor S. Gerasimov" Date: Wed, 27 Nov 2024 05:24:35 +0100 Subject: [PATCH] Compute once those that does not depend on k Signed-off-by: Igor S. Gerasimov --- src/constrain_pot.f90 | 33 ++++++++++++++++----------------- 1 file changed, 16 insertions(+), 17 deletions(-) diff --git a/src/constrain_pot.f90 b/src/constrain_pot.f90 index 335fb8264..71fad3460 100644 --- a/src/constrain_pot.f90 +++ b/src/constrain_pot.f90 @@ -106,14 +106,14 @@ subroutine qpothess2(fix,n,at,xyz,h) ii = (fix%atoms(i)-1)*3 do j = 1, fix%n !inner loop for sum over all atoms if (i.ne.j) then + rij = xyz(:,fix%atoms(j))-xyz(:,fix%atoms(i)) + r = norm2(rij) + r2 = r*r + r3 = r2*r + ij = lin(i-1,j-1)+1 + r0 = fix%val(ij) + dr = r-r0 do k = 1, 3 !loop diagonal elements (xyz components) - rij = xyz(:,fix%atoms(j))-xyz(:,fix%atoms(i)) - r = norm2(rij) - r2 = r*r - r3 = r2*r - ij = lin(i-1,j-1)+1 - r0 = fix%val(ij) - dr = r-r0 dx = xyz(k,fix%atoms(i))-xyz(k,fix%atoms(j)) dx2 = dx*dx iik = lin(ii+k,ii+k) @@ -127,23 +127,22 @@ subroutine qpothess2(fix,n,at,xyz,h) endif enddo !end inner loop for sum over all atoms do j = i+1, fix%n !loop over the rest (mixed atoms) + rij = xyz(:,fix%atoms(j))-xyz(:,fix%atoms(i)) + r = norm2(rij) + r2 = r*r + r3 = r2*r + ij = lin(i-1,j-1)+1 + r0 = fix%val(ij) + dr = r-r0 + jj = (fix%atoms(j)-1)*3 do k = 1, 3 !loop diagonal elements (xyz components) - rij = xyz(:,fix%atoms(j))-xyz(:,fix%atoms(i)) - r = norm2(rij) - r2 = r*r - r3 = r2*r - ij = lin(i-1,j-1)+1 - r0 = fix%val(ij) - dr = r-r0 - jj = (fix%atoms(j)-1)*3 + dx = xyz(k,fix%atoms(i))-xyz(k,fix%atoms(j)) do m = 1, 3 if (k.eq.m) then !same component case - dx = xyz(k,fix%atoms(i))-xyz(k,fix%atoms(j)) dx2 = dx*dx ijk = lin(ii+k,jj+k) h(ijk) = h(ijk) - 2.0_wp*fix%fc*(1.0_wp+dx2/r2-dx2*dr/r3-r0/r) else !different component case - dx = xyz(k,fix%atoms(i))-xyz(k,fix%atoms(j)) dy = xyz(m,fix%atoms(i))-xyz(m,fix%atoms(j)) ikjm = lin(ii+k,jj+m) h(ikjm) = h(ikjm) - 2.0_wp*fix%fc*r0*dx*dy/r3