1 #include "spglib.h"
2 #include <stdio.h>
3 #include <stdlib.h>
4 #include <math.h>
5
6 static void test_spg_get_symmetry(void);
7 static void test_spg_get_symmetry_with_collinear_spin(void);
8 static void test_spg_get_multiplicity(void);
9 static void test_spg_find_primitive_BCC(void);
10 static void test_spg_find_primitive_corundum(void);
11 static void test_spg_standardize_cell_BCC(void);
12 static void test_spg_standardize_cell_BCC_prim(void);
13 static void test_spg_standardize_cell_corundum(void);
14 static int sub_spg_standardize_cell(double lattice[3][3],
15 double position[][3],
16 int types[],
17 const int num_atom,
18 const double symprec,
19 const int to_primitive,
20 const int no_idealize);
21 static void test_spg_get_international(void);
22 static void test_spg_get_schoenflies(void);
23 static void test_spg_get_spacegroup_type(void);
24 static void test_spg_get_symmetry_from_database(void);
25 static void test_spg_refine_cell_BCC(void);
26 static void test_spg_get_dataset(void);
27 static void test_spg_get_ir_reciprocal_mesh(void);
28 static void test_spg_get_stabilized_reciprocal_mesh(void);
29 static void test_spg_relocate_BZ_grid_address(void);
30 static void test_spg_get_error_message(void);
31 static void show_spg_dataset(double lattice[3][3],
32 const double origin_shift[3],
33 double position[][3],
34 const int num_atom,
35 const int types[]);
36 static void show_cell(double lattice[3][3],
37 double position[][3],
38 const int types[],
39 const int num_atom);
40
main(void)41 int main(void)
42 {
43 test_spg_find_primitive_BCC();
44 test_spg_find_primitive_corundum();
45 test_spg_standardize_cell_BCC();
46 test_spg_standardize_cell_BCC_prim();
47 test_spg_standardize_cell_corundum();
48 test_spg_get_multiplicity();
49 test_spg_get_symmetry();
50 test_spg_get_symmetry_with_collinear_spin();
51 test_spg_get_international();
52 test_spg_get_schoenflies();
53 test_spg_get_spacegroup_type();
54 test_spg_get_symmetry_from_database();
55 test_spg_refine_cell_BCC();
56 test_spg_get_dataset();
57 test_spg_get_ir_reciprocal_mesh();
58 test_spg_get_stabilized_reciprocal_mesh();
59 test_spg_relocate_BZ_grid_address();
60 test_spg_get_error_message();
61
62 return 0;
63 }
64
test_spg_find_primitive_BCC(void)65 static void test_spg_find_primitive_BCC(void)
66 {
67 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 4}};
68 double position[][3] = {
69 {0, 0, 0},
70 {0.5, 0.5, 0.5}
71 };
72 int types[] = {1, 1};
73 int i, num_atom = 2, num_primitive_atom;
74 double symprec = 1e-5;
75
76 /* lattice, position, and types are overwirtten. */
77 printf("*** Example of spg_find_primitive (BCC unitcell --> primitive) ***:\n");
78 num_primitive_atom = spg_find_primitive(lattice, position, types, num_atom, symprec);
79 if (num_primitive_atom == 0) {
80 printf("Primitive cell was not found.\n");
81 } else {
82 show_cell(lattice, position, types, num_primitive_atom);
83 }
84 }
85
test_spg_find_primitive_corundum(void)86 static void test_spg_find_primitive_corundum(void)
87 {
88 double lattice[3][3] = {{4.8076344022756095, -2.4038172011378047, 0},
89 {0, 4.1635335244786962, 0},
90 {0, 0, 13.1172699198127543}};
91 double position[][3] = {
92 {0.0000000000000000, 0.0000000000000000, 0.3521850942289043},
93 {0.6666666666666643, 0.3333333333333357, 0.6855184275622400},
94 {0.3333333333333357, 0.6666666666666643, 0.0188517608955686},
95 {0.0000000000000000, 0.0000000000000000, 0.6478149057711028},
96 {0.6666666666666643, 0.3333333333333357, 0.9811482391044314},
97 {0.3333333333333357, 0.6666666666666643, 0.3144815724377600},
98 {0.0000000000000000, 0.0000000000000000, 0.1478149057710957},
99 {0.6666666666666643, 0.3333333333333357, 0.4811482391044314},
100 {0.3333333333333357, 0.6666666666666643, 0.8144815724377600},
101 {0.0000000000000000, 0.0000000000000000, 0.8521850942288972},
102 {0.6666666666666643, 0.3333333333333357, 0.1855184275622400},
103 {0.3333333333333357, 0.6666666666666643, 0.5188517608955686},
104 {0.3061673906454899, 0.0000000000000000, 0.2500000000000000},
105 {0.9728340573121541, 0.3333333333333357, 0.5833333333333357},
106 {0.6395007239788255, 0.6666666666666643, 0.9166666666666643},
107 {0.6938326093545102, 0.0000000000000000, 0.7500000000000000},
108 {0.3604992760211744, 0.3333333333333357, 0.0833333333333357},
109 {0.0271659426878458, 0.6666666666666643, 0.4166666666666643},
110 {0.0000000000000000, 0.3061673906454899, 0.2500000000000000},
111 {0.6666666666666643, 0.6395007239788255, 0.5833333333333357},
112 {0.3333333333333357, 0.9728340573121541, 0.9166666666666643},
113 {0.0000000000000000, 0.6938326093545102, 0.7500000000000000},
114 {0.6666666666666643, 0.0271659426878458, 0.0833333333333357},
115 {0.3333333333333357, 0.3604992760211744, 0.4166666666666643},
116 {0.6938326093545102, 0.6938326093545102, 0.2500000000000000},
117 {0.3604992760211744, 0.0271659426878458, 0.5833333333333357},
118 {0.0271659426878458, 0.3604992760211744, 0.9166666666666643},
119 {0.3061673906454899, 0.3061673906454899, 0.7500000000000000},
120 {0.9728340573121541, 0.6395007239788255, 0.0833333333333357},
121 {0.6395007239788255, 0.9728340573121541, 0.4166666666666643},
122 };
123 int types[30];
124 int i, num_atom = 30, num_primitive_atom;
125 double symprec = 1e-5;
126
127 for (i = 0; i < 12; i++) {
128 types[i] = 1;
129 }
130 for (i = 12; i < 30; i++) {
131 types[i] = 2;
132 }
133
134 /* lattice, position, and types are overwirtten. */
135 printf("*** Example of spg_find_primitive (Corundum) ***:\n");
136 num_primitive_atom = spg_find_primitive(lattice, position, types, num_atom, symprec);
137 if (num_primitive_atom == 0) {
138 printf("Primitive cell was not found.\n");
139 } else {
140 show_cell(lattice, position, types, num_primitive_atom);
141 }
142 }
143
test_spg_refine_cell_BCC(void)144 static void test_spg_refine_cell_BCC(void)
145 {
146 double lattice[3][3] = {{0, 2, 2}, {2, 0, 2}, {2, 2, 0}};
147
148 /* 4 times larger memory space must be prepared. */
149 double position[4][3];
150 int types[4];
151
152 int i, num_atom_bravais, num_atom = 1;
153 double symprec = 1e-5;
154
155 position[0][0] = 0;
156 position[0][1] = 0;
157 position[0][2] = 0;
158 types[0] = 1;
159
160 /* lattice, position, and types are overwirtten. */
161 printf("*** Example of spg_refine_cell ***:\n");
162 num_atom_bravais = spg_refine_cell( lattice,
163 position,
164 types,
165 num_atom,
166 symprec );
167 show_cell(lattice, position, types, num_atom_bravais);
168 }
169
test_spg_standardize_cell_BCC(void)170 static void test_spg_standardize_cell_BCC(void)
171 {
172 double lattice[3][3] = {{3.97, 0, 0}, {0, 4.03, 0}, {0, 0, 4.0}};
173 double position[][3] = {
174 {0.002, 0, 0},
175 {0.5, 0.5001, 0.5}
176 };
177 int types[] = {1, 1};
178 int i, j, k, num_atom = 2, num_primitive_atom;
179 double symprec = 1e-1;
180
181 /* lattice, position, and types are overwirtten. */
182 printf("*** Example of spg_standardize_cell (BCC unitcell) ***:\n");
183 printf("------------------------------------------------------\n");
184 for (j = 0; j < 2; j++) {
185 for (k = 0; k < 2; k++) {
186 sub_spg_standardize_cell(lattice,
187 position,
188 types,
189 num_atom,
190 symprec,
191 j,
192 k);
193 printf("------------------------------------------------------\n");
194 }
195 }
196 }
197
test_spg_standardize_cell_BCC_prim(void)198 static void test_spg_standardize_cell_BCC_prim(void)
199 {
200 double lattice[3][3] = {{-2.01, 2, 2}, {2, -2.02, 2}, {2, 2, -2.03}};
201 double position[][3] = {
202 {0.002, 0, 0},
203 };
204 int types[] = {1};
205 int i, j, k, num_atom = 1, num_primitive_atom;
206 double symprec = 1e-1;
207
208 /* lattice, position, and types are overwirtten. */
209 printf("*** Example of spg_standardize_cell (BCC primitive) ***:\n");
210 printf("------------------------------------------------------\n");
211 for (j = 0; j < 2; j++) {
212 for (k = 0; k < 2; k++) {
213 sub_spg_standardize_cell(lattice,
214 position,
215 types,
216 num_atom,
217 symprec,
218 j,
219 k);
220 printf("------------------------------------------------------\n");
221 }
222 }
223 }
224
test_spg_standardize_cell_corundum(void)225 static void test_spg_standardize_cell_corundum(void)
226 {
227 double lattice[3][3] = {{4.8076344022756095, -2.4038172011378047, 0},
228 {0, 4.1635335244786962, 0},
229 {0, 0, 13.1172699198127543}};
230 double position[][3] = {
231 {0.0000000000000000, 0.0000000000000000, 0.3521850942289043},
232 {0.6666666666666643, 0.3333333333333357, 0.6855184275622400},
233 {0.3333333333333357, 0.6666666666666643, 0.0188517608955686},
234 {0.0000000000000000, 0.0000000000000000, 0.6478149057711028},
235 {0.6666666666666643, 0.3333333333333357, 0.9811482391044314},
236 {0.3333333333333357, 0.6666666666666643, 0.3144815724377600},
237 {0.0000000000000000, 0.0000000000000000, 0.1478149057710957},
238 {0.6666666666666643, 0.3333333333333357, 0.4811482391044314},
239 {0.3333333333333357, 0.6666666666666643, 0.8144815724377600},
240 {0.0000000000000000, 0.0000000000000000, 0.8521850942288972},
241 {0.6666666666666643, 0.3333333333333357, 0.1855184275622400},
242 {0.3333333333333357, 0.6666666666666643, 0.5188517608955686},
243 {0.3061673906454899, 0.0000000000000000, 0.2500000000000000},
244 {0.9728340573121541, 0.3333333333333357, 0.5833333333333357},
245 {0.6395007239788255, 0.6666666666666643, 0.9166666666666643},
246 {0.6938326093545102, 0.0000000000000000, 0.7500000000000000},
247 {0.3604992760211744, 0.3333333333333357, 0.0833333333333357},
248 {0.0271659426878458, 0.6666666666666643, 0.4166666666666643},
249 {0.0000000000000000, 0.3061673906454899, 0.2500000000000000},
250 {0.6666666666666643, 0.6395007239788255, 0.5833333333333357},
251 {0.3333333333333357, 0.9728340573121541, 0.9166666666666643},
252 {0.0000000000000000, 0.6938326093545102, 0.7500000000000000},
253 {0.6666666666666643, 0.0271659426878458, 0.0833333333333357},
254 {0.3333333333333357, 0.3604992760211744, 0.4166666666666643},
255 {0.6938326093545102, 0.6938326093545102, 0.2500000000000000},
256 {0.3604992760211744, 0.0271659426878458, 0.5833333333333357},
257 {0.0271659426878458, 0.3604992760211744, 0.9166666666666643},
258 {0.3061673906454899, 0.3061673906454899, 0.7500000000000000},
259 {0.9728340573121541, 0.6395007239788255, 0.0833333333333357},
260 {0.6395007239788255, 0.9728340573121541, 0.4166666666666643},
261 };
262 int types[30];
263 int i, j, k, num_atom = 30;
264 double symprec = 1e-5;
265
266 for (i = 0; i < 12; i++) {
267 types[i] = 1;
268 }
269 for (i = 12; i < 30; i++) {
270 types[i] = 2;
271 }
272
273 /* lattice, position, and types are overwirtten. */
274 printf("*** Example of spg_standardize_cell (Corundum) ***:\n");
275 printf("------------------------------------------------------\n");
276 for (j = 0; j < 2; j++) {
277 for (k = 0; k < 2; k++) {
278 sub_spg_standardize_cell(lattice,
279 position,
280 types,
281 num_atom,
282 symprec,
283 j,
284 k);
285 printf("------------------------------------------------------\n");
286 }
287 }
288 }
289
sub_spg_standardize_cell(double lattice[3][3],double position[][3],int types[],const int num_atom,const double symprec,const int to_primitive,const int no_idealize)290 static int sub_spg_standardize_cell(double lattice[3][3],
291 double position[][3],
292 int types[],
293 const int num_atom,
294 const double symprec,
295 const int to_primitive,
296 const int no_idealize)
297 {
298 int i, num_primitive_atom;
299 double lat[3][3], pos[num_atom][3];
300 int typ[num_atom];
301
302 for (i = 0; i < 3; i++) {
303 lat[i][0] = lattice[i][0];
304 lat[i][1] = lattice[i][1];
305 lat[i][2] = lattice[i][2];
306 }
307
308 for (i = 0; i < num_atom; i++) {
309 pos[i][0] = position[i][0];
310 pos[i][1] = position[i][1];
311 pos[i][2] = position[i][2];
312 typ[i] = types[i];
313 }
314
315 /* lattice, position, and types are overwirtten. */
316 num_primitive_atom = spg_standardize_cell(lat,
317 pos,
318 typ,
319 num_atom,
320 to_primitive,
321 no_idealize,
322 symprec);
323 printf("VASP POSCAR format: ");
324 if (to_primitive == 0) {
325 printf("to_primitive=0 and ");
326 } else {
327 printf("to_primitive=1 and ");
328 }
329
330 if (no_idealize == 0) {
331 printf("no_idealize=0\n");
332 } else {
333 printf("no_idealize=1\n");
334 }
335 printf("1.0\n");
336 for (i = 0; i < 3; i++) {
337 printf("%f %f %f\n", lat[0][i], lat[1][i], lat[2][i]);
338 }
339 printf("%d\n", num_primitive_atom);
340 printf("Direct\n");
341 for (i = 0; i < num_primitive_atom; i++) {
342 printf("%f %f %f\n", pos[i][0], pos[i][1], pos[i][2]);
343 }
344 }
345
test_spg_get_international(void)346 static void test_spg_get_international(void)
347 {
348 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 3}};
349 double position[][3] =
350 {
351 {0, 0, 0},
352 {0.5, 0.5, 0.5},
353 {0.3, 0.3, 0},
354 {0.7, 0.7, 0},
355 {0.2, 0.8, 0.5},
356 {0.8, 0.2, 0.5},
357 };
358 int types[] = {1, 1, 2, 2, 2, 2};
359 int num_spg, num_atom = 6;
360 char symbol[21];
361
362 num_spg = spg_get_international(symbol, lattice, position, types, num_atom, 1e-5);
363 printf("*** Example of spg_get_international ***:\n");
364 if ( num_spg > 0 ) {
365 printf("%s (%d)\n", symbol, num_spg);
366 } else {
367 printf("Space group could not be found.\n");
368 }
369 }
370
test_spg_get_schoenflies(void)371 static void test_spg_get_schoenflies(void)
372 {
373 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 3}};
374 double position[][3] =
375 {
376 {0, 0, 0},
377 {0.5, 0.5, 0.5},
378 {0.3, 0.3, 0},
379 {0.7, 0.7, 0},
380 {0.2, 0.8, 0.5},
381 {0.8, 0.2, 0.5},
382 };
383 int types[] = {1, 1, 2, 2, 2, 2};
384 int num_atom = 6;
385 char symbol[7];
386
387 spg_get_schoenflies(symbol, lattice, position, types, num_atom, 1e-5);
388 printf("*** Example of spg_get_schoenflies ***:\n");
389 printf("Schoenflies: %s\n", symbol);
390 }
391
test_spg_get_spacegroup_type(void)392 static void test_spg_get_spacegroup_type(void)
393 {
394 SpglibSpacegroupType spgtype;
395 spgtype = spg_get_spacegroup_type(446);
396
397 printf("*** Example of spg_get_spacegroup_type ***:\n");
398 printf("Number: %d\n", spgtype.number);
399 printf("Schoenflies: %s\n", spgtype.schoenflies);
400 printf("International: %s\n", spgtype.international);
401 printf("International: %s\n", spgtype.international_full);
402 printf("International: %s\n", spgtype.international_short);
403 printf("Hall symbol: %s\n", spgtype.hall_symbol);
404
405 }
406
test_spg_get_symmetry_from_database(void)407 static void test_spg_get_symmetry_from_database(void)
408 {
409 int rotations[192][3][3];
410 double translations[192][3];
411 int max_size, i, j, size;
412
413 size = spg_get_symmetry_from_database(rotations,
414 translations,
415 460);
416
417 printf("*** Example of spg_get_symmetry_from_database ***:\n");
418 for (i = 0; i < size; i++) {
419 printf("--- %d ---\n", i + 1);
420 for (j = 0; j < 3; j++) {
421 printf("%2d %2d %2d\n",
422 rotations[i][j][0], rotations[i][j][1], rotations[i][j][2]);
423 }
424 printf("%f %f %f\n",
425 translations[i][0], translations[i][1], translations[i][2]);
426 }
427 }
428
test_spg_get_multiplicity(void)429 static void test_spg_get_multiplicity(void)
430 {
431 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 4}};
432 double position[][3] = {
433 {0, 0, 0},
434 {0.5, 0.5, 0.5}
435 };
436 int types[] = {1, 2};
437 int num_atom = 2;
438 int size;
439
440 size = spg_get_multiplicity(lattice, position, types, num_atom, 1e-5);
441
442 printf("*** Example of spg_get_multiplicity ***:\n");
443 printf("Number of symmetry operations: %d\n", size);
444 }
445
446
test_spg_get_symmetry(void)447 static void test_spg_get_symmetry(void)
448 {
449 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 3}};
450 double position[][3] =
451 {
452 {0, 0, 0},
453 {0.5, 0.5, 0.25},
454 {0.3, 0.3, 0},
455 {0.7, 0.7, 0},
456 {0.2, 0.8, 0.25},
457 {0.8, 0.2, 0.25},
458 {0, 0, 0.5},
459 {0.5, 0.5, 0.75},
460 {0.3, 0.3, 0.5},
461 {0.7, 0.7, 0.5},
462 {0.2, 0.8, 0.75},
463 {0.8, 0.2, 0.75}
464 };
465 int types[] = {1, 1, 2, 2, 2, 2, 1, 1, 2, 2, 2, 2};
466 int num_atom = 12;
467 int max_size = 50;
468 int i, j, size;
469 int rotation[max_size][3][3];
470 double translation[max_size][3];
471
472 double origin_shift[3] = {0.1, 0.1, 0};
473 for (i = 0; i < num_atom; i++) {
474 for (j = 0; j < 3; j++) {
475 position[i][j] += origin_shift[j];
476 }
477 }
478
479 printf("*** Example of spg_get_symmetry (Rutile two unit cells) ***:\n");
480 size = spg_get_symmetry(rotation,
481 translation,
482 max_size,
483 lattice,
484 position,
485 types,
486 num_atom,
487 1e-5);
488 for (i = 0; i < size; i++) {
489 printf("--- %d ---\n", i + 1);
490 for (j = 0; j < 3; j++)
491 printf("%2d %2d %2d\n", rotation[i][j][0], rotation[i][j][1], rotation[i][j][2]);
492 printf("%f %f %f\n",
493 translation[i][0], translation[i][1], translation[i][2]);
494 }
495
496 }
497
test_spg_get_symmetry_with_collinear_spin(void)498 static void test_spg_get_symmetry_with_collinear_spin(void) {
499 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 4}};
500 double position[][3] =
501 {
502 {0, 0, 0},
503 {0.5, 0.5, 0.5}
504 };
505 int types[] = {1, 1};
506 int equivalent_atoms[2];
507 double spins[2];
508 int num_atom = 2;
509 int max_size = 300;
510 int i, j, size;
511 int rotation[max_size][3][3];
512 double translation[max_size][3];
513
514 printf("*** Example of spg_get_symmetry_with_spin (BCC ferro) ***:\n");
515 spins[0] = 1;
516 spins[1] = 1;
517 size = spg_get_symmetry_with_collinear_spin(rotation,
518 translation,
519 equivalent_atoms,
520 max_size,
521 lattice,
522 position,
523 types,
524 spins,
525 num_atom,
526 1e-5);
527 for (i = 0; i < size; i++) {
528 printf("--- %d ---\n", i + 1);
529 for (j = 0; j < 3; j++) {
530 printf("%2d %2d %2d\n",
531 rotation[i][j][0], rotation[i][j][1], rotation[i][j][2]);
532 }
533 printf("%f %f %f\n", translation[i][0], translation[i][1],
534 translation[i][2]);
535 }
536
537 printf("*** Example of spg_get_symmetry_with_spin (BCC antiferro) ***:\n");
538 spins[0] = 1;
539 spins[1] = -1;
540 size = spg_get_symmetry_with_collinear_spin(rotation,
541 translation,
542 equivalent_atoms,
543 max_size,
544 lattice,
545 position,
546 types,
547 spins,
548 num_atom,
549 1e-5);
550 for (i = 0; i < size; i++) {
551 printf("--- %d ---\n", i + 1);
552 for (j = 0; j < 3; j++) {
553 printf("%2d %2d %2d\n",
554 rotation[i][j][0], rotation[i][j][1], rotation[i][j][2]);
555 }
556 printf("%f %f %f\n", translation[i][0], translation[i][1],
557 translation[i][2]);
558 }
559
560 printf("*** Example of spg_get_symmetry_with_spin (BCC broken spin) ***:\n");
561 spins[0] = 1;
562 spins[1] = 2;
563 size = spg_get_symmetry_with_collinear_spin(rotation,
564 translation,
565 equivalent_atoms,
566 max_size,
567 lattice,
568 position,
569 types,
570 spins,
571 num_atom,
572 1e-5);
573 for (i = 0; i < size; i++) {
574 printf("--- %d ---\n", i + 1);
575 for (j = 0; j < 3; j++) {
576 printf("%2d %2d %2d\n",
577 rotation[i][j][0], rotation[i][j][1], rotation[i][j][2]);
578 }
579 printf("%f %f %f\n", translation[i][0], translation[i][1],
580 translation[i][2]);
581 }
582 }
583
584
test_spg_get_dataset(void)585 static void test_spg_get_dataset(void)
586 {
587 printf("*** Example of spg_get_dataset (Rutile two unit cells) ***:\n");
588 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 3}};
589 double origin_shift[3] = {0.1, 0.1, 0};
590 double position[][3] =
591 {
592 {0, 0, 0},
593 {0.5, 0.5, 0.25},
594 {0.3, 0.3, 0},
595 {0.7, 0.7, 0},
596 {0.2, 0.8, 0.25},
597 {0.8, 0.2, 0.25},
598 {0, 0, 0.5},
599 {0.5, 0.5, 0.75},
600 {0.3, 0.3, 0.5},
601 {0.7, 0.7, 0.5},
602 {0.2, 0.8, 0.75},
603 {0.8, 0.2, 0.75}
604 };
605 int types[] = {1, 1, 2, 2, 2, 2, 1, 1, 2, 2, 2, 2};
606 int num_atom = 12;
607
608 show_spg_dataset(lattice, origin_shift, position, num_atom, types);
609
610 double lattice_2[3][3] = {{3.7332982433264039, -1.8666491216632011, 0},
611 {0, 3.2331311186244847, 0},
612 {0, 0, 6.0979971306362799}};
613 double origin_shift_2[3] = {0.1, 0.1, 0};
614 double position_2[][3] =
615 {
616 {0, 0, 0},
617 {1.0 / 3, 2.0 / 3, 0.4126},
618 {1.0 / 3, 2.0 / 3, 0.776},
619 {2.0 / 3, 1.0 / 3, 0.2542},
620 };
621 int types_2[] = {1, 2, 3, 3};
622 int num_atom_2 = 4;
623
624 show_spg_dataset(lattice_2, origin_shift_2, position_2, num_atom_2, types_2);
625 }
626
test_spg_get_ir_reciprocal_mesh(void)627 static void test_spg_get_ir_reciprocal_mesh(void)
628 {
629 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 3}};
630 double position[][3] =
631 {
632 {0, 0, 0},
633 {0.5, 0.5, 0.5},
634 {0.3, 0.3, 0},
635 {0.7, 0.7, 0},
636 {0.2, 0.8, 0.5},
637 {0.8, 0.2, 0.5},
638 };
639 int types[] = {1, 1, 2, 2, 2, 2};
640 int num_atom = 6;
641 int m = 40;
642 int mesh[] = {m, m, m};
643 int is_shift[] = {1, 1, 1};
644 int grid_address[m * m * m][3];
645 int grid_mapping_table[m * m * m];
646
647 printf("*** Example of spg_get_ir_reciprocal_mesh of Rutile structure ***:\n");
648
649 int num_ir = spg_get_ir_reciprocal_mesh(grid_address,
650 grid_mapping_table,
651 mesh,
652 is_shift,
653 1,
654 lattice,
655 position,
656 types,
657 num_atom,
658 1e-5);
659
660 printf("Number of irreducible k-points of Rutile with\n");
661 printf("40x40x40 Monkhorst-Pack mesh is %d (4200).\n", num_ir);
662 }
663
test_spg_get_stabilized_reciprocal_mesh(void)664 static void test_spg_get_stabilized_reciprocal_mesh(void)
665 {
666 SpglibDataset *dataset;
667 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 3}};
668 double position[][3] =
669 {
670 {0, 0, 0},
671 {0.5, 0.5, 0.5},
672 {0.3, 0.3, 0},
673 {0.7, 0.7, 0},
674 {0.2, 0.8, 0.5},
675 {0.8, 0.2, 0.5},
676 };
677 int types[] = {1, 1, 2, 2, 2, 2};
678 int num_atom = 6;
679
680 dataset = spg_get_dataset(lattice,
681 position,
682 types,
683 num_atom,
684 1e-5);
685
686
687 int m = 40;
688 int mesh[] = {m, m, m};
689 int is_shift[] = {1, 1, 1};
690 int grid_address[m * m * m][3];
691 int grid_mapping_table[m * m * m];
692 double q[] = {0, 0.5, 0.5};
693
694 printf("*** Example of spg_get_stabilized_reciprocal_mesh of Rutile structure ***:\n");
695
696 int num_ir = spg_get_stabilized_reciprocal_mesh(grid_address,
697 grid_mapping_table,
698 mesh,
699 is_shift,
700 1,
701 dataset->n_operations,
702 dataset->rotations,
703 1,
704 (double(*)[3])q);
705 spg_free_dataset(dataset);
706
707 printf("Number of irreducible k-points stabilized by q=(0, 1/2, 1/2) of Rutile with\n");
708 printf("40x40x40 Monkhorst-Pack mesh is %d (8000).\n", num_ir);
709 }
710
test_spg_relocate_BZ_grid_address(void)711 static void test_spg_relocate_BZ_grid_address(void)
712 {
713 double rec_lattice[3][3] = {{-0.17573761, 0.17573761, 0.17573761},
714 { 0.17573761, -0.17573761, 0.17573761},
715 { 0.17573761, 0.17573761, -0.17573761}};
716 int rotations[][3][3] = {{{1, 0, 0},
717 {0, 1, 0},
718 {0, 0, 1}}};
719
720
721 int m = 40;
722 int mesh[] = {m, m, m};
723 int is_shift[] = {0, 0, 0};
724 int bz_grid_address[(m + 1) * (m + 1) * (m + 1)][3];
725 int bz_map[m * m * m * 8];
726 int grid_address[m * m * m][3];
727 int grid_mapping_table[m * m * m];
728 double q[] = {0, 0, 0};
729
730 int num_ir = spg_get_stabilized_reciprocal_mesh(grid_address,
731 grid_mapping_table,
732 mesh,
733 is_shift,
734 1,
735 1,
736 rotations,
737 1,
738 (double(*)[3])q);
739
740
741 printf("*** Example of spg_relocate_BZ_grid_address of NaCl structure ***:\n");
742
743 int num_q = spg_relocate_BZ_grid_address(bz_grid_address,
744 bz_map,
745 grid_address,
746 mesh,
747 rec_lattice,
748 is_shift);
749
750 printf("Number of k-points of NaCl Brillouin zone\n");
751 printf("with Gamma-centered 40x40x40 Monkhorst-Pack mesh is %d (65861).\n", num_q);
752 }
753
show_spg_dataset(double lattice[3][3],const double origin_shift[3],double position[][3],const int num_atom,const int types[])754 static void show_spg_dataset(double lattice[3][3],
755 const double origin_shift[3],
756 double position[][3],
757 const int num_atom,
758 const int types[])
759 {
760 SpglibDataset *dataset;
761 char ptsymbol[6];
762 int pt_trans_mat[3][3];
763
764 int i, j, size;
765 const char *wl = "abcdefghijklmnopqrstuvwxyz";
766
767 for ( i = 0; i < num_atom; i++ ) {
768 for ( j = 0; j < 3; j++ ) {
769 position[i][j] += origin_shift[j];
770 }
771 }
772
773 dataset = spg_get_dataset(lattice,
774 position,
775 types,
776 num_atom,
777 1e-5);
778
779 printf("International: %s (%d)\n", dataset->international_symbol, dataset->spacegroup_number );
780 printf("Hall symbol: %s\n", dataset->hall_symbol );
781 spg_get_pointgroup(ptsymbol,
782 pt_trans_mat,
783 dataset->rotations,
784 dataset->n_operations);
785 printf("Point group: %s\n", ptsymbol);
786 printf("Transformation matrix:\n");
787 for ( i = 0; i < 3; i++ ) {
788 printf("%f %f %f\n",
789 dataset->transformation_matrix[i][0],
790 dataset->transformation_matrix[i][1],
791 dataset->transformation_matrix[i][2]);
792 }
793 printf("Wyckoff letters:\n");
794 for ( i = 0; i < dataset->n_atoms; i++ ) {
795 printf("%c ", wl[dataset->wyckoffs[i]]);
796 }
797 printf("\n");
798 printf("Equivalent atoms:\n");
799 for (i = 0; i < dataset->n_atoms; i++) {
800 printf("%d ", dataset->equivalent_atoms[i]);
801 }
802 printf("\n");
803
804 for (i = 0; i < dataset->n_operations; i++) {
805 printf("--- %d ---\n", i + 1);
806 for (j = 0; j < 3; j++) {
807 printf("%2d %2d %2d\n",
808 dataset->rotations[i][j][0],
809 dataset->rotations[i][j][1],
810 dataset->rotations[i][j][2]);
811 }
812 printf("%f %f %f\n",
813 dataset->translations[i][0],
814 dataset->translations[i][1],
815 dataset->translations[i][2]);
816 }
817
818 spg_free_dataset(dataset);
819
820 }
821
show_cell(double lattice[3][3],double position[][3],const int types[],const int num_atom)822 static void show_cell(double lattice[3][3],
823 double position[][3],
824 const int types[],
825 const int num_atom)
826 {
827 int i;
828
829 printf("Lattice parameter:\n");
830 for (i = 0; i < 3; i++) {
831 printf("%f %f %f\n", lattice[0][i], lattice[1][i], lattice[2][i]);
832 }
833 printf("Atomic positions:\n");
834 for (i = 0; i < num_atom; i++) {
835 printf("%d: %f %f %f\n",
836 types[i], position[i][0], position[i][1], position[i][2]);
837 }
838 }
839
test_spg_get_error_message(void)840 static void test_spg_get_error_message(void)
841 {
842 double lattice[3][3] = {{4, 0, 0}, {0, 4, 0}, {0, 0, 4}};
843 double position[][3] = {
844 {0, 0, 0},
845 {0.5, 0.5, 0.5},
846 {0.5, 0.5, 0.5}
847 };
848 int types[] = {1, 1, 1};
849 int i, num_atom = 3, num_primitive_atom;
850 double symprec = 1e-5;
851 SpglibError error;
852
853
854 /* lattice, position, and types are overwirtten. */
855 printf("*** Example of spg_get_error_message ***:\n");
856 num_primitive_atom = spg_find_primitive(lattice, position, types, num_atom, symprec);
857 if (num_primitive_atom == 0) {
858 printf("Primitive cell was not found.\n");
859 error = spg_get_error_code();
860 printf("%s\n", spg_get_error_message(error));
861 }
862 }
863