1 /* This code is part of the tng binary trajectory format.
2  *
3  * Written by Magnus Lundborg
4  * Copyright (c) 2012-2014, The GROMACS development team.
5  * Check out http://www.gromacs.org for more information.
6  *
7  *
8  * This program is free software; you can redistribute it and/or
9  * modify it under the terms of the Revised BSD License.
10  */
11 
12 #include "tng/tng_io.h"
13 
14 #ifdef USE_STD_INTTYPES_H
15 #include <inttypes.h>
16 #endif
17 
18 #include <stdlib.h>
19 #include <string.h>
20 #include <math.h>
21 #include "tng/version.h"
22 
23 #define BOX_SHAPE_X 150.0
24 #define BOX_SHAPE_Y 145.5
25 #define BOX_SHAPE_Z 155.5
26 #define N_FRAME_SETS 100
27 #define MEDIUM_STRIDE_LEN 5
28 #define LONG_STRIDE_LEN 25
29 #define TIME_PER_FRAME 2e-15
30 #define COMPRESSION_PRECISION 1000
31 #define USER_NAME "USER 1"
32 #define PROGRAM_NAME "tng_testing"
33 #define COMPUTER_NAME "Unknown computer"
34 #define FORCEFIELD_NAME "No forcefield"
35 
tng_test_setup_molecules(tng_trajectory_t traj)36 static tng_function_status tng_test_setup_molecules(tng_trajectory_t traj)
37 {
38     tng_molecule_t molecule;
39     tng_chain_t chain;
40     tng_residue_t residue;
41     tng_atom_t atom;
42     tng_bond_t bond;
43     int64_t cnt;
44 
45     tng_molecule_add(traj, "water", &molecule);
46     tng_molecule_chain_add(traj, molecule, "W", &chain);
47     tng_chain_residue_add(traj, chain, "WAT", &residue);
48     if(tng_residue_atom_add(traj, residue, "O", "O", &atom) == TNG_CRITICAL)
49     {
50         return(TNG_CRITICAL);
51     }
52     if(tng_residue_atom_add(traj, residue, "HO1", "H", &atom) == TNG_CRITICAL)
53     {
54         return(TNG_CRITICAL);
55     }
56     if(tng_residue_atom_add(traj, residue, "HO2", "H", &atom) == TNG_CRITICAL)
57     {
58         return(TNG_CRITICAL);
59     }
60     tng_molecule_bond_add(traj, molecule, 0, 1, &bond);
61     tng_molecule_bond_add(traj, molecule, 0, 2, &bond);
62     tng_molecule_cnt_set(traj, molecule, 200);
63     tng_molecule_cnt_get(traj, molecule, &cnt);
64 
65     if(cnt != 200)
66     {
67         return(TNG_CRITICAL);
68     }
69 /*     printf("Created %"PRId64" %s molecules.\n", cnt, molecule->name); */
70 
71 /*     traj->molecule_cnt_list[traj->n_molecules-1] = 5;
72 //     tng_molecule_name_set(traj, &traj->molecules[1], "ligand");
73 //     tng_molecule_name_set(traj, &traj->molecules[2], "water");
74 //     tng_molecule_name_set(traj, &traj->molecules[3], "dummy");
75 //     traj->molecules[0].id = 0;
76 //     traj->molecules[1].id = 1;
77 //     traj->molecules[2].id = 2;
78 //     traj->molecules[3].id = 3;
79 
80 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom1", "type1") == TNG_CRITICAL)
81 //     {
82 //         return(TNG_CRITICAL);
83 //     }
84 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom2", "type1") == TNG_CRITICAL)
85 //     {
86 //         return(TNG_CRITICAL);
87 //     }
88 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom3", "type1") == TNG_CRITICAL)
89 //     {
90 //         return(TNG_CRITICAL);
91 //     }
92 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom4", "type2") == TNG_CRITICAL)
93 //     {
94 //         return(TNG_CRITICAL);
95 //     }
96 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom5", "type2") == TNG_CRITICAL)
97 //     {
98 //         return(TNG_CRITICAL);
99 //     }
100 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom6", "type2") == TNG_CRITICAL)
101 //     {
102 //         return(TNG_CRITICAL);
103 //     }
104 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[0], "atom7", "type3") == TNG_CRITICAL)
105 //     {
106 //         return(TNG_CRITICAL);
107 //     }
108 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[1], "C1", "C") == TNG_CRITICAL)
109 //     {
110 //         return(TNG_CRITICAL);
111 //     }
112 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[1], "O1", "O") == TNG_CRITICAL)
113 //     {
114 //         return(TNG_CRITICAL);
115 //     }
116 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[1], "H11", "H") == TNG_CRITICAL)
117 //     {
118 //         return(TNG_CRITICAL);
119 //     }
120 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[1], "H12", "H") == TNG_CRITICAL)
121 //     {
122 //         return(TNG_CRITICAL);
123 //     }
124 //     if(tng_add_atom_to_molecule(traj, &traj->molecules[1], "H13", "H") == TNG_CRITICAL)
125 //     {
126 //         return(TNG_CRITICAL);
127 //     }
128 */
129     return(TNG_SUCCESS);
130 }
131 
tng_test_molecules(tng_trajectory_t traj)132 static tng_function_status tng_test_molecules(tng_trajectory_t traj)
133 {
134     tng_molecule_t molecule, molecule_new;
135     tng_chain_t chain;
136     tng_residue_t residue;
137     tng_atom_t atom;
138     int64_t cnt, *bonds_to, *bonds_from;
139     char var_atoms, str[TNG_MAX_STR_LEN];
140     tng_function_status stat;
141 
142     stat = tng_num_molecule_types_get(traj, &cnt);
143     if(stat != TNG_SUCCESS || cnt != 1)
144     {
145         printf("Molecule reading error. %s: %d\n",
146                __FILE__, __LINE__);
147         return(TNG_FAILURE);
148     }
149 
150     stat = tng_num_molecules_get(traj, &cnt);
151     if(stat != TNG_SUCCESS || cnt != 200)
152     {
153         printf("Molecule reading error. %s: %d\n",
154                __FILE__, __LINE__);
155         return(TNG_FAILURE);
156     }
157 
158     stat = tng_num_particles_variable_get(traj, &var_atoms);
159     if(stat != TNG_SUCCESS || var_atoms)
160     {
161         printf("Molecule reading error. %s: %d\n",
162                __FILE__, __LINE__);
163         return(TNG_FAILURE);
164     }
165 
166     stat = tng_molecule_of_index_get(traj, 0, &molecule);
167     if(stat != TNG_SUCCESS)
168     {
169         printf("Cannot find molecule. %s: %d\n",
170                __FILE__, __LINE__);
171         return(stat);
172     }
173     stat = tng_molecule_find(traj, "water", -1, &molecule);
174     if(stat != TNG_SUCCESS)
175     {
176         printf("Cannot find molecule. %s: %d\n",
177                __FILE__, __LINE__);
178         return(stat);
179     }
180     stat =tng_molecule_name_get(traj, molecule, str, TNG_MAX_STR_LEN);
181     if(stat != TNG_SUCCESS)
182     {
183         printf("Cannot get molecule name. %s: %d\n",
184                __FILE__, __LINE__);
185         return(stat);
186     }
187     stat = tng_molecule_num_chains_get(traj, molecule, &cnt);
188     if(stat != TNG_SUCCESS || cnt != 1)
189     {
190         printf("Cannot get number of chains in molecule. %s: %d\n",
191                __FILE__, __LINE__);
192         return(TNG_FAILURE);
193     }
194     stat = tng_molecule_chain_of_index_get(traj, molecule, 0, &chain);
195     if(stat != TNG_SUCCESS)
196     {
197         printf("Cannot get chain in molecule. %s: %d\n",
198                __FILE__, __LINE__);
199         return(stat);
200     }
201     stat = tng_molecule_chain_find(traj, molecule, "W", -1, &chain);
202     if(stat != TNG_SUCCESS)
203     {
204         printf("Cannot get chain in molecule. %s: %d\n",
205                __FILE__, __LINE__);
206         return(stat);
207     }
208     stat = tng_molecule_num_residues_get(traj, molecule, &cnt);
209     if(stat != TNG_SUCCESS || cnt != 1)
210     {
211         printf("Cannot get number of residues in molecule. %s: %d\n",
212                __FILE__, __LINE__);
213         return(TNG_FAILURE);
214     }
215     stat = tng_molecule_residue_of_index_get(traj, molecule, 0, &residue);
216     if(stat != TNG_SUCCESS)
217     {
218         printf("Cannot get residue in molecule. %s: %d\n",
219                __FILE__, __LINE__);
220         return(stat);
221     }
222     stat = tng_molecule_num_atoms_get(traj, molecule, &cnt);
223     if(stat != TNG_SUCCESS || cnt != 3)
224     {
225         printf("Cannot get number of atoms in molecule. %s: %d\n",
226                __FILE__, __LINE__);
227         return(TNG_FAILURE);
228     }
229     stat = tng_molecule_atom_of_index_get(traj, molecule, 0, &atom);
230     if(stat != TNG_SUCCESS)
231     {
232         printf("Cannot get atom in molecule. %s: %d\n",
233                __FILE__, __LINE__);
234         return(stat);
235     }
236     stat = tng_molecule_atom_find(traj, molecule, "O", -1, &atom);
237     if(stat != TNG_SUCCESS)
238     {
239         printf("Cannot get atom in molecule. %s: %d\n",
240                __FILE__, __LINE__);
241         return(stat);
242     }
243 
244     stat =tng_chain_name_get(traj, chain, str, TNG_MAX_STR_LEN);
245     if(stat != TNG_SUCCESS)
246     {
247         printf("Cannot get name of chain. %s: %d\n",
248                __FILE__, __LINE__);
249         return(stat);
250     }
251     stat = tng_chain_num_residues_get(traj, chain, &cnt);
252     if(stat != TNG_SUCCESS || cnt != 1)
253     {
254         printf("Cannot get number of residues in chain. %s: %d\n",
255                __FILE__, __LINE__);
256         return(TNG_FAILURE);
257     }
258     stat = tng_chain_residue_of_index_get(traj, chain, 0, &residue);
259     if(stat != TNG_SUCCESS)
260     {
261         printf("Cannot get residue in chain. %s: %d\n",
262                __FILE__, __LINE__);
263         return(stat);
264     }
265     stat = tng_chain_residue_find(traj, chain, "WAT", -1, &residue);
266     if(stat != TNG_SUCCESS)
267     {
268         printf("Cannot get residue in chain. %s: %d\n",
269                __FILE__, __LINE__);
270         return(stat);
271     }
272 
273     stat = tng_residue_name_get(traj, residue, str, TNG_MAX_STR_LEN);
274     if(stat != TNG_SUCCESS)
275     {
276         printf("Cannot get name of residue. %s: %d\n",
277                __FILE__, __LINE__);
278         return(stat);
279     }
280     stat = tng_residue_num_atoms_get(traj, residue, &cnt);
281     if(stat != TNG_SUCCESS || cnt != 3)
282     {
283         printf("Cannot get number of atoms in residue. %s: %d\n",
284                __FILE__, __LINE__);
285         return(TNG_FAILURE);
286     }
287     stat = tng_residue_atom_of_index_get(traj, residue, 0, &atom);
288     if(stat != TNG_SUCCESS)
289     {
290         printf("Cannot get residue of atom. %s: %d\n",
291                __FILE__, __LINE__);
292         return(stat);
293     }
294 
295     stat = tng_atom_name_get(traj, atom, str, TNG_MAX_STR_LEN);
296     if(stat != TNG_SUCCESS)
297     {
298         printf("Cannot get name of atom. %s: %d\n",
299                __FILE__, __LINE__);
300         return(stat);
301     }
302     stat = tng_atom_type_get(traj, atom, str, TNG_MAX_STR_LEN);
303     if(stat != TNG_SUCCESS)
304     {
305         printf("Cannot get atom type of atom. %s: %d\n",
306                __FILE__, __LINE__);
307         return(stat);
308     }
309 
310     stat = tng_molecule_id_of_particle_nr_get(traj, 0, &cnt);
311     if(stat != TNG_SUCCESS || cnt != 1)
312     {
313         printf("Cannot get molecule id of atom. %s: %d\n",
314                __FILE__, __LINE__);
315         return(TNG_FAILURE);
316     }
317     stat = tng_residue_id_of_particle_nr_get(traj, 0, &cnt);
318     if(stat != TNG_SUCCESS || cnt != 0)
319     {
320         printf("Cannot get residue id of atom. %s: %d\n",
321                __FILE__, __LINE__);
322         return(TNG_FAILURE);
323     }
324     stat = tng_global_residue_id_of_particle_nr_get(traj, 599, &cnt);
325     if(stat != TNG_SUCCESS || cnt != 199)
326     {
327         printf("Cannot get global residue id of atom. %s: %d\n",
328                __FILE__, __LINE__);
329         return(TNG_FAILURE);
330     }
331     stat = tng_molecule_name_of_particle_nr_get(traj, 0, str, TNG_MAX_STR_LEN);
332     if(stat != TNG_SUCCESS)
333     {
334         printf("Cannot get molecule name of atom. %s: %d\n",
335                __FILE__, __LINE__);
336         return(stat);
337     }
338     stat = tng_chain_name_of_particle_nr_get(traj, 0, str, TNG_MAX_STR_LEN);
339     if(stat != TNG_SUCCESS)
340     {
341         printf("Cannot get chain name of atom. %s: %d\n",
342                __FILE__, __LINE__);
343         return(stat);
344     }
345     stat = tng_residue_name_of_particle_nr_get(traj, 0, str, TNG_MAX_STR_LEN);
346     if(stat != TNG_SUCCESS)
347     {
348         printf("Cannot get residue name of atom. %s: %d\n",
349                __FILE__, __LINE__);
350         return(stat);
351     }
352     stat = tng_atom_name_of_particle_nr_get(traj, 0, str, TNG_MAX_STR_LEN);
353     if(stat != TNG_SUCCESS)
354     {
355         printf("Cannot get atom name of atom. %s: %d\n",
356                __FILE__, __LINE__);
357         return(stat);
358     }
359 
360     stat = tng_molecule_alloc(traj, &molecule_new);
361     if(stat != TNG_SUCCESS)
362     {
363         printf("Cannot setup new molecule. %s: %d\n",
364                __FILE__, __LINE__);
365         return(stat);
366     }
367     stat = tng_molecule_name_set(traj, molecule_new, "TEST");
368     if(stat != TNG_SUCCESS)
369     {
370         printf("Cannot set name of new molecule. %s: %d\n",
371                __FILE__, __LINE__);
372         return(stat);
373     }
374     stat = tng_molecule_existing_add(traj, &molecule_new);
375     if(stat != TNG_SUCCESS)
376     {
377         printf("Cannot add new molecule to molecule system. %s: %d\n",
378                __FILE__, __LINE__);
379         return(stat);
380     }
381 
382     stat = tng_molsystem_bonds_get(traj, &cnt, &bonds_from, &bonds_to);
383     if(stat != TNG_SUCCESS || cnt != 400)
384     {
385         printf("Cannot get bonds in molecule system. %s: %d\n",
386                __FILE__, __LINE__);
387         return(stat);
388     }
389 
390     free(bonds_from);
391     free(bonds_to);
392 
393     return(TNG_SUCCESS);
394 }
395 
tng_test_read_and_write_file(tng_trajectory_t traj,const char hash_mode)396 static tng_function_status tng_test_read_and_write_file
397                 (tng_trajectory_t traj, const char hash_mode)
398 {
399     char file_name[TNG_MAX_STR_LEN];
400     tng_function_status stat;
401 
402     stat = tng_input_file_get(traj, file_name, TNG_MAX_STR_LEN);
403     if(stat != TNG_SUCCESS)
404     {
405         printf("Could not get name of input file. %s: %d\n",
406                __FILE__, __LINE__);
407         return(stat);
408     }
409     stat = tng_output_file_get(traj, file_name, TNG_MAX_STR_LEN);
410     if(stat != TNG_SUCCESS)
411     {
412         printf("Could not get name of output file. %s: %d\n",
413                __FILE__, __LINE__);
414         return(stat);
415     }
416 
417     stat = tng_file_headers_read(traj, hash_mode);
418     if(stat != TNG_SUCCESS)
419     {
420         printf("Could not read headers. %s: %d\n",
421                __FILE__, __LINE__);
422         return(stat);
423     }
424     stat = tng_file_headers_write(traj, hash_mode);
425     if(stat != TNG_SUCCESS)
426     {
427         printf("Could not write headers. %s: %d\n",
428                __FILE__, __LINE__);
429         return(stat);
430     }
431 
432     while(stat == TNG_SUCCESS)
433     {
434         stat = tng_frame_set_read_next(traj, hash_mode);
435         if(stat == TNG_CRITICAL)
436         {
437             printf("Could not read frame set. %s: %d\n",
438                    __FILE__, __LINE__);
439             return(stat);
440         }
441         if(stat == TNG_FAILURE)
442         {
443             return(TNG_SUCCESS);
444         }
445         stat = tng_frame_set_write(traj, hash_mode);
446     }
447 
448     return(stat);
449 }
450 
tng_test_write_and_read_traj(tng_trajectory_t * traj,const char hash_mode)451 static tng_function_status tng_test_write_and_read_traj(tng_trajectory_t *traj,
452                                                         const char hash_mode)
453 {
454     int i, j, k, nr, cnt, dependency;
455     float *data, *molpos, *charges, *masses;
456     int64_t mapping[300], n_particles;
457     int64_t read_n_frames, read_stride_length, read_n_particles, read_n_values_per_frame;
458     int64_t n_frames_per_frame_set, tot_n_mols;
459     int64_t codec_id;
460     int64_t dist_exp = -9, temp_int, temp_int2;
461 //     int64_t frame_nr;
462     double box_shape[9], temp_double;
463     char atom_type[16], annotation[128];
464     char temp_str[TNG_MAX_STR_LEN];
465     char read_data_type;
466     tng_trajectory_frame_set_t frame_set;
467     tng_function_status stat = TNG_SUCCESS;
468 
469     tng_medium_stride_length_set(*traj, MEDIUM_STRIDE_LEN);
470     tng_long_stride_length_set(*traj, LONG_STRIDE_LEN);
471 
472     tng_first_user_name_set(*traj, USER_NAME);
473     tng_first_program_name_set(*traj, PROGRAM_NAME);
474     tng_first_computer_name_set(*traj, COMPUTER_NAME);
475     tng_forcefield_name_set(*traj, FORCEFIELD_NAME);
476 
477     tng_compression_precision_set(*traj, COMPRESSION_PRECISION);
478 
479     tng_distance_unit_exponential_set(*traj, dist_exp);
480 
481     tng_time_per_frame_set(*traj, TIME_PER_FRAME);
482 
483     /* Create molecules */
484     if(tng_test_setup_molecules(*traj) == TNG_CRITICAL)
485     {
486         return(TNG_CRITICAL);
487     }
488 
489     /* Set the box shape */
490     box_shape[1] = box_shape[2] = box_shape[3] = box_shape[5] = box_shape[6] =
491     box_shape[7] = 0;
492     box_shape[0] = BOX_SHAPE_X;
493     box_shape[4] = BOX_SHAPE_Y;
494     box_shape[8] = BOX_SHAPE_Z;
495     if(tng_data_block_add(*traj, TNG_TRAJ_BOX_SHAPE, "BOX SHAPE", TNG_DOUBLE_DATA,
496                        TNG_NON_TRAJECTORY_BLOCK, 1, 9, 1, TNG_UNCOMPRESSED,
497                        box_shape) == TNG_CRITICAL)
498     {
499         tng_trajectory_destroy(traj);
500         printf("Cannot write trajectory box shape.\n");
501         exit(1);
502     }
503 
504     /* Set partial charges (treat the water as TIP3P). */
505     tng_num_particles_get(*traj, &n_particles);
506     charges = malloc(sizeof(float) * n_particles);
507     for(i = 0; i < n_particles; i++)
508     {
509         stat = tng_atom_type_of_particle_nr_get(*traj, i, atom_type,
510                                                 sizeof(atom_type));
511         if(stat == TNG_CRITICAL)
512         {
513             break;
514         }
515         /* We only have water in the system. If the atom is oxygen set its
516          * partial charge to -0.834, if it's a hydrogen set its partial charge to
517          * 0.417. */
518         switch(atom_type[0])
519         {
520             case 'O':
521                 charges[i] = -0.834;
522                 break;
523             case 'H':
524                 charges[i] = 0.417;
525                 break;
526             default:
527                 printf("Failed setting partial charges. %s: %d\n",
528                        __FILE__, __LINE__);
529                 return(TNG_CRITICAL);
530         }
531     }
532     if(stat == TNG_CRITICAL)
533     {
534         free(charges);
535         printf("Failed setting partial charges. %s: %d\n",
536                __FILE__, __LINE__);
537         return(TNG_CRITICAL);
538     }
539 
540     stat = tng_particle_data_block_add(*traj, TNG_TRAJ_PARTIAL_CHARGES, "PARTIAL CHARGES",
541                                        TNG_FLOAT_DATA, TNG_NON_TRAJECTORY_BLOCK,
542                                        1, 1, 1, 0, n_particles,
543                                        TNG_UNCOMPRESSED, charges);
544     free(charges);
545     charges = 0;
546     if(stat != TNG_SUCCESS)
547     {
548         printf("Failed adding partial charges. %s: %d\n",
549                __FILE__, __LINE__);
550         return(TNG_CRITICAL);
551     }
552 
553     /* Set atom masses. */
554     masses = malloc(sizeof(float) * n_particles);
555     for(i = 0; i < n_particles; i++)
556     {
557         stat = tng_atom_type_of_particle_nr_get(*traj, i, atom_type,
558                                                 sizeof(atom_type));
559         if(stat == TNG_CRITICAL)
560         {
561             break;
562         }
563         /* We only have water in the system. If the atom is oxygen set its
564          * mass to 16.00000, if it's a hydrogen set its mass to
565          * 1.00800. */
566         switch(atom_type[0])
567         {
568             case 'O':
569                 masses[i] = 16.00000;
570                 break;
571             case 'H':
572                 masses[i] = 1.00800;
573                 break;
574             default:
575                 printf("Failed setting atom masses. %s: %d\n",
576                        __FILE__, __LINE__);
577                 return(TNG_CRITICAL);
578         }
579     }
580     if(stat == TNG_CRITICAL)
581     {
582         free(masses);
583         printf("Failed setting atom masses. %s: %d\n",
584                __FILE__, __LINE__);
585         return(TNG_CRITICAL);
586     }
587 
588     stat = tng_particle_data_block_add(*traj, TNG_TRAJ_MASSES, "ATOM MASSES",
589                                        TNG_FLOAT_DATA, TNG_NON_TRAJECTORY_BLOCK,
590                                        1, 1, 1, 0, n_particles,
591                                        TNG_GZIP_COMPRESSION, masses);
592     free(masses);
593     masses = 0;
594     if(stat != TNG_SUCCESS)
595     {
596         printf("Failed adding atom masses. %s: %d\n",
597                __FILE__, __LINE__);
598         return(TNG_CRITICAL);
599     }
600 
601     /* Generate a custom annotation data block */
602     strcpy(annotation, "This trajectory was generated from tng_io_testing. "
603                        "It is not a real MD trajectory.");
604     if(tng_data_block_add(*traj, TNG_TRAJ_GENERAL_COMMENTS, "COMMENTS", TNG_CHAR_DATA,
605                           TNG_NON_TRAJECTORY_BLOCK, 1, 1, 1, TNG_UNCOMPRESSED,
606                           annotation) != TNG_SUCCESS)
607     {
608         printf("Failed adding details annotation data block. %s: %d\n",
609                __FILE__, __LINE__);
610         return(TNG_CRITICAL);
611     }
612 
613     /* Write file headers (includes non trajectory data blocks */
614     if(tng_file_headers_write(*traj, hash_mode) == TNG_CRITICAL)
615     {
616         printf("Cannot write file headers. %s: %d\n",
617                __FILE__, __LINE__);
618     }
619 
620 
621     tng_num_frames_per_frame_set_get(*traj, &n_frames_per_frame_set);
622 
623     data = malloc(sizeof(float) * n_particles *
624                   n_frames_per_frame_set * 3);
625     if(!data)
626     {
627         printf("Cannot allocate memory. %s: %d\n", __FILE__, __LINE__);
628         return(TNG_CRITICAL);
629     }
630 
631     tng_num_molecules_get(*traj, &tot_n_mols);
632     molpos = malloc(sizeof(float) * tot_n_mols * 3);
633 
634     /* Set initial coordinates */
635     for(i = 0; i < tot_n_mols; i++)
636     {
637         nr = i * 3;
638         /* Somewhat random coordinates (between 0 and 100),
639          * but not specifying a random seed */
640         molpos[nr] = 100.0 * rand() / (RAND_MAX + 1.0);
641         molpos[nr+1] = 100.0 * rand() / (RAND_MAX + 1.0);
642         molpos[nr+2] = 100.0 * rand() / (RAND_MAX + 1.0);
643     }
644 
645     /* Generate frame sets - each with 100 frames (by default) */
646     for(i = 0; i < N_FRAME_SETS; i++)
647     {
648         cnt = 0;
649         if(i < N_FRAME_SETS/2)
650         {
651             codec_id = TNG_GZIP_COMPRESSION;
652         }
653         else
654         {
655             codec_id = TNG_TNG_COMPRESSION;
656         }
657         for(j = 0; j < n_frames_per_frame_set; j++)
658         {
659             for(k = 0; k < tot_n_mols; k++)
660             {
661                 nr = k * 3;
662                 /* Move -1 to 1 */
663                 molpos[nr] += 2 * (rand() / (RAND_MAX + 1.0)) - 1;
664                 molpos[nr+1] += 2 * (rand() / (RAND_MAX + 1.0)) - 1;
665                 molpos[nr+2] += 2 * (rand() / (RAND_MAX + 1.0)) - 1;
666 
667                 data[cnt++] = molpos[nr];
668                 data[cnt++] = molpos[nr + 1];
669                 data[cnt++] = molpos[nr + 2];
670                 data[cnt++] = molpos[nr] + 1;
671                 data[cnt++] = molpos[nr + 1] + 1;
672                 data[cnt++] = molpos[nr + 2] + 1;
673                 data[cnt++] = molpos[nr] - 1;
674                 data[cnt++] = molpos[nr + 1] - 1;
675                 data[cnt++] = molpos[nr + 2] - 1;
676             }
677         }
678         if(tng_frame_set_with_time_new(*traj, i * n_frames_per_frame_set,
679                                        n_frames_per_frame_set, 2e-15 * (i*n_frames_per_frame_set)) != TNG_SUCCESS)
680         {
681             printf("Error creating frame set %d. %s: %d\n",
682                    i, __FILE__, __LINE__);
683             free(molpos);
684             free(data);
685             return(TNG_CRITICAL);
686         }
687 
688         tng_frame_set_particle_mapping_free(*traj);
689 
690         /* Setup particle mapping. Use 4 different mapping blocks with arbitrary
691          * mappings. */
692         for(k=0; k<150; k++)
693         {
694             mapping[k]=k;
695         }
696         if(tng_particle_mapping_add(*traj, 0, 150, mapping) != TNG_SUCCESS)
697         {
698             printf("Error creating particle mapping. %s: %d\n",
699                    __FILE__, __LINE__);
700             free(molpos);
701             free(data);
702             return(TNG_CRITICAL);
703         }
704         for(k=0; k<150; k++)
705         {
706             mapping[k]=599-k;
707         }
708         if(tng_particle_mapping_add(*traj, 150, 150, mapping) != TNG_SUCCESS)
709         {
710             printf("Error creating particle mapping. %s: %d\n",
711                    __FILE__, __LINE__);
712             free(molpos);
713             free(data);
714             return(TNG_CRITICAL);
715         }
716         for(k=0; k<150; k++)
717         {
718             mapping[k]=k+150;
719         }
720         if(tng_particle_mapping_add(*traj, 300, 150, mapping) != TNG_SUCCESS)
721         {
722             printf("Error creating particle mapping. %s: %d\n",
723                    __FILE__, __LINE__);
724             free(molpos);
725             free(data);
726             return(TNG_CRITICAL);
727         }
728         for(k=0; k<150; k++)
729         {
730             mapping[k]=449-k;
731         }
732         if(tng_particle_mapping_add(*traj, 450, 150, mapping) != TNG_SUCCESS)
733         {
734             printf("Error creating particle mapping. %s: %d\n",
735                    __FILE__, __LINE__);
736             free(molpos);
737             free(data);
738             return(TNG_CRITICAL);
739         }
740 
741         /* Add the positions in a data block */
742         if(tng_particle_data_block_add(*traj, TNG_TRAJ_POSITIONS,
743                                        "POSITIONS",
744                                        TNG_FLOAT_DATA,
745                                        TNG_TRAJECTORY_BLOCK,
746                                        n_frames_per_frame_set, 3,
747                                        1, 0, n_particles,
748 /*                                        TNG_UNCOMPRESSED, */
749                                        codec_id,
750                                        data) != TNG_SUCCESS)
751         {
752             printf("Error adding data. %s: %d\n", __FILE__, __LINE__);
753             free(molpos);
754             free(data);
755             return(TNG_CRITICAL);
756         }
757         /* Write the frame set */
758         if(tng_frame_set_write(*traj, hash_mode) != TNG_SUCCESS)
759         {
760             printf("Error writing frame set. %s: %d\n", __FILE__, __LINE__);
761             free(molpos);
762             free(data);
763             return(TNG_CRITICAL);
764         }
765     }
766 
767     free(molpos);
768     free(data);
769 
770     tng_trajectory_destroy(traj);
771     tng_trajectory_init(traj);
772     tng_input_file_set(*traj, TNG_EXAMPLE_FILES_DIR "tng_test.tng");
773 
774     stat = tng_file_headers_read(*traj, hash_mode);
775 
776     tng_first_user_name_get(*traj, temp_str, TNG_MAX_STR_LEN);
777     if(strcmp(USER_NAME, temp_str) != 0)
778     {
779         printf("User name does not match when reading written file. %s: %d\n",
780                __FILE__, __LINE__);
781         return(TNG_FAILURE);
782     }
783 
784     tng_first_program_name_get(*traj, temp_str, TNG_MAX_STR_LEN);
785     if(strcmp(PROGRAM_NAME, temp_str) != 0)
786     {
787         printf("Program name does not match when reading written file. %s: %d\n",
788                __FILE__, __LINE__);
789         return(TNG_FAILURE);
790     }
791 
792     tng_first_computer_name_get(*traj, temp_str, TNG_MAX_STR_LEN);
793     if(strcmp(COMPUTER_NAME, temp_str) != 0)
794     {
795         printf("Computer name does not match when reading written file. %s: %d\n",
796                __FILE__, __LINE__);
797         return(TNG_FAILURE);
798     }
799 
800     tng_forcefield_name_get(*traj, temp_str, TNG_MAX_STR_LEN);
801     if(strcmp(FORCEFIELD_NAME, temp_str) != 0)
802     {
803         printf("Forcefield name does not match when reading written file. %s: %d\n",
804                __FILE__, __LINE__);
805         return(TNG_FAILURE);
806     }
807 
808     tng_medium_stride_length_get(*traj, &temp_int);
809     if(temp_int != MEDIUM_STRIDE_LEN)
810     {
811         printf("Stride length does not match when reading written file. %s: %d\n",
812                __FILE__, __LINE__);
813         return(TNG_FAILURE);
814     }
815 
816     tng_long_stride_length_get(*traj, &temp_int);
817     if(temp_int != LONG_STRIDE_LEN)
818     {
819         printf("Stride length does not match when reading written file. %s: %d\n",
820                __FILE__, __LINE__);
821         return(TNG_FAILURE);
822     }
823 
824     tng_compression_precision_get(*traj, &temp_double);
825     if(temp_double != COMPRESSION_PRECISION)
826     {
827         printf("Compression precision does not match when reading written file. %s: %d\n",
828                __FILE__, __LINE__);
829         return(TNG_FAILURE);
830     }
831 
832     tng_distance_unit_exponential_get(*traj, &temp_int);
833     if(temp_int != dist_exp)
834     {
835         printf("Distance unit exponential does not match when reading written file. %s: %d\n",
836                __FILE__, __LINE__);
837         return(TNG_FAILURE);
838     }
839 
840     stat = tng_test_molecules(*traj);
841     if(stat != TNG_SUCCESS)
842     {
843         return(stat);
844     }
845 
846     stat = tng_particle_data_vector_get(*traj, TNG_TRAJ_MASSES, (void **)&masses, &read_n_frames,
847                                         &read_stride_length, &read_n_particles,
848                                         &read_n_values_per_frame, &read_data_type);
849     if(stat != TNG_SUCCESS)
850     {
851         free(masses);
852         return(stat);
853     }
854     if(read_n_particles != n_particles)
855     {
856         printf("Number of particles does not match when reading atom masses. %s: %d\n",
857                __FILE__, __LINE__);
858         return(TNG_FAILURE);
859     }
860 
861     /* Above we have written only water molecules (in the order oxygen, hydrogen, hydrogen ...).
862      * Test that the first and second as well as the very last atoms (oxygen, hydrogen and hydrogen)
863      * have the correct atom masses. */
864     if(fabs(masses[0] - 16.00000) > 0.0001 || fabs(masses[1] - 1.00800) > 0.0001 ||
865        fabs(masses[read_n_particles-1] - 1.00800) > 0.0001)
866     {
867         printf("Atom masses do not match when reading written file. %s: %d\n",
868                __FILE__, __LINE__);
869         return(TNG_FAILURE);
870     }
871     free(masses);
872 
873     i = 0;
874     while(stat == TNG_SUCCESS)
875     {
876         stat = tng_frame_set_read_next(*traj, hash_mode);
877         tng_current_frame_set_get(*traj, &frame_set);
878         tng_frame_set_prev_frame_set_file_pos_get(*traj, frame_set, &temp_int);
879         tng_frame_set_next_frame_set_file_pos_get(*traj, frame_set, &temp_int2);
880         if(i > 0)
881         {
882             if(temp_int == -1)
883             {
884                 printf("File position of previous frame set not correct. %s: %d\n",
885                        __FILE__, __LINE__);
886                 return(TNG_FAILURE);
887             }
888         }
889         else if(temp_int != -1)
890         {
891             printf("File position of previous frame set not correct. %s: %d\n",
892                     __FILE__, __LINE__);
893             return(TNG_FAILURE);
894         }
895         if(i < N_FRAME_SETS -1)
896         {
897             if(temp_int2 == -1)
898             {
899                 printf("File position of next frame set not correct. %s: %d\n",
900                        __FILE__, __LINE__);
901                 return(TNG_FAILURE);
902             }
903         }
904         else if(temp_int2 != -1)
905         {
906             printf("File position of previous next set not correct. %s: %d\n",
907                     __FILE__, __LINE__);
908             return(TNG_FAILURE);
909         }
910         i++;
911     }
912     if(stat == TNG_CRITICAL)
913     {
914         return(stat);
915     }
916 
917     tng_time_per_frame_get(*traj, &temp_double);
918     if(fabs(TIME_PER_FRAME - temp_double) > 0.000001)
919     {
920         printf("Time per frame does not match when reading written file. %s: %d\n",
921                __FILE__, __LINE__);
922         printf("Value: %e, expected value: %e\n", temp_double, TIME_PER_FRAME);
923         return(TNG_FAILURE);
924     }
925 
926     stat = tng_frame_set_nr_find(*traj, (int64_t)(0.30*N_FRAME_SETS));
927     if(stat != TNG_SUCCESS)
928     {
929         printf("Could not find frame set %"PRId64". %s: %d\n", (int64_t)0.30*N_FRAME_SETS,
930                __FILE__, __LINE__);
931         return(stat);
932     }
933 
934     stat = tng_frame_set_nr_find(*traj, (int64_t)(0.75*N_FRAME_SETS));
935     if(stat != TNG_SUCCESS)
936     {
937         printf("Could not find frame set %"PRId64". %s: %d\n", (int64_t)0.75*N_FRAME_SETS,
938                __FILE__, __LINE__);
939         return(stat);
940     }
941 
942     tng_current_frame_set_get(*traj, &frame_set);
943     tng_frame_set_frame_range_get(*traj, frame_set, &temp_int, &temp_int2);
944     if(temp_int !=  75 * n_frames_per_frame_set)
945     {
946         printf("Unexpected first frame in frame set. %s: %d\n",
947                __FILE__, __LINE__);
948         return(TNG_FAILURE);
949     }
950 
951     stat = tng_frame_set_read_current_only_data_from_block_id(*traj, hash_mode, TNG_TRAJ_POSITIONS);
952     if(stat != TNG_SUCCESS)
953     {
954         printf("Cannot read positions in current frame set. %s: %d\n",
955                __FILE__, __LINE__);
956         return(TNG_FAILURE);
957     }
958     stat = tng_frame_set_read_next_only_data_from_block_id(*traj, hash_mode, TNG_TRAJ_POSITIONS);
959     if(stat != TNG_SUCCESS)
960     {
961         printf("Cannot read positions in next frame set. %s: %d\n",
962                __FILE__, __LINE__);
963         return(TNG_FAILURE);
964     }
965     stat = tng_data_block_name_get(*traj, TNG_TRAJ_POSITIONS, temp_str, TNG_MAX_STR_LEN);
966     if(stat != TNG_SUCCESS || strcmp("POSITIONS", temp_str) != 0)
967     {
968         printf("Cannot get name of data block or unexpected name. %s: %d\n",
969                __FILE__, __LINE__);
970         return(TNG_FAILURE);
971     }
972     stat = tng_data_block_name_get(*traj, TNG_TRAJ_FORCES, temp_str, TNG_MAX_STR_LEN);
973     if(stat != TNG_FAILURE)
974     {
975         printf("Trying to retrieve name of non-existent data block did not return failure. %s: %d\n",
976                __FILE__, __LINE__);
977         return(TNG_FAILURE);
978     }
979     stat = tng_data_block_dependency_get(*traj, TNG_TRAJ_POSITIONS, &dependency);
980     if(stat != TNG_SUCCESS || dependency != TNG_FRAME_DEPENDENT + TNG_PARTICLE_DEPENDENT)
981     {
982         printf("Cannot get dependency of data block or unexpected dependency. %s: %d\n",
983                __FILE__, __LINE__);
984         return(TNG_FAILURE);
985     }
986     stat = tng_data_block_num_values_per_frame_get(*traj, TNG_TRAJ_POSITIONS, &temp_int);
987     if(stat != TNG_SUCCESS || temp_int != 3)
988     {
989         printf("Cannot get number of values per frame of data block or unexpected value. %s: %d\n",
990                __FILE__, __LINE__);
991         return(TNG_FAILURE);
992     }
993     stat = tng_data_get_stride_length(*traj, TNG_TRAJ_POSITIONS, 100, &temp_int);
994     if(stat != TNG_SUCCESS || temp_int != 1)
995     {
996         printf("Cannot get stride length of data block or unexpected value. %s: %d\n",
997                __FILE__, __LINE__);
998         return(TNG_FAILURE);
999     }
1000 
1001     return(TNG_SUCCESS);
1002 }
1003 
1004 /* This test relies on knowing that the box shape is stored as double */
tng_test_get_box_data(tng_trajectory_t traj)1005 tng_function_status tng_test_get_box_data(tng_trajectory_t traj)
1006 {
1007     int64_t n_frames, n_values_per_frame;
1008     union data_values **values = 0;
1009     char type;
1010 
1011     if(tng_data_get(traj, TNG_TRAJ_BOX_SHAPE, &values, &n_frames,
1012                     &n_values_per_frame, &type) != TNG_SUCCESS)
1013     {
1014         printf("Failed getting box shape. %s: %d\n", __FILE__, __LINE__);
1015         return(TNG_CRITICAL);
1016     }
1017 
1018     /* The X dimension in the example file is 50 */
1019     if(fabs(values[0][0].d - 50) > 0.000001)
1020     {
1021         printf("Unexpected value in box shape. %s: %d\n", __FILE__, __LINE__);
1022         return(TNG_FAILURE);
1023     }
1024 
1025     tng_data_values_free(traj, values, n_frames, n_values_per_frame, type);
1026 
1027     return(TNG_SUCCESS);
1028 }
1029 
1030 /* This test relies on knowing that the positions are stored as float
1031  * and that the data is not sparse (i.e. as many frames in the data
1032  * as in the frame set */
tng_test_get_positions_data(tng_trajectory_t traj,const char hash_mode)1033 tng_function_status tng_test_get_positions_data(tng_trajectory_t traj,
1034                                                 const char hash_mode)
1035 {
1036     int64_t i, j, k, n_frames, n_particles, n_values_per_frame;
1037     union data_values ***values = 0;
1038     char type;
1039 
1040     if(tng_particle_data_get(traj, TNG_TRAJ_POSITIONS, &values, &n_frames,
1041                              &n_particles, &n_values_per_frame, &type) !=
1042        TNG_SUCCESS)
1043     {
1044         printf("Failed getting particle positions. %s: %d\n", __FILE__, __LINE__);
1045         return(TNG_CRITICAL);
1046     }
1047 
1048     if(n_values_per_frame != 3)
1049     {
1050         printf("Number of values per frame does not match expected value. %s: %d\n",
1051                __FILE__, __LINE__);
1052         tng_particle_data_values_free(traj, values, n_frames, n_particles,
1053                                       n_values_per_frame, type);
1054         return(TNG_FAILURE);
1055     }
1056 
1057     for(i = 0; i < n_frames; i++)
1058     {
1059 //         printf("%"PRId64"\n", i);
1060         for(j = 0; j < n_particles; j++)
1061         {
1062             for(k = 0; k < n_values_per_frame; k++)
1063             {
1064 //                 printf("%f ", values[i][j][k].f);
1065                 if(values[i][j][k].f < -500 || values[i][j][k].f > 500)
1066                 {
1067                     printf("Coordinates not in range. %s: %d\n",
1068                            __FILE__, __LINE__);
1069                     tng_particle_data_values_free(traj, values, n_frames, n_particles,
1070                                                   n_values_per_frame, type);
1071                     return(TNG_FAILURE);
1072                 }
1073             }
1074 //             printf("\n");
1075         }
1076     }
1077 
1078     if(tng_particle_data_interval_get(traj, TNG_TRAJ_POSITIONS, 111000, 111499,
1079                                       hash_mode, &values, &n_particles,
1080                                       &n_values_per_frame, &type) == TNG_SUCCESS)
1081     {
1082         printf("Getting particle positions succeeded when it should have failed. %s: %d\n",
1083                __FILE__, __LINE__);
1084         return(TNG_CRITICAL);
1085     }
1086 
1087     if(tng_particle_data_interval_get(traj, TNG_TRAJ_POSITIONS, 1000, 1050,
1088                                       hash_mode, &values, &n_particles,
1089                                       &n_values_per_frame, &type) != TNG_SUCCESS)
1090     {
1091         printf("Failed getting particle positions. %s: %d\n", __FILE__, __LINE__);
1092         return(TNG_CRITICAL);
1093     }
1094 
1095     for(i = 0; i < 50; i++)
1096     {
1097 //         printf("%"PRId64"\n", i);
1098         for(j = 0; j < n_particles; j++)
1099         {
1100             for(k = 0; k < n_values_per_frame; k++)
1101             {
1102 //                 printf("%f ", values[i][j][k].f);
1103                 if(values[i][j][k].f < -500 || values[i][j][k].f > 500)
1104                 {
1105                     printf("Coordinates not in range. %s: %d\n",
1106                            __FILE__, __LINE__);
1107                     tng_particle_data_values_free(traj, values, n_frames, n_particles,
1108                                                   n_values_per_frame, type);
1109                     return(TNG_FAILURE);
1110                 }
1111             }
1112 //             printf("\n");
1113         }
1114     }
1115 
1116     tng_particle_data_values_free(traj, values, n_frames, n_particles,
1117                                   n_values_per_frame, type);
1118 
1119     return(TNG_SUCCESS);
1120 }
1121 
tng_test_utility_functions(tng_trajectory_t traj,const char hash_mode)1122 tng_function_status tng_test_utility_functions(tng_trajectory_t traj, const char hash_mode)
1123 {
1124     tng_function_status stat;
1125     int64_t n_particles, i, j, k, codec_id, n_frames, n_frames_per_frame_set;
1126     int64_t n_frames_to_read=30, stride_len, next_frame, n_blocks, *block_ids = 0;
1127     double time, multiplier;
1128     float *positions = 0;
1129 
1130     stat = tng_util_trajectory_open(TNG_EXAMPLE_FILES_DIR "tng_test.tng", 'r', &traj);
1131     if(stat != TNG_SUCCESS)
1132     {
1133         return(stat);
1134     }
1135 
1136     stat = tng_util_time_of_frame_get(traj, 50, &time);
1137     if(stat != TNG_SUCCESS || fabs(time - 100e-13) > 0.000001)
1138     {
1139         printf("Unexpected time at frame 50. %s: %d\n", __FILE__, __LINE__);
1140         printf("Value: %e, expected value: %e\n", time, 100e-13);
1141         return(stat);
1142     }
1143     stat = tng_util_time_of_frame_get(traj, 100, &time);
1144     if(stat != TNG_SUCCESS || fabs(time - 200e-13) > 0.000001)
1145     {
1146         printf("Unexpected time at frame 100. %s: %d\n", __FILE__, __LINE__);
1147         printf("Value: %e, expected value: %e\n", time, 100e-13);
1148         return(stat);
1149     }
1150 
1151     tng_num_frames_per_frame_set_get(traj, &n_frames_per_frame_set);
1152 
1153     stat = tng_util_num_frames_with_data_of_block_id_get(traj, TNG_TRAJ_POSITIONS, &n_frames);
1154     if(stat != TNG_SUCCESS || n_frames != n_frames_per_frame_set * N_FRAME_SETS)
1155     {
1156         printf("Unexpected number of frames with positions data. %s: %d\n",
1157                __FILE__, __LINE__);
1158         printf("Value: %"PRId64", expected value: %"PRId64"\n", n_frames,
1159                n_frames_per_frame_set * N_FRAME_SETS);
1160         return(stat);
1161     }
1162 
1163     tng_num_frames_per_frame_set_get(traj, &n_frames_per_frame_set);
1164 
1165     stat = tng_util_num_frames_with_data_of_block_id_get(traj, TNG_TRAJ_POSITIONS, &n_frames);
1166     if(stat != TNG_SUCCESS || n_frames != n_frames_per_frame_set * N_FRAME_SETS)
1167     {
1168         return(stat);
1169     }
1170 
1171     tng_num_particles_get(traj, &n_particles);
1172 
1173     stat = tng_util_pos_read_range(traj, 1, n_frames_to_read, &positions, &stride_len);
1174     if(stat != TNG_SUCCESS)
1175     {
1176         if(positions)
1177         {
1178             free(positions);
1179         }
1180         return(stat);
1181     }
1182 
1183     for(i = 0; i < n_frames_to_read / stride_len; i++)
1184     {
1185         for(j = 0; j < n_particles; j++)
1186         {
1187             for(k = 0; k < 3; k++)
1188             {
1189                 if(positions[i*n_particles + j*3 + k] < -500 || positions[i*n_particles + j*3 + k] > 500)
1190                 {
1191                     printf("Coordinates not in range. %s: %d\n",
1192                            __FILE__, __LINE__);
1193                     free(positions);
1194                     return(TNG_FAILURE);
1195                 }
1196             }
1197         }
1198     }
1199 
1200     free(positions);
1201 
1202     stat=tng_util_trajectory_next_frame_present_data_blocks_find(traj, n_frames_to_read,
1203                                                                  0, 0, &next_frame,
1204                                                                  &n_blocks, &block_ids);
1205     if(block_ids)
1206     {
1207         free(block_ids);
1208     }
1209     if(stat != TNG_SUCCESS || n_blocks != 1 || next_frame != n_frames_to_read + stride_len)
1210     {
1211         printf("Unexpected data blocks in next frame. %s: %d\n",
1212                __FILE__, __LINE__);
1213         return(TNG_FAILURE);
1214     }
1215 
1216     stat = tng_util_frame_current_compression_get(traj, TNG_TRAJ_POSITIONS, &codec_id, &multiplier);
1217     if(stat != TNG_SUCCESS || codec_id != TNG_GZIP_COMPRESSION)
1218     {
1219         printf("Could not get compression. %s: %d\n",
1220                __FILE__, __LINE__);
1221         return(TNG_FAILURE);
1222     }
1223     stat = tng_util_trajectory_close(&traj);
1224     if(stat != TNG_SUCCESS)
1225     {
1226         return(stat);
1227     }
1228 
1229     return(TNG_SUCCESS);
1230 }
1231 
1232 
tng_test_append(tng_trajectory_t traj,const char hash_mode)1233 tng_function_status tng_test_append(tng_trajectory_t traj, const char hash_mode)
1234 {
1235     char str[TNG_MAX_STR_LEN];
1236     int64_t n_frames, n_particles, i;
1237     double time, *velocities;
1238     tng_function_status stat;
1239 
1240     stat = tng_util_trajectory_open(TNG_EXAMPLE_FILES_DIR "tng_test.tng", 'a', &traj);
1241     if(stat != TNG_SUCCESS)
1242     {
1243         printf("Cannot open trajectory. %s: %d\n",
1244                __FILE__, __LINE__);
1245         return(stat);
1246     }
1247 
1248     stat = tng_last_user_name_set(traj, USER_NAME);
1249     if(stat != TNG_SUCCESS)
1250     {
1251         printf("Cannot set last user name. %s: %d\n",
1252                __FILE__, __LINE__);
1253         return(stat);
1254     }
1255     stat = tng_last_user_name_get(traj, str, TNG_MAX_STR_LEN);
1256     if(stat != TNG_SUCCESS)
1257     {
1258         printf("Cannot get last user name. %s: %d\n",
1259                __FILE__, __LINE__);
1260         return(stat);
1261     }
1262     stat = tng_last_program_name_set(traj, PROGRAM_NAME);
1263     if(stat != TNG_SUCCESS)
1264     {
1265         printf("Cannot set last program name. %s: %d\n",
1266                __FILE__, __LINE__);
1267         return(stat);
1268     }
1269     stat = tng_last_program_name_get(traj, str, TNG_MAX_STR_LEN);
1270     if(stat != TNG_SUCCESS)
1271     {
1272         printf("Cannot get last program name. %s: %d\n",
1273                __FILE__, __LINE__);
1274         return(stat);
1275     }
1276     stat = tng_last_computer_name_set(traj, "Still " COMPUTER_NAME);
1277     if(stat != TNG_SUCCESS)
1278     {
1279         printf("Cannot set last computer name. %s: %d\n",
1280                __FILE__, __LINE__);
1281         return(stat);
1282     }
1283     stat = tng_last_computer_name_get(traj, str, TNG_MAX_STR_LEN);
1284     if(stat != TNG_SUCCESS)
1285     {
1286         printf("Cannot get last computer name. %s: %d\n",
1287                __FILE__, __LINE__);
1288         return(stat);
1289     }
1290 
1291     stat = tng_file_headers_write(traj, hash_mode);
1292     if(stat != TNG_SUCCESS)
1293     {
1294         printf("Cannot write file headers. %s: %d\n",
1295                __FILE__, __LINE__);
1296         return(stat);
1297     }
1298 
1299     tng_num_frames_get(traj, &n_frames);
1300     tng_frame_set_of_frame_find(traj, n_frames - 1);
1301     tng_util_time_of_frame_get(traj, n_frames - 1, &time);
1302     time += TIME_PER_FRAME;
1303     tng_num_particles_get(traj, &n_particles);
1304 
1305     velocities = malloc(sizeof(double) * n_particles * 3);
1306     if(!velocities)
1307     {
1308         printf("Cannot allocate memory. %s: %d\n",
1309                __FILE__, __LINE__);
1310         return(TNG_CRITICAL);
1311     }
1312 
1313     for(i = 0; i < n_particles * 3; i++)
1314     {
1315         velocities[i] = i;
1316     }
1317 
1318     stat = tng_util_vel_with_time_double_write(traj, n_frames, time, velocities);
1319 
1320     free(velocities);
1321 
1322     stat = tng_util_trajectory_close(&traj);
1323 
1324     return(stat);
1325 }
1326 
tng_test_copy_container(tng_trajectory_t traj,const char hash_mode)1327 tng_function_status tng_test_copy_container(tng_trajectory_t traj, const char hash_mode)
1328 {
1329     tng_trajectory_t dest;
1330     tng_function_status stat;
1331 
1332     stat = tng_util_trajectory_open(TNG_EXAMPLE_FILES_DIR "tng_test.tng", 'r', &traj);
1333     if(stat != TNG_SUCCESS)
1334     {
1335         printf("Cannot open trajectory. %s: %d\n",
1336                __FILE__, __LINE__);
1337         return(stat);
1338     }
1339 
1340     stat = tng_trajectory_init_from_src(traj, &dest);
1341     if(stat != TNG_SUCCESS)
1342     {
1343         return(stat);
1344     }
1345 
1346     stat = tng_molecule_system_copy(traj, dest);
1347     if(stat != TNG_SUCCESS)
1348     {
1349         return(stat);
1350     }
1351 
1352     stat = tng_util_trajectory_close(&traj);
1353     if(stat != TNG_SUCCESS)
1354     {
1355         return(stat);
1356     }
1357     stat = tng_util_trajectory_close(&dest);
1358 
1359     return(stat);
1360 }
1361 
main()1362 int main()
1363 {
1364     tng_trajectory_t traj = 0;
1365     char time_str[TNG_MAX_DATE_STR_LEN];
1366     char version_str[TNG_MAX_STR_LEN];
1367     char hash_mode = TNG_USE_HASH;
1368 
1369     tng_version(traj, version_str, TNG_MAX_STR_LEN);
1370     printf("Test version control:\t\t\t\t");
1371     if(strncmp(TNG_VERSION, version_str, TNG_MAX_STR_LEN) == 0)
1372     {
1373         printf("Succeeded.\n");
1374     }
1375     else
1376     {
1377         printf("Failed.\n");
1378     }
1379 
1380     printf("Test Init trajectory:\t\t\t\t");
1381     if(tng_trajectory_init(&traj) != TNG_SUCCESS)
1382     {
1383         tng_trajectory_destroy(&traj);
1384         printf("Failed. %s: %d.\n", __FILE__, __LINE__);
1385         exit(1);
1386     }
1387     printf("Succeeded.\n");
1388 
1389     tng_time_get_str(traj, time_str);
1390 
1391     printf("Creation time: %s\n", time_str);
1392 
1393     tng_input_file_set(traj, TNG_EXAMPLE_FILES_DIR "tng_example.tng");
1394     tng_output_file_set(traj, TNG_EXAMPLE_FILES_DIR "tng_example_out.tng");
1395 
1396     printf("Test Read and write file:\t\t\t");
1397     if(tng_test_read_and_write_file(traj, hash_mode) != TNG_SUCCESS)
1398     {
1399         printf("Failed. %s: %d\n", __FILE__, __LINE__);
1400     }
1401     else
1402     {
1403         printf("Succeeded.\n");
1404     }
1405 
1406     printf("Test Get data:\t\t\t\t\t");
1407     if(tng_test_get_box_data(traj) != TNG_SUCCESS)
1408     {
1409         printf("Failed. %s: %d\n", __FILE__, __LINE__);
1410     }
1411     else
1412     {
1413         printf("Succeeded.\n");
1414     }
1415 
1416     printf("Test Destroy and init trajectory:\t\t");
1417     if(tng_trajectory_destroy(&traj) != TNG_SUCCESS ||
1418        tng_trajectory_init(&traj) != TNG_SUCCESS)
1419     {
1420         printf("Failed. %s: %d\n", __FILE__, __LINE__);
1421     }
1422     else
1423     {
1424         printf("Succeeded.\n");
1425     }
1426 
1427 
1428     tng_output_file_set(traj, TNG_EXAMPLE_FILES_DIR "tng_test.tng");
1429 
1430     printf("Test Write and read file:\t\t\t");
1431     if(tng_test_write_and_read_traj(&traj, hash_mode) != TNG_SUCCESS)
1432     {
1433         printf("Failed. %s: %d\n", __FILE__, __LINE__);
1434     }
1435     else
1436     {
1437         printf("Succeeded.\n");
1438     }
1439 
1440     printf("Test Get particle data:\t\t\t\t");
1441     if(tng_test_get_positions_data(traj, hash_mode) != TNG_SUCCESS)
1442     {
1443         printf("Failed. %s: %d\n",
1444                __FILE__, __LINE__);
1445     }
1446     else
1447     {
1448         printf("Succeeded.\n");
1449     }
1450 
1451     printf("Test Destroy trajectory:\t\t\t");
1452     if(tng_trajectory_destroy(&traj) != TNG_SUCCESS)
1453     {
1454         printf("Failed. %s: %d.\n", __FILE__, __LINE__);
1455         exit(1);
1456     }
1457     else
1458     {
1459         printf("Succeeded.\n");
1460     }
1461 
1462     printf("Test Utility functions:\t\t\t\t");
1463     if(tng_test_utility_functions(traj, hash_mode) != TNG_SUCCESS)
1464     {
1465         printf("Failed. %s: %d.\n", __FILE__, __LINE__);
1466         exit(1);
1467     }
1468     else
1469     {
1470         printf("Succeeded.\n");
1471     }
1472 
1473     printf("Test Append:\t\t\t\t\t");
1474     if(tng_test_append(traj, hash_mode) != TNG_SUCCESS)
1475     {
1476         printf("Failed. %s: %d.\n", __FILE__, __LINE__);
1477     }
1478     else
1479     {
1480         printf("Succeeded.\n");
1481     }
1482 
1483     printf("Test Copy trajectory container:\t\t\t");
1484     if(tng_test_copy_container(traj, hash_mode) != TNG_SUCCESS)
1485     {
1486         printf("Failed. %s: %d.\n", __FILE__, __LINE__);
1487     }
1488     else
1489     {
1490         printf("Succeeded.\n");
1491     }
1492 
1493     printf("Tests finished\n");
1494 
1495     exit(0);
1496 }
1497