Merge pull request #10 from aselimov/bg--fix-to-pycac-format

Fixes to element definitions

src/io.f90

@@ -363,7 +363,7 @@ module io
         !NOTE: This code doesn't work for arbitrary number of basis atoms per node. It assumes that the 
         !each element has only 1 atom type at the node.
         character(len=100), intent(in) :: file
-        integer :: interp_max, i, j,  lat_size, inod, ibasis, ip
+        integer :: interp_max, i, j, lat_size, inod, ibasis, ip, unique_index(10), unique_num
         real(kind=dp) :: box_vec(3)
 
 1  format('time' / i16, f23.15)
@@ -442,14 +442,21 @@ module io
         !write the element information
         if(ele_num > 0) then
             write(11,12)
-            do i = 1, lattice_types
+            !First figure out all of the unique element types
-                do j = 1, ele_num
+            unique_num = 0
-                    if (lat_ele(j) == i) then 
+            unique_index(:) = 0
-                        lat_size = size_ele(j)-1
+            eleloop:do i = 1, ele_num
-                        exit
+                do j =1 , unique_num
+                    if ( ( size_ele(i) == size_ele( unique_index(j) ) ).and. &
+                        ( lat_ele(i) == lat_ele(unique_index(j)) ) ) then
+                        cycle eleloop
                     end if
                 end do
-                write(11,'(3i16)') i, lat_size, basis_type(1,i)
+                unique_num = unique_num + 1
+                unique_index(unique_num) = i
+            end do eleloop
+            do i = 1, unique_num
+                write(11,'(3i16)') i, size_ele(i)-1, basis_type(1,i)
             end do
             ip = 0
             write(11,13)
@@ -600,8 +607,8 @@ module io
         real(kind=dp), dimension(3), intent(in) :: displace
         real(kind = dp), dimension(6), intent(out) :: temp_box_bd
 
-        integer :: i, j, k, n, inod, ibasis, type, size, in_atoms, in_eles, new_atom_types, &
+        integer :: i, j, k, l, n, inod, ibasis, type, size, in_atoms, in_eles, new_atom_types, &
-                   new_type_to_type(10), new_lattice_types, sbox
+                   new_type_to_type(10), new_lattice_types, sbox, new_lattice_map(10)
         character(len=100) :: etype
         real(kind=dp) ::  r(3), newdisplace(3)
         real(kind=dp), allocatable :: r_innode(:,:,:)
@@ -665,6 +672,30 @@ module io
         end do
         !Read the lattice parameters for every lattice type
         read(11,*) (lapa(i), i = lattice_types+1, lattice_types+new_lattice_types)
+
+        !Now we loop over all new lattice types and check to see if they are exactly the same as any old lattice types
+        !If they are then we map the new type to the old type.
+        k = lattice_types + 1
+        new_lattice_map(:) = 0
+        new_loop:do i = 1, new_lattice_types
+            old_loop:do j = 1, lattice_types
+                !First check all the lattice level variables
+                if ((basisnum(i) == basisnum(j)).and. &
+                    (ng_node(i) == ng_node(j))) then 
+                    !Now check the basis  level variables
+                    do ibasis =1, basisnum(i)
+                        if(basis_type(ibasis,i) /= basis_type(ibasis,j)) then 
+                            cycle old_loop
+                        end if
+                    end do
+                    new_lattice_map(i) = j
+                    cycle new_loop
+                end if
+            end do old_loop
+            new_lattice_map(i) = k 
+            k = k+1
+        end do new_loop
+
         !Read number of elements and atoms and allocate arrays
         read(11, *) in_atoms, in_eles
         call grow_ele_arrays(in_eles, in_atoms)
@@ -678,15 +709,14 @@ module io
 
         !Read the elements
         do i = 1, in_eles
-            read(11, *) n, type, size, sbox, etype
+            read(11, *) l, type, size, sbox, etype
             do inod = 1, ng_node(type)
                 do ibasis =1, basisnum(type)
                     read(11,*) j, k, r_innode(:, ibasis, inod)
                     r_innode(:,ibasis,inod) = r_innode(:, ibasis, inod) + newdisplace
                 end do
             end do
-            type = type + lattice_types
+            call add_element(etype, size, new_lattice_map(type), sbox+n, r_innode)
-            call add_element(etype, size, type, sbox+n, r_innode)
         end do
 
         !Close the file being read
@@ -694,7 +724,7 @@ module io
 
         !Only increment the lattice types if there are elements, if there are no elements then we 
         !just overwrite the arrays
-        if(in_eles > 0) lattice_types = lattice_types + new_lattice_types
+        if(in_eles > 0) lattice_types =  maxval(new_lattice_map)
 
 
         sub_box_num = sub_box_num + n