#include "precompilerdefinitions" program crystal_structure_info !!{!src/crystal_structure_info/manual.md!} ! Prints general info about the lattice, high symmetry points and so on. use konstanter, only: flyt,lo_tol,lo_sqtol,lo_bohr_to_A use helpers, only: lo_frobnorm,lo_mpi_helper,lo_chop,open_file,tochar,lo_clean_fractional_coordinates,lo_determ,& lo_fetch_tolerance,walltime,lo_return_unique,lo_invert3x3matrix,qsort use options, only: lo_opts use type_crystalstructure, only: lo_crystalstructure use type_qpointmesh, only: lo_bandstructure ! implicit none ! for normal phonon things type(lo_opts) :: opts type(lo_crystalstructure) :: uc type(lo_bandstructure) :: bs type(lo_mpi_helper) :: mw ! real(flyt), dimension(3) :: v0,v1 real(flyt) :: tstart,f0 integer :: i,j,u character(len=1000) :: opf ! ! Get the command line arguments ! call mw%init() tstart=walltime() call opts%parse() ! Read unitcell call uc%readfromfile('infile.ucposcar',verbosity=opts%verbosity) call uc%classify('wedge',timereversal=opts%timereversal) write(*,*) 'Info about the crystal structure:' write(*,*) '... I believe the basis is ',trim(uc%info%bravaislatticelongname),' (',trim(uc%info%bravaislattice),') with the basis vectors' write(*,"(1X,'a1:',3(2X,F18.12))") uc%latticevectors(:,1)*lo_bohr_to_A write(*,"(1X,'a2:',3(2X,F18.12))") uc%latticevectors(:,2)*lo_bohr_to_A write(*,"(1X,'a3:',3(2X,F18.12))") uc%latticevectors(:,3)*lo_bohr_to_A write(*,*) '... and the atoms at these positions (in fractional coordinates)' do i=1,uc%na write(*,'(1X,A8,3(2x,F15.13))') trim(uc%atomic_symbol(uc%species(i))),uc%r(:,i) enddo ! Print what the conventional basis is f0=lo_frobnorm(uc%latticevectors-uc%info%unique_conventional_basis) if ( f0 .gt. lo_tol ) then write(*,*) 'I believe this is the conventional lattice:' write(*,"(1X,'a1:',3(2X,F18.12))") uc%info%unique_conventional_basis(:,1)*lo_bohr_to_A write(*,"(1X,'a2:',3(2X,F18.12))") uc%info%unique_conventional_basis(:,2)*lo_bohr_to_A write(*,"(1X,'a3:',3(2X,F18.12))") uc%info%unique_conventional_basis(:,3)*lo_bohr_to_A endif !! It might not agree with my canonical form. !f0=lo_frobnorm(uc%latticevectors-uc%info%unique_primitive_basis) !if ( f0 .gt. lo_tol ) then ! write(*,*) ' ' ! write(*,*) 'I like this primitive basis better:' ! write(*,"(1X,'a1:',3(2X,F18.12))") uc%info%unique_primitive_basis(:,1) ! write(*,"(1X,'a2:',3(2X,F18.12))") uc%info%unique_primitive_basis(:,2) ! write(*,"(1X,'a3:',3(2X,F18.12))") uc%info%unique_primitive_basis(:,3) ! write(*,*) '... with the atoms at these positions (in fractional coordinates)' ! ! Convert the coordinates ! allocate(dum(3,uc%na)) ! do i=1,uc%na ! ! get the current coordinate ! v0=uc%r(:,i) ! v0=matmul(uc%info%permutation_to_unique,v0) ! v0=lo_chop(lo_clean_fractional_coordinates(v0),lo_sqtol) ! dum(:,i)=v0 ! write(*,'(1X,A8,3(2x,F15.13))') trim(uc%atomic_symbol(uc%species(i))),dum(:,i) ! enddo !endif ! Get the high symmetry points if ( uc%info%havewedge .eqv. .false. ) call uc%classify('wedge',timereversal=opts%timereversal) write(*,*) ' ' write(*,*) 'Available high symmetry points:' write(*,*) '(the ones called NP are in fact high symmetry points, but noone has bothered to name them)' write(*,*) ' Cartesian coordinates Fractional coordinates' do i=1,uc%irrw%nnodes v0=lo_chop(uc%irrw%r(:,i),lo_tol) v1=lo_chop(matmul(uc%inv_reciprocal_latticevectors,v0),lo_tol) write(*,"(1X,I2,1X,A4,3(1X,F13.10),3X,3(1X,F13.10))") i,uc%irrw%label(i),v0/lo_bohr_to_A,v1 enddo call uc%write_bz_to_hdf5('outfile.brillouin_zone.hdf5') call bs%standardpath(uc) bs%npts=100 write(*,*) ' ' write(*,*) 'Per default, the following path will be used:' do i=1,bs%npath write(*,*) trim(bs%symb_q_start(i)),' -> ',trim(bs%symb_q_end(i)) enddo write(*,*) 'It is written in "outfile.qpoints_dispersion", modify and' write(*,*) 'copy it to "infile.qpoints_dispersion" if you want' ! Dump a path file u=open_file('out','outfile.qpoints_dispersion') opf="(A5,22X,' ! Bravais lattice type')" write(u,opf) uc%info%bravaislattice opf="(I5,22X,' ! Number of points on each path')" write(u,opf) 100 opf="(I5,22X,' ! Number paths between special points')" write(u,opf) bs%npath opf="(A3,1X,A3,20X,' ! Starting and ending special point')" write(u,opf) bs%symb_q_start(1),bs%symb_q_end(1) opf="(A3,1X,A3,20X,' !')" do i=2,bs%npath write(u,opf) bs%symb_q_start(i),bs%symb_q_end(i) enddo close(u) if ( opts%printsymmetry ) then write(*,*) ' Symmetry operations:' do i=1,uc%sym%n write(*,*) 'Operation ',tochar(i),' ',uc%sym%op(i)%whatkind write(*,*) ' Axis:',lo_chop(uc%sym%op(i)%axis,lo_sqtol) if ( uc%sym%op(i)%fold .gt. 0 ) write(*,*) ' Angle: ',tochar(int(anint(360.0_flyt/uc%sym%op(i)%fold))) opf=' FMAP:' do j=1,size(uc%sym%op(i)%fmap) opf=trim(opf)//' '//tochar(j)//'->'//tochar(uc%sym%op(i)%fmap(j)) enddo write(*,*) trim(opf) !write(*,*) ' FMAP:',uc%sym%op(i)%fmap write(*,'(1X,A,3(1X,F9.6))') ' translation:',uc%sym%op(i)%ftr write(*,*) ' ... in Cartesian coordinates' do j=1,3 write(*,"(6(3X,F15.11))") uc%sym%op(i)%m(j,:) !write(*,"(3(3X,A10))") lo_fancy_real2char(uc%sym%op(i)%m(j,1)),lo_fancy_real2char(uc%sym%op(i)%m(j,2)),lo_fancy_real2char(uc%sym%op(i)%m(j,3)) enddo write(*,*) ' ... in fractional coordinates' do j=1,3 write(*,"(6(3X,F15.11))") uc%sym%op(i)%fm(j,:) !write(*,"(6(3X,A3))") lo_fancy_real2char(uc%sym%op(i)%fm(j,1)),lo_fancy_real2char(uc%sym%op(i)%fm(j,2)),lo_fancy_real2char(uc%sym%op(i)%fm(j,3)) enddo write(*,*) ' ' enddo ! Dump symmetry operations to Mathematica-format, in case you need to play with it u=open_file('out','outfile.symops_mathematica') write(u,*) 'ops=ConstantArray[0,{'//tochar(uc%sym%n)//',3,3}];' do i=1,uc%sym%n do j=1,3 write(u,*) 'ops[['//tochar(i)//','//tochar(j)//& ']]={ '//lo_fancy_real2char(uc%sym%op(i)%m(1,j)) //','//lo_fancy_real2char(uc%sym%op(i)%m(2,j))//','//lo_fancy_real2char(uc%sym%op(i)%m(3,j))//'};' enddo enddo close(u) endif write(*,*) 'all done' write(*,*) ' ' call mw%destroy() contains !!> Find the primitive cell !subroutine find_primitive_cell(basis,r,species,flavor,primbasis,tolerance) ! !> cell to be reduced ! real(flyt), dimension(3,3), intent(in) :: basis ! !> atoms in the cell ! real(flyt), dimension(:,:), intent(in) :: r ! !> species classification ! integer, dimension(:), intent(in) :: species ! !> flavor classification ! integer, dimension(:), intent(in) :: flavor ! !> primitive lattice vectors ! real(flyt), dimension(3,3), intent(out) :: primbasis ! !> tolerance ! real(flyt), intent(in) :: tolerance ! ! ! real(flyt), dimension(:,:), allocatable :: r0,r1 ! real(flyt), dimension(3,3) :: basis0,basis1 ! integer, dimension(:), allocatable :: species0,flavor0,species1,flavor1 ! integer :: i,na ! ! ! na=size(r,2) ! allocate(r0(3,na)) ! allocate(species0(na)) ! allocate(flavor0(na)) ! basis0=basis ! r0=r ! species0=species ! flavor0=flavor ! ! ! Recursively find a smaller cell ! redloop: do ! ! Try to find a smaller cell ! call reduce_cell(basis0,r0,species0,flavor0,basis1,r1,species1,flavor1,tolerance) ! if ( abs(abs(lo_determ(basis0))-abs(lo_determ(basis1))) .lt. tolerance ) then ! ! this means no new cell, we are converged. ! exit redloop ! else ! ! this means a new cell, update and try again! ! if ( allocated(r0) ) deallocate(r0) ! if ( allocated(species0) ) deallocate(species0) ! if ( allocated(flavor0) ) deallocate(flavor0) ! i=size(r1,2) ! allocate(r0(3,i)) ! allocate(species0(i)) ! allocate(flavor0(i)) ! basis0=basis1 ! r0=r1 ! species0=species1 ! flavor0=flavor1 ! if ( allocated(r1) ) deallocate(r1) ! if ( allocated(species1) ) deallocate(species1) ! if ( allocated(flavor1) ) deallocate(flavor1) ! endif ! enddo redloop ! ! ! Here we should be done. ! primbasis=basis0 ! ! if ( allocated(r0) ) deallocate(r0) ! if ( allocated(r1) ) deallocate(r1) ! if ( allocated(species0) ) deallocate(species0) ! if ( allocated(species1) ) deallocate(species1) ! if ( allocated(flavor0) ) deallocate(flavor0) ! if ( allocated(flavor1) ) deallocate(flavor1) !end subroutine ! !> single iteration of the cell reduction ! subroutine reduce_cell(latticevectors,r,species,flavor,newlatticevectors,newr,newspecies,newflavor,tolerance) ! !> original lattice vectors ! real(flyt), dimension(3,3), intent(in) :: latticevectors ! !> original positions ! real(flyt), dimension(:,:), intent(in) :: r ! !> original species ! integer, dimension(:), intent(in) :: species ! !> additional flavor for symmetry thing ! integer, dimension(:), intent(in) :: flavor ! !> tolerance ! real(flyt), intent(in) :: tolerance ! ! real(flyt), dimension(3,3), intent(out) :: newlatticevectors ! real(flyt), dimension(:,:), allocatable, intent(out) :: newr ! integer, dimension(:), allocatable, intent(out) :: newspecies ! integer, dimension(:), allocatable,intent(out) :: newflavor ! ! real(flyt), dimension(:,:), allocatable :: vecs,dum ! real(flyt), dimension(3,3) :: m0,m1 ! real(flyt) :: vol,f0,tol,reltol,frtol ! integer :: i,j,k,l,nvecs,na ! logical :: cellvalid ! ! ! set the tolerance ! call lo_fetch_tolerance(tolerance,latticevectors,realspace_cart_tol=tol,realspace_fractional_tol=frtol,relative_tol=reltol) ! ! na=size(r,2) ! ! Get a list of possible lattice vectors: All the vectors in the cell + old latticevectors ! l=na+3 ! allocate(dum(3,na+3)) ! dum(:,1:na)=lo_clean_fractional_coordinates(r,frtol**2) ! do i=1,na ! dum(:,i)=matmul(latticevectors,dum(:,i)) ! enddo ! do i=1,3 ! dum(:,i+na)=latticevectors(:,i) ! enddo ! ! Reduce these to the unique. Should not be necessary, but never hurts. ! call lo_return_unique(dum,vecs,frtol) ! nvecs=size(vecs,2) ! ! ! volume of original cell ! vol=abs(lo_determ(latticevectors)) ! ! Loop over all vectors to try to find a new cell. ! vl1: do i=1,nvecs ! vl2: do j=i+1,nvecs ! vl3: do k=j+1,nvecs ! ! possible new set of lattice vectors ! m0(:,1)=vecs(:,i) ! m0(:,2)=vecs(:,j) ! m0(:,3)=vecs(:,k) ! ! some quick tests to see if there is any point in continuing: ! f0=abs(lo_determ(m0)) ! if ( f0 .lt. tol ) cycle vl3 ! stupid cell if volume is zero ! if ( abs(mod(vol/f0,1.0_flyt)) .gt. reltol ) cycle vl3 ! has to be an integer division ! if ( int(anint(vol/f0)) .eq. 1 ) cycle vl3 ! has to make the cell smaller, not the same ! ! I can get really retarded cells from this, do the Delaunay thing on it ! call reduce_basis_to_something_decent(m0,m1) ! ! Now do an actual test ! call test_cell(latticevectors,r,species,flavor,m1,cellvalid,newr,newspecies,newflavor,tolerance) ! if ( cellvalid ) then ! ! If I find a new cell, we will exit this routine and start over. Makes things way way faster. ! newlatticevectors=m1 ! deallocate(vecs) ! return ! endif ! enddo vl3 ! enddo vl2 ! enddo vl1 ! ! If we made it here, the cell is already the smallest. Return the same basis again. ! newlatticevectors=latticevectors ! end subroutine ! ! !> test if a new cell is a smaller unit cell ! subroutine test_cell(oldcell,oldr,oldspecies,oldflavor,cell,valid,newr,newspecies,newflavor,tolerance) ! !> original cell ! real(flyt), dimension(3,3), intent(in) :: oldcell ! !> original positions ! real(flyt), dimension(:,:), intent(in) :: oldr ! !> original species ! integer, dimension(:), intent(in) :: oldspecies ! !> original flavor ! integer, dimension(:), intent(in) :: oldflavor ! !> new cell ! real(flyt), dimension(3,3), intent(in) :: cell ! !> is it a new cell? ! logical, intent(out) :: valid ! !> new positions in case of a valid cell ! real(flyt), dimension(:,:), allocatable, intent(out) :: newr ! !> new species in case of a valid cell ! integer, dimension(:), allocatable, intent(out) :: newspecies ! !> new flavor in case of a valid cell ! integer, dimension(:), allocatable, intent(out) :: newflavor ! !> tolerance ! real(flyt), intent(in) :: tolerance ! ! ! real(flyt), dimension(3,3) :: icell,m0 ! real(flyt), dimension(:,:), allocatable :: dumr1,dumr2 ! real(flyt) :: tol,sqtol ! integer :: i,j,l,mult,oldna,newna ! ! ! set the tolerance ! call lo_fetch_tolerance(tolerance,oldcell,realspace_fractional_tol=tol) ! sqtol=tol**2 ! ! valid=.false. ! ! How many atoms should I get of each? ! mult=int(anint(abs(lo_determ(oldcell)/lo_determ(cell)))) ! ! ! m0 is a matrix that converts to fractional coordinates in the new cell ! icell=lo_invert3x3matrix( cell ) ! m0=matmul(icell,oldcell) ! ! ! convert coordinates to fractional with respect to the new cell. Also ! ! I decorate the positions with the species and flavor, to keep track of that. ! oldna=size(oldr,2) ! allocate(dumr1(5,oldna)) ! do i=1,oldna ! dumr1(1:3,i)=lo_clean_fractional_coordinates(matmul(m0,oldr(:,i)),sqtol) ! dumr1(4:5,i)=[oldspecies(i),oldflavor(i)] ! enddo ! ! Reduce this to the union of the list, the unique ! call lo_return_unique(dumr1,dumr2,tol) ! ! ! Did it become the correct number of atoms? ! if ( abs(size(dumr2,2)*mult-oldna) .ne. 0 ) then ! valid=.false. ! return ! endif ! ! ! Slightly better check: Each of the atoms in the reduced cell must appear ! ! exactly mult times in the old cell. ! do i=1,size(dumr2,2) ! l=0 ! do j=1,oldna ! if ( sum(abs(dumr1(:,j)-dumr2(:,i))) .lt. tol ) l=l+1 ! enddo ! ! kill if not true ! if ( l .ne. mult ) then ! valid=.false. ! return ! endif ! enddo ! ! ! ok, now I say this is valid, return the new cell ! newna=size(dumr2,2) ! allocate(newr(3,newna)) ! allocate(newspecies(newna)) ! allocate(newflavor(newna)) ! newr=dumr2(1:3,:) ! newspecies=int(anint(dumr2(4,:))) ! newflavor=int(anint(dumr2(5,:))) ! valid=.true. ! end subroutine !!> reduce a set of lattice vectors to canonical form !subroutine reduce_basis_to_something_decent(basis,canonical_basis) ! !> original lattice vectors ! real(flyt), dimension(3,3), intent(in) :: basis ! !> canonical lattice vectors ! real(flyt), dimension(3,3), intent(out) :: canonical_basis ! ! real(flyt), dimension(3,4) :: ext_basis ! real(flyt), dimension(3,7) :: dum ! real(flyt), dimension(7) :: d ! real(flyt) :: f0 ! integer, dimension(7) :: ind ! integer :: i,j,k ! ! ! First do the Delaunay reduction thing ! ext_basis=0.0_flyt ! ext_basis(:,1:3)=basis ! do i=1,3 ! ext_basis(:,4)=ext_basis(:,4)-basis(:,i) ! enddo ! ! Iterate and stuff ! iterloop: do ! do i=1,4 ! do j=i+1,4 ! ! check scalar product ! f0=dot_product(ext_basis(:,i),ext_basis(:,j)) ! if ( f0 .gt. lo_sqtol ) then ! ! subtract i from all not i and j ! do k=1,4 ! if ( k .ne. i .and. k .ne. j ) then ! ext_basis(:,k)=ext_basis(:,k)+ext_basis(:,i) ! endif ! enddo ! ! flip sign of i ! ext_basis(:,i)=-ext_basis(:,i) ! ! start over! ! cycle iterloop ! endif ! enddo ! enddo ! ! If we make it here, we are done! ! exit iterloop ! enddo iterloop ! ! ! Look through all variants to get shortest vectors ! dum=0.0_flyt ! dum(:,1:4)=ext_basis ! dum(:,5)=ext_basis(:,1)+ext_basis(:,2) ! dum(:,6)=ext_basis(:,2)+ext_basis(:,3) ! dum(:,7)=ext_basis(:,1)+ext_basis(:,3) ! do i=1,7 ! d(i)=norm2(dum(:,i)) ! enddo ! call qsort(d,ind) ! ! ! do i=3,7 ! canonical_basis(:,1)=dum(:,ind(1)) ! canonical_basis(:,2)=dum(:,ind(2)) ! canonical_basis(:,3)=dum(:,ind(i)) ! if ( abs(lo_determ(canonical_basis)) .gt. lo_tol ) exit ! enddo ! if ( lo_determ(canonical_basis) .lt. 0.0_flyt ) canonical_basis=-canonical_basis !end subroutine ! !> Check if an array has repeated values ! logical function lo_is_vector_degenerate_brute(v,thres) ! !> array to test ! real(flyt), dimension(:), intent(in) :: v ! !> tolerance ! real(flyt), intent(in) :: thres ! ! ! integer :: i,j ! real(flyt) :: f0,f1 ! ! ! f0=lo_huge ! do i=1,size(v,1) ! do j=i+1,size(v,1) ! f1=abs(v(i)-v(j)) ! if ( f1 .lt. f0 ) f0=f1 ! enddo ! enddo ! ! ! if ( f0 .lt. thres ) then ! lo_is_vector_degenerate_brute=.true. ! else ! lo_is_vector_degenerate_brute=.false. ! endif ! end function ! ! !> Check if an array has repeated values ! logical function lo_is_vector_degenerate_qsort(v,thres) ! !> array to test ! real(flyt), dimension(:), intent(in) :: v ! !> tolerance ! real(flyt), intent(in) :: thres ! ! ! integer :: i,j ! real(flyt) :: f0,f1 ! real(flyt), dimension(:), allocatable :: dum ! ! ! lo_allocate(dum(size(v,1))) ! dum=v ! call qsort(dum) ! ! ! lo_is_vector_degenerate_qsort=.false. ! do i=1,size(dum,1)-1 ! if ( abs(dum(i)-dum(i+1)) .lt. thres ) then ! lo_is_vector_degenerate_qsort=.true. ! return ! endif ! enddo ! end function !> converts a double to a character, and tries to find well defined numbers function lo_fancy_real2char(x,tolerance) result(s) !> number to be converted real(flyt), intent(in) :: x !> optional tolerance real(flyt), intent(in), optional :: tolerance !> the string character(len=:), allocatable :: s real(flyt) :: tol real(flyt), dimension(67), parameter :: welldefined=[& 0.0000000000000000_flyt,& 0.1000000000000000_flyt,& -0.1000000000000000_flyt,& 0.1111111111111111_flyt,& -0.1111111111111111_flyt,& 0.1250000000000000_flyt,& -0.1250000000000000_flyt,& 0.1428571428571428_flyt,& -0.1428571428571428_flyt,& 0.1666666666666667_flyt,& -0.1666666666666667_flyt,& 0.2000000000000000_flyt,& -0.2000000000000000_flyt,& 0.2222222222222222_flyt,& -0.2222222222222222_flyt,& 0.2500000000000000_flyt,& -0.2500000000000000_flyt,& 0.2857142857142857_flyt,& -0.2857142857142857_flyt,& 0.3000000000000000_flyt,& -0.3000000000000000_flyt,& 0.3333333333333333_flyt,& -0.3333333333333333_flyt,& 0.3750000000000000_flyt,& -0.3750000000000000_flyt,& 0.4000000000000000_flyt,& -0.4000000000000000_flyt,& 0.4285714285714285_flyt,& -0.4285714285714285_flyt,& 0.4444444444444444_flyt,& -0.4444444444444444_flyt,& 0.5000000000000000_flyt,& -0.5000000000000000_flyt,& 0.5555555555555556_flyt,& -0.5555555555555556_flyt,& 0.5714285714285714_flyt,& -0.5714285714285714_flyt,& 0.6000000000000000_flyt,& -0.6000000000000000_flyt,& 0.6250000000000000_flyt,& -0.6250000000000000_flyt,& 0.6666666666666667_flyt,& -0.6666666666666667_flyt,& 0.7000000000000000_flyt,& -0.7000000000000000_flyt,& 0.7142857142857143_flyt,& -0.7142857142857143_flyt,& 0.7500000000000000_flyt,& -0.7500000000000000_flyt,& 0.7777777777777778_flyt,& -0.7777777777777778_flyt,& 0.8000000000000000_flyt,& -0.8000000000000000_flyt,& 0.8333333333333334_flyt,& -0.8333333333333334_flyt,& 0.8571428571428571_flyt,& -0.8571428571428571_flyt,& 0.8660254037844386_flyt,& -0.8660254037844386_flyt,& 0.8750000000000000_flyt,& -0.8750000000000000_flyt,& 0.8888888888888888_flyt,& -0.8888888888888888_flyt,& 0.9000000000000000_flyt,& -0.9000000000000000_flyt,& 1.0000000000000000_flyt,& -1.0000000000000000_flyt] character(len=12), dimension(67), parameter :: welldefined_string=[& '0 ',& ! 0.0000000000000000_flyt,& '1/10 ',& ! 0.1000000000000000_flyt,& '-1/10 ',& ! -0.1000000000000000_flyt,& 'x ',& ! 0.1111111111111111_flyt,& 'x ',& ! -0.1111111111111111_flyt,& '1/8 ',& ! 0.1250000000000000_flyt,& '/1/8 ',& ! -0.1250000000000000_flyt,& 'x ',& ! 0.1428571428571428_flyt,& 'x ',& ! -0.1428571428571428_flyt,& '1/6 ',& ! 0.1666666666666667_flyt,& '-1/6 ',& ! -0.1666666666666667_flyt,& '1/5 ',& ! 0.2000000000000000_flyt,& '-1/5 ',& ! -0.2000000000000000_flyt,& 'x ',& ! 0.2222222222222222_flyt,& 'x ',& ! -0.2222222222222222_flyt,& '1/4 ',& ! 0.2500000000000000_flyt,& '-1/4 ',& ! -0.2500000000000000_flyt,& 'x ',& ! 0.2857142857142857_flyt,& 'x ',& ! -0.2857142857142857_flyt,& 'x ',& ! 0.3000000000000000_flyt,& 'x ',& ! -0.3000000000000000_flyt,& 'x ',& ! 0.3333333333333333_flyt,& 'x ',& ! -0.3333333333333333_flyt,& 'x ',& ! 0.3750000000000000_flyt,& 'x ',& ! -0.3750000000000000_flyt,& 'x ',& ! 0.4000000000000000_flyt,& 'x ',& ! -0.4000000000000000_flyt,& 'x ',& ! 0.4285714285714285_flyt,& 'x ',& ! -0.4285714285714285_flyt,& 'x ',& ! 0.4444444444444444_flyt,& 'x ',& ! -0.4444444444444444_flyt,& '1/2 ',& ! 0.5000000000000000_flyt,& '-1/2 ',& ! -0.5000000000000000_flyt,& 'x ',& ! 0.5555555555555556_flyt,& 'x ',& ! -0.5555555555555556_flyt,& 'x ',& ! 0.5714285714285714_flyt,& 'x ',& ! -0.5714285714285714_flyt,& 'x ',& ! 0.6000000000000000_flyt,& 'x ',& ! -0.6000000000000000_flyt,& 'x ',& ! 0.6250000000000000_flyt,& 'x ',& ! -0.6250000000000000_flyt,& 'x ',& ! 0.6666666666666667_flyt,& 'x ',& ! -0.6666666666666667_flyt,& 'x ',& ! 0.7000000000000000_flyt,& 'x ',& ! -0.7000000000000000_flyt,& 'x ',& ! 0.7142857142857143_flyt,& 'x ',& ! -0.7142857142857143_flyt,& 'x ',& ! 0.7500000000000000_flyt,& 'x ',& ! -0.7500000000000000_flyt,& 'x ',& ! 0.7777777777777778_flyt,& 'x ',& ! -0.7777777777777778_flyt,& 'x ',& ! 0.8000000000000000_flyt,& 'x ',& ! -0.8000000000000000_flyt,& 'x ',& ! 0.8333333333333334_flyt,& 'x ',& ! -0.8333333333333334_flyt,& 'x ',& ! 0.8571428571428571_flyt,& 'x ',& ! -0.8571428571428571_flyt,& 'Sqrt[3]/2 ',& ! 0.8660254037844386_flyt,& '-Sqrt[3]/2 ',& ! -0.8660254037844386_flyt,& 'x ',& ! 0.8750000000000000_flyt,& 'x ',& ! -0.8750000000000000_flyt,& 'x ',& ! 0.8888888888888888_flyt,& 'x ',& ! -0.8888888888888888_flyt,& '9/10 ',& ! 0.9000000000000000_flyt,& '-9/10 ',& ! -0.9000000000000000_flyt,& '1 ',& ! 1.0000000000000000_flyt,& '-1 '] ! -1.0000000000000000_flyt] integer :: i if ( present(tolerance) ) then tol=tolerance else tol=lo_sqtol endif do i=1,67 if ( abs(x-welldefined(i)) .lt. tol ) then s=trim(welldefined_string(i)) return endif enddo s='dunno' end function end program