From 84a84e578a28f7bdbf652cf51b747d04b9c55a72 Mon Sep 17 00:00:00 2001
From: Alex Selimov <alex@alexselimov.xyz>
Date: Tue, 27 Oct 2020 16:44:03 -0400
Subject: [PATCH] Working version of slip_plane ooption with effill

---
 src/opt_slip_plane.f90 | 18 +++++++++++++-----
 1 file changed, 13 insertions(+), 5 deletions(-)

diff --git a/src/opt_slip_plane.f90 b/src/opt_slip_plane.f90
index c6c8468..e607ae8 100644
--- a/src/opt_slip_plane.f90
+++ b/src/opt_slip_plane.f90
@@ -20,7 +20,8 @@ module opt_slip_plane
         integer :: ie, ia, slip_enum, old_atom_num, esize, new_ele_num, n, m, o, ele(3,8), nump_ele, inod, vlat(3), ibasis
                 
         integer, allocatable :: slip_eles(:), temp_int(:)
-        real(kind=dp) :: r_interp(3, max_basisnum*max_esize**3), rfill(3,max_basisnum, max_ng_node), ratom(3,max_basisnum)
+        real(kind=dp) :: r_interp(3, max_basisnum*max_esize**3), rfill(3,max_basisnum, max_ng_node), ratom(3,max_basisnum), &
+                         maxp, minp
 
         integer :: type_interp(max_basisnum*max_esize**3)
         logical :: lat_points(max_esize,max_esize, max_esize)
@@ -77,23 +78,29 @@ module opt_slip_plane
                             ele = cubic_cell*(esize -1) 
                             do o = 1, size_ele(ie) - esize
                                 do n = 1, size_ele(ie) - esize
-                                    do m = 1, size_ele(ie) - esize
+                                    latloop:do m = 1, size_ele(ie) - esize
                                         do inod = 1, ng_node(lat_ele(ie))
                                             vlat = ele(:,inod) + (/ m, n, o /)
+                                            if (.not.lat_points(vlat(1), vlat(2),vlat(3))) cycle latloop
                                             call get_interp_pos(vlat(1), vlat(2), vlat(3), ie, rfill(:,:,inod))
                                         end do
+
+                                        !Check to make sure all lattice points exist for the current element
+                                        if(any(lat_points(m:m+esize-1, n:n+esize-1, o:o+esize-1) == 0)) cycle latloop 
                                         !Check to see if the plane intersects this element if not then add it 
-                                        if((spos < maxval(rfill(sdim,1:basisnum(lat_ele(ie)),1:ng_node(lat_ele(ie))))).and. &
-                                            (spos > minval(rfill(sdim,1:basisnum(lat_ele(ie)),1:ng_node(lat_ele(ie)))))) then 
+                                        maxp = maxval(rfill(sdim,1:basisnum(lat_ele(ie)),1:ng_node(lat_ele(ie))))
+                                        minp = minval(rfill(sdim,1:basisnum(lat_ele(ie)),1:ng_node(lat_ele(ie))))
+                                        if(.not.(spos < maxp).and.(spos > minp))then
                                             nump_ele = nump_ele - esize**3 
                                             lat_points(m:m+esize-1, n:n+esize-1, o:o+esize-1) = .false.
                                             call add_element(0, type_ele(ie), esize, lat_ele(ie), sbox_ele(ie), rfill) 
                                             new_ele_num = new_ele_num + 1
                                         end if
-                                    end do
+                                    end do latloop
                                 end do
                             end do
                         end if
+                        esize= esize-2
                     end do
                     ! Now add the leftover lattice points as atoms
                     do o = 1, size_ele(ie)
@@ -102,6 +109,7 @@ module opt_slip_plane
                                 if(lat_points(m,n,o)) then 
                                     call get_interp_pos(m,n,o, ie, ratom(:,:))
                                     do ibasis = 1, basisnum(lat_ele(ie))
+                                        call apply_periodic(r_atom(:,ibasis))
                                         call add_atom(0, basis_type(ibasis,lat_ele(ie)), sbox_ele(ie), ratom(:,ibasis))
                                     end do
                                 end if