diff --git a/src/opt_group.f90 b/src/opt_group.f90 index 73ea62c..9a01f1f 100644 --- a/src/opt_group.f90 +++ b/src/opt_group.f90 @@ -302,12 +302,12 @@ module opt_group !This command is used to remesh the group to a desired element size integer :: i, j, k, ix, iy, iz, inod, ibasis, ie, type_interp(max_basisnum*max_esize**3), add_atom_num, orig_atom_num, & - current_esize, dof, max_lat(3), r_lat(3), ele(3,8), vlat(3), bd_in_lat(6), bd_in_array(3), old_ele, old_atom, & - max_loops, working_esize, group_lat_num, lat_list(10), group_sbox_num, sbox_list(100), is, ilat + current_esize, dof, max_lat(3), r_lat(3), ele(3,8), vlat(3), bd_in_lat(6), bd_in_array(3), new_ele, new_atom, & + max_loops, working_esize, group_lat_num, lat_list(10), group_sbox_num, sbox_list(100), is, ilat, tot_dof real(kind=dp) :: r_interp(3, max_basisnum*max_esize**3), ori_inv(3,3), r(3), & r_new_node(3,max_basisnum, max_ng_node), orient(3,3), group_in_lat(3,8), group_bd(6) - logical, allocatable :: lat_points(:,:,:) + logical, allocatable :: lat_points(:,:,:), new_sbox, new_lat character(len=100) :: remesh_ele_type @@ -318,35 +318,67 @@ module opt_group group_sbox_num = 0 group_lat_num = 0 do i = 1, group_atom_num + new_sbox=.true. + new_lat=.true. do j = 1, group_sbox_num - if (sbox_list(j) == sbox_atom(atom_index(i))) exit + if (sbox_list(j) == sbox_atom(atom_index(i))) then + new_sbox=.false. + exit + end if + end do + + if(new_sbox) then group_sbox_num = group_sbox_num + 1 sbox_list(group_sbox_num) = sbox_atom(atom_index(i)) - end do + end if do j = 1, group_lat_num - if (basis_type(1,lat_list(j)) == type_atom(atom_index(i))) exit + if (basis_type(1,lat_list(j)) == type_atom(atom_index(i))) then + new_lat = .false. + exit + end if + end do + + if (new_lat) then group_lat_num = group_lat_num + 1 do k = 1, lattice_types - if (basis_type(1,k) == type_atom(atom_index(i))) lat_list(group_lat_num) = k + if (basis_type(1,k) == type_atom(atom_index(i))) lat_list(group_lat_num) = k end do - end do + end if + end do do i = 1, group_ele_num + new_sbox=.true. + new_lat = .true. do j = 1, group_sbox_num - if (sbox_list(j) == sbox_ele(element_index(i))) exit - group_sbox_num = group_sbox_num + 1 - sbox_list(group_sbox_num) = sbox_ele(element_index(i)) + if (sbox_list(j) == sbox_ele(element_index(i))) then + new_sbox = .false. + exit + end if end do + + if (new_sbox) then + group_sbox_num = group_sbox_num + 1 + sbox_list(group_sbox_num) = sbox_ele(element_index(i)) + end if do j = 1, group_lat_num - if (lat_list(group_lat_num) == lat_ele(element_index(i))) exit + if (lat_list(group_lat_num) == lat_ele(element_index(i))) then + new_lat=.false. + exit + end if + end do + + if (new_lat) then group_lat_num = group_lat_num + 1 lat_list(group_lat_num) = lat_ele(element_index(i)) - end do + end if end do + new_atom = 0 + new_ele=0 + tot_dof=0 do is = 1, group_sbox_num orient = sub_box_ori(:, :, sbox_list(is)) @@ -388,6 +420,7 @@ module opt_group !If for some reason there are no dof in this loop then cycle out if(dof == 0) cycle + tot_dof = tot_dof+dof group_in_lat = reshape((/ group_bd(1),group_bd(3),group_bd(5), & group_bd(2),group_bd(3),group_bd(5), & @@ -410,6 +443,7 @@ module opt_group !Now place all group atoms and group interpolated atoms into lat_points do i = 1, group_atom_num + if (.not.((sbox_atom(atom_index(i)) == is).and.(basis_type(1,ilat) == type_atom(atom_index(i))))) cycle r = r_atom(:,atom_index(i))/lapa(ilat) r = matmul(fcc_inv,matmul(ori_inv,r)) do j = 1, 3 @@ -427,6 +461,7 @@ module opt_group !Now place interpolated atoms within lat_points array do i =1, group_ele_num + if (.not.((sbox_ele(element_index(i)) == is).and.( lat_ele(element_index(i)) == ilat))) cycle ie = element_index(i) call interpolate_atoms(type_ele(ie), size_ele(ie), lat_ele(ie), r_node(:,:,:,ie), type_interp, r_interp) do j = 1, size_ele(ie)**3 * basisnum(lat_ele(ie)) @@ -446,17 +481,6 @@ module opt_group end do end do - print *, "Group has ", dof, " degrees of freedom to remesh" - - !Delete all elements and atoms to make space for new elements and atoms - call delete_atoms(group_atom_num, atom_index) - call delete_elements(group_ele_num, element_index) - - old_atom = atom_num - old_ele = ele_num - - - !Now run remeshing algorithm, not the most optimized or efficient but gets the job done !Figure out new looping boundaries bd_in_array(1) = bd_in_lat(2) - bd_in_lat(1) + 10 @@ -499,6 +523,7 @@ module opt_group lat_points(ix:ix+working_esize-1,iy:iy+working_esize-1,iz:iz+working_esize-1) = .false. !Add the element, for the sbox we just set it to the same sbox that we get the orientation !from. In this case it is from the sbox of the first atom in the group. + new_ele = new_ele+1 call add_element(remesh_ele_type, working_esize, ilat, sbox_atom(atom_index(1)),r_new_node) end if @@ -519,15 +544,22 @@ module opt_group end do lat_points(ix,iy,iz) = .false. r = matmul(orient, matmul(fcc_mat, vlat))*lapa(ilat) - call add_atom(basis_type(1,ilat), sbox_atom(atom_index(1)), r) + new_atom=new_atom+1 + call add_atom(basis_type(1,ilat), is, r) end if end do end do end do + deallocate(lat_points) end do end do - print *, "Remeshing has created ", ele_num-old_ele, " elements and ", atom_num-old_atom, " atoms." + !Delete all elements and atoms to make space for new elements and atoms + call delete_atoms(group_atom_num, atom_index) + call delete_elements(group_ele_num, element_index) + + print *, tot_dof, " degrees of freedom in group" + print *, "remesh_group creates ", new_atom, " atoms and ", new_ele, " elements." end subroutine remesh_group subroutine delete_group