From fd143783ec889798da0377052c0ddf8001f1794a Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Tue, 28 Jan 2020 09:36:00 -0500 Subject: [PATCH] Added remesh option to group --- src/elements.f90 | 36 +++++++++++++++++++++++++++++++++- src/opt_group.f90 | 49 ++++++++++++++++++++++++++++++++++++++++++++++- 2 files changed, 83 insertions(+), 2 deletions(-) diff --git a/src/elements.f90 b/src/elements.f90 index b4e59be..9aa5f54 100644 --- a/src/elements.f90 +++ b/src/elements.f90 @@ -223,6 +223,8 @@ module elements real(kind=dp), intent(in) :: r(3, max_basisnum, max_ng_node) ele_num = ele_num + 1 + !Check to see if we need to grow the arrays + call grow_ele_arrays(1,0) type_ele(ele_num) = type size_ele(ele_num) = size lat_ele(ele_num) = lat @@ -239,6 +241,8 @@ module elements real(kind=dp), intent(in), dimension(3) :: r atom_num = atom_num+1 + !Check to see if we need to grow the arrays + call grow_ele_arrays(0,1) type_atom(atom_num) = type r_atom(:,atom_num) = r(:) @@ -356,7 +360,7 @@ module elements end do end select - if (ia /= esize**3) then + if (ia /= bnum*esize**3) then print *, "Incorrect interpolation" stop 3 end if @@ -403,4 +407,34 @@ module elements atom_num = atom_num - 1 end do end subroutine + + subroutine delete_elements(num, index) + !This subroutine deletes elements from the element array + integer, intent(in) :: num + integer, intent(inout), dimension(num) :: index + + integer :: i, j + + call heapsort(index) + + !We go from largest index to smallest index just to make sure that we don't miss + !accidentally overwrite values which need to be deleted + do i = num, 1, -1 + if(index(i) == ele_num) then + r_node(:,:,:,index(i)) = 0.0_dp + type_ele(index(i)) ='' + size_ele(index(i)) = 0 + lat_ele(index(i)) = 0 + sbox_ele(index(i)) = 0 + else + node_num = node_num - ng_node(lat_ele(index(i))) + r_node(:,:,:,index(i)) = r_node(:,:,:,ele_num) + type_ele(index(i)) = type_ele(ele_num) + size_ele(index(i)) = size_ele(ele_num) + lat_ele(index(i)) = lat_ele(ele_num) + sbox_ele(index(i)) = sbox_ele(ele_num) + end if + ele_num = ele_num - 1 + end do + end subroutine delete_elements end module elements \ No newline at end of file diff --git a/src/opt_group.f90 b/src/opt_group.f90 index 13b5b03..25d8671 100644 --- a/src/opt_group.f90 +++ b/src/opt_group.f90 @@ -8,7 +8,7 @@ module opt_group use box implicit none - integer :: group_ele_num, group_atom_num + integer :: group_ele_num, group_atom_num, remesh_size character(len=15) :: type, shape !Type indicates what element type is selected and shape is the group shape real(kind=dp) :: block_bd(6), disp_vec(3) logical :: displace, wrap @@ -24,6 +24,7 @@ module opt_group group_ele_num = 0 group_atom_num = 0 + remesh_size=0 if(allocated(element_index)) deallocate(element_index) if(allocated(atom_index)) deallocate(atom_index) @@ -34,6 +35,7 @@ module opt_group !Now call the transformation functions for the group if(displace) call displace_group + if(remesh_size > 0) call remesh_group end subroutine group subroutine parse_group(arg_pos) @@ -94,6 +96,11 @@ module opt_group end do case('wrap') wrap = .true. + case('remesh') + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder, arglen) + if (arglen==0) stop "Missing remesh element size in group command" + read(textholder, *) remesh_size case default !If it isn't an available option to opt_disl then we just exit exit @@ -197,4 +204,44 @@ module opt_group end subroutine displace_group + subroutine remesh_group + !This command is used to remesh the group to a desired element size + + integer :: i, j, ie, type_interp(max_basisnum*max_esize**3), add_atom_num + real(kind=dp) :: r_interp(3, max_basisnum*max_esize**3) + + !Refining to atoms and remeshing to elements are different processes so check which code we need to run + select case(remesh_size) + + !Refining to atoms + case(2) + if(group_ele_num > 0) then + !Estimate number of atoms we are adding, this doesn't have to be exact + add_atom_num = group_ele_num*basisnum(lat_ele(element_index(1)))*size_ele(element_index(1))**3 + call grow_ele_arrays(0,add_atom_num) + + do i = 1, group_ele_num + ie = element_index(i) + !Get the interpolated atom positions + call interpolate_atoms(type_ele(ie), size_ele(ie), lat_ele(ie), r_node(:,:,:,ie), type_interp, r_interp) + + !Loop over all interpolated atoms and add them to the system, we apply periodic boundaries here as well to make sure + !they are in the box + do j = 1, basisnum(lat_ele(ie))*size_ele(ie)**3 + call apply_periodic(r_interp(:,j)) + call add_atom(type_interp(j), r_interp(:,j)) + end do + end do + + !Once all atoms are added we delete all of the elements + call delete_elements(group_ele_num, element_index) + + end if + !Remeshing to elements, currently not available + case default + print *, "Remeshing to elements is currently not available. Please refine to atoms by passing a remsh size of 2" + end select + + + end subroutine remesh_group end module opt_group \ No newline at end of file