1! How to make
2!
3! gcc -c spglib_f.c
4! gfortran -c example.f90
5! gfortran -o example example.o spglib_f.o ~/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, is_shift, is_time_reversal )
15
16  use defs_basis
17
18  implicit none
19
20  ! Arguments ------------------------------------
21  ! scalars
22  integer, intent(in) :: num_atom, max_num_sym, is_time_reversal
23  real(dp), intent(in) :: symprec
24  ! arrays
25  integer, intent(in), dimension(num_atom) :: atom_types
26  integer, intent(in), dimension(3) :: mesh, is_shift
27  real(dp), intent(in), dimension(3, 3) :: lattice
28  real(dp), intent(in), dimension(3, num_atom) :: positions
29  ! Local variables-------------------------------
30  ! scalars
31  integer :: i, j, counter, weight, sum_weight, space_group, num_sym, indent
32  integer :: num_ir_grid
33  character(len=11) :: international
34  character(len=7) :: schoenflies
35  character(len=30) :: space
36  ! arrays
37  integer, dimension(3, 3, max_num_sym) :: rotations
38  integer, dimension(3, mesh(1)*mesh(2)*mesh(3)) :: grid_point
39  integer, dimension(mesh(1)*mesh(2)*mesh(3)) :: map
40  real(dp), dimension(3, max_num_sym) :: translations
41  real(dp), dimension(3, 3) :: lattice_t
42  !**************************************************************************
43
44  space = "                              "
45
46
47  ! transpose due to array order difference between C and fortran
48  lattice_t = transpose( lattice )
49
50  call spg_get_symmetry( num_sym, rotations, translations, max_num_sym, &
51       & lattice_t, positions, atom_types, num_atom, symprec)
52
53  ! transpose due to array order difference between C and fortran
54  do i = 1, num_sym
55     rotations(:,:,i) = transpose(rotations(:,:,i))
56  end do
57
58
59  indent = 1
60  call spg_get_international( space_group, international, &
61       & lattice_t, positions, atom_types, num_atom, symprec );
62
63  if (space_group /= 0) then
64     call spg_get_schoenflies( space_group, schoenflies, &
65          & lattice_t, positions, atom_types, num_atom, symprec );
66
67     print('(a, "space_group: ", i3)'), space(1:indent*2), space_group
68     print('(a, "international: ", a, a)' ), space(1:indent*2), trim(international)
69     print('(a, "schoenflies: ", a)'), space(1:indent*2), trim(schoenflies)
70  else
71     print '("Space group could not be found,")'
72  end if
73
74  print ('(a, "atom-type:")'), space(1:indent*2)
75  do i = 1, num_atom
76     print('(a, "- { type: ", i3, " }")'), space(1:indent*2), atom_types(i)
77  end do
78  print('(a, "real-basis:")'), space(1:indent*2)
79  do i = 1, 3
80     print('(a, "- [", f19.14, ", ", f19.14, ", ", f19.14, " ]")'), space(1:indent*2), lattice(:, i)
81  end do
82  print('(a, "position:")'), space(1:indent*2)
83  do i = 1, num_atom
84     print('(a, "- [", f17.14, ", ", f17.14, ", ", f17.14, " ]")'), space(1:indent*2), positions(:, i)
85  end do
86  print('(a, "operation:")'), space(1:indent*2)
87  do i = 1, num_sym
88     print('(a, "- rotation: #", i4)'), space(1:indent*2), i
89     do j = 1, 3
90        print('(a, "  - [", i3,",", i3,",", i3," ]")'), space(1:indent*2), rotations(j,:, i)
91     end do
92     print('(a, "  translation: [ ", f10.7,", ", f10.7,", ", f10.7," ]")'), space(1:indent*2), translations(:,i)
93  end do
94
95  print('(a, "reciprocal-mesh: [", i3, ", ", i3, ", ", i3, " ]")'), space(1:indent*2), mesh(:)
96  print('(a, "- is_shift: [", i3, ", ", i3, ", ", i3, " ]")'), space(1:indent*2), is_shift(:)
97  print('(a, "- is_time_reversal: ", i3)'), space(1:indent*2), is_time_reversal
98  call spg_get_ir_reciprocal_mesh( num_ir_grid, grid_point, map, &
99       mesh, is_shift, 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  sum_weight = 0
106  do i = 1, mesh(1)*mesh(2)*mesh(3)
107     if ( i == map(i) ) then
108        ! Ad-hoc and intuitive implementation of weight
109        weight = 0
110        do j = 1, mesh(1)*mesh(2)*mesh(3)
111           if ( i == map(j) ) then
112              weight = weight + 1
113           end if
114        end do
115
116        counter = counter + 1
117        sum_weight = sum_weight + weight
118
119        print('(a, "  - #", i4)'), space(1:indent*2), counter
120        print('(a, "    address: [", i3, ", ", i3, ", ", i3, " ]")'), space(1:indent*2), grid_point(:, map(i))
121        print('(a, "    weight: ", i4)'), space(1:indent*2), weight
122     end if
123  end do
124  ! print('(a, "  sum_weight: ", i4)'), space(1:indent*2), sum_weight
125
126end subroutine write_syminfo
127
128program spglib_example
129
130  use defs_basis
131
132  implicit none
133
134  ! max_sym is the expected maximum number of symmetry operations.
135  ! This can be very big if it is supercell.
136  integer :: max_num_sym=192
137
138  integer :: num_atom
139  real(dp) :: symprec
140  integer, dimension(3, 3, 192) :: rotations
141  real(dp), dimension(3, 192) :: translations
142  real(dp), dimension(3, 3) :: lattice
143  real(dp), dimension(3, 1000) :: positions
144  integer, dimension(1000) :: atom_types
145
146  integer, dimension(3) :: mesh=(/ 20, 20, 20 /)
147  integer, dimension(3) :: is_shift=(/ 0, 0, 0/)
148  integer, dimension(3, 8000) :: grid_point
149  integer, dimension(8000) :: map
150  integer :: is_time_reversal=1
151
152  symprec = 1e-5
153
154  ! Rutile
155
156  print '("Rutile:")'
157  num_atom = 6
158  lattice(1:3, 1) = (/ 4.0, 0.0, 0.0 /)
159  lattice(1:3, 2) = (/ 0.0, 4.0, 0.0 /)
160  lattice(1:3, 3) = (/ 0.0, 0.0, 2.6 /)
161  positions(1:3, 1) = (/ 0.0, 0.0, 0.0 /)
162  positions(1:3, 2) = (/ 0.5, 0.5, 0.5 /)
163  positions(1:3, 3) = (/ 0.3, 0.3, 0.0 /)
164  positions(1:3, 4) = (/ 0.7, 0.7, 0.0 /)
165  positions(1:3, 5) = (/ 0.2, 0.8, 0.5 /)
166  positions(1:3, 6) = (/ 0.8, 0.2, 0.5 /)
167  atom_types(1:6) = (/ 1, 1, 2, 2, 2, 2 /)
168
169  call write_syminfo( max_num_sym, num_atom, &
170       lattice, symprec, atom_types, positions, &
171       mesh, is_shift, is_time_reversal )
172
173  ! FCC
174  print '("")'
175  print '("FCC:")'
176  num_atom = 1
177  lattice(1:3, 1) = (/ 0.0, 2.0, 2.0 /)
178  lattice(1:3, 2) = (/ 2.0, 0.0, 2.0 /)
179  lattice(1:3, 3) = (/ 2.0, 2.0, 0.0 /)
180  positions(1:3, 1) = (/ 0.0, 0.0, 0.0 /)
181  atom_types(1) = 1
182
183  call write_syminfo( max_num_sym, num_atom, &
184       lattice, symprec, atom_types, positions, &
185       mesh, is_shift, is_time_reversal )
186
187end program spglib_example
188
189