From 724e5732877d2dd42d1d3189184fad55fe211935 Mon Sep 17 00:00:00 2001 From: Alex Selimov Date: Mon, 10 Aug 2020 15:52:19 -0400 Subject: [PATCH] Added deform box option and allowed for reading time information from restart files --- src/Makefile | 2 +- src/box.f90 | 3 +- src/call_option.f90 | 3 ++ src/elements.f90 | 2 + src/io.f90 | 20 +++++---- src/opt_deform.f90 | 98 +++++++++++++++++++++++++++++++++++++++++++++ 6 files changed, 117 insertions(+), 11 deletions(-) create mode 100644 src/opt_deform.f90 diff --git a/src/Makefile b/src/Makefile index 87f9306..fde4484 100644 --- a/src/Makefile +++ b/src/Makefile @@ -2,7 +2,7 @@ FC=ifort FFLAGS=-mcmodel=large -g -O0 -stand f08 -fpe0 -traceback -check bounds,uninit -warn all -implicitnone -no-wrap-margin -heap-arrays #FFLAGS=-mcmodel=large -Ofast -no-wrap-margin -heap-arrays MODES=mode_create.o mode_merge.o mode_convert.o -OPTIONS=opt_disl.o opt_group.o opt_orient.o opt_delete.o +OPTIONS=opt_disl.o opt_group.o opt_orient.o opt_delete.o opt_deform.o OBJECTS=main.o elements.o io.o subroutines.o functions.o atoms.o call_mode.o box.o $(MODES) $(OPTIONS) call_option.o sorts.o .SUFFIXES: diff --git a/src/box.f90 b/src/box.f90 index c6e7948..faf268c 100644 --- a/src/box.f90 +++ b/src/box.f90 @@ -14,7 +14,6 @@ module 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 - !Below are some simulation parameters which may be adjusted if reading in restart files integer :: timestep=0 real(kind=dp) :: total_time=0.0_dp @@ -98,4 +97,4 @@ module box end do return end subroutine in_sub_box -end module box \ No newline at end of file +end module box diff --git a/src/call_option.f90 b/src/call_option.f90 index 368e542..1c3759e 100644 --- a/src/call_option.f90 +++ b/src/call_option.f90 @@ -3,6 +3,7 @@ subroutine call_option(option, arg_pos) use opt_disl use opt_group use opt_orient + use opt_deform use opt_delete use box implicit none @@ -31,6 +32,8 @@ subroutine call_option(option, arg_pos) bound_called = .true. case('-sbox_ori') call sbox_ori(arg_pos) + case('-deform') + call deform(arg_pos) case('-delete') call run_delete(arg_pos) case('-set_cac') diff --git a/src/elements.f90 b/src/elements.f90 index 7f8cef4..6df859b 100644 --- a/src/elements.f90 +++ b/src/elements.f90 @@ -503,6 +503,7 @@ module elements size_ele(sorted_index(i)) = 0 lat_ele(sorted_index(i)) = 0 sbox_ele(sorted_index(i)) = 0 + tag_ele(sorted_index(i)) = 0 else node_num = node_num - ng_node(lat_ele(sorted_index(i))) r_node(:,:,:,sorted_index(i)) = r_node(:,:,:,ele_num) @@ -510,6 +511,7 @@ module elements size_ele(sorted_index(i)) = size_ele(ele_num) lat_ele(sorted_index(i)) = lat_ele(ele_num) sbox_ele(sorted_index(i)) = sbox_ele(ele_num) + tag_ele(sorted_index(i)) = tag_ele(ele_num) end if ele_num = ele_num - 1 end do diff --git a/src/io.f90 b/src/io.f90 index ba3a389..aba4dd9 100644 --- a/src/io.f90 +++ b/src/io.f90 @@ -491,7 +491,7 @@ module io exit endif end do - write(11, '(4i16)') i, etype, 1, basis_type(1,lat_ele(i)) + write(11, '(4i16)') tag_ele(i), etype, 1, basis_type(1,lat_ele(i)) do inod = 1, ng_node(lat_ele(i)) do ibasis = 1, basisnum(lat_ele(i)) ip = ip + 1 @@ -505,7 +505,7 @@ module io if(atom_num /= 0) then write(11,14) do i = 1, atom_num - write(11, '(3i16, 3f23.15)') i, 1, type_atom(i), r_atom(:,i) + write(11, '(3i16, 3f23.15)') tag_atom(i), 1, type_atom(i), r_atom(:,i) end do end if @@ -781,19 +781,19 @@ module io integer :: i, inod, ibasis, j, k, l, in_eles, in_atoms, ele_types, in_lat_num, in_atom_types, & atom_type_map(100), etype_map(100), etype, lat_type, new_lattice_map(100), & - atom_type + atom_type, stat real(kind=dp) :: newdisplace(3), r_in(3,1,8), r_in_atom(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 the timestep information + read(11,*) textholder + read(11,*) timestep, total_time !Read element number + read(11,*) textholder read(11,*) in_eles !Discard info and read ng_max_node @@ -941,7 +941,11 @@ module io end if do i = 1, in_atoms - read(11,*) j, k, atom_type, r_in_atom(:) + read(11,*, iostat=stat) j, k, atom_type, r_in_atom(:) + if(stat > 0) then + print *, j + stop + end if r_in_atom = r_in_atom + newdisplace call add_atom(j,atom_type_map(atom_type), sub_box_num + 1, r_in_atom) end do diff --git a/src/opt_deform.f90 b/src/opt_deform.f90 new file mode 100644 index 0000000..c9be654 --- /dev/null +++ b/src/opt_deform.f90 @@ -0,0 +1,98 @@ +module opt_deform + !This module constains the deform option which applies a uniaxial strain to the system + use parameters + use subroutines + use elements + use box + + implicit none + + + real(kind=dp), save :: applied_strain + integer, save :: sdim + + public + contains + + subroutine deform(arg_pos) + !This subroutine applies the simulation box deformation to the system + + integer, intent(inout) :: arg_pos + + character(len=1) :: dims(3) + integer :: i, j, k + real(kind=dp) :: frac_atom(atom_num), frac_node(max_basisnum, max_ng_node, ele_num) + + !initialize some variables + dims(1) = 'x' + dims(2) = 'y' + dims(3) = 'z' + + !First parse the input command + call parse_deform(arg_pos) + print *, '-----------------------Option Deform------------------------' + !Now convert all positions in the specified dimension to fractional coordinates + do i = 1, atom_num + frac_atom(i) = (r_atom(sdim, i) - box_bd(2*sdim-1))/(box_bd(2*sdim)-box_bd(2*sdim-1)) + end do + do i = 1, ele_num + do j = 1, ng_node(lat_ele(i)) + do k = 1, basisnum(lat_ele(i)) + frac_node(k,j,i) = (r_node(sdim,k,j,i) - box_bd(2*sdim-1))/(box_bd(2*sdim)-box_bd(2*sdim-1)) + end do + end do + end do + + print *, "Original box bounds in ", dims(sdim), " are ", box_bd(2*sdim-1:2*sdim) + box_bd(2*sdim) = box_bd(2*sdim) + applied_strain + print *, "New box bounds are ", box_bd(2*sdim-1:2*sdim) + + !Now reassign the positions + do i = 1, atom_num + r_atom(sdim,i) = frac_atom(i)*(box_bd(2*sdim)-box_bd(2*sdim-1)) + box_bd(2*sdim-1) + end do + do i = 1, ele_num + do j = 1, ng_node(lat_ele(i)) + do k = 1, basisnum(lat_ele(i)) + r_node(sdim,k,j,i) = frac_node(k,j,i)*(box_bd(2*sdim)-box_bd(2*sdim-1)) + box_bd(2*sdim-1) + end do + end do + end do + + end subroutine deform + + subroutine parse_deform(arg_pos) + + integer, intent(inout) :: arg_pos + integer :: arg_len + character(len=100) :: textholder + + !Pull out the dimension to be strained + arg_pos=arg_pos+1 + call get_command_argument(arg_pos, textholder, arg_len) + if (arg_len == 0) stop "Missing dim in deform command" + + select case(trim(adjustl(textholder))) + case('x','X') + sdim = 1 + case('y','Y') + sdim = 2 + case('z','Z') + sdim = 3 + case default + print *, "Dimension ", trim(adjustl(textholder)), " is not accepted. Please select either x, y, or z" + end select + + !Now pull out the strain vector, currently only accepts a real number by which to reduce the simulation cell size by in + !that dimension + arg_pos = arg_pos + 1 + call get_command_argument(arg_pos, textholder, arg_len) + if (arg_len == 0) stop "Missing strain in deform command" + read(textholder, *) applied_strain + + arg_pos = arg_pos + 1 + + end subroutine parse_deform + + +end module opt_deform