From 14dd8ad755685add1d7db37614b3d42b3df89d92 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Wed, 15 Jan 2020 16:03:18 -0500 Subject: [PATCH 01/18] Fixes to mode_merge --- src/mode_convert.f90 | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/mode_convert.f90 b/src/mode_convert.f90 index 1f8a6bd..32ded60 100644 --- a/src/mode_convert.f90 +++ b/src/mode_convert.f90 @@ -13,8 +13,6 @@ module mode_convert integer, intent(out) :: arg_pos character(len=100) :: infile, outfile real(kind = dp) :: temp_box_bd(6) - !We have to allocate the element and atom arrays with a size of 1 for the read in code to work - call alloc_ele_arrays(1,1) !First read in the file call get_command_argument(2, infile) call get_in_file(infile) From cde96402e29578a89bbfadea2297292f4d6d3144 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Wed, 15 Jan 2020 16:44:13 -0500 Subject: [PATCH 02/18] Working shift command for merge --- src/io.f90 | 1 + src/mode_merge.f90 | 76 +++++++++++++++++++++++++++++++++++++++++++-- src/subroutines.f90 | 19 ++++++++++++ 3 files changed, 94 insertions(+), 2 deletions(-) diff --git a/src/io.f90 b/src/io.f90 index e5b098d..99f9e96 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -194,6 +194,7 @@ module io select case(trim(adjustl(type_ele(i)))) case('fcc') do iatom = 1, basisnum(lat_ele(i))*size_ele(i)**3 + call apply_periodic(r_interp(:,iatom)) write(11, '(2i16, 3f23.15)') atom_num+iatom, type_interp(iatom), r_interp(:,iatom) end do end select diff --git a/src/mode_merge.f90 b/src/mode_merge.f90 index cfb4273..567b19b 100644 --- a/src/mode_merge.f90 +++ b/src/mode_merge.f90 @@ -8,7 +8,9 @@ module mode_merge use elements character(len=4) :: dim - integer :: in_num + integer :: in_num, new_starts(2) + real(kind=dp) :: shift_vec(3) + logical :: wrap, shift_flag public contains @@ -18,12 +20,23 @@ module mode_merge integer :: i real(kind=dp) :: displace(3), temp_box_bd(6) + wrap = .false. + shift_flag = .false. + + shift_vec(:) = 0.0_dp !First we parse the merge command call parse_command(arg_pos) !Now loop over all files and stack them do i = 1, in_num displace(:) = 0.0_dp + + !The new starts variable dictate where in the atom and element array each new + !file starts. This is used for additional options that can be applied to solely + !these new atoms/elements that are read in. + new_starts(1) = atom_num + 1 + new_starts(2) = ele_num + 1 + if ((i==1).or.(trim(adjustl(dim)) == 'none')) then call read_in(i, displace, temp_box_bd) call grow_box(temp_box_bd) @@ -40,6 +53,8 @@ module mode_merge call read_in(i, displace, temp_box_bd) call grow_box(temp_box_bd) end if + + if(shift_flag) call shift(new_starts, i) end do return @@ -89,7 +104,16 @@ module mode_merge !Choose what to based on what the option string is select case(trim(textholder)) - + case('shift') + shift_flag = .true. + do i = 1,3 + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder, arglen) + if (arglen==0) stop "Missing vector component for shift command" + read(textholder, *) shift_vec(i) + end do + case('wrap') + wrap = .true. case default !If it isn't an available option to mode merge then we just exit exit @@ -97,4 +121,52 @@ module mode_merge end do end subroutine parse_command + + subroutine shift(array_start, filenum) + !This subroutine applies a shift to newly added atoms and elements. It also wraps the atoms + !if the user provides the wrap flag + + integer, dimension(2), intent(in) :: array_start + integer, intent(in) :: filenum + + integer :: i, j, ibasis, inod + real(kind=dp), dimension(3) :: current_shift + + + !Calculate the current shift which is the filenum-1 multiplied by the user specified shift + current_shift = (filenum-1)*shift_vec + + !First shift all the atoms + do i = array_start(1), atom_num + r_atom(:,i) = r_atom(:,i) + current_shift + end do + + !Now shift all the elements + do i = array_start(2), ele_num + do inod = 1, ng_node(lat_ele(i)) + do ibasis = 1, basisnum(lat_ele(i)) + r_node(:,ibasis, inod, i) = r_node(:,ibasis, inod, i) + current_shift + end do + end do + end do + + !Now we check if we have to wrap the atoms, nodes are not wrapped. For elements the periodic + !boundary conditions are applied in the actual CAC codes + if(wrap) then + do i = array_start(1), atom_num + call apply_periodic(r_atom(:,i)) + end do + + !If we don't include the wrap command then we have to increase the size of the box + else + do i = 1,3 + if (current_shift(i) < -lim_zero) then + box_bd(2*i-1) = box_bd(2*i-1) - current_shift(i) + else if (current_shift(i) > lim_zero) then + box_bd(2*i) = box_bd(2*i) + current_shift(i) + end if + end do + end if + + end subroutine shift end module mode_merge \ No newline at end of file diff --git a/src/subroutines.f90 b/src/subroutines.f90 index 330fbca..7085df1 100644 --- a/src/subroutines.f90 +++ b/src/subroutines.f90 @@ -249,4 +249,23 @@ module subroutines end do end subroutine siftdown + + + subroutine apply_periodic(r) + !This function checks if an atom is outside the box and wraps it back in. This is generally used when some + !kind of displacement is applied but the simulation cell is desired to be maintained as the same size. + + real(kind=dp), dimension(3), intent(inout) :: r + + integer :: j + + do j = 1, 3 + if (r(j) > box_bd(2*j)) then + r(j) = r(j) - box_bd(2*j) + else if (r(j) < box_bd(2*j-1)) then + r(j) = r(j) + box_bd(2*j-1) + end if + end do + end subroutine + end module subroutines \ No newline at end of file From f0665ce3ef75d9cfcf34ea191af680eb7dae5dd5 Mon Sep 17 00:00:00 2001 From: Alex Date: Thu, 16 Jan 2020 12:47:09 -0500 Subject: [PATCH 03/18] Add installation instructions --- INSTALL.md | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) create mode 100644 INSTALL.md diff --git a/INSTALL.md b/INSTALL.md new file mode 100644 index 0000000..e760b93 --- /dev/null +++ b/INSTALL.md @@ -0,0 +1,16 @@ +# Install CACmb + +In order to build CACmb switch to the src repository and run: + +``` +make +``` + +Once this command is complete you can run: + +``` +sudo make install +``` + +in order to copy the installation to `/usr/local/bin`. + From 12aa13b94bc3bc1497e4da5647400b614ae45681 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Fri, 17 Jan 2020 13:53:21 -0500 Subject: [PATCH 04/18] Fixes to element definitions --- src/io.f90 | 60 ++++++++++++++++++++++++++++++++++++++++-------------- 1 file changed, 45 insertions(+), 15 deletions(-) diff --git a/src/io.f90 b/src/io.f90 index 99f9e96..27d81f6 100644 --- a/src/io.f90 +++ b/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 - do j = 1, ele_num - if (lat_ele(j) == i) then - lat_size = size_ele(j)-1 - exit + !First figure out all of the unique element types + unique_num = 0 + unique_index(:) = 0 + eleloop:do i = 1, ele_num + 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) + end do + 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, & - new_type_to_type(10), new_lattice_types, sbox + 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_lattice_map(10) character(len=100) :: etype real(kind=dp) :: r(3), newdisplace(3) real(kind=dp), allocatable :: r_innode(:,:,:) @@ -662,9 +669,33 @@ module io do i = 1, basisnum(j) basis_type(i,j) = new_type_to_type(basis_type(i,j)) end do - end do + 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, type, sbox+n, r_innode) + call add_element(etype, size, new_lattice_map(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 From f80045058fdacca564908e97becb84a7014846aa Mon Sep 17 00:00:00 2001 From: Alex Date: Tue, 21 Jan 2020 12:01:45 -0500 Subject: [PATCH 05/18] quick fix to bug to write pycac format --- src/io.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/io.f90 b/src/io.f90 index 27d81f6..e4535a5 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -684,7 +684,7 @@ module io (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 + if(basis_type(ibasis,lattice_types+i) /= basis_type(ibasis,j)) then cycle old_loop end if end do From 347b73054b8dd044c6e3edb0b9c45e099352e64c Mon Sep 17 00:00:00 2001 From: Alex Date: Mon, 27 Jan 2020 10:01:22 -0500 Subject: [PATCH 06/18] Added group option and shift command --- README.md | 21 ++++- src/Makefile | 8 +- src/io.f90 | 2 + src/opt_group.f90 | 200 ++++++++++++++++++++++++++++++++++++++++++++ src/subroutines.f90 | 9 +- 5 files changed, 232 insertions(+), 8 deletions(-) create mode 100644 src/opt_group.f90 diff --git a/README.md b/README.md index cf7fcfb..47f7045 100644 --- a/README.md +++ b/README.md @@ -153,4 +153,23 @@ This option deletes vacancies on a plane which when minimized should result in a `bx by bz` - The burgers vector for the dislocation -`poisson` - Poisson ratio for continuum solution \ No newline at end of file +`poisson` - Poisson ratio for continuum solution + +### Option Group + +`-group select_type group_shape shape_arguments additional keywords` + +This option selects a group of either elements, nodes, atoms and applies some transformation to them. + +`select_type` - Either `nodes`, `atoms`, `elements`, `nodes/atoms`, `all`. When using the option `nodes` only nodes which are within the group are selected, `elements` selects elements based on whether the element center is within the group, `nodes/atoms` selects both nodes and atoms for the group. `all` selects elements based on the element center and atoms based on their position. + +`group_shape` - Specifies what shape the group takes and dictates which options must be passed. Each shape requires different arguments and these arguments are represented by the placeholder `shape_arguments`. The accepted group shapes and arguments are below: + +*Block:* + +`-group nodes block xlo xhi ylo yhi zlo zhi` + +This selects a group residing in a block with edges perpendicular to the simulation cell. The block boundaries are given by `xlo xhi ylo yhi zlo zhi`. + +`additional keywords`- Represents the various transformations which can be performed on a group. These additional keywords are given below. + diff --git a/src/Makefile b/src/Makefile index 0cdd008..16cac59 100644 --- a/src/Makefile +++ b/src/Makefile @@ -2,8 +2,8 @@ FC=ifort FFLAGS=-mcmodel=large -g -O0 -stand f08 -fpe0 -traceback -check bounds,uninit -warn all -implicitnone -no-wrap-margin #FFLAGS=-mcmodel=large -Ofast -no-wrap-margin MODES=mode_create.o mode_merge.o mode_convert.o -OPTIONS=opt_disl.o -OBJECTS=main.o elements.o io.o subroutines.o functions.o atoms.o call_mode.o box.o call_option.o $(MODES) $(OPTIONS) +OPTIONS=opt_disl.o opt_group.o +OBJECTS=main.o elements.o io.o subroutines.o functions.o atoms.o call_mode.o box.o $(MODES) $(OPTIONS) call_option.o .SUFFIXES: .SUFFIXES: .c .f .f90 .F90 .o @@ -38,6 +38,6 @@ atoms.o subroutines.o testfuncs.o box.o : functions.o main.o io.o $(MODES) $(OPTIONS) : elements.o call_mode.o : $(MODES) call_option.o : $(OPTIONS) -$(MODES) io.o: atoms.o box.o +$(MODES) $(OPTIONS) io.o : atoms.o box.o $(MODES) main.o : io.o -testfuncs.o elements.o mode_create.o opt_disl.o: subroutines.o +testfuncs.o elements.o mode_create.o $(OPTIONS): subroutines.o diff --git a/src/io.f90 b/src/io.f90 index e4535a5..9602a77 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -399,6 +399,8 @@ module io select case(max_ng_node) case(8) interp_max = (max_esize)**3 + case default + interp_max = 0 end select write(11,20) interp_max write(11,3) node_num diff --git a/src/opt_group.f90 b/src/opt_group.f90 new file mode 100644 index 0000000..13b5b03 --- /dev/null +++ b/src/opt_group.f90 @@ -0,0 +1,200 @@ +module opt_group + + !This module contains all code associated with dislocations + + use parameters + use elements + use subroutines + use box + implicit none + + integer :: group_ele_num, group_atom_num + 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 + + integer, allocatable :: element_index(:), atom_index(:) + + public + contains + + subroutine group(arg_pos) + !Main calling function for the group option + integer, intent(inout) :: arg_pos + + group_ele_num = 0 + group_atom_num = 0 + if(allocated(element_index)) deallocate(element_index) + if(allocated(atom_index)) deallocate(atom_index) + + call parse_group(arg_pos) + + call get_group + + !Now call the transformation functions for the group + if(displace) call displace_group + + end subroutine group + + subroutine parse_group(arg_pos) + !Parse the group command + integer, intent(inout) :: arg_pos + + integer :: i,arglen + character(len=100) :: textholder + + !Parse type and shape command + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, type, arglen) + if (arglen==0) STOP "Missing select_type in group command" + select case(trim(adjustl(type))) + case('atoms', 'elements', 'both') + continue + case default + print *, "Select_type ", trim(adjustl(type)), " is not an accept group selection criteria. ", & + "Please select from atoms, nodes, or both." + end select + + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, shape, arglen) + if (arglen==0) STOP "Missing group_shape in group command" + + !Now parse the arguments required by the user selected shape + select case(trim(adjustl(shape))) + case('block') + do i= 1, 6 + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder, arglen) + if (arglen==0) STOP "Missing block boundary in dislgen command" + call parse_pos(int((i+1)/2), textholder, block_bd(i)) + end do + + case default + print *, "Group shape ", trim(adjustl(shape)), " is not currently accepted. Please check documentation ", & + "for accepted group shapes." + end select + + !Now parse the additional options which may be present + do while(.true.) + if(arg_pos > command_argument_count()) exit + !Pull out the next argument which should either be a keyword or an option + arg_pos=arg_pos+1 + call get_command_argument(arg_pos, textholder) + textholder=adjustl(textholder) + + !Choose what to based on what the option string is + select case(trim(textholder)) + case('displace') + displace = .true. + do i = 1,3 + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder, arglen) + if (arglen==0) stop "Missing vector component for shift command" + read(textholder, *) disp_vec(i) + end do + case('wrap') + wrap = .true. + case default + !If it isn't an available option to opt_disl then we just exit + exit + end select + end do + end subroutine parse_group + + subroutine get_group + !This subroutine finds all elements and/or atoms within the group boundaries + !specified by the user. + integer :: i, inod, ibasis + integer, allocatable :: resize_array(:) + real(kind=dp) :: r_center(3) + + select case(trim(adjustl(shape))) + case('block') + !Allocate variables to arbitrary size + allocate(element_index(1024), atom_index(1024)) + !Check the type to see whether we need to find the elements within the group + select case(trim(adjustl(type))) + case('elements', 'both') + do i = 1, ele_num + r_center(:) = 0.0_dp + do inod = 1, ng_node(lat_ele(i)) + do ibasis = 1, basisnum(lat_ele(i)) + r_center = r_center + r_node(:,ibasis,inod,i)/(basisnum(lat_ele(i))*ng_node(lat_ele(i))) + end do + end do + + if (in_block_bd(r_center, block_bd)) then + group_ele_num = group_ele_num + 1 + if(group_ele_num > size(element_index)) then + allocate(resize_array(size(element_index) + 1024)) + resize_array(1:group_ele_num-1) = element_index + resize_array(group_ele_num:) = 0 + call move_alloc(resize_array, element_index) + end if + + element_index(group_ele_num) = i + end if + end do + end select + + !Check the type to see if we need to find the atoms within the group + select case(trim(adjustl(type))) + case('atoms','both') + do i = 1, atom_num + if(in_block_bd(r_atom(:,i),block_bd)) then + group_atom_num = group_atom_num + 1 + if (group_atom_num > size(atom_index)) then + allocate(resize_array(size(atom_index) + 1024)) + resize_array(1:group_atom_num -1) = atom_index + resize_array(group_atom_num:) = 0 + call move_alloc(resize_array, atom_index) + end if + + atom_index(group_atom_num) = i + end if + end do + end select + end select + end subroutine get_group + + subroutine displace_group + !This subroutine applies a displacement to elements/atoms in the groups + + integer :: i, inod, ibasis + + !Displace atoms + do i = 1, group_atom_num + r_atom(:,atom_index(i)) = r_atom(:,atom_index(i)) + disp_vec + end do + + !Displace elements + do i = 1, group_ele_num + do inod = 1, ng_node(lat_ele(element_index(i))) + do ibasis = 1, basisnum(lat_ele(element_index(i))) + r_node(:,ibasis,inod,element_index(i)) = r_node(:,ibasis,inod,element_index(i)) + disp_vec + end do + end do + end do + + !Now either apply periodic boundaries if wrap command was passed or adjust box dimensions + !Now we check if we have to wrap the atoms, nodes are not wrapped. For elements the periodic + !boundary conditions are applied in the actual CAC codes + if(wrap) then + do i = 1, atom_num + call apply_periodic(r_atom(:,i)) + end do + + !If we don't include the wrap command then we have to increase the size of the box + else + do i = 1,3 + if (disp_vec(i) < -lim_zero) then + box_bd(2*i-1) = box_bd(2*i-1) - disp_vec(i) + else if (disp_vec(i) > lim_zero) then + box_bd(2*i) = box_bd(2*i) + disp_vec(i) + end if + end do + end if + + end subroutine displace_group + +end module opt_group \ No newline at end of file diff --git a/src/subroutines.f90 b/src/subroutines.f90 index 7085df1..6efb97e 100644 --- a/src/subroutines.f90 +++ b/src/subroutines.f90 @@ -178,6 +178,8 @@ module subroutines real(kind=dp), intent(out) :: pos !The output parsed position value integer :: iospara + + iospara = 0 if(trim(adjustl(pos_string)) == 'inf') then pos=box_bd(2*i) else if(trim(adjustl(pos_string)) == '-inf') then @@ -258,12 +260,13 @@ module subroutines real(kind=dp), dimension(3), intent(inout) :: r integer :: j - + real(kind=dp) ::box_len do j = 1, 3 + box_len = box_bd(2*j) - box_bd(2*j-1) if (r(j) > box_bd(2*j)) then - r(j) = r(j) - box_bd(2*j) + r(j) = r(j) - box_len else if (r(j) < box_bd(2*j-1)) then - r(j) = r(j) + box_bd(2*j-1) + r(j) = r(j) + box_len end if end do end subroutine From 8a1bbcbc4e3c657eaba83644761dfb3209a4e4c5 Mon Sep 17 00:00:00 2001 From: Alex Date: Mon, 27 Jan 2020 10:06:20 -0500 Subject: [PATCH 07/18] Added call_option which has been needed for a while but was missing --- src/call_option.f90 | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) create mode 100644 src/call_option.f90 diff --git a/src/call_option.f90 b/src/call_option.f90 new file mode 100644 index 0000000..befb46c --- /dev/null +++ b/src/call_option.f90 @@ -0,0 +1,18 @@ +subroutine call_option(option, arg_pos) + use parameters + use opt_disl + use opt_group + implicit none + + integer, intent(inout) :: arg_pos + character(len=100), intent(in) :: option + + select case(trim(adjustl(option))) + case('-dislgen', '-disloop') + call dislocation(option, arg_pos) + case('-group') + call group(arg_pos) + case default + print *, 'Option ', trim(adjustl(option)), ' is not currently accepted.' + end select +end subroutine call_option \ No newline at end of file From e76187bd3c177ce8ecc536efc68210fa852f137d Mon Sep 17 00:00:00 2001 From: Alex Date: Mon, 27 Jan 2020 10:16:07 -0500 Subject: [PATCH 08/18] Square loop option --- src/opt_disl.f90 | 75 +++++++++++++++++++++++++++++------------------- 1 file changed, 45 insertions(+), 30 deletions(-) diff --git a/src/opt_disl.f90 b/src/opt_disl.f90 index 4510ede..00886e5 100644 --- a/src/opt_disl.f90 +++ b/src/opt_disl.f90 @@ -291,36 +291,51 @@ module opt_disl a3 = 3 end select - !Calculate loop perimeter - perimeter = 2.0_dp*pi*loop_radius - - !Define the number of points forming the loop - ! The following criteria are used as a trade-off between - ! good accuracy and computational efficiency: - ! - each dislocation segment should have a length of 5 angströms; - ! - the loop should contain at least 3 points - ! (for very small loops, this will result in segments shorter than 5 A); - ! - there should not be more than 100 points - ! (for very large loops, this will result in segments longer than 5 A). - Npoints = MAX( 3 , MIN( NINT(perimeter/5.d0) , 100 ) ) - - !angle between two consecutive points - theta = 2.0_dp*pi / dble(Npoints) - - !allocate xLoop - allocate(xLoop(Npoints,3)) - xLoop(:,:) = 0.0_dp - - !Calculate the position of each point in the loop - angle = 0.0_dp - do i = 1, size(xLoop,1) - xLoop(i,a1) = centroid(a1) + loop_radius*dcos(angle) - xLoop(i,a2) = centroid(a2) + loop_radius*dsin(angle) - xLoop(i,a3) = centroid(a3) - ! Increment angle for next point - angle = angle + theta - end do - + if(loop_radius < 0.0_dp) then + xLoop(:,:) = 0.d0 + xLoop(1,a1) = centroid(1) - loop_radius + xLoop(1,a2) = centroid(2) - loop_radius + xLoop(1,a3) = centroid(3) + xLoop(2,a1) = centroid(1) + loop_radius + xLoop(2,a2) = centroid(2) - loop_radius + xLoop(2,a3) = centroid(3) + xLoop(3,a1) = centroid(1) + loop_radius + xLoop(3,a2) = centroid(2) + loop_radius + xLoop(3,a3) = centroid(3) + xLoop(4,a1) = centroid(1) - loop_radius + xLoop(4,a2) = centroid(2) + loop_radius + xLoop(4,a3) = centroid(3) + else + !Calculate loop perimeter + perimeter = 2.0_dp*pi*loop_radius + + !Define the number of points forming the loop + ! The following criteria are used as a trade-off between + ! good accuracy and computational efficiency: + ! - each dislocation segment should have a length of 5 angströms; + ! - the loop should contain at least 3 points + ! (for very small loops, this will result in segments shorter than 5 A); + ! - there should not be more than 100 points + ! (for very large loops, this will result in segments longer than 5 A). + Npoints = MAX( 3 , MIN( NINT(perimeter/5.d0) , 100 ) ) + + !angle between two consecutive points + theta = 2.0_dp*pi / dble(Npoints) + + !allocate xLoop + allocate(xLoop(Npoints,3)) + xLoop(:,:) = 0.0_dp + + !Calculate the position of each point in the loop + angle = 0.0_dp + do i = 1, size(xLoop,1) + xLoop(i,a1) = centroid(a1) + loop_radius*dcos(angle) + xLoop(i,a2) = centroid(a2) + loop_radius*dsin(angle) + xLoop(i,a3) = centroid(a3) + ! Increment angle for next point + angle = angle + theta + end do + end if !Now actually calculate the displacement created by a loop for every atom do i = 1, atom_num u = 0.0_dp From 94e1c9fd7bc5d7fedbbec5b01436bdd3f44a3fbc Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 27 Jan 2020 11:56:31 -0500 Subject: [PATCH 09/18] Missing allocation in disloop code --- src/opt_disl.f90 | 1 + 1 file changed, 1 insertion(+) diff --git a/src/opt_disl.f90 b/src/opt_disl.f90 index 00886e5..79c45bf 100644 --- a/src/opt_disl.f90 +++ b/src/opt_disl.f90 @@ -292,6 +292,7 @@ module opt_disl end select if(loop_radius < 0.0_dp) then + ALLOCATE(xLoop(4,3)) xLoop(:,:) = 0.d0 xLoop(1,a1) = centroid(1) - loop_radius xLoop(1,a2) = centroid(2) - loop_radius From 1b099a17d5adc42b24b59d8868a86b91c9a00a2e Mon Sep 17 00:00:00 2001 From: Kevin Chu Date: Mon, 27 Jan 2020 12:00:47 -0500 Subject: [PATCH 10/18] Fix to makefile build order --- src/Makefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Makefile b/src/Makefile index 16cac59..a3f926f 100644 --- a/src/Makefile +++ b/src/Makefile @@ -38,6 +38,6 @@ atoms.o subroutines.o testfuncs.o box.o : functions.o main.o io.o $(MODES) $(OPTIONS) : elements.o call_mode.o : $(MODES) call_option.o : $(OPTIONS) -$(MODES) $(OPTIONS) io.o : atoms.o box.o +$(MODES) $(OPTIONS) subroutines.o io.o : atoms.o box.o $(MODES) main.o : io.o testfuncs.o elements.o mode_create.o $(OPTIONS): subroutines.o From 38872d8d9f1d9ba6c3019e6eb42788b0e870b69c Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 27 Jan 2020 20:28:47 -0500 Subject: [PATCH 11/18] Update to readme with correct format for dislocation command --- README.md | 47 ++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 46 insertions(+), 1 deletion(-) diff --git a/README.md b/README.md index 47f7045..1851f36 100644 --- a/README.md +++ b/README.md @@ -113,6 +113,26 @@ This mode merges multiple data files and creates one big simulation cell. The pa `dim` - the dimension they are to be stacked along, can be either `x`, `y`, or `z`. If the argument `none` is passed then the cells are just overlaid. Future options will include a delete overlap command. +**Additional options:** + +**Shift** + +``` +shift x y z +``` + +If the shift command is passed to mode merge then each file after the first file in the merge command is displaced by the vector `[x, y, z]`. This is additive so if you are merging three files and this command is passed then the second file is shifted by `[x,y,z]` and the third file is shifted by `2*[x,y,z]`. + +Example: `cacmb --merge z 2 Cu.mb Cu2.mb Cu3.mb Cumerged.mb shift 2 0 0` will shift the atomic and element positions in the `Cu2.mb` file by `[2,0,0]` and the positions in `Cu3.mb` by `[4,0,0]`. + +**Wrap** + +``` +wrap +``` + +This will wrap atomic positions back inside the box. Effectively as if periodic boundary conditions are applied so that atoms which exit from one side of the simulation cell enter back in through the other. + ## Options Options are additional portions of code which have additional functionality. Options are performed in the order that they appear in the argument list and can be added to any mode. If wanting to use strictly options use `--convert` to specify input and output files. @@ -120,7 +140,7 @@ Options are additional portions of code which have additional functionality. Opt ### Option dislgen ``` --dislgen [ijk] [hkl] x y z char_angle poisson +-dislgen [ijk] [hkl] x y z burgers char_angle poisson ``` @@ -131,6 +151,8 @@ This options adds an arbitrarily oriented dislocation into your model based on u `[hkl]` - The vector for the slip plane +`burgers` - The magnitude of the burgers vector for the dislocation to be inserted + `x y z` - The position of the dislocation centroid `char_angle` - Character angle of the dislocation (0 is screw and 90 is edge) @@ -173,3 +195,26 @@ This selects a group residing in a block with edges perpendicular to the simulat `additional keywords`- Represents the various transformations which can be performed on a group. These additional keywords are given below. +**Displace** + +``` +displace x y z +``` + +This is similar to the mode merge `-shift` argument. This simply shift atoms and elements within the group by a vector `[x,y,z]`. + +**Wrap** + +``` +wrap +``` + +This command wraps atoms back into the simulation cell as though periodic boundary conditions are being used. + +**Remesh** + +``` +remesh esize +``` + +This command remeshes the atoms/elements within the group to the new element size `esize`. Currently only accepts an `esize` of 2 which refines it to full atomistics. \ No newline at end of file From 0d2a351cd576750a827ade223abf83060afd7976 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 27 Jan 2020 20:30:40 -0500 Subject: [PATCH 12/18] Fix to disl_gen. Problems with sub_box_bd and element dislocation displacement field --- src/io.f90 | 6 +++++- src/opt_disl.f90 | 28 ++++++---------------------- 2 files changed, 11 insertions(+), 23 deletions(-) diff --git a/src/io.f90 b/src/io.f90 index 9602a77..fdde73d 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -642,7 +642,11 @@ module io read(11,*) ((sub_box_ori(j, k, sub_box_num+i), j = 1, 3), k = 1, 3) !Read in subbox boundaries read(11,*) sub_box_bd(:,sub_box_num+i) - sub_box_bd(:,sub_box_num+i) = sub_box_bd(:, sub_box_num+i) + displace(:) + + do j = 1, 3 + sub_box_bd(2*j-1,sub_box_num+i) = sub_box_bd(2*j-1, sub_box_num+i) + displace(j) + sub_box_bd(2*j,sub_box_num+i) = sub_box_bd(2*j, sub_box_num+i) + displace(j) + end do !Read in sub_box_array_bd read(11,*) ((sub_box_array_bd(j, k, sub_box_num+i), j = 1, 2), k = 1, 2) end do diff --git a/src/opt_disl.f90 b/src/opt_disl.f90 index 79c45bf..62ceae6 100644 --- a/src/opt_disl.f90 +++ b/src/opt_disl.f90 @@ -70,7 +70,11 @@ module opt_disl call parse_pos(i, textholder, centroid(i)) end do - print *, centroid + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder, arglen) + if (arglen==0) STOP "Missing parameter in dislgen command" + read(textholder, *) b + arg_pos = arg_pos + 1 call get_command_argument(arg_pos, textholder, arglen) @@ -84,26 +88,6 @@ module opt_disl read(textholder, *) poisson arg_pos = arg_pos + 1 - call get_command_argument(arg_pos, lattice, arglen) - if (arglen==0) STOP "Missing lattice in dislgen command" - - arg_pos = arg_pos + 1 - call get_command_argument(arg_pos, textholder, arglen) - if (arglen==0) STOP "Missing lattice parameter in dislgen command" - read(textholder, *) lattice_parameter - - arg_pos = arg_pos + 1 - !Now set the vurgers vector based on the lattice paarameter and lattice type - select case(lattice) - case('fcc') - b = lattice_parameter / sqrt(2.0_dp) - ! case('bcc') - ! b = lattice_parameter * sqrt(3.0_dp) / 2.0_dp - case default - print *, 'Error: Lattice structure', lattice, ' is not accepted for dislgen option' - STOP - end select - end subroutine parse_dislgen subroutine dislgen @@ -168,7 +152,7 @@ module opt_disl do i = 1, ele_num do inod=1, ng_node(lat_ele(i)) do ibasis = 1, basisnum(lat_ele(i)) - r = r_node(:,ibasis,inod,i) + r = r_node(:,ibasis,inod,i)-centroid r = matmul(inv_transform, r) if (r(1) == 0) then actan = pi/2 From 99a680542d9f95a5168f2943f29f9f4bea55e4d3 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 27 Jan 2020 20:32:18 -0500 Subject: [PATCH 13/18] Extended inf commands to be position independent of the inf command --- src/subroutines.f90 | 19 +++++++++++++++---- 1 file changed, 15 insertions(+), 4 deletions(-) diff --git a/src/subroutines.f90 b/src/subroutines.f90 index 6efb97e..e8c365b 100644 --- a/src/subroutines.f90 +++ b/src/subroutines.f90 @@ -186,15 +186,26 @@ module subroutines pos=box_bd(2*i-1) else if ((index(pos_string,'-') > 0).and.(index(pos_string,'inf')>0)) then !Now extract the number we are reducing from infinity - read(pos_string(index(pos_string,'-')+1:), *, iostat=iospara) pos + if(index(pos_string,'inf') < index(pos_string,'-')) then + read(pos_string(index(pos_string,'-')+1:), *, iostat=iospara) pos + else + read(pos_string(1:index(pos_string,'-')-1), *, iostat=iospara) pos + end if pos = box_bd(2*i) - pos else if ((index(pos_string,'+') > 0).and.(index(pos_string,'inf')>0)) then !Now extract the number we are reducing from infinity - read(pos_string(index(pos_string,'+')+1:), *, iostat=iospara) pos - pos = box_bd(2*i-1) + pos + if(index(pos_string,'inf') < index(pos_string,'-')) then + read(pos_string(index(pos_string,'+')+1:), *, iostat=iospara) pos + else + read(pos_string(1:index(pos_string,'+')-1), *, iostat=iospara) pos + end if else if ((index(pos_string,'*') > 0).and.(index(pos_string,'inf')>0)) then !Now extract the number we are reducing from infinity - read(pos_string(index(pos_string,'*')+1:), *, iostat=iospara) pos + if(index(pos_string,'inf') < index(pos_string,'-')) then + read(pos_string(index(pos_string,'*')+1:), *, iostat=iospara) pos + else + read(pos_string(1:index(pos_string,'*')-1), *, iostat=iospara) pos + end if pos = (box_bd(2*i)-box_bd(2*i-1))*pos else read(pos_string, *, iostat=iospara) pos From 69ab46b6e36ac374f4cc82d449193883d70e699e Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 27 Jan 2020 21:06:22 -0500 Subject: [PATCH 14/18] Fix to inf*val format for positions --- src/subroutines.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/subroutines.f90 b/src/subroutines.f90 index e8c365b..ae15a65 100644 --- a/src/subroutines.f90 +++ b/src/subroutines.f90 @@ -201,7 +201,7 @@ module subroutines end if else if ((index(pos_string,'*') > 0).and.(index(pos_string,'inf')>0)) then !Now extract the number we are reducing from infinity - if(index(pos_string,'inf') < index(pos_string,'-')) then + if(index(pos_string,'inf') < index(pos_string,'*')) then read(pos_string(index(pos_string,'*')+1:), *, iostat=iospara) pos else read(pos_string(1:index(pos_string,'*')-1), *, iostat=iospara) pos From 1db2b24499ca2447392b0657b148870cbc60b818 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 27 Jan 2020 21:10:31 -0500 Subject: [PATCH 15/18] Fixes to parse_pos subroutine --- src/subroutines.f90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/subroutines.f90 b/src/subroutines.f90 index ae15a65..9fd159c 100644 --- a/src/subroutines.f90 +++ b/src/subroutines.f90 @@ -194,11 +194,12 @@ module subroutines pos = box_bd(2*i) - pos else if ((index(pos_string,'+') > 0).and.(index(pos_string,'inf')>0)) then !Now extract the number we are reducing from infinity - if(index(pos_string,'inf') < index(pos_string,'-')) then + if(index(pos_string,'inf') < index(pos_string,'+')) then read(pos_string(index(pos_string,'+')+1:), *, iostat=iospara) pos else read(pos_string(1:index(pos_string,'+')-1), *, iostat=iospara) pos end if + pos = box_bd(2*i-1) + pos else if ((index(pos_string,'*') > 0).and.(index(pos_string,'inf')>0)) then !Now extract the number we are reducing from infinity if(index(pos_string,'inf') < index(pos_string,'*')) then @@ -206,7 +207,7 @@ module subroutines else read(pos_string(1:index(pos_string,'*')-1), *, iostat=iospara) pos end if - pos = (box_bd(2*i)-box_bd(2*i-1))*pos + pos = (box_bd(2*i)-box_bd(2*i-1))*pos + box_bd(2*i-1) else read(pos_string, *, iostat=iospara) pos end if From fd143783ec889798da0377052c0ddf8001f1794a Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Tue, 28 Jan 2020 09:36:00 -0500 Subject: [PATCH 16/18] 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 From 9ebdfff0a137179cb30384e1f20016f4f642ceed Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Tue, 28 Jan 2020 10:42:30 -0500 Subject: [PATCH 17/18] Added print messages to let user know whats going on --- src/io.f90 | 4 ++++ src/main.f90 | 40 +++++++++++++++++++++++++++------------- src/mode_create.f90 | 7 +++++++ src/mode_merge.f90 | 5 +++++ src/opt_disl.f90 | 7 ++++++- src/opt_group.f90 | 14 +++++++++++++- 6 files changed, 62 insertions(+), 15 deletions(-) diff --git a/src/io.f90 b/src/io.f90 index fdde73d..8fd139b 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -627,6 +627,7 @@ module io temp_box_bd(2*i) = temp_box_bd(2*i) + newdisplace(i) end do + call grow_box(temp_box_bd) !Read in the number of sub_boxes and allocate the variables read(11, *) n @@ -707,6 +708,9 @@ module io call grow_ele_arrays(in_eles, in_atoms) allocate(r_innode(3,max_basisnum, max_ng_node)) + print *, "Read in ", in_eles, " elements and ", in_atoms, " atoms from ", trim(adjustl(file)) + print *, "New box dimensions are: ", box_bd + !Read the atoms do i = 1, in_atoms read(11,*) j, type, r(:) diff --git a/src/main.f90 b/src/main.f90 index ed34c53..7d24d6d 100644 --- a/src/main.f90 +++ b/src/main.f90 @@ -1,25 +1,39 @@ program main - !**************************** CACmb ******************************* - !* CAC model building toolkit * - ! ____________ * - ! / / * - ! / / * - ! /___________/ * - ! _|_ _|_ _|____________ * - ! / / * - ! / / * - ! /___________/ * - ! * - !******************************************************************* + !**************************** CACmb ******************** + !* CAC model building toolkit * + !* ____________ * + !* / / * + !* / / * + !* /___________/ * + !* _|_ _|_ _|____________ * + !* / / * + !* / / * + !* /___________/ * + !* * + !******************************************************** use parameters use elements use io - + integer :: i, end_mode_arg, arg_num, arg_pos character(len=100) :: argument + !Print introduction text + print *, '*********************** CACmb *********************' + print *, '* CAC model building toolkit *' + print *, '* _______ *' + print *, '* / / *' + print *, '* / / *' + print *, '* /______ / *' + print *, '* _|_ _|_ _|_______ *' + print *, '* / / *' + print *, '* / / *' + print *, '* /______ / *' + print *, '* *' + print *, '****************************************************' + !Call initialization functions call lattice_init call box_init diff --git a/src/mode_create.f90 b/src/mode_create.f90 index 9463e92..c3e6bcc 100644 --- a/src/mode_create.f90 +++ b/src/mode_create.f90 @@ -30,6 +30,8 @@ module mode_create integer :: i, ibasis, inod real(kind=dp), allocatable :: r_node_temp(:,:,:) + print *, '-----------------------Mode Create---------------------------' + !Initialize default parameters orient = reshape((/ 1.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 1.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 1.0_dp /), shape(orient)) cell_mat(:,:)=0.0_dp @@ -91,6 +93,8 @@ module mode_create end do else + print *, "Creating 1 element" + call cell_init(lattice_parameter, esize, element_type, orient, cell_mat) !If the user doesn't pass any build instructions than we just put the cell mat into the element_array call alloc_ele_arrays(1,0) @@ -110,6 +114,7 @@ module mode_create !If we passed the dup_flag or dim_flag then we have to convert the lattice points and add them to the atom/element arrays if(dup_flag.or.dim_flag) then + !Call the build function with the correct transformation matrix select case(trim(adjustl(element_type))) case('fcc') @@ -123,6 +128,8 @@ module mode_create !Now that it is built multiply by the lattice parameter box_bd = box_bd*lattice_parameter + print *, "Using mode create, ", lat_ele_num, " elements are created and ", lat_atom_num*basisnum(1), " atoms are created." + !Allocate variables call alloc_ele_arrays(lat_ele_num, lat_atom_num*basisnum(1)) if(lat_atom_num > 0) then diff --git a/src/mode_merge.f90 b/src/mode_merge.f90 index 567b19b..14a015d 100644 --- a/src/mode_merge.f90 +++ b/src/mode_merge.f90 @@ -20,6 +20,8 @@ module mode_merge integer :: i real(kind=dp) :: displace(3), temp_box_bd(6) + print *, '-----------------------Mode Merge---------------------------' + wrap = .false. shift_flag = .false. @@ -34,6 +36,7 @@ module mode_merge !The new starts variable dictate where in the atom and element array each new !file starts. This is used for additional options that can be applied to solely !these new atoms/elements that are read in. + new_starts(1) = atom_num + 1 new_starts(2) = ele_num + 1 @@ -136,6 +139,8 @@ module mode_merge !Calculate the current shift which is the filenum-1 multiplied by the user specified shift current_shift = (filenum-1)*shift_vec + print *, "Atoms/elements from file ", trim(adjustl(infiles(filenum))), " are shifted by ", current_shift + !First shift all the atoms do i = array_start(1), atom_num r_atom(:,i) = r_atom(:,i) + current_shift diff --git a/src/opt_disl.f90 b/src/opt_disl.f90 index 62ceae6..025d680 100644 --- a/src/opt_disl.f90 +++ b/src/opt_disl.f90 @@ -23,12 +23,12 @@ module opt_disl contains subroutine dislocation(option, arg_pos) - !Main calling function for all codes related to dislocations character(len=100), intent(in) :: option integer, intent(inout) :: arg_pos + print *, '--------------------Option Dislocation-----------------------' select case(trim(adjustl(option))) case('-dislgen') @@ -97,6 +97,8 @@ module opt_disl real(kind=dp) :: ss_ori(3,3), ss_inv(3,3), be, bs, slipx(3), disp_transform(3,3), inv_transform(3,3), & actan, r(3), disp(3) + print *, "Dislocation with centroid ", centroid, " is inserted" + !Calculate screw and edge burgers vectors be = sin(char_angle*pi/180.0_dp)*b bs = cos(char_angle*pi/180.0_dp)*b @@ -254,6 +256,9 @@ module opt_disl real(kind = dp) :: perimeter, angle, theta, omega, xA(3), xB(3), xC(3), u(3) real(kind=dp), dimension(:,:), allocatable :: xloop !coordinate of points forming loop + + print *, "Dislocation loop with centroid ", centroid, " is inserted" + if(allocated(xLoop)) deallocate(xLoop) !Define new directions diff --git a/src/opt_group.f90 b/src/opt_group.f90 index 25d8671..6d129f1 100644 --- a/src/opt_group.f90 +++ b/src/opt_group.f90 @@ -22,6 +22,8 @@ module opt_group !Main calling function for the group option integer, intent(inout) :: arg_pos + print *, '-----------------------Option Group-------------------------' + group_ele_num = 0 group_atom_num = 0 remesh_size=0 @@ -117,6 +119,9 @@ module opt_group select case(trim(adjustl(shape))) case('block') + + print *, "Group has block shape with boundaries: ", block_bd + !Allocate variables to arbitrary size allocate(element_index(1024), atom_index(1024)) !Check the type to see whether we need to find the elements within the group @@ -162,6 +167,8 @@ module opt_group end do end select end select + + print *, 'Group contains ', group_ele_num, " elements and ", group_atom_num, " atoms." end subroutine get_group subroutine displace_group @@ -169,6 +176,8 @@ module opt_group integer :: i, inod, ibasis + print *, "Elements/atoms in group displaced by ", disp_vec + !Displace atoms do i = 1, group_atom_num r_atom(:,atom_index(i)) = r_atom(:,atom_index(i)) + disp_vec @@ -207,7 +216,7 @@ module opt_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 + integer :: i, j, ie, type_interp(max_basisnum*max_esize**3), add_atom_num, orig_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 @@ -216,6 +225,7 @@ module opt_group !Refining to atoms case(2) if(group_ele_num > 0) then + orig_atom_num = atom_num !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) @@ -236,6 +246,8 @@ module opt_group !Once all atoms are added we delete all of the elements call delete_elements(group_ele_num, element_index) + print *, group_ele_num, " elements of group are refined to ", atom_num -orig_atom_num, " atoms." + end if !Remeshing to elements, currently not available case default From a942b8e2a1ccad6cc883cc4355d4afd6c8515357 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Tue, 28 Jan 2020 10:59:49 -0500 Subject: [PATCH 18/18] Added -ow option which allows for automatic overwriting --- README.md | 10 +++++++++- src/call_option.f90 | 6 +++++- src/io.f90 | 4 ++-- src/main.f90 | 12 ++++++++++++ 4 files changed, 28 insertions(+), 4 deletions(-) diff --git a/README.md b/README.md index 1851f36..5d45ef3 100644 --- a/README.md +++ b/README.md @@ -217,4 +217,12 @@ This command wraps atoms back into the simulation cell as though periodic bounda remesh esize ``` -This command remeshes the atoms/elements within the group to the new element size `esize`. Currently only accepts an `esize` of 2 which refines it to full atomistics. \ No newline at end of file +This command remeshes the atoms/elements within the group to the new element size `esize`. Currently only accepts an `esize` of 2 which refines it to full atomistics. + +### Option overwrite + +``` +-ow +``` + +If this option is passed then all files are automatically overwritten without asking the user. \ No newline at end of file diff --git a/src/call_option.f90 b/src/call_option.f90 index befb46c..60551af 100644 --- a/src/call_option.f90 +++ b/src/call_option.f90 @@ -12,7 +12,11 @@ subroutine call_option(option, arg_pos) call dislocation(option, arg_pos) case('-group') call group(arg_pos) + case('-ow') + arg_pos = arg_pos + 1 + continue case default - print *, 'Option ', trim(adjustl(option)), ' is not currently accepted.' + print *, 'Option ', trim(adjustl(option)), ' is not currently accepted. Skipping to next argument' + arg_pos = arg_pos + 1 end select end subroutine call_option \ No newline at end of file diff --git a/src/io.f90 b/src/io.f90 index 8fd139b..26d4250 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -9,7 +9,7 @@ module io integer :: outfilenum = 0, infilenum = 0 character(len=100) :: outfiles(10), infiles(10) - + logical :: force_overwrite public contains @@ -38,7 +38,7 @@ module io !Check to see if file exists, if it does then ask user if they would like to overwrite the file inquire(file=trim(temp_outfile), exist=file_exists) - if (file_exists) then + if (file_exists.and.(.not.(force_overwrite))) then if (overwrite == 'r') print *, "File ", trim(temp_outfile), " already exists. Would you like to overwrite? (Y/N)" read(*,*) overwrite if((scan(overwrite, "n") > 0).or.(scan(overwrite, "N") > 0)) then diff --git a/src/main.f90 b/src/main.f90 index 7d24d6d..31f5f58 100644 --- a/src/main.f90 +++ b/src/main.f90 @@ -37,11 +37,23 @@ program main !Call initialization functions call lattice_init call box_init + force_overwrite=.false. + end_mode_arg = 0 ! Command line parsing arg_num = command_argument_count() + !Check to see if overwrite flag is passed + do i = 1, arg_num + call get_command_argument(i,argument) + select case(trim(adjustl(argument))) + + case('-ow') + force_overwrite = .true. + print *, "Overwrite flag passed, output files will be overwritten" + end select + end do !Determine if a mode is being used and what it is. The first argument has to be the mode !if a mode is being used call get_command_argument(1, argument)