1! How to make 2! 3! gfortran -c spglib_f08.f90 4! gfortran -c example_f08.f90 5! gfortran -o example_f08 example_f08.o spglib_f08.o ~/code/spglib/src/.libs/libsymspg.a 6 7module defs_basis 8 implicit none 9 integer, parameter :: dp=kind(1.0d0) 10end module defs_basis 11 12subroutine write_syminfo( max_num_sym, num_atom, & 13 & lattice, symprec, atom_types, positions, & 14 & mesh, shifted, is_time_reversal ) 15 16 use defs_basis 17 18 use spglib_f08 19 20 implicit none 21 22 ! Arguments ------------------------------------ 23 ! scalars 24 integer, intent(in) :: num_atom, max_num_sym, is_time_reversal 25 real(dp), intent(in) :: symprec 26 ! arrays 27 integer, intent(in), dimension(num_atom) :: atom_types 28 integer, intent(in), dimension(3) :: mesh, shifted 29 real(dp), intent(in), dimension(3, 3) :: lattice 30 real(dp), intent(in), dimension(3, num_atom) :: positions 31 ! Local variables------------------------------- 32 ! scalars 33 integer :: i, j, counter, weight, space_group, num_sym, indent 34 integer :: num_ir_grid 35 character(len=21) :: international 36 character(len=7) :: schoenflies 37 character(len=30) :: space 38 character(len=128) :: fmt 39 ! arrays 40 integer, dimension(3, 3, max_num_sym) :: rotations 41 integer, dimension(3, mesh(1)*mesh(2)*mesh(3)) :: grid_point 42 integer, dimension(mesh(1)*mesh(2)*mesh(3)) :: map 43 real(dp), dimension(3, max_num_sym) :: translations 44 45 type(SpglibDataset) :: dset 46 type(SpglibSpacegroupType) :: spgtype 47 48 !************************************************************************** 49 50 space = " " 51 52! The allocatable components of dset get deallocated on scope exit 53 dset = spg_get_dataset(lattice, positions, atom_types, num_atom, symprec) 54 55 num_sym = dset % n_operations 56 57 indent = 1 58 if (dset % spacegroup_number /= 0) then 59 spgtype = spg_get_spacegroup_type(dset % hall_number) 60 61 print('(a, "space_group: ", i3)'), space(1:indent*2), dset % spacegroup_number 62 print('(a, "international: ", a, a)' ), space(1:indent*2), trim(dset % international_symbol) 63 print('(a, "schoenflies: ", a)'), space(1:indent*2), trim(spgtype % schoenflies) 64 else 65 print '("Space group could not be found,")' 66 end if 67 68 print ('(a, "atom-type:")'), space(1:indent*2) 69 do i = 1, num_atom 70 print('(a, "- { type: ", i3, "}")'), space(1:indent*2), atom_types(i) 71 end do 72 print('(a, "real-basis:")'), space(1:indent*2) 73 do i = 1, 3 74 print('(a, "- [", f19.14, ", ", f19.14, ", ", f19.14, "]")'), space(1:indent*2), lattice(i, :) 75 end do 76 print('(a, "position:")'), space(1:indent*2) 77 do i = 1, num_atom 78 print('(a, "- [", f17.14, ", ", f17.14, ", ", f17.14, "]")'), space(1:indent*2), positions(:, i) 79 end do 80 print('(a, "operation:")'), space(1:indent*2) 81 do i = 1, num_sym 82 print('(a, "- rotation: #", i4)'), space(1:indent*2), i 83 do j = 1, 3 84 print('(a, " - [", i3,",", i3,",", i3,"]")'), space(1:indent*2), dset % rotations(:,j, i) 85 end do 86 print('(a, " translation: [ ", f10.7,", ", f10.7,", ", f10.7,"]")'), space(1:indent*2), dset % translations(:,i) 87 end do 88 89 90 write(fmt, '(i0,a)') num_atom, "(x,i0,':',i0))" 91 print "(a,'wyckoffs: ', "//fmt,space(1:indent*2),(i, dset % wyckoffs(i), i = 1, dset % n_atoms) 92 print "(a,'equivalent_atoms: ', "//fmt,space(1:indent*2),(i, dset % equivalent_atoms(i) + 1, i = 1, dset % n_atoms) 93 94 print('(a, "reciprocal-mesh: [", i3, ", ", i3, ", ", i3, "]")'), space(1:indent*2), mesh(:) 95 print('(a, "- shifted: [", i3, ", ", i3, ", ", i3, "]")'), space(1:indent*2), shifted(:) 96 print('(a, "- is_time_reversal: ", i3)'), space(1:indent*2), is_time_reversal 97 98 num_ir_grid = spg_get_ir_reciprocal_mesh( grid_point, map, & 99 & mesh, shifted, is_time_reversal, lattice, positions, & 100 & atom_types, num_atom, symprec ) 101 102 print('(a, "- num-ir-grid-point:", i4)'), space(1:indent*2), num_ir_grid 103 print('(a, "- ir-grid-point:")'), space(1:indent*2) 104 counter = 0 105 106 do i = 1, product(mesh) 107 if ( i-1 == map(i) ) then 108 ! Ad-hoc and intuitive implementation of weight. Not optimal for very large size 109 weight = count(map == i-1) 110 counter = counter + 1 111 112 print '(i5, 3(x,f9.6), 2x,i4,2x, f9.6, 3(x,i4))', counter, & 113 & real(shifted + 2*grid_point(:, i))/real(2*mesh), map(i)+1, & 114 & real(weight)/real(product(mesh)), grid_point(:, i) 115 end if 116 end do 117 118end subroutine write_syminfo 119 120program spglib_example 121 122 use defs_basis 123 124 implicit none 125 126 ! max_sym is the expected maximum number of symmetry operations. 127 ! This can be very big in a supercell. 128 integer :: max_num_sym=192 129 130 integer :: num_atom 131 real(dp) :: symprec 132 integer, dimension(3, 3, 192) :: rotations 133 real(dp), dimension(3, 192) :: translations 134 real(dp), dimension(3, 3) :: lattice 135 real(dp), dimension(3, 1000) :: positions 136 integer, dimension(1000) :: atom_types 137 138 integer, dimension(3) :: mesh = [20,20,20] 139 integer, dimension(3) :: shifted 140 integer, dimension(3, 8000) :: grid_point 141 integer, dimension(8000) :: map 142 integer :: is_time_reversal = 1 143 144 real(dp) :: a 145 146! All numbers are improperly given in default (single) precision to avoid clutter. It should be fine as long as symprec > the error of sp 147 148 symprec = 1e-6_dp 149 150 151 ! Rutile 152 153 print '("Rutile:")' 154 num_atom = 6 155 lattice(1,:) = (/ 4.0_dp, 0.0_dp, 0.0_dp /) 156 lattice(2,:) = (/ 0.0_dp, 4.0_dp, 0.0_dp /) 157 lattice(3,:) = (/ 0.0_dp, 0.0_dp, 2.6_dp /) 158 positions(1:3, 1) = (/ 0.0_dp, 0.0_dp, 0.0_dp /) 159 positions(1:3, 2) = (/ 0.5_dp, 0.5_dp, 0.5_dp /) 160 positions(1:3, 3) = (/ 0.3_dp, 0.3_dp, 0.0_dp /) 161 positions(1:3, 4) = (/ 0.7_dp, 0.7_dp, 0.0_dp /) 162 positions(1:3, 5) = (/ 0.2_dp, 0.8_dp, 0.5_dp /) 163 positions(1:3, 6) = (/ 0.8_dp, 0.2_dp, 0.5_dp /) 164 atom_types(1:6) = (/ 1, 1, 2, 2, 2, 2 /) 165 shifted = [0,0,0] 166 167 call write_syminfo( max_num_sym, num_atom, & 168 lattice, symprec, atom_types, positions, & 169 mesh, shifted, is_time_reversal ) 170 171 ! FCC 172 print '("")' 173 print '("FCC:")' 174 num_atom = 1 175 lattice(1,:) = (/ 0.0_dp, 2.0_dp, 2.0_dp /) 176 lattice(2,:) = (/ 2.0_dp, 0.0_dp, 2.0_dp /) 177 lattice(3,:) = (/ 2.0_dp, 2.0_dp, 0.0_dp /) 178 positions(1:3, 1) = (/ 0.0_dp, 0.0_dp, 0.0_dp /) 179 atom_types(1) = 1 180 shifted = [0,0,0] 181 182 call write_syminfo( max_num_sym, num_atom, & 183 lattice, symprec, atom_types, positions, & 184 mesh, shifted, is_time_reversal ) 185 186 187 188 189 190 print '(/,"D0_{24}:")' 191 a = 2.89_dp 192 193 lattice(1,:) = [ sqrt(3.0_dp), -1.0_dp, 0.0_dp] 194 lattice(2,:) = [ sqrt(3.0_dp), 1.0_dp, 0.0_dp] 195 lattice(3,:) = [ 0.0_dp , 0.0_dp, sqrt(32/3.0_dp)] 196 197 198 num_atom = 16 199 positions(:, 1) = [ 0.0_dp, 0.0_dp, 0.0_dp] 200 positions(:, 2) = [ 3.0_dp, 0.0_dp, 0.0_dp] 201 positions(:, 3) = [ 0.0_dp, 3.0_dp, 0.0_dp] 202 positions(:, 4) = [ 3.0_dp, 3.0_dp, 0.0_dp] 203 positions(:, 5) = [ 2.0_dp, 2.0_dp, 1.0_dp] 204 positions(:, 6) = [ 5.0_dp, 2.0_dp, 1.0_dp] 205 positions(:, 7) = [ 2.0_dp, 5.0_dp, 1.0_dp] 206 positions(:, 8) = [ 5.0_dp, 5.0_dp, 1.0_dp] 207 positions(:, 9) = [ 0.0_dp, 0.0_dp, 2.0_dp] 208 positions(:,10) = [ 3.0_dp, 0.0_dp, 2.0_dp] 209 positions(:,11) = [ 0.0_dp, 3.0_dp, 2.0_dp] 210 positions(:,12) = [ 3.0_dp, 3.0_dp, 2.0_dp] 211 positions(:,13) = [ 1.0_dp, 1.0_dp, 3.0_dp] 212 positions(:,14) = [ 4.0_dp, 1.0_dp, 3.0_dp] 213 positions(:,15) = [ 1.0_dp, 4.0_dp, 3.0_dp] 214 positions(:,16) = [ 4.0_dp, 4.0_dp, 3.0_dp] 215 216 positions(1:2,:num_atom) = positions(1:2,:num_atom)/6 217 positions(3,:num_atom) = positions(3,:num_atom)/4 218 219 220 atom_types(:num_atom) = [ 1, 2, 2, 2, 1, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 1] 221 mesh = [ 12, 12, 6] 222 shifted = [1,1,1] 223 224 call write_syminfo( max_num_sym, num_atom, & 225 & lattice, symprec, atom_types, positions, & 226 & mesh, shifted, is_time_reversal ) 227 228end program spglib_example 229