diff --git a/README.md b/README.md index 346fa58..fe71ed0 100644 --- a/README.md +++ b/README.md @@ -277,6 +277,14 @@ Delete requires the usage of additional keywords to specify which delete action This command will delete all overlapping atoms within a specific cutoff radius `rc_off`. This currently does not affect elements. +### Option sbox_ori + +``` +-sbox_ori sbox [hkl] [hkl] [hkl] +``` + +This option is primarily used when reading data from non .mb formats. This code simply sets the orientation variable for the specified sub box `sbox`. + **** ## Position Specification diff --git a/src/box.f90 b/src/box.f90 index f98ed88..c6e7948 100644 --- a/src/box.f90 +++ b/src/box.f90 @@ -11,7 +11,6 @@ module box !command. Currently only mode_merge will require sub_boxes, for mode_create it will always !allocate to only 1 sub_box integer :: sub_box_num = 0 - integer, allocatable :: sub_box_array_bd(:,:,:)!Boundaries in the atom and element arrays for each sub_box real(kind=dp), allocatable :: sub_box_ori(:,:,:)!Orientations for each of the subboxes real(kind=dp), allocatable :: sub_box_bd(:,:)!Boundaries for each of the sub_boxes @@ -38,11 +37,10 @@ module box integer :: i - allocate(sub_box_ori(3,3,n), sub_box_bd(6,n), sub_box_array_bd(2,2,n)) + allocate(sub_box_ori(3,3,n), sub_box_bd(6,n)) do i = 1, n sub_box_ori(:,:,i) = identity_mat(3) sub_box_bd(:,i) = 0.0_dp - sub_box_array_bd(:,:,i) = 1 end do end subroutine alloc_sub_box @@ -65,10 +63,6 @@ module box temp_bd(:, sub_box_num+1:) = 0.0_dp call move_alloc(temp_bd, sub_box_bd) - temp_array_bd(:,:,1:sub_box_num) = sub_box_array_bd - temp_array_bd(:,:,sub_box_num+1:) = 1 - call move_alloc(temp_array_bd, sub_box_array_bd) - return end subroutine grow_sub_box diff --git a/src/call_option.f90 b/src/call_option.f90 index 73cd1de..ead7af0 100644 --- a/src/call_option.f90 +++ b/src/call_option.f90 @@ -29,6 +29,8 @@ subroutine call_option(option, arg_pos) call get_command_argument(arg_pos, box_bc) arg_pos=arg_pos+1 bound_called = .true. + case('-sbox_ori') + call sbox_ori(arg_pos) case('-delete') call run_delete(arg_pos) case default diff --git a/src/io.f90 b/src/io.f90 index 4aab723..90efe0f 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -389,7 +389,7 @@ module io 10 format('box matrix') 11 format(3f23.15) 12 format('coarse-grained domain') -13 format('ie ele_type grain_ele lat_type_ele'/ 'ip ibasis type x y z') +13 format('ie ele_type grain_ele lat_type_ele'/ 'ip ibasis x y z') 14 format('atomistic domain' / 'ia grain_atom type_atom x y z') 15 format('maximum lattice periodicity length' / 3f23.15) 16 format('Number of lattice types and atom types '/ 2i16) @@ -521,11 +521,10 @@ module io write(11,*) box_bd(:) !Write the number of sub_boxes in the system write(11,*) sub_box_num - !For every subbox write the orientation, sub box boundary, and sub_box_array_bds + !For every subbox write the orientation, sub box boundar do i = 1, sub_box_num write(11,*) sub_box_ori(:,:,i) write(11,*) sub_box_bd(:,i) - write(11,*) ((sub_box_array_bd(j,k,i), j = 1, 2), k = 1, 2) end do !Write the number of atom types in the current model and all of their names @@ -590,13 +589,13 @@ module io end if select case(temp_infile(scan(temp_infile,'.',.true.)+1:)) - case('xyz', 'lmp', 'vtk', 'mb') + case('restart', 'mb') infilenum=infilenum+1 infiles(infilenum) = temp_infile exit case default print *, "File type: ", trim(temp_infile(scan(temp_infile,'.',.true.):)), "not currently accepted. ", & - "please input a filename with extension from following list: mb." + "please input a filename with extension from following list: mb, restart." read(*,*) temp_infile end select @@ -615,6 +614,8 @@ module io select case(trim(adjustl(infiles(i)(scan(infiles(i),'.',.true.)+1:)))) case('mb') call read_mb(infiles(i), displace, temp_box_bd) + case('restart') + call read_pycac(infiles(i), displace, temp_box_bd) case default print *, "The extension ", trim(adjustl(outfiles(i)(scan(outfiles(i),'.',.true.)+1:))), & " is not accepted for writing. Please select from: mb and try again" @@ -643,7 +644,6 @@ module io !Read in the box boundary and grow the current active box bd read(11, *) temp_box_bd(:) - print *, "displace", displace do i = 1, 3 if (abs(displace(i)) > lim_zero) then newdisplace(i) = displace(i) - temp_box_bd(2*i-1) @@ -654,7 +654,6 @@ module io temp_box_bd(2*i) = temp_box_bd(2*i) + newdisplace(i) end do - print *, "newdisplace", newdisplace call grow_box(temp_box_bd) !Read in the number of sub_boxes and allocate the variables read(11, *) n @@ -676,14 +675,8 @@ module io 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 - !Add the existing element boundaries - sub_box_array_bd(:,1,sub_box_num+1:) = sub_box_array_bd(:,1,sub_box_num+1:) + atom_num - sub_box_array_bd(:,2,sub_box_num+1:) = sub_box_array_bd(:,2,sub_box_num+1:) + ele_num - !Read in the number of atom types and all their names read(11, *) new_atom_types, (new_type_to_name(i), i = 1, new_atom_types) !Now fit these into the global list of atom types, after this new_type_to_type is the actual global @@ -715,8 +708,8 @@ module io 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 + if ((basisnum(lattice_types + i) == basisnum(j)).and. & + (ng_node(lattice_types + i) == ng_node(j))) then !Now check the basis level variables do ibasis =1, basisnum(i) if(basis_type(ibasis,lattice_types+i) /= basis_type(ibasis,j)) then @@ -770,4 +763,181 @@ module io call set_max_esize end subroutine read_mb + + subroutine read_pycac(file, displace, temp_box_bd) + !Add subroutine for reading in restart files from PyCAC. This code currently only + !works for 8 node elements. + character(len=100), intent(in) :: file + real(kind=dp), dimension(3), intent(in) :: displace + real(kind = dp), dimension(6), intent(out) :: temp_box_bd + + integer :: i, inod, ibasis, j, k, in_eles, in_atoms, ele_types, in_lat_num, in_atom_types, & + atom_type_map(10), etype_map(10), etype, lat_type, new_lattice_map(10), & + atom_type + real(kind=dp) :: newdisplace(3), r_in(3,1,8), r_in_atom(3), new_displace(3) + character(len=100) :: textholder, in_lattype_map(10) + character(len=2) :: atomic_element + !First open the file + open(unit=11, file=trim(adjustl(file)), action='read',position='rewind') + + !Disregard unneeded information + do i = 1, 3 + read(11,*) textholder + end do + + !Read element number + read(11,*) in_eles + + !Discard info and read ng_max_node + do i = 1, 5 + read(11,*) textholder + end do + read(11,*) max_ng_node + + !Read element types (only needed inside this subroutine) + read(11,*) textholder + read(11,*) ele_types + + !Read in atom num + read(11,*) textholder + read(11,*) in_atoms + + !read in lattice_types and atom types + do i = 1,3 + read(11,*) textholder + end do + read(11,*) in_lat_num, in_atom_types + + + !Read define atom_types by name + read(11,*) textholder + do i = 1, in_atom_types + read(11,*) j, atomic_element + call add_atom_type(atomic_element, atom_type_map(i)) + end do + + !Read in the boundary + read(11,*) textholder + read(11,*) box_bc + + !Disregard useless info + do i = 1, 3 + read(11,*) textholder + end do + + !Read in lat_type map + do i = 1, in_lat_num + read(11,*) j, in_lattype_map(i) + ng_node(lattice_types+i) = 8 !Only cubic elements are currently supported in pyCAC + end do + + !Disregard useless info + do i =1 , 3 + read(11,*) textholder + end do + + !Read box boundaries and displace them if necessary + read(11,*) temp_box_bd(:) + do i = 1, 3 + if (abs(displace(i)) > lim_zero) then + newdisplace(i) = displace(i) - temp_box_bd(2*i-1) + else + newdisplace(i)=displace(i) + end if + temp_box_bd(2*i-1) = temp_box_bd(2*i-1) + newdisplace(i) + temp_box_bd(2*i) = temp_box_bd(2*i) + newdisplace(i) + end do + + call grow_box(temp_box_bd) + + !Allocate sub_box + if (sub_box_num == 0) then + call alloc_sub_box(1) + else + call grow_sub_box(1) + end if + + !Because orientations and other needed sub_box information isn't really + !present within the restart file we just default a lot of it. + sub_box_ori(:,:,sub_box_num+1) = identity_mat(3) + sub_box_bd(:, sub_box_num+1) = temp_box_bd + + !Read in more useless info + do i = 1, 10 + read(11,*) textholder + end do + + !Now read the ele_type to size and lat map + do i = 1, ele_types + read(11,*) j, etype_map(i) + etype_map(i) = etype_map(i) + 1 + end do + + + !Now set up the lattice types. In this code it assumes that lattice_type 1 maps to + !atom type 1 because it only allows 1 atom per basis in pyCAC at the moment. + do i = 1, in_lat_num + basisnum(lattice_types+i) = 1 + basis_type(1,lattice_types+i) = atom_type_map(i) + end do + + !Figure out the lattice type maps in case we have repeated lattice_types + k = lattice_types + 1 + new_lattice_map(:) = 0 + new_loop:do i = 1, in_lat_num + old_loop:do j = 1, lattice_types + !First check all the lattice level variables + if ((basisnum(lattice_types+i) == basisnum(j)).and. & + (ng_node(lattice_types+i) == ng_node(j))) then + !Now check the basis level variables + do ibasis =1, basisnum(i) + if(basis_type(ibasis,lattice_types+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 more useless data + read(11,*) textholder + read(11,*) textholder + + !set max values and allocate variables + max_basisnum = maxval(basisnum) + max_ng_node = maxval(ng_node) + call grow_ele_arrays(in_eles, in_atoms) + + !Now start reading the elements + do i = 1, in_eles + read(11,*) j, etype, k, lat_type + do inod = 1, 8 + read(11, *) j, k, r_in(:,1,inod) + r_in(:,1,inod) = r_in(:,1,inod) + newdisplace + end do + call add_element(in_lattype_map(lat_type), etype_map(etype), new_lattice_map(lat_type), sub_box_num + 1, r_in) + end do + + !Read useless data + read(11,*) textholder + read(11,*) textholder + + do i = 1, in_atoms + read(11,*) j, k, atom_type, r_in_atom(:) + r_in_atom = r_in_atom + newdisplace + call add_atom(atom_type_map(atom_type), sub_box_num + 1, r_in_atom) + end do + !Close file + close(11) + + lattice_types = maxval(new_lattice_map) + + sub_box_num = sub_box_num + 1 + + call set_max_esize + end subroutine read_pycac end module io diff --git a/src/mode_create.f90 b/src/mode_create.f90 index 3cbafa2..ffdfc05 100644 --- a/src/mode_create.f90 +++ b/src/mode_create.f90 @@ -158,9 +158,6 @@ module mode_create sub_box_num = 1 sub_box_ori(:,:,1) = orient sub_box_bd(:,1) = box_bd - sub_box_array_bd(1,:,1) = 1 - sub_box_array_bd(2,1,1) = atom_num - sub_box_array_bd(2,2,1) = ele_num end subroutine create !This subroutine parses the command and pulls out information needed for mode_create subroutine parse_command(arg_pos) diff --git a/src/opt_group.f90 b/src/opt_group.f90 index 03b3063..be87611 100644 --- a/src/opt_group.f90 +++ b/src/opt_group.f90 @@ -364,7 +364,7 @@ module opt_group if (max_remesh) then - max_loops = (remesh_size-2)/2 + max_loops = (remesh_size-3)/2 else max_loops = 1 end if diff --git a/src/opt_orient.f90 b/src/opt_orient.f90 index 812754f..0f0dea9 100644 --- a/src/opt_orient.f90 +++ b/src/opt_orient.f90 @@ -133,4 +133,29 @@ module opt_orient box_bd = orig_box_bd end subroutine unorient + subroutine sbox_ori(arg_pos) + integer, intent(inout) :: arg_pos + + integer :: i, sbox_in, arg_len + real(kind = dp) :: new_orient(3,3) + character(len=100) :: textholder + character(len=8) :: ori_string + + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder,arg_len) + if (arg_len== 0) stop 'Missing sbox in sbox_ori command' + read(textholder,*) sbox_in + + do i = 1, 3 + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, ori_string, arg_len) + if (arg_len == 0) print *, "Missing orientation vector in -orient option" + call parse_ori_vec(ori_string, new_orient(i,:)) + new_orient(i,:) = new_orient(i,:) / norm2(new_orient(i,:)) + end do + + sub_box_ori(:,:,sbox_in) = new_orient + arg_pos = arg_pos + 1 + + end subroutine sbox_ori end module opt_orient \ No newline at end of file