atomes 1.2.1
atomes: an atomic scale modeling tool box
Loading...
Searching...
No Matches
dlp_init.c
Go to the documentation of this file.
1/* This file is part of the 'atomes' software
2
3'atomes' is free software: you can redistribute it and/or modify it under the terms
4of the GNU Affero General Public License as published by the Free Software Foundation,
5either version 3 of the License, or (at your option) any later version.
6
7'atomes' is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
8without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
9See the GNU General Public License for more details.
10
11You should have received a copy of the GNU Affero General Public License along with 'atomes'.
12If not, see <https://www.gnu.org/licenses/>
13
14Copyright (C) 2022-2025 by CNRS and University of Strasbourg */
15
22/*
23* This file: 'dlp_init.c'
24*
25* Contains:
26*
27
28 - The functions to initialize the creation of a classical force field
29
30*
31* List of functions:
32
33 int get_position_in_field_atom_from_model_id (int fat, int at);
34 int get_field_atom_id_from_model_id (field_molecule * fmol, int at);
35 int get_fragment_atom_id_from_model_id (field_molecule * fmol, int at);
36 int get_atom_id_in_fragment_from_model_id (int frag, int at);
37 int set_atom_id (field_atom* at, int c, int p, int id);
38 int get_struct_id_from_atom_id (int ids, int * aid);
39 int prepare_field_atom (int i, int j, int k, int l, int m);
40 int test_for_bonds (field_atom* at, field_atom* bt);
41 int prepare_field_struct (int ids, int sid, int yes_no_num, int * aid);
42 int bonds_between_atoms (int n, field_atom* at, field_atom* bt, int a, int b);
43 int test_for_angles (field_atom* at,
44 field_atom* bt,
45 field_atom* ct);
46 int angles_from_bonds (int n,
47 field_atom* at,
48 field_atom* bt,
49 field_atom* ct);
50 int test_for_dihedrals (field_atom* at,
51 field_atom* bt,
52 field_atom* ct,
53 field_atom* dt);
54 int dihedrals_from_angles (int n,
55 field_atom* at,
56 field_atom* bt,
57 field_atom* ct,
58 field_atom* dt);
59 int impropers_inversion (int n, int stru,
60 int at, int bt, int ct, int dt,
61 int a, int b, int c, int d);
62 int coord_sphere_multiplicity (atom * at, int id, gboolean set_atom_id, gboolean in_field_atom);
63 int find_neighbor_loop (int frag, int aid, int fai,
64 field_atom* fat,
65 field_neighbor * ngmb,
66 int sid, gboolean save_it);
67 int setup_atomic_weight (int seq);
68 int init_vdw (gboolean init);
69
70 gboolean was_not_created_struct (int ids, int num, int * aid);
71 gboolean in_bond (int at, int bd[2]);
72 gboolean are_neighbors (field_neighbor * ngb, int at);
73 gboolean are_in_bond (atom ato, int at);
74 gboolean is_numbering_possible (int frag);
75 gboolean is_this_numbering_possible_for_this_atom (int frag, field_neighbor * ngma, int atom);
76 gboolean id_n_fold_atoms_in_fragment (int frag, int limit, int num_ngb, int search_type,
77 int init, field_neighbor * ngma_init);
78 gboolean id_atoms_in_fragment (int frag, int seed);
79
80 gchar * set_field_atom_name (field_atom* ato, field_molecule * mol);
81
82 void init_all_atoms (int i);
83 void init_all_bonds ();
84 void init_all_angles ();
85 void init_all_dihedrals ();
86 void init_all_impropers_inversions (int stru);
87 void get_weight (int seq);
88 void find_atom_id_in_field_molecule ();
89 void clean_field_struct_list (field_struct * stru);
90 void init_all_field_struct (gboolean init);
91 void init_field_molecule (int i);
92 void init_dlpoly_field (classical_field * new_field);
93 void init_lammps_field (classical_field * new_field);
94
95 field_atom* init_field_atom (int id, int type, int at, int nat, int coo, int * list);
96 field_shell * init_field_shell (int id, int ia, int ib);
97 field_constraint * init_field_constraint (int id, int ia, int ib);
98 field_pmf * init_field_pmf (int id, int num[2], int * list[2], float * w[2]);
99 field_rigid * init_field_rigid (int id, int num, int * list);
100 field_tethered * init_field_tethered (int id, int num);
101 field_prop * init_field_prop (int ti, int key, gboolean show, gboolean use);
102 field_struct * init_field_struct (int st, int ai, int an, int * aid);
103 field_nth_body * init_field_nth_body (int bi, int bd, int * na, int ** ma, int ** ba);
104 field_external * init_field_external (int bi);
105 field_neighbor * get_init_neighbor (int a);
106
107 classical_field * create_force_field_data_structure (int ai);
108
109*/
110
111#include "dlp_field.h"
112#include "interface.h"
113
114extern void print_all_field_struct (field_molecule * mol, int str);
115
116int * astr;
117float val;
118extern int * duplicate_int (int num, int * old_val);
119extern float * duplicate_float (int num, float * old_val);
120int id_atom;
121int multi;
122int a_multi;
123
132int get_position_in_field_atom_from_model_id (int fat, int at)
133{
134 int i;
135 field_atom* tmp_ft = get_active_atom (tmp_fmol -> id, fat);
136 for (i=0; i<tmp_ft -> num; i++)
137 {
138 if (tmp_ft -> list[i] == at) return i;
139 }
140 return -1;
141}
142
151int get_field_atom_id_from_model_id (field_molecule * fmol, int at)
152{
153 int i;
154 field_atom* fat = fmol -> first_atom;
155 while (fat)
156 {
157 for (i=0; i< fat -> num; i++)
158 {
159 if (at == fat -> list[i]) return fat -> id;
160 }
161 fat = fat -> next;
162 }
163 return -1;
164}
165
174int get_fragment_atom_id_from_model_id (field_molecule * fmol, int at)
175{
176 int i;
177 field_atom* fat = fmol -> first_atom;
178 while (fat)
179 {
180 for (i=0; i< fat -> num; i++)
181 {
182 if (at == fat -> list[i]) return fat -> list_id[i];
183 }
184 fat = fat -> next;
185 }
186 return -1;
187}
188
197int get_atom_id_in_fragment_from_model_id (int frag, int at)
198{
199 int i, j, k;
200 field_atom* tmp_ft;
201 for (i=0; i< tmp_fmol -> mol -> natoms; i++)
202 {
203 if (tmp_fmol -> atoms_id[i][frag].a > -1 && tmp_fmol -> atoms_id[i][frag].b > -1)
204 {
205 tmp_ft = get_active_atom (tmp_fmol -> id, tmp_fmol -> atoms_id[i][frag].a);
206 j = tmp_fmol -> atoms_id[i][frag].b;
207 k = tmp_ft -> list[j];
208 if (k == at) return i;
209 }
210 }
211 return -1;
212}
213
224int set_atom_id (field_atom* at, int c, int p, int id)
225{
226 if (c == tmp_fmol -> fragments[0])
227 {
228 id_atom ++;
229 tmp_fmol -> atoms_id[id_atom][0].a = id;
230 tmp_fmol -> atoms_id[id_atom][0].b = p;
231 return id_atom;
232 }
233 else
234 {
235 return -1;
236 }
237}
238
247gchar * set_field_atom_name (field_atom* ato, field_molecule * mol)
248{
249 int i, j;
250 gchar * str;
251 field_atom* tmpatf;
252 switch (ato -> type)
253 {
254 case SPEC_ONLY:
255 str = exact_name(tmp_proj -> chemistry -> label[ato -> sp]);
256 break;
257 case OTHER:
258 tmpatf = mol -> first_atom;
259 j = 1;
260 while (tmpatf)
261 {
262 if (tmpatf -> sp == ato -> sp) j ++;
263 tmpatf = tmpatf -> next;
264 }
265 str = exact_name(g_strdup_printf ("%s-%d-%d", tmp_proj -> chemistry -> label[ato -> sp], j, mol -> id+1));
266 break;
267 default:
268 i = ato -> list[0];
269 j = tmp_proj -> atoms[0][i].coord[ato -> type-1];
270 switch (ato -> type)
271 {
272 case TOTAL_COORD_AND_SPEC:
273 str = g_strdup_printf ("%s-%d-%d",
274 exact_name(tmp_proj -> chemistry -> label[ato -> sp]),
275 tmp_coord -> geolist[0][ato -> sp][j],
276 mol -> id+1);
277 break;
278 case PARTIAL_COORD_AND_SPEC:
279 str = g_strdup_printf ("%s-%d-%d",
280 exact_name(tmp_proj -> chemistry -> label[ato -> sp]),
281 j+1,
282 mol -> id+1);
283 /*str = g_strdup_printf ("%s-%d",
284 exact_name (env_name (tmp_proj, j, ato -> sp, ato -> type, NULL)),
285 mol -> id+1);*/
286 break;
287 }
288 }
289 return str;
290}
291
304field_atom* init_field_atom (int id, int type, int at, int nat, int coo, int * list)
305{
306 int h, i, j, k, l, m, n;
307 field_atom* ato;
308 ato = g_malloc0 (sizeof*ato);
309 ato -> id = id;
310 ato -> type = type;
311 ato -> num = nat;
312 ato -> list = allocint (ato -> num);
313 ato -> list_id = allocint (ato -> num);
314 ato -> frozen_id = allocbool (ato -> num);
315 h = 0;
316 switch (type)
317 {
318 case OTHER:
319 for (i=0; i<tmp_fmol -> multi; i++)
320 {
321 for (j=0; j < nat/tmp_fmol -> multi; j++)
322 {
323 k = list[j];
324 l = tmp_fmol -> atoms_id[k][i].a;
325 m = tmp_fmol -> atoms_id[k][i].b;
326 n = get_active_atom (tmp_fmol -> id, l) -> list[m];
327 ato -> list[h] = n;
328 ato -> list_id[h] = k;
329 ato -> frozen_id[h] = FALSE;
330 tmp_fmol -> atoms_id[k][i].a = ato -> id;
331 tmp_fmol -> atoms_id[k][i].b = h;
332 tmp_proj -> atoms[0][n].faid = ato -> id;
333 h ++;
334 }
335 }
336 break;
337 default:
338 for (i=0; i < tmp_fmol -> multi; i++)
339 {
340 j = tmp_fmol -> fragments[i];
341 for (k=0; k < tmp_proj -> natomes; k++)
342 {
343 if (tmp_proj -> atoms[0][k].coord[2] == j && tmp_proj -> atoms[0][k].sp == at)
344 {
345 if (type == SPEC_ONLY)
346 {
347 ato -> list[h] = k;
348 ato -> frozen_id[h] = FALSE;
349 ato -> list_id[h] = set_atom_id (ato, tmp_proj -> atoms[0][k].coord[2], h, id);
350 h ++;
351 }
352 else
353 {
354 if (tmp_proj -> atoms[0][k].coord[type-1] == coo)
355 {
356 ato -> list[h] = k;
357 ato -> frozen_id[h] = FALSE;
358 ato -> list_id[h] = set_atom_id (ato, tmp_proj -> atoms[0][k].coord[2], h, id);
359 h ++;
360 }
361 }
362 }
363 }
364 }
365 break;
366 }
367#ifdef DEBUG
368 if (h != ato -> num) g_debug ("Something is wrong in with this field atom: id= %d, sp= %d, co= %d", id, at, coo);
369#endif
370 ato -> sp = at;
371
372 ato -> mass = 0.0;//get_force_field_atom_mass (at, 0);
373 ato -> fid = (ato -> mass == 0.0) ? -1 : tmp_field -> type;
374 ato -> afid = (ato -> mass == 0.0) ? -1 : 0;
375 ato -> mass = (ato -> mass == 0.0) ? tmp_proj -> chemistry -> chem_prop[CHEM_M][at] : ato -> mass;
376 ato -> name = g_strdup_printf ("%s", set_field_atom_name (ato, tmp_fmol));
377 ato -> charge = 0.0;
378 ato -> frozen = 0;
379 ato -> show = FALSE;
380 ato -> next = NULL;
381 ato -> prev = NULL;
382 return ato;
383}
384
394field_shell * init_field_shell (int id, int ia, int ib)
395{
396 field_shell * shell;
397 shell = g_malloc0 (sizeof*shell);
398 shell -> id = id;
399 shell -> ia[0] = ia;
400 shell -> ia[1] = ib;
401 shell -> k2 = 0.0;
402 shell -> k4 = 0.0;
403 shell -> show = FALSE;
404 shell -> use = TRUE;
405 shell -> next = NULL;
406 shell -> prev = NULL;
407 return shell;
408}
409
419field_constraint * init_field_constraint (int id, int ia, int ib)
420{
421 field_constraint * cons;
422 cons = g_malloc0 (sizeof*cons);
423 cons -> id = id;
424 cons -> ia[0] = ia;
425 cons -> ia[1] = ib;
426 cons -> av = -1.0;
427 cons -> length = 0.0;
428 cons -> show = FALSE;
429 cons -> use = TRUE;
430 cons -> next = NULL;
431 cons -> prev = NULL;
432 return cons;
433}
434
445field_pmf * init_field_pmf (int id, int num[2], int * list[2], float * w[2])
446{
447 field_pmf * pmf;
448 pmf = g_malloc0 (sizeof*pmf);
449 pmf -> id = id;
450 pmf -> av = -1.0;
451 pmf -> length = 0.0;
452 int i;
453 if (num != NULL)
454 {
455 for (i=0; i<2; i++)
456 {
457 pmf -> num[i] = num[i];
458 pmf -> list[i] = duplicate_int (num[i], list[i]);
459 pmf -> weight[i] = duplicate_float (num[i], w[i]);
460 }
461 }
462 else
463 {
464 for (i=0; i<2; i++)
465 {
466 pmf -> num[i] = 0;
467 pmf -> list[i] = NULL;
468 pmf -> weight[i] = NULL;
469 }
470 }
471 pmf -> show = FALSE;
472 pmf -> use = TRUE;
473 pmf -> next = NULL;
474 pmf -> prev = NULL;
475 return pmf;
476}
477
487field_rigid * init_field_rigid (int id, int num, int * list)
488{
489 field_rigid * rig;
490 rig = g_malloc0 (sizeof*rig);
491 rig -> id = id;
492 rig -> num = num;
493 rig -> list = NULL;
494 if (list != NULL) rig -> list = duplicate_int (num, list);
495 rig -> show = FALSE;
496 rig -> use = TRUE;
497 rig -> next = NULL;
498 rig -> prev = NULL;
499 return rig;
500}
501
510field_tethered * init_field_tethered (int id, int num)
511{
512 field_tethered * tet;
513 tet = g_malloc0 (sizeof*tet);
514 tet -> id = id;
515 tet -> num = num;
516 tet -> key = 0;
517 tet -> val = allocfloat (fvalues[activef][0][0]);
518 tet -> show = FALSE;
519 tet -> use = TRUE;
520 tet -> next = NULL;
521 tet -> prev = NULL;
522 return tet;
523}
524
535field_prop * init_field_prop (int ti, int key, gboolean show, gboolean use)
536{
537 field_prop * prop;
538 prop = g_malloc0 (sizeof*prop);
539 prop -> aid = allocint (struct_id(ti+7));
540 prop -> key = key;
541 prop -> pid = -1;
542 prop -> fpid = -1;
543 prop -> val = allocfloat (fvalues[activef][ti+1][key]);
544 prop -> show = show;
545 prop -> use = use;
546 prop -> next = NULL;
547 prop -> prev = NULL;
548 return prop;
549}
550
559int get_struct_id_from_atom_id (int ids, int * aid)
560{
561 int i, j;
562 field_struct * stru = tmp_fmol -> first_struct[ids];
563 j = struct_id (ids+7);
564 int res = 0;
565 while (stru)
566 {
567 res = stru -> id + 1;
568 for (i=0; i<j; i++)
569 {
570 if (stru -> aid[i] != aid[i])
571 {
572 res = 0;
573 break;
574 }
575 }
576 if (! res)
577 {
578 res = -(stru -> id + 1);
579 for (i=0; i<j; i++)
580 {
581 if (stru -> aid[i] != aid[j-1-i])
582 {
583 res = 0;
584 break;
585 }
586 }
587 }
588 if (res != 0) break;
589 stru = stru -> next;
590 }
591 return res;
592}
593
603gboolean was_not_created_struct (int ids, int num, int * aid)
604{
605 if (num == 0) return TRUE;
606 int j, k;
607 k = struct_id(ids+7);
608 tmp_fstr = tmp_fmol -> first_struct[ids];
609 gboolean res;
610 while (tmp_fstr)
611 {
612 res = FALSE;
613 if (ids < 6)
614 {
615 for (j=0; j<k; j++)
616 {
617 if (tmp_fstr -> aid[j] != aid[j])
618 {
619 res = TRUE;
620 break;
621 }
622 }
623 if (res)
624 {
625 res = FALSE;
626 for (j=0; j<k; j++)
627 {
628 if (tmp_fstr -> aid[k-1-j] != aid[j])
629 {
630 res = TRUE;
631 break;
632 }
633 }
634 }
635 }
636 else if (ids == 6)
637 {
638 for (j=0; j<k; j++)
639 {
640 if (tmp_fstr -> aid[j] != aid[j])
641 {
642 res = TRUE;
643 break;
644 }
645 }
646 if (res)
647 {
648 if (tmp_fstr -> aid[0] != aid[0])
649 {
650 res = TRUE;
651 }
652 else if (tmp_fstr -> aid[3] != aid[3])
653 {
654 res = TRUE;
655 }
656 else if (tmp_fstr -> aid[1] == aid[2] && tmp_fstr -> aid[2] == aid[1])
657 {
658 res = FALSE;
659 }
660 }
661 }
662 else
663 {
664 if (tmp_fstr -> aid[0] != aid[0])
665 {
666 res = TRUE;
667 }
668 else
669 {
670 if (tmp_fstr -> aid[1] != aid[1])
671 {
672 res = TRUE;
673 }
674 else if (tmp_fstr -> aid[2] == aid[2] && tmp_fstr -> aid[3] == aid[3])
675 {
676 res = FALSE;
677 }
678 else if (tmp_fstr -> aid[2] == aid[3] && tmp_fstr -> aid[3] == aid[2])
679 {
680 res = FALSE;
681 }
682 else
683 {
684 res = TRUE;
685 }
686 }
687 }
688 if (! res) return FALSE;
689 tmp_fstr = tmp_fstr -> next;
690 }
691 return TRUE;
692}
693
704field_struct * init_field_struct (int st, int ai, int an, int * aid)
705{
706 field_struct * str;
707 str = g_malloc0 (sizeof*str);
708 str -> st = st;
709 str -> id = ai;
710 if (st == 6 || st == 7)
711 {
712 str -> num = 0;
713 str -> av = 0.0;
714 }
715 else
716 {
717 str -> num = an;
718 str -> av = val;
719 }
720 str -> aid = NULL;
721 if (aid != NULL) str -> aid = duplicate_int (struct_id(st+7), aid);
722 str -> def = init_field_prop (st, 0, FALSE, TRUE);
723 str -> other = NULL;
724 str -> next = NULL;
725 str -> prev = NULL;
726 return str;
727}
728
740field_nth_body * init_field_nth_body (int bi, int bd, int * na, int ** ma, int ** ba)
741{
742 int i, j;
743 field_nth_body * nthbd;
744 nthbd = g_malloc0 (sizeof*nthbd);
745 nthbd -> id = bi;
746 nthbd -> bd = bd;
747 if (! bd)
748 {
749 nthbd -> fpid = allocint(2);
750 nthbd -> fpid[0] = nthbd -> fpid[1] = -1;
751 }
752 j = body_at (bd);
753 // Id in Vdw list
754 if (na == NULL)
755 {
756 nthbd -> na = allocint (j);
757 for (i=0; i<j; i++) nthbd -> na[i] = -1;
758 }
759 else
760 {
761 nthbd -> na = duplicate_int (j, na);
762 }
763 // Molecule id
764 nthbd -> ma = NULL;
765 nthbd -> ma = g_malloc0 (j*sizeof*nthbd -> ma);
766 // Field atom id in molecule
767 nthbd -> a = NULL;
768 nthbd -> a = g_malloc0 (j*sizeof*nthbd -> a);
769 for (i=0; i < j; i++)
770 {
771 if (ma != NULL) nthbd -> ma[i] = duplicate_int (na[i], ma[i]);
772 if (ba != NULL) nthbd -> a[i] = duplicate_int (na[i], ba[i]);
773 }
774 nthbd -> key = 0;
775 nthbd -> val = allocfloat (fvalues[activef][9+bd][nthbd -> key]);
776 nthbd -> show = FALSE;
777 nthbd -> use = TRUE;
778 nthbd -> next = NULL;
779 nthbd -> prev = NULL;
780 return nthbd;
781}
782
790field_external * init_field_external (int bi)
791{
792 field_external * nfext;
793 nfext = g_malloc0 (sizeof*nfext);
794 nfext -> id = bi;
795 nfext -> key = -1;
796 nfext -> val = NULL;
797 nfext -> use = FALSE;
798 nfext -> next = NULL;
799 nfext -> prev = NULL;
800 return nfext;
801}
802
814int prepare_field_atom (int i, int j, int k, int l, int m)
815{
816 if (tmp_fmol -> first_atom == NULL)
817 {
818 tmp_fmol -> first_atom = init_field_atom (i, j, k, l, m, NULL);
819 tmp_fat = tmp_fmol -> first_atom;
820 }
821 else
822 {
823 tmp_fat -> next = init_field_atom (i, j, k, l, m, NULL);
824 tmp_fat -> next -> prev = g_malloc0 (sizeof*tmp_fat -> next -> prev);
825 tmp_fat -> next -> prev = tmp_fat;
826 tmp_fat = tmp_fat -> next;
827 }
828 return 1;
829}
830
838void init_all_atoms (int i)
839{
840 int j, k, l, m, n;
841
842 tmp_fmol -> first_atom = NULL;
843 k = 0;
844 switch (tmp_field -> atom_init)
845 {
846 case SPEC_ONLY:
847 for (j=0; j<tmp_proj -> nspec; j++)
848 {
849 if (tmp_mol -> species[j] > 0)
850 {
851 k += prepare_field_atom (k, tmp_field -> atom_init, j, tmp_mol -> species[j]*tmp_fmol -> multi, 0);
852 }
853 }
854 break;
855 case OTHER:
856
857 break;
858 default:
859 for (j=0; j<tmp_proj -> nspec; j++)
860 {
861 if (tmp_mol -> species[j] > 0)
862 {
863 for (l=0; l < tmp_coord -> totcoord[tmp_field -> atom_init-1]; l++)
864 {
865 n = 0;
866 for (m=0; m < tmp_proj -> natomes; m++)
867 {
868 if (tmp_proj -> atoms[0][m].sp == j
869 && tmp_proj -> atoms[0][m].coord[3] == i
870 && tmp_proj -> atoms[0][m].coord[tmp_field -> atom_init -1] == l) n ++;
871 }
872 if (n > 0)
873 {
874 k += prepare_field_atom (k, tmp_field -> atom_init, j, n, l);
875 }
876 }
877 }
878 }
879 break;
880 }
881 tmp_fmol -> atoms = k;
882}
883
892gboolean in_bond (int at, int bd[2])
893{
894 if (at == bd[0] || at == bd[1])
895 {
896 return TRUE;
897 }
898 else
899 {
900 return FALSE;
901 }
902}
903
912gboolean are_neighbors (field_neighbor * ngb, int at)
913{
914 int i;
915 for (i=0; i<ngb -> num; i++)
916 {
917 if (ngb -> vois[i] == at) return TRUE;
918 }
919 return FALSE;
920}
921
930gboolean are_in_bond (atom ato, int at)
931{
932 int i;
933 for (i=0; i<ato.numv; i++)
934 {
935 if (ato.vois[i] == at) return TRUE;
936 }
937 return FALSE;
938}
939
948int test_for_bonds (field_atom* at, field_atom* bt)
949{
950 int i, j, k, l, m;
951 m = 0;
952 val = 0.0;
953 for (i=0; i < at -> num; i++)
954 {
955 j = at -> list[i];
956 for (k= 0; k < tmp_proj -> atoms[0][j].numv; k++)
957 {
958 l = tmp_proj -> atoms[0][j].vois[k];
959 if (tmp_proj -> atoms[0][l].faid == bt -> id)
960 {
961 m ++;
962 val += distance_3d (& tmp_proj -> cell, 0, & tmp_proj -> atoms[0][j], & tmp_proj -> atoms[0][l]).length;
963 }
964 }
965 }
966 if (m > 0) val /= m;
967 if (at -> id == bt -> id) m /= 2;
968 return m / tmp_fmol -> multi;
969}
970
981int prepare_field_struct (int ids, int sid, int yes_no_num, int * aid)
982{
983 if (yes_no_num > 0)
984 {
985 if (tmp_fmol -> first_struct[ids] == NULL)
986 {
987 tmp_fmol -> first_struct[ids] = init_field_struct (ids, sid, yes_no_num, aid);
988 tmp_fstr = tmp_fmol -> first_struct[ids];
989 return 1;
990 }
991 else if (was_not_created_struct (ids, sid, aid))
992 {
993 tmp_fstr = tmp_fmol -> first_struct[ids];
994 while (tmp_fstr -> next) tmp_fstr = tmp_fstr -> next;
995 tmp_fstr -> next = init_field_struct (ids, sid, yes_no_num, aid);
996 tmp_fstr -> next -> prev = tmp_fstr;
997 tmp_fstr = tmp_fstr -> next;
998 return 1;
999 }
1000 }
1001 return 0;
1002}
1003
1015int bonds_between_atoms (int n, field_atom* at, field_atom* bt, int a, int b)
1016{
1017 if ((at -> sp == a && bt -> sp == b) || (at -> sp == b && bt -> sp == a))
1018 {
1019 astr[0] = at -> id;
1020 astr[1] = bt -> id;
1021 n += prepare_field_struct (0, n, test_for_bonds (at, bt), astr);
1022 }
1023 return n;
1024}
1025
1031void init_all_bonds ()
1032{
1033 int j, k, l;
1034 tmp_fmol -> first_struct[0] = NULL;
1035 k = 0;
1036 for (j=0; j< tmp_proj -> nspec; j++)
1037 {
1038 tmp_fat = tmp_fmol -> first_atom;
1039 while (tmp_fat)
1040 {
1041 tmp_fbt = tmp_fat;
1042 while (tmp_fbt)
1043 {
1044 k = bonds_between_atoms (k, tmp_fat, tmp_fbt, j, j);
1045 tmp_fbt = tmp_fbt -> next;
1046 }
1047 tmp_fat = tmp_fat -> next;
1048 }
1049 }
1050 for (j=0; j<tmp_proj -> nspec-1; j++)
1051 {
1052 for (l=j+1; l<tmp_proj -> nspec; l++)
1053 {
1054 tmp_fat = tmp_fmol -> first_atom;
1055 while (tmp_fat)
1056 {
1057 tmp_fbt = tmp_fat -> next;
1058 while (tmp_fbt)
1059 {
1060 k = bonds_between_atoms (k, tmp_fat, tmp_fbt, j, l);
1061 tmp_fbt = tmp_fbt -> next;
1062 }
1063 tmp_fat = tmp_fat -> next;
1064 }
1065 }
1066 }
1067 tmp_fmol -> nstruct[0] = k;
1068}
1069
1079int test_for_angles (field_atom* at,
1080 field_atom* bt,
1081 field_atom* ct)
1082{
1083 int i, j, k, l, m, n, o;
1084 val = 0.0;
1085 o = 0;
1086 for (i=0; i<at -> num; i++)
1087 {
1088 j = at -> list[i];
1089 for (k=0; k<tmp_proj -> atoms[0][j].numv; k++)
1090 {
1091 l = tmp_proj -> atoms[0][j].vois[k];
1092 if (tmp_proj -> atoms[0][l].numv >= 2 && tmp_proj -> atoms[0][l].faid == bt -> id)
1093 {
1094 for (m=0; m<tmp_proj -> atoms[0][l].numv; m++)
1095 {
1096 n = tmp_proj -> atoms[0][l].vois[m];
1097 if (n != j && tmp_proj -> atoms[0][n].faid == ct -> id)
1098 {
1099 o ++;
1100 val += angle_3d (& tmp_proj -> cell, 0,
1101 & tmp_proj -> atoms[0][j],
1102 & tmp_proj -> atoms[0][l],
1103 & tmp_proj -> atoms[0][n]).angle;
1104 }
1105 }
1106 }
1107 }
1108 }
1109 if (o > 0) val /= o;
1110 if (at -> id == ct -> id) o /= 2;
1111 return o / tmp_fmol -> multi;
1112}
1113
1124int angles_from_bonds (int n, field_atom* at,
1125 field_atom* bt,
1126 field_atom* ct)
1127{
1128 astr[0] = at -> id;
1129 astr[1] = bt -> id;
1130 astr[2] = ct -> id;
1131 n += prepare_field_struct (2, n, test_for_angles (at, bt, ct), astr);
1132 return n;
1133}
1134
1140void init_all_angles ()
1141{
1142 int m, p;
1143
1144 tmp_fmol -> first_struct[2] = NULL;
1145 p = 0;
1146 field_struct * tmp_fst = tmp_fmol -> first_struct[0];
1147 for (m=0; m < tmp_fmol -> nstruct[0]; m++)
1148 {
1149 tmp_fat = get_active_atom (tmp_fmol -> id, tmp_fst -> aid[0]);
1150 tmp_fbt = get_active_atom (tmp_fmol -> id, tmp_fst -> aid[1]);
1151 tmp_fct = tmp_fmol -> first_atom;
1152 while (tmp_fct)
1153 {
1154 p = angles_from_bonds (p, tmp_fat, tmp_fbt, tmp_fct);
1155 p = angles_from_bonds (p, tmp_fbt, tmp_fat, tmp_fct);
1156 tmp_fct = tmp_fct -> next;
1157 }
1158 if (tmp_fst -> next != NULL) tmp_fst = tmp_fst -> next;
1159 }
1160 tmp_fmol -> nstruct[2] = p;
1161}
1162
1173int test_for_dihedrals (field_atom* at,
1174 field_atom* bt,
1175 field_atom* ct,
1176 field_atom* dt)
1177{
1178 int i, j, k, l, m, n, o, p, q;
1179 q = 0;
1180 val = 0.0;
1181 for (i=0; i<at -> num; i++)
1182 {
1183 j = at -> list[i];
1184 for (k=0; k<tmp_proj -> atoms[0][j].numv; k++)
1185 {
1186 l = tmp_proj -> atoms[0][j].vois[k];
1187 if (tmp_proj -> atoms[0][l].faid == bt -> id)
1188 {
1189 for (m=0; m<tmp_proj -> atoms[0][l].numv; m++)
1190 {
1191 n = tmp_proj -> atoms[0][l].vois[m];
1192 if (n != j && tmp_proj -> atoms[0][n].faid == ct -> id)
1193 {
1194 for (o=0; o<tmp_proj -> atoms[0][n].numv; o++)
1195 {
1196 p = tmp_proj -> atoms[0][n].vois[o];
1197 if (p != j && p != l && tmp_proj -> atoms[0][p].faid == dt -> id)
1198 {
1199 q ++;
1200 val += dihedral_3d (& tmp_proj -> cell, 0,
1201 & tmp_proj -> atoms[0][j],
1202 & tmp_proj -> atoms[0][l],
1203 & tmp_proj -> atoms[0][n],
1204 & tmp_proj -> atoms[0][p]).angle;
1205 }
1206 }
1207 }
1208 }
1209 }
1210 }
1211 }
1212 if (q > 0) val /= q;
1213 if (at -> id == dt -> id && bt -> id == ct -> id) q /= 2;
1214 return q / tmp_fmol -> multi;
1215}
1216
1228int dihedrals_from_angles (int n,
1229 field_atom* at,
1230 field_atom* bt,
1231 field_atom* ct,
1232 field_atom* dt)
1233{
1234 astr[0] = at -> id;
1235 astr[1] = bt -> id;
1236 astr[2] = ct -> id;
1237 astr[3] = dt -> id;
1238 n += prepare_field_struct (4, n, test_for_dihedrals (at, bt, ct, dt), astr);
1239 return n;
1240}
1241
1247void init_all_dihedrals ()
1248{
1249 int n, p;
1250
1251 tmp_fmol -> first_struct[4] = NULL;
1252 p = 0;
1253 field_struct * tmp_fst = tmp_fmol -> first_struct[2];
1254 for (n=0; n< tmp_fmol -> nstruct[2]; n++)
1255 {
1256 tmp_fat = get_active_atom (tmp_fmol -> id, tmp_fst -> aid[0]);
1257 tmp_fbt = get_active_atom (tmp_fmol -> id, tmp_fst -> aid[1]);
1258 tmp_fct = get_active_atom (tmp_fmol -> id, tmp_fst -> aid[2]);
1259 tmp_fdt = tmp_fmol -> first_atom;
1260 while (tmp_fdt)
1261 {
1262 p = dihedrals_from_angles (p, tmp_fat, tmp_fbt, tmp_fct, tmp_fdt);
1263 p = dihedrals_from_angles (p, tmp_fct, tmp_fbt, tmp_fat, tmp_fdt);
1264 tmp_fdt = tmp_fdt -> next;
1265 }
1266 if (tmp_fst -> next != NULL) tmp_fst = tmp_fst -> next;
1267 }
1268 tmp_fmol -> nstruct[4] = p;
1269}
1270
1289int impropers_inversion (int n, int stru,
1290 int at, int bt, int ct, int dt,
1291 int a, int b, int c, int d)
1292{
1293 astr[0] = at;
1294 astr[1] = bt;
1295 astr[2] = ct;
1296 astr[3] = dt;
1297 double v;
1298 //g_debug ("Improp:: at= %d, bt= %d, ct= %d, dt= %d", at, bt, ct, dt);
1299 //g_debug ("Improp:: a= %d, b= %d, c= %d, d= %d", a, b, c, d);
1300 n += prepare_field_struct (stru, n, 1, astr);
1301 //g_debug ("tmp_fstr -> id= %d", tmp_fstr -> id);
1302 if (stru == 6)
1303 {
1304 v = dihedral_3d (& tmp_proj -> cell, 0,
1305 & tmp_proj -> atoms[0][b],
1306 & tmp_proj -> atoms[0][c],
1307 & tmp_proj -> atoms[0][a],
1308 & tmp_proj -> atoms[0][d]).angle;
1309 }
1310 else
1311 {
1312 v = inversion_3d (& tmp_proj -> cell, 0,
1313 & tmp_proj -> atoms[0][a],
1314 & tmp_proj -> atoms[0][b],
1315 & tmp_proj -> atoms[0][c],
1316 & tmp_proj -> atoms[0][d]).angle;
1317 }
1318 tmp_fstr -> av += v;
1319 tmp_fstr -> num ++;
1320 //g_debug ("tmp_fstr -> id = %d, tmp_fstr -> num = %d, v= %f, av= %f\n", tmp_fstr -> id, tmp_fstr -> num, v, tmp_fstr -> av);
1321 return n;
1322}
1323
1331void init_all_impropers_inversions (int stru)
1332{
1333 int i, j, k, l, m, n, p;
1334 int * atid, * matid;
1335 atid = allocint (tmp_proj -> coord -> cmax+1);
1336 matid = allocint (tmp_proj -> coord -> cmax+1);
1337 tmp_fmol -> first_struct[stru] = NULL;
1338 p = 0;
1339 for (i=0; i<tmp_fmol -> mol -> natoms; i++)
1340 {
1341 matid[0] = tmp_fmol -> atoms_id[i][0].a;
1342 j = tmp_fmol -> atoms_id[i][0].b;
1343 k = atid[0] = get_active_atom (tmp_fmol -> id, matid[0]) -> list[j];
1344 if ((tmp_proj -> atoms[0][k].numv > 2 && stru == 6) || (tmp_proj -> atoms[0][k].numv == 3 && stru == 7))
1345 {
1346 for (m=0; m<tmp_proj -> atoms[0][k].numv; m++)
1347 {
1348 atid[m+1] = tmp_proj -> atoms[0][k].vois[m];
1349 l = tmp_proj -> atoms[0][k].vois[m];
1350 matid[m+1] = tmp_proj -> atoms[0][l].faid;
1351 }
1352 if (stru == 6)
1353 {
1354 for (l=0; l<tmp_proj -> atoms[0][k].numv-2; l++)
1355 {
1356 for (m=l+1; m<tmp_proj -> atoms[0][k].numv-1; m++)
1357 {
1358 for (n=m+1; n<tmp_proj -> atoms[0][k].numv; n++)
1359 {
1360 p = impropers_inversion (p, stru,
1361 matid[0], matid[l+1], matid[m+1], matid[n+1],
1362 atid[0], atid[l+1], atid[m+1], atid[n+1]);
1363 p = impropers_inversion (p, stru,
1364 matid[0], matid[l+1], matid[n+1], matid[m+1],
1365 atid[0], atid[l+1], atid[n+1], atid[m+1]);
1366 p = impropers_inversion (p, stru,
1367 matid[0], matid[m+1], matid[n+1], matid[l+1],
1368 atid[0], atid[m+1], atid[n+1], atid[l+1]);
1369 }
1370 }
1371 }
1372 }
1373 else
1374 {
1375 p = impropers_inversion (p, stru,
1376 matid[0], matid[1], matid[2], matid[3],
1377 atid[0], atid[1], atid[2], atid[3]);
1378 p = impropers_inversion (p, stru,
1379 matid[0], matid[2], matid[1], matid[3],
1380 atid[0], atid[2], atid[1], atid[3]);
1381 p = impropers_inversion (p, stru,
1382 matid[0], matid[3], matid[1], matid[2],
1383 atid[0], atid[3], atid[1], atid[2]);
1384 }
1385 }
1386 }
1387 tmp_fmol -> nstruct[stru] = p;
1388 tmp_fstr = tmp_fmol -> first_struct[stru];
1389 for (i=0; i<tmp_fmol -> nstruct[stru]; i++)
1390 {
1391 tmp_fstr -> def -> use = FALSE;
1392 tmp_fstr -> av /= tmp_fstr -> num;
1393 if (tmp_fstr -> next != NULL) tmp_fstr = tmp_fstr -> next;
1394 }
1395 g_free (atid);
1396 g_free (matid);
1397}
1398
1408int coord_sphere_multiplicity (atom * at, int id, gboolean set_atom_id)
1409{
1410 int j, k, l, m;
1411 gboolean not_alone;
1412 multi = 0;
1413 a_multi = 0;
1414 tmp_fct = tmp_fmol -> first_atom;
1415 while (tmp_fct)
1416 {
1417 if (tmp_fct -> sp == at -> sp)
1418 {
1419 for (j=0; j<tmp_fct -> num; j++)
1420 {
1421 k = tmp_fct -> list[j];
1422 not_alone = TRUE;
1423 if (k != at -> id)
1424 {
1425 for (l=0; l<2; l++)
1426 {
1427 if (tmp_proj -> atoms[0][k].coord[l] != at -> coord[l])
1428 {
1429 not_alone = FALSE;
1430 break;
1431 }
1432 }
1433 if (not_alone)
1434 {
1435 multi ++;
1436 if (tmp_fct -> id == tmp_fat -> id) a_multi ++;
1437 if (set_atom_id)
1438 {
1439 l = tmp_proj -> atoms[0][k].coord[2];
1440 for (m=1; m<tmp_fmol -> multi; m++)
1441 {
1442 if (l == tmp_fmol -> fragments[m])
1443 {
1444 tmp_fmol -> atoms_id[id][m].a = tmp_fct -> id;
1445 tmp_fmol -> atoms_id[id][m].b = j;
1446 tmp_fct -> list_id[j] = id;
1447 }
1448 }
1449 }
1450 }
1451 }
1452 }
1453 }
1454 tmp_fct = tmp_fct -> next;
1455 }
1456 return multi;
1457}
1458
1459int assigned;
1460int * i_weight; // For the molecule, indice = atom id for the molecule
1461int * n_weight; // For the analyzed multi, indice = atom id for the project
1462field_neighbor * init_ngb = NULL;
1463
1471field_neighbor * get_init_neighbor (int a)
1472{
1473 field_neighbor * ngb = NULL;
1474 ngb = init_ngb;
1475 while (ngb)
1476 {
1477 if (ngb -> id == a) break;
1478 ngb = ngb -> next;
1479 }
1480 return ngb;
1481}
1482
1490gboolean is_numbering_possible (int frag)
1491{
1492 int i, j, k, l, m, n, o , p;
1493 field_neighbor * ngma, * ngmb;
1494 field_atom* fat, * fbt;
1495
1496 for (i=0; i<tmp_fmol -> mol -> natoms-1; i++)
1497 {
1498 if (tmp_fmol -> atoms_id[i][frag].a > -1 && tmp_fmol -> atoms_id[i][frag].b > -1)
1499 {
1500 ngma = get_init_neighbor (i);
1501 j = tmp_fmol -> atoms_id[i][frag].a;
1502 k = tmp_fmol -> atoms_id[i][frag].b;
1503 fat = get_active_atom(tmp_fmol -> id, j);
1504 l = fat -> list[k];
1505 for (m=i+1; m<tmp_fmol -> mol -> natoms; m++)
1506 {
1507 if (tmp_fmol -> atoms_id[m][frag].a > -1 && tmp_fmol -> atoms_id[m][frag].b > -1)
1508 {
1509 ngmb = get_init_neighbor (m);
1510 n = tmp_fmol -> atoms_id[m][frag].a;
1511 o = tmp_fmol -> atoms_id[m][frag].b;
1512 fbt = get_active_atom(tmp_fmol -> id, n);
1513 p = fbt -> list[o];
1514 if (are_in_bond (tmp_proj -> atoms[0][l], p) != are_neighbors (ngma, m))
1515 {
1516 return FALSE;
1517 }
1518 if (are_in_bond (tmp_proj -> atoms[0][p], l) != are_neighbors (ngmb, i))
1519 {
1520 return FALSE;
1521 }
1522 }
1523 }
1524 }
1525 }
1526 return TRUE;
1527}
1528
1538gboolean is_this_numbering_possible_for_this_atom (int frag, field_neighbor * ngma, int atom)
1539{
1540 int m, n, o, p;
1541 field_atom* fbt;
1542
1543 for (m=0; m<tmp_fmol -> mol -> natoms; m++)
1544 {
1545 if (m != ngma -> id)
1546 {
1547 if (tmp_fmol -> atoms_id[m][frag].a > -1 && tmp_fmol -> atoms_id[m][frag].b > -1)
1548 {
1549 n = tmp_fmol -> atoms_id[m][frag].a;
1550 o = tmp_fmol -> atoms_id[m][frag].b;
1551 fbt = get_active_atom(tmp_fmol -> id, n);
1552 p = fbt -> list[o];
1553 if (are_in_bond (tmp_proj -> atoms[0][atom], p) != are_neighbors (ngma, m))
1554 {
1555 return FALSE;
1556 }
1557 }
1558 }
1559 }
1560 return TRUE;
1561}
1562
1579int find_neighbor_loop (int frag, int aid, int fai,
1580 field_atom* fat,
1581 field_neighbor * ngmb,
1582 int sid, gboolean save_it)
1583{
1584 int i, j, k;
1585 int iter = 0;
1586
1587 i = fat -> num / tmp_fmol -> multi;
1588 for (j=frag*i; j < (frag+1)*i; j++)
1589 {
1590 if (fat -> list_id[j] < 0)
1591 {
1592 k = fat -> list[j];
1593 if (are_in_bond (tmp_proj -> atoms[0][k], aid))
1594 {
1595 //g_debug ("Are linked, n_weight[%d]= %d, i_weight[%d]= %d", k+1, n_weight[k], j+1, i_weight[j]);
1596 if (n_weight[k] == i_weight[ngmb -> id])
1597 {
1598 tmp_fmol -> atoms_id[fai][frag].a = fat -> id;
1599 tmp_fmol -> atoms_id[fai][frag].b = j;
1600 fat -> list_id[j] = fai;
1601 if (save_it && iter == sid) return k;
1602 if (is_this_numbering_possible_for_this_atom (frag, ngmb, k)) iter ++;
1603 tmp_fmol -> atoms_id[fai][frag].a = -1;
1604 tmp_fmol -> atoms_id[fai][frag].b = -1;
1605 fat -> list_id[j] = -1;
1606 }
1607 }
1608 }
1609 }
1610 return iter;
1611}
1612
1626gboolean id_n_fold_atoms_in_fragment (int frag, int limit, int num_ngb, int search_type,
1627 int init, field_neighbor * ngma_init)
1628{
1629 int i, j, k, l, m, n, o, p, q, r, s, t;
1630 field_neighbor * ngma, * ngmb, * ngmc;
1631 field_atom* fat, * fbt;
1632 if (assigned < limit)
1633 {
1634 if (search_type)
1635 {
1636 i = ngma_init -> id;
1637 j = tmp_fmol -> atoms_id[i][frag].a;
1638 for (k=0; k<ngma_init -> num; k++)
1639 {
1640 l = ngma_init -> vois[k];
1641 ngma = get_init_neighbor (l);
1642 if (tmp_fmol -> atoms_id[l][frag].a < 0 && tmp_fmol -> atoms_id[l][frag].b < 0)
1643 {
1644 if (ngma -> num > num_ngb)
1645 {
1646 m = tmp_fmol -> atoms_id[l][0].a;
1647 fat = get_active_atom (tmp_fmol -> id, m);
1648
1649 n = find_neighbor_loop (frag, init, l, fat, ngma, 0, FALSE);
1650 for (o=0; o<n; o++)
1651 {
1652 p = find_neighbor_loop (frag, init, l, fat, ngma, o, TRUE);
1653 assigned ++;
1654 if (assigned == limit)
1655 {
1656 return TRUE;
1657 }
1658 else if (id_n_fold_atoms_in_fragment(frag, limit, num_ngb, 1, p, ngma))
1659 {
1660 return TRUE;
1661 }
1662 else
1663 {
1664 for (q=0; q<tmp_fmol -> mol -> natoms; q++)
1665 {
1666 if (q != l && tmp_fmol -> atoms_id[q][frag].a > -1 && tmp_fmol -> atoms_id[q][frag].b > -1)
1667 {
1668 ngmb = get_init_neighbor(q);
1669 t = 0;
1670 for (r=0; r<ngmb -> num; r++)
1671 {
1672 s = ngmb -> vois[r];
1673 ngmc = get_init_neighbor(s);
1674 if (tmp_fmol -> atoms_id[s][frag].a < 0 && tmp_fmol -> atoms_id[s][frag].b < 0 && ngmc -> num > num_ngb) t ++;
1675 }
1676 if (t)
1677 {
1678 r = tmp_fmol -> atoms_id[q][frag].a;
1679 fbt = get_active_atom (tmp_fmol -> id, r);
1680 s = tmp_fmol -> atoms_id[q][frag].b;
1681 t = fbt -> list[s];
1682 if (id_n_fold_atoms_in_fragment(frag, limit, num_ngb, 1, t, ngmb))
1683 {
1684 return TRUE;
1685 }
1686 }
1687 }
1688 }
1689 if (id_n_fold_atoms_in_fragment(frag, limit, num_ngb, 0, -1, NULL))
1690 {
1691 return TRUE;
1692 }
1693 else
1694 {
1695 assigned --;
1696 tmp_fmol -> atoms_id[l][frag].a = -1;
1697 p = tmp_fmol -> atoms_id[l][frag].b;
1698 fat -> list_id[p] = -1;
1699 tmp_fmol -> atoms_id[l][frag].b = -1;
1700 }
1701 }
1702 }
1703 }
1704 }
1705 }
1706 return FALSE;
1707 }
1708 else
1709 {
1710 for (i=0; i<tmp_fmol -> mol -> natoms; i++)
1711 {
1712 if (tmp_fmol -> atoms_id[i][frag].a < 0 && tmp_fmol -> atoms_id[i][frag].b <0)
1713 {
1714 j = 0;
1715 ngma = get_init_neighbor (i);
1716 if (ngma -> num > num_ngb)
1717 {
1718 for (k=0; k<tmp_proj -> natomes; k++)
1719 {
1720 if (tmp_proj -> atoms[0][k].coord[2] == tmp_fmol -> fragments[frag])
1721 {
1722 j ++;
1723 if ((n_weight[k] == i_weight[i]) && (ngma -> num == tmp_proj -> atoms[0][k].numv))
1724 {
1725 if (get_atom_id_in_fragment_from_model_id (frag, k) < 0)
1726 {
1727 assigned ++;
1728 l = get_field_atom_id_from_model_id (tmp_fmol, k);
1729 tmp_fmol -> atoms_id[i][frag].a = l;
1730 m = get_position_in_field_atom_from_model_id (l, k);
1731 tmp_fmol -> atoms_id[i][frag].b = m;
1732 get_active_atom(tmp_fmol -> id, l) -> list_id[m] = i;
1733 if (is_this_numbering_possible_for_this_atom (frag, ngma, k))
1734 {
1735 if (id_n_fold_atoms_in_fragment(frag, limit, num_ngb, 1, k, ngma))
1736 {
1737 return TRUE;
1738 }
1739 else if (id_n_fold_atoms_in_fragment(frag, limit, num_ngb, 0, -1, NULL))
1740 {
1741 return TRUE;
1742 }
1743 }
1744 tmp_fmol -> atoms_id[i][frag].a = -1;
1745 tmp_fmol -> atoms_id[i][frag].b = -1;
1746 get_active_atom(tmp_fmol -> id, l) -> list_id[m] = -1;
1747 assigned --;
1748 }
1749 }
1750 if (j == tmp_fmol -> mol -> natoms)
1751 {
1752 break;
1753 }
1754 }
1755 }
1756 return FALSE;
1757 }
1758 }
1759 }
1760 return FALSE;
1761 }
1762 }
1763 else
1764 {
1765 return TRUE;
1766 }
1767}
1768
1777gboolean id_atoms_in_fragment (int frag, int seed)
1778{
1779 int i, j;
1780 i = 0;
1781 for (j=0; j<tmp_proj -> natomes; j++)
1782 {
1783 if (tmp_proj -> atoms[0][j].coord[2] == tmp_fmol -> fragments[frag])
1784 {
1785 if (tmp_proj -> atoms[0][j].numv > 1) i++;
1786 }
1787 }
1788 if (seed > -1)
1789 {
1790 field_neighbor * ngb = get_init_neighbor (seed);
1791 if (id_n_fold_atoms_in_fragment (frag, i, 1, 0, seed, ngb))
1792 {
1793 return id_n_fold_atoms_in_fragment (frag, tmp_fmol -> mol -> natoms, 0, 0, seed, ngb);
1794 }
1795 else
1796 {
1797 return FALSE;
1798 }
1799 }
1800 else
1801 {
1802 if (id_n_fold_atoms_in_fragment (frag, i, 1, 0, -1, NULL))
1803 {
1804 return id_n_fold_atoms_in_fragment (frag, tmp_fmol -> mol -> natoms, 0, 0, -1, NULL);
1805 }
1806 else
1807 {
1808 return FALSE;
1809 }
1810 }
1811}
1812
1820void get_weight (int seq)
1821{
1822 int h, i, j, k;
1823 for (i=0; i<tmp_proj -> natomes; i++)
1824 {
1825 for (j=0; j<tmp_proj -> atoms[0][i].numv; j++)
1826 {
1827 k = tmp_proj -> atoms[0][i].vois[j];
1828 n_weight[i] += tmp_proj -> atoms[0][k].numv;
1829 }
1830 }
1831 int * weight = duplicate_int (tmp_proj -> natomes, n_weight);
1832 for (h=1; h<seq; h++)
1833 {
1834 for (i=0; i<tmp_proj -> natomes; i++)
1835 {
1836 for (j=0; j<tmp_proj -> atoms[0][i].numv; j++)
1837 {
1838 k = tmp_proj -> atoms[0][i].vois[j];
1839 weight[i] += n_weight[k];
1840 }
1841 }
1842 for (i=0; i<tmp_proj -> natomes; i++)
1843 {
1844 n_weight[i] = weight[i];
1845 }
1846 }
1847 g_free (weight);
1848
1849#ifdef DEBUG
1850 /*for (i=0; i<tmp_proj -> natomes; i++)
1851 {
1852 j = tmp_proj -> atoms[0][i].sp;
1853 g_debug ("N_WEIGHT:: %s i= %d, n_weight[i] = %d", tmp_proj -> chemistry -> label[j], i+1, n_weight[i]);
1854 }*/
1855#endif
1856}
1857
1865int setup_atomic_weight (int seq)
1866{
1867 int h, i, j, k, l;
1868 field_neighbor * ngm;
1869 i_weight = allocint (tmp_fmol -> mol -> natoms);
1870 init_ngb = NULL;
1871 for (i=0; i<tmp_proj -> natomes; i++)
1872 {
1873 if (tmp_proj -> atoms[0][i].coord[2] == tmp_fmol -> fragments[0])
1874 {
1875 if (init_ngb)
1876 {
1877 ngm -> next = g_malloc0(sizeof*ngm);
1878 ngm -> next -> prev = ngm;
1879 ngm = ngm -> next;
1880 }
1881 else
1882 {
1883 init_ngb = g_malloc0(sizeof*init_ngb);
1884 ngm = init_ngb;
1885 }
1886 ngm -> num = tmp_proj -> atoms[0][i].numv;
1887 ngm -> vois = duplicate_int (ngm -> num, tmp_proj -> atoms[0][i].vois);
1888 j = get_atom_id_in_fragment_from_model_id (0, i);
1889 ngm -> id = j;
1890 for (k=0; k<ngm -> num; k++)
1891 {
1892 l = ngm -> vois[k];
1893 ngm -> vois[k] = get_atom_id_in_fragment_from_model_id (0, l);
1894 i_weight[j] += tmp_proj -> atoms[0][l].numv;
1895 }
1896 }
1897 }
1898
1899 int * weight = duplicate_int (tmp_fmol -> mol -> natoms, i_weight);
1900 for (h=1; h<seq; h++)
1901 {
1902 for (i=0; i<tmp_fmol -> mol -> natoms; i++)
1903 {
1904 ngm = get_init_neighbor (i);
1905 if (ngm)
1906 {
1907 for (j=0; j<ngm -> num; j++)
1908 {
1909 k = ngm -> vois[j];
1910 weight[i] += i_weight[k];
1911 }
1912 }
1913 }
1914 for (i=0; i<tmp_fmol -> mol -> natoms; i++)
1915 {
1916 i_weight[i] = weight[i];
1917 }
1918 }
1919
1920#ifdef DEBUG
1921 /*field_atom* fat;
1922 for (i=0; i<tmp_fmol -> mol -> natoms; i++)
1923 {
1924 fat = get_active_atom (tmp_fmol -> id, tmp_fmol -> atoms_id[i][0].a);
1925 j = tmp_fmol -> atoms_id[i][0].b;
1926 k = fat -> list[j];
1927 l = tmp_proj -> atoms[0][k].sp;
1928 g_debug ("I_WEIGHT:: %s i= %d, i_weight[i] = %d", tmp_proj -> chemistry -> label[l], i, i_weight[i]);
1929 }*/
1930#endif
1931
1932 g_free (weight);
1933
1934 int seed = -1;
1935 for (i=0; i<tmp_fmol -> mol -> natoms; i++)
1936 {
1937 seed = i;
1938 for (j=0; j<tmp_fmol -> mol -> natoms; j++)
1939 {
1940 if (j!=i)
1941 {
1942 if (i_weight[i] == i_weight[j])
1943 {
1944 seed = -1;
1945 break;
1946 }
1947 }
1948 }
1949 if (seed > -1) break;
1950 }
1951 get_weight (seq);
1952 return seed;
1953}
1954
1960void find_atom_id_in_field_molecule ()
1961{
1962 int i, j, k, l;
1963 gboolean sym = FALSE;
1964 gboolean done;
1965 int ** coordnum;
1966 int start = 1;
1967 int seed = -1;
1968 assigned = 0;
1969
1970 coordnum = allocdint(tmp_proj -> nspec, tmp_coord -> totcoord[1]);
1971 tmp_fat = tmp_fmol -> first_atom;
1972 while (tmp_fat)
1973 {
1974 if (tmp_fat -> num == tmp_fmol -> multi)
1975 {
1976 i = tmp_fat -> list_id[0];
1977 // Multiplicity = num of atoms ...
1978 // Unic atom, ie. as much as identical molecular fragments
1979 for (j=1; j<tmp_fmol -> multi; j++)
1980 {
1981 for (k=1; k<tmp_fat -> num; k++)
1982 {
1983 l = tmp_fat -> list[k];
1984 if (tmp_proj -> atoms[0][l].coord[2] == tmp_fmol -> fragments[j])
1985 {
1986 tmp_fmol -> atoms_id[i][j].a = tmp_fmol -> atoms_id[i][0].a;
1987 tmp_fmol -> atoms_id[i][j].b = k;
1988 tmp_fat -> list_id[j] = i;
1989 break;
1990 }
1991 }
1992 }
1993 assigned ++;
1994 }
1995 tmp_fat = tmp_fat -> next;
1996 }
1997
1998 for (i=0; i< tmp_fmol -> mol -> natoms; i++)
1999 {
2000 tmp_fat = get_active_atom (tmp_fmol -> id, tmp_fmol -> atoms_id[i][0].a);
2001 done = TRUE;
2002 for (j=0; j<tmp_fat -> num; j++)
2003 {
2004 if (tmp_fat -> list_id[j] < 0)
2005 {
2006 done = FALSE;
2007 break;
2008 }
2009 }
2010 if (! done)
2011 {
2012 j = tmp_fmol -> atoms_id[i][0].b;
2013 k = tmp_fat -> list[j];
2014 l = coord_sphere_multiplicity (& tmp_proj -> atoms[0][k], i, FALSE);
2015 coordnum[tmp_proj -> atoms[0][k].sp][tmp_proj -> atoms[0][k].coord[1]] = l+1;
2016 if (l == tmp_fmol -> multi - 1 || a_multi == tmp_fmol -> multi - 1)
2017 {
2018 coord_sphere_multiplicity (& tmp_proj -> atoms[0][k], i, TRUE);
2019 assigned ++;
2020 }
2021 }
2022 }
2023 for (i=0; i<tmp_proj -> nspec; i++)
2024 {
2025 for (j=0; j<tmp_coord -> totcoord[1]; j++)
2026 {
2027 if (coordnum[i][j] > 1 && coordnum[i][j] % tmp_fmol -> multi != 0) g_debug ("This should not happen");
2028 coordnum[i][j] /= tmp_fmol -> multi;
2029 }
2030 }
2031 for (k=6; k>0; k--)
2032 {
2033 sym = TRUE;
2034 for (i=0; i<tmp_proj -> nspec; i++)
2035 {
2036 for (j=0; j<tmp_coord -> totcoord[1]; j++)
2037 {
2038 if (coordnum[i][j] > 1 && coordnum[i][j] >= k && coordnum[i][j] % k != 0) sym = FALSE;
2039 }
2040 }
2041 if (sym) break;
2042 }
2043 g_free (coordnum);
2044
2045 if (assigned < tmp_fmol -> mol -> natoms)
2046 {
2047#ifdef DEBUG
2048 if ((sym || assigned == 0))
2049 {
2050 g_debug ("Molecule seems to be symmetric, k= %d", k);
2051 g_debug ("Mol -> id: %d", tmp_fmol -> id);
2052 g_debug ("Mol -> multi= %d", tmp_fmol -> multi);
2053 g_debug ("Mol -> field atoms= %d", tmp_fmol -> atoms);
2054 g_debug ("Mol -> atomes= %d", tmp_fmol -> mol -> natoms);
2055 g_debug ("Assigned so far= %d", assigned);
2056 }
2057#endif
2058
2059 // k is a weight factor
2060 // The values of the total coord of each atom's neighbors
2061 // is sum k times over, this is a safety for systems with
2062 // a peculiar symetry, otherwise it could be very very long
2063 k = 3;
2064 j = assigned;
2065 gboolean done = FALSE;
2066 if (n_weight) g_free (n_weight);
2067 n_weight = NULL;
2068 n_weight = allocint (tmp_proj -> natomes);
2069 while (! done && (seed < 0 || (tmp_fmol -> mol -> natoms && k < tmp_fmol -> mol -> natoms)))
2070 {
2071 seed = setup_atomic_weight (k);
2072 if (seed < 0) k ++;
2073 done = TRUE;
2074 }
2075#ifdef DEBUG
2076 g_debug ("Assigning using: seed= %d, w_factor= %d", seed, k);
2077#endif
2078
2079 for (i=start; i<tmp_fmol -> multi; i++)
2080 {
2081#ifdef DEBUG
2082 g_debug (" -> fragment id: %d", i+1);
2083#endif
2084 assigned = j;
2085 while (! id_atoms_in_fragment (i, seed) && k < tmp_fmol -> mol -> natoms)
2086 {
2087#ifdef DEBUG
2088 g_debug ("**********************************");
2089 g_debug ("*** NOT ABLE TO ASSIGN NUMBERS ***");
2090 g_debug ("*** Mol -> Id= %3d ***", tmp_fmol -> id);
2091 g_debug ("*** Mol -> Multi= %3d ***", tmp_fmol -> multi);
2092 g_debug ("*** Fragment= %3d ***", i+1);
2093 g_debug ("*** k= %3d ***", k);
2094 g_debug ("*** assigned= %10d ***", assigned);
2095 g_debug ("**********************************");
2096#endif
2097 k ++;
2098 if (n_weight) g_free (n_weight);
2099 n_weight = NULL;
2100 n_weight = allocint (tmp_proj -> natomes);
2101 seed = setup_atomic_weight (k);
2102 }
2103 for (l=0; l<tmp_proj -> natomes; l++)
2104 {
2105 if (tmp_proj -> atoms[0][l].coord[2] == tmp_fmol -> fragments[i])
2106 {
2107 tmp_proj -> atoms[0][l].faid = get_field_atom_id_from_model_id (tmp_fmol, l);
2108 tmp_proj -> atoms[0][l].fid = get_atom_id_in_fragment_from_model_id (i, l);
2109 }
2110 }
2111#ifdef DEBUG
2112 g_debug (" -> done !");
2113#endif
2114 }
2115 }
2116 if (i_weight) g_free (i_weight);
2117 i_weight = NULL;
2118 if (n_weight) g_free (n_weight);
2119 n_weight = NULL;
2120 if (init_ngb) g_free (init_ngb);
2121 init_ngb = NULL;
2122#ifdef DEBUG
2123 g_debug ("Assigning completed !");
2124#endif
2125}
2126
2134void clean_field_struct_list (field_struct * stru)
2135{
2136 while (stru)
2137 {
2138 if (stru -> next)
2139 {
2140 stru = stru -> next;
2141 g_free (stru -> prev);
2142 }
2143 else
2144 {
2145 g_free (stru);
2146 stru = NULL;
2147 }
2148 }
2149}
2150
2158void init_all_field_struct (gboolean init)
2159{
2160 int i;
2161 for (i=0; i<tmp_proj -> natomes; i++)
2162 {
2163 if (tmp_proj -> atoms[0][i].coord[2] == tmp_fmol -> fragments[0])
2164 {
2165 tmp_proj -> atoms[0][i].faid = get_field_atom_id_from_model_id (tmp_fmol, i);
2166 tmp_proj -> atoms[0][i].fid = get_atom_id_in_fragment_from_model_id (0, i);
2167 }
2168 }
2169 if (tmp_fmol -> multi > 1 && init) find_atom_id_in_field_molecule ();
2170
2171 for (i=0; i<8; i++)
2172 {
2173 tmp_fmol -> nstruct[i] = 0;
2174 if (tmp_fmol -> first_struct[i] && tmp_fmol -> nstruct[i])
2175 {
2176 clean_field_struct_list (tmp_fmol -> first_struct[i]);
2177 }
2178 tmp_fmol -> first_struct[i] = NULL;
2179 }
2180
2181 // Bonds start
2182 astr = g_malloc0 (2*sizeof*astr);
2183 init_all_bonds ();
2184 g_free (astr);
2185
2186 if (tmp_fmol -> nstruct[0] > 0)
2187 {
2188 tmp_fmol -> first_struct[1] = duplicate_field_struct_list (tmp_fmol -> first_struct[0], TRUE);
2189 tmp_fmol -> nstruct[1] = tmp_fmol -> nstruct[0];
2190 }
2191 // Bonds end
2192
2193 // Angle start
2194 astr = g_malloc0 (3*sizeof*astr);
2195 init_all_angles ();
2196 g_free (astr);
2197 if (tmp_fmol -> nstruct[2] > 0)
2198 {
2199 tmp_fmol -> nstruct[3] = tmp_fmol -> nstruct[2];
2200 tmp_fmol -> first_struct[3] = duplicate_field_struct_list (tmp_fmol -> first_struct[2], TRUE);
2201 }
2202 // Angle end
2203
2204 // Dihedral start
2205 astr = g_malloc0 (4*sizeof*astr);
2206 init_all_dihedrals ();
2207
2208 // Torsional restraints
2209 if (tmp_fmol -> nstruct[4] > 0)
2210 {
2211 tmp_fmol -> nstruct[5] = tmp_fmol -> nstruct[4];
2212 tmp_fmol -> first_struct[5] = duplicate_field_struct_list (tmp_fmol -> first_struct[4], TRUE);
2213 }
2214
2215// Improper and inversion angles
2216 for (i=6; i<8; i++) init_all_impropers_inversions (i);
2217 g_free (astr);
2218}
2219
2227void init_field_molecule (int i)
2228{
2229 int j, k;
2230 tmp_fmol -> id = i;
2231 tmp_fmol -> multi = tmp_mol -> multiplicity;
2232 tmp_fmol -> name = g_strdup_printf("MOL-%d", i+1);
2233 tmp_fmol -> fragments = NULL;
2234 tmp_fmol -> fragments = allocint (tmp_fmol -> multi);
2235 for (j=0; j<tmp_fmol -> multi; j++)
2236 {
2237 tmp_fmol -> fragments[j] = tmp_mol -> fragments[j];
2238 }
2239 tmp_fmol -> mol = tmp_mol;
2240
2241 // Atoms start
2242 id_atom = -1;
2243 tmp_fmol -> atoms_id = g_malloc (tmp_mol -> natoms*sizeof*tmp_fmol -> atoms_id);
2244 for (j=0; j<tmp_mol -> natoms; j++)
2245 {
2246 tmp_fmol -> atoms_id[j] = g_malloc0 (tmp_fmol -> multi*sizeof*tmp_fmol -> atoms_id[j]);
2247 for (k=0; k<tmp_fmol -> multi; k++)
2248 {
2249 tmp_fmol -> atoms_id[j][k].a = -1;
2250 tmp_fmol -> atoms_id[j][k].b = -1;
2251 }
2252 }
2253
2254 init_all_atoms (i);
2255 // Atoms end
2256
2257 // setup_field_molecule_neighbors (i, tmp_proj);
2258
2259 tmp_fmol -> shells = 0;
2260 // Field struct
2261 init_all_field_struct (TRUE);
2262
2263 tmp_fmol -> shells = 0;
2264 tmp_fmol -> first_shell = NULL;
2265 tmp_fmol -> constraints = 0;
2266 tmp_fmol -> first_constraint = NULL;
2267 tmp_fmol -> pmfs = 0;
2268 tmp_fmol -> first_pmf = NULL;
2269 tmp_fmol -> rigids = 0;
2270 tmp_fmol -> first_rigid = NULL;
2271 tmp_fmol -> tethered = 0;
2272 tmp_fmol -> first_tethered = NULL;
2273
2274 tmp_fmol -> show = FALSE;
2275 tmp_fmol -> show_id = FALSE;
2276 tmp_fmol -> next = NULL;
2277 tmp_fmol -> prev = NULL;
2278}
2279
2287int init_vdw (gboolean init)
2288{
2289 int nvdw;
2290 int i, j, k, l;
2291 tmp_fmol = tmp_field -> first_molecule;
2292 i = 0;
2293 while (tmp_fmol)
2294 {
2295 i += tmp_fmol -> atoms;
2296 tmp_fmol = tmp_fmol -> next;
2297 }
2298 nvdw = i;
2299 gchar ** to_be_vdw = g_malloc (nvdw*sizeof*to_be_vdw);
2300 int * vdw_mlist = allocint (nvdw);
2301 int ** vdw_aids = allocdint (nvdw,nvdw);
2302 int ** vdw_mids = allocdint (nvdw,nvdw);
2303 int vdw_na[2];
2304 int * vdw_a[2], * vdw_ma[2];
2305 gboolean add_vdw;
2306 i = 0;
2307 tmp_fmol = tmp_field -> first_molecule;
2308 while (tmp_fmol)
2309 {
2310 tmp_fat = tmp_fmol -> first_atom;
2311 while (tmp_fat)
2312 {
2313 add_vdw = TRUE;
2314 for (j=0; j<i; j++)
2315 {
2316 if (g_strcmp0 (to_be_vdw[j], tmp_fat -> name) == 0)
2317 {
2318 add_vdw = FALSE;
2319 vdw_mlist[j] ++;
2320 vdw_mids[j][vdw_mlist[j]] = tmp_fmol -> id;
2321 vdw_aids[j][vdw_mlist[j]] = tmp_fat -> id;
2322 break;
2323 }
2324 }
2325 if (add_vdw)
2326 {
2327 to_be_vdw[j] = g_strdup_printf ("%s", tmp_fat -> name);
2328 vdw_mids[j][0] = tmp_fmol -> id;
2329 vdw_aids[j][0] = tmp_fat -> id;
2330 i ++;
2331 }
2332 tmp_fat = tmp_fat -> next;
2333 }
2334 tmp_fmol = tmp_fmol -> next;
2335 }
2336
2337 j = 0 ;
2338 tmp_fbody = NULL;
2339 for (k=0; k<i; k++)
2340 {
2341 vdw_na[0] = vdw_na[1] = vdw_mlist[k]+1;
2342 vdw_a[0] = vdw_a[1] = duplicate_int (vdw_na[0], vdw_aids[k]);
2343 vdw_ma[0] = vdw_ma[1] = duplicate_int (vdw_na[0], vdw_mids[k]);
2344 if (tmp_fbody == NULL)
2345 {
2346 if (init)
2347 {
2348 tmp_field -> first_body[0] = init_field_nth_body (j, 0, vdw_na, vdw_ma, vdw_a);
2349 tmp_fbody = tmp_field -> first_body[0];
2350 }
2351 else
2352 {
2353 comp_fbody = init_field_nth_body (j, 0, vdw_na, vdw_ma, vdw_a);
2354 tmp_fbody = comp_fbody;
2355 }
2356 }
2357 else
2358 {
2359 tmp_fbody -> next = init_field_nth_body (j, 0, vdw_na, vdw_ma, vdw_a);
2360 tmp_fbody -> next -> prev = tmp_fbody;
2361 tmp_fbody = tmp_fbody -> next;
2362 }
2363 j ++;
2364 }
2365 for (k=0; k<i-1; k++)
2366 {
2367 for (l=k+1; l<i; l++)
2368 {
2369 vdw_na[0] = vdw_mlist[k]+1;
2370 vdw_na[1] = vdw_mlist[l]+1;
2371 vdw_a[0] = duplicate_int (vdw_na[0], vdw_aids[k]);
2372 vdw_a[1] = duplicate_int (vdw_na[1], vdw_aids[l]);
2373 vdw_ma[0] = duplicate_int (vdw_na[0], vdw_mids[k]);
2374 vdw_ma[1] = duplicate_int (vdw_na[1], vdw_mids[l]);
2375 tmp_fbody -> next = init_field_nth_body (j, 0, vdw_na, vdw_ma, vdw_a);
2376 tmp_fbody -> next -> prev = tmp_fbody;
2377 tmp_fbody = tmp_fbody -> next;
2378 j ++;
2379 }
2380 }
2381 g_free (to_be_vdw);
2382 g_free (vdw_aids);
2383 g_free (vdw_mids);
2384 g_free (vdw_mlist);
2385 return j;
2386}
2387
2395void init_dlpoly_field (classical_field * new_field)
2396{
2397 // Preparing DLPOLY Control file
2398 new_field -> sys_opts = allocdouble (17);
2399 new_field -> io_opts = allocdouble (23);
2400 new_field -> ana_opts = allocdouble (17);
2401 new_field -> elec_opts = allocdouble (11);
2402 new_field -> vdw_opts = allocdouble (6);
2403 new_field -> met_opts = allocdouble (2);
2404 new_field -> equi_opts = allocdouble (17);
2405 new_field -> thermo_opts= allocdouble (10);
2406 new_field -> md_opts = allocdouble (20);
2407 new_field -> out_opts = allocdouble (31);
2408/*
2409 new_field -> io_opts[]
2410 */
2411 // Set job time control (yes/no)
2412 new_field -> io_opts[0] = 1.0;
2413 // Job time control (s):
2414 new_field -> io_opts[1] = 10000.0;
2415 // Set job closure time (yes/no)
2416 new_field -> io_opts[2] = 1.0;
2417 // Job closure time (s):
2418 new_field -> io_opts[3] = 100.0;
2419 // 4 (1/0) General I/O read interace
2420 // 5 General I/O method (0= MPI, 1= FORTRAN, 2= Traditional, 3= netCDF)
2421 // Reader count:
2422 new_field -> io_opts[6] = 12;
2423 // Batch size:
2424 new_field -> io_opts[7] = 2000000;
2425 // Buffer size:
2426 new_field -> io_opts[8] = 20000;
2427 // 9 (1/0) Error check on general read I/O
2428 // 10 (1/0) General I/O write interace
2429 // 11 General I/O method (0= MPI, 1= FORTRAN, 2= Traditional, 3= netCDF)
2430 // Real precision (0 = 32 bits, 1= 64 bits)
2431 new_field -> io_opts[12] = 1;
2432 // 13 Type (0= Sorted, 1= Unsorted)
2433 // Reader count:
2434 new_field -> io_opts[14] = 12;
2435 // Batch size:
2436 new_field -> io_opts[15] = 2000000;
2437 // Buffer size:
2438 new_field -> io_opts[16] = 20000;
2439 // 17 (1/0) Error check on general write I/O
2440 // 18 (1/0) Define random number generator seeds
2441 // n1 seed:
2442 new_field -> io_opts[19] = 1;
2443 // n2 seed:
2444 new_field -> io_opts[20] = 20;
2445 // n3 seed:
2446 new_field -> io_opts[21] = 300;
2447 // 22 (1/0) limit the number of processors in the z direction
2448
2449 /*
2450 new_field -> sys_opts[]
2451 */
2452 // Relative dieletric constant:
2453 new_field -> sys_opts[0] = 1.0;
2454 // 1 = allowed local variation of the system density (%)
2455 // 2 (1/0) Ingore the particle indices in CONFIG file
2456 // 3 (1/0) Ignore strict checks
2457 // 4 (1/0) Skip detailed topology reporting
2458 // 5 (1/0) Ignore center of mass removal
2459 // 6 (1/0) Set tolerance of the relaxed shell model
2460 // Tolerance for the relaxed shell model:
2461 new_field -> sys_opts[7] = 1.0;
2462 // 8 (1/0) Set the subscelling threshold density of particles per link cell
2463 // Subscelling threshold density of particles per link cell:
2464 new_field -> sys_opts[9] = 50.0;
2465 // 10 (1/0) Create an expanded version of the current model ...
2466 // ... on x:
2467 new_field -> sys_opts[11] = 1.0;
2468 // ... on y:
2469 new_field -> sys_opts[12] = 1.0;
2470 // .. on z:
2471 new_field -> sys_opts[13] = 1.0;
2472 // 14 (1/0) Restart a calculation
2473 // 15 (0= Continue current, 1= Start new without T reset, 2= Start new with T reset)
2474
2475 /*
2476 new_field -> equi_opts[]
2477 */
2478 new_field -> equi_opts[0] = 1.0;
2479 // Equilibrate for the first N steps, N=
2480 new_field -> equi_opts[1] = 20.0;
2481 // Scale temperature every N steps, N=
2482 new_field -> equi_opts[3] = 1.0;
2483 // Cap forces fmax=
2484 new_field -> equi_opts[5] = 1000.0;
2485 // Resample the momenta distribution during N steps N=
2486 new_field -> equi_opts[7] = 1.0;
2487 // 8 (1/0) = Minimize system configuration
2488 // 9 Minimize target (0=Force, 1=Energy, 2= Distance)
2489 // Minimization threshold:
2490 new_field -> equi_opts[10] = 50.0;
2491 // Minimization to be done every N steps, N=
2492 new_field -> equi_opts[11] = 1.0;
2493 // 12 (1/0) = Optimize system configuration from start
2494 // 13 Minimize target (0=Force, 1=Energy, 2= Distance)
2495 // Minimization threshold:
2496 new_field -> equi_opts[14] = 50.0;
2497 // 15 (1/0) = Perform zero temperature (10K) optimization
2498 // 16 (1/0) = Include equilibration data in statistics
2499
2500 /*
2501 new_field -> vdw_opts[]
2502 */
2503 // Vdw cutoff
2504 // 0 (1/0) Use VdW interactions
2505 new_field -> vdw_opts[0] = 1.0;
2506 // Vdw Cutoff radius
2507 new_field -> vdw_opts[1] = 12.0;
2508 // 2 (1/0) Enforce direct VdW calculation
2509 // 3 (1/0) Apply force shifting
2510 // 4 (1/0) Apply mixing rule
2511 // 5 Mixing rule (0= LB, 1= FH, 2= H, 3= HHG, 4= TT, 5= F)
2512
2513 /*
2514 new_field -> elec_opts[]
2515 */
2516 // 0 (1/0) Evaluate electrostatics
2517 new_field -> elec_opts[0] = 1.0;
2518 // Long range cutoff:
2519 new_field -> elec_opts[1] = 12.0;
2520 // 2 (1/0) Add padding to rcut
2521 // 3 Padding in angstrom
2522 // 4 (1/0) Use extended coulombic exclusion
2523 // 5 Electrostatic calculation method (0= Direct, ...)
2524 new_field -> elec_opts[5] = 2.0;
2525 // Ewald precision:
2526 new_field -> elec_opts[6] = 1e-20;
2527 // 7 to 9 method options
2528 // 10 Evaluation every N steps, N=
2529 tmp_field -> elec_opts[10] = 1.0;
2530
2531 /*
2532 new_field -> met_opts[]
2533 */
2534 // 0 (1/0) Enable direct calculation of metal interactions
2535 // 1 (1/0) Switch the default embeding fuctions
2536
2537 // Ensemble (0= NVE, 1= NVT, 2= NPT, 4= NST)
2538 // new_field -> ensemble
2539 // Thermostat = f(ensemble):
2540 // new_field -> thermostat
2541 /*
2542 new_field -> thermo_opts[]
2543 */
2544 // 0-5 ensemble + thermostat related options
2545 // 6 (1/0) Use pseudo thermal bath
2546 // 7 Thermostat type (0= L+D, 1= L, 2= G, 3= D)
2547 // Thickness:
2548 new_field -> thermo_opts[8] = 2.0;
2549 // 9 = Target temperature
2550
2551 /*
2552 new_field -> md_opts[]
2553 */
2554 // MD target temperature (K)
2555 new_field -> md_opts[0] = 300.0;
2556 // 1 Verlet integrator (0 = velocity, 1 = leapfrog)
2557 // Number of MD steps:
2558 new_field -> md_opts[2] = 100.0;
2559 // 3 Time step (0 = fixed, 1 = variable)
2560 // Time step delta(t)= (ps)
2561 new_field -> md_opts[4] = 0.001;
2562 // Target pressure (katms):
2563 new_field -> md_opts[5] = 0.001;
2564 // Variable time step, maximum time step:
2565 new_field -> md_opts[6] = 0.001;
2566 // Variable time step, maximum distance allowed:
2567 new_field -> md_opts[7] = 0.5;
2568 // Variable time step, minimun distance allowed:
2569 new_field -> md_opts[8] = 0.01;
2570 // Shake iterations limit, in cycles:
2571 new_field -> md_opts[9] = 250.0;
2572 // Shake tolerance:
2573 new_field -> md_opts[10] = 1e-6;
2574 // FIQA iterations limit, in cycles:
2575 new_field -> md_opts[11] = 100;
2576 // FIQA quaternions tolerance:
2577 new_field -> md_opts[12] = 1e-8;
2578 // 13 (1/0) Initiate impact on particles
2579 // Index if the particle to impact:
2580 new_field -> md_opts[14] = 1.0;
2581 // Time step of the impact:
2582 new_field -> md_opts[15] = 1.0;
2583 // Energy of the imapct (keV):
2584 new_field -> md_opts[16] = 1.0;
2585 // Direction of the impact (from center of mass) ...
2586 // ... on x:
2587 new_field -> md_opts[17] = 1.0;
2588 // ... on y:
2589 new_field -> md_opts[18] = 1.0;
2590 // ... on z:
2591 new_field -> md_opts[19] = 1.0;
2592
2593 /*
2594 new_field -> out_opts[]
2595 */
2596 // 0 (1/0) Print Defect file ...
2597 // 1 ... starting from step
2598 // ... every N steps, N=
2599 new_field -> out_opts[2] = 1.0;
2600 // ... with a site inteerstitial cutoff:
2601 new_field -> out_opts[3] = 0.75;
2602 // 4 (1/0) Print Displacement file ..
2603 // 5 ... starting from step
2604 // ... every N steps, N=
2605 new_field -> out_opts[6] = 1.0;
2606 // ... with a qualifying cutoff:
2607 new_field -> out_opts[7] = 0.15;
2608 // 8 (1/0) Print HISTORY file ...
2609 // 9 ... starting from step
2610 // ... every N steps, N=
2611 new_field -> out_opts[10] = 1.0;
2612 // 11 Data level (0, 1, 2)
2613 // 12 (1/0) Print Atom's MSD ...
2614 // 13 ... starting from step
2615 // ... every N steps, N=
2616 new_field -> out_opts[14] = 1.0;
2617 // 15 (1/0) Print system data ...
2618 // ... every N steps, N=
2619 new_field -> out_opts[16] = 1.0;
2620 // 17 (1/0) Accumulate statistics ...
2621 // ... every N steps, N=
2622 new_field -> out_opts[18] = 1.0;
2623 // 19 (1/0) Rolling average stack ...
2624 // ... to N steps, N=
2625 new_field -> out_opts[20] = 1.0;
2626 // 21 (1/0) RDF ...
2627 // ... every N steps, N=
2628 new_field -> out_opts[22] = 1.0;
2629 // ... with a bin size=
2630 new_field -> out_opts[23] = 0.05;
2631 // 24 (1/0) Overall VAF analysis ...
2632 // ... every N steps, N=
2633 new_field -> out_opts[25] = 1.0;
2634 // ... with a bin size=
2635 new_field -> out_opts[26] = 50.0;
2636 // 27 (1/0) Z density profile ...
2637 // ... every N steps, N=
2638 new_field -> out_opts[28] = 1.0;
2639 // 29 (1/0) ignoring time-av. VAFs
2640 // Dump restart information every N steps, N=
2641 new_field -> out_opts[30] = 1000.0;
2642
2643 /*
2644 new_field -> ana_opts[]
2645 */
2646 // 0 (1/0) Activate all intra-molecular PDF analysis ...
2647 // ... every N steps, N=
2648 new_field -> ana_opts[1] = 1.0;
2649 // ... number of delta r in:
2650 new_field -> ana_opts[2] = 100.0;
2651 // ... Cutoff:
2652 new_field -> ana_opts[3] = 2.0;
2653 // 4 (1/0) Activate bonds PDF analysis ...
2654 // ... every N steps, N=
2655 new_field -> ana_opts[5] = 1.0;
2656 // ... number of delta r in:
2657 new_field -> ana_opts[6] = 100.0;
2658 // ... Cutoff:
2659 new_field -> ana_opts[7] = 2.0;
2660 // 8 (1/0) Activate angles analysis ...
2661 // ... every N steps, N=
2662 new_field -> ana_opts[9] = 1.0;
2663 // ... number of delta angle in:
2664 new_field -> ana_opts[10] = 90.0;
2665 // 11 (1/0) Activate dihedrals analysis ...
2666 // ... every N steps, N=
2667 new_field -> ana_opts[12] = 1.0;
2668 // ... number of delta angle in:
2669 new_field -> ana_opts[13] = 90.0;
2670 // 14 (1/0) Activate inversions analysis ...
2671 // ... every N steps, N=
2672 new_field -> ana_opts[15] = 1.0;
2673 // ... number of delta angle in:
2674 new_field -> ana_opts[16] = 90.0;
2675
2676}
2677
2685void init_lammps_field (classical_field * new_field)
2686{
2687 // Preparing LAMMPS Control file
2688 new_field -> sys_opts = allocdouble (17);
2689 new_field -> io_opts = allocdouble (23);
2690 new_field -> ana_opts = allocdouble (17);
2691 new_field -> elec_opts = allocdouble (11);
2692 new_field -> vdw_opts = allocdouble (6);
2693 new_field -> met_opts = allocdouble (2);
2694 new_field -> equi_opts = allocdouble (17);
2695 new_field -> thermo_opts= allocdouble (10);
2696 new_field -> md_opts = allocdouble (20);
2697 new_field -> out_opts = allocdouble (31);
2698/*
2699 new_field -> io_opts[]
2700 */
2701 // Set job time control (yes/no)
2702 new_field -> io_opts[0] = 1.0;
2703 // Job time control (s):
2704 new_field -> io_opts[1] = 10000.0;
2705 // Set job closure time (yes/no)
2706 new_field -> io_opts[2] = 1.0;
2707 // Job closure time (s):
2708 new_field -> io_opts[3] = 100.0;
2709 // 4 (1/0) General I/O read interace
2710 // 5 General I/O method (0= MPI, 1= FORTRAN, 2= Traditional, 3= netCDF)
2711 // Reader count:
2712 new_field -> io_opts[6] = 12;
2713 // Batch size:
2714 new_field -> io_opts[7] = 2000000;
2715 // Buffer size:
2716 new_field -> io_opts[8] = 20000;
2717 // 9 (1/0) Error check on general read I/O
2718 // 10 (1/0) General I/O write interace
2719 // 11 General I/O method (0= MPI, 1= FORTRAN, 2= Traditional, 3= netCDF)
2720 // Real precision (0 = 32 bits, 1= 64 bits)
2721 new_field -> io_opts[12] = 1;
2722 // 13 Type (0= Sorted, 1= Unsorted)
2723 // Reader count:
2724 new_field -> io_opts[14] = 12;
2725 // Batch size:
2726 new_field -> io_opts[15] = 2000000;
2727 // Buffer size:
2728 new_field -> io_opts[16] = 20000;
2729 // 17 (1/0) Error check on general write I/O
2730 // 18 (1/0) Define random number generator seeds
2731 // n1 seed:
2732 new_field -> io_opts[19] = 1;
2733 // n2 seed:
2734 new_field -> io_opts[20] = 20;
2735 // n3 seed:
2736 new_field -> io_opts[21] = 300;
2737 // 22 (1/0) limit the number of processors in the z direction
2738
2739 /*
2740 new_field -> sys_opts[]
2741 */
2742 // Relative dieletric constant:
2743 new_field -> sys_opts[0] = 1.0;
2744 // 1 = allowed local variation of the system density (%)
2745 // 2 (1/0) Ingore the particle indices in CONFIG file
2746 // 3 (1/0) Ignore strict checks
2747 // 4 (1/0) Skip detailed topology reporting
2748 // 5 (1/0) Ignore center of mass removal
2749 // 6 (1/0) Set tolerance of the relaxed shell model
2750 // Tolerance for the relaxed shell model:
2751 new_field -> sys_opts[7] = 1.0;
2752 // 8 (1/0) Set the subscelling threshold density of particles per link cell
2753 // Subscelling threshold density of particles per link cell:
2754 new_field -> sys_opts[9] = 50.0;
2755 // 10 (1/0) Create an expanded version of the current model ...
2756 // ... on x:
2757 new_field -> sys_opts[11] = 1.0;
2758 // ... on y:
2759 new_field -> sys_opts[12] = 1.0;
2760 // .. on z:
2761 new_field -> sys_opts[13] = 1.0;
2762 // 14 (1/0) Restart a calculation
2763 // 15 (0= Continue current, 1= Start new without T reset, 2= Start new with T reset)
2764
2765 /*
2766 new_field -> equi_opts[]
2767 */
2768 new_field -> equi_opts[0] = 1.0;
2769 // Equilibrate for the first N steps, N=
2770 new_field -> equi_opts[1] = 20.0;
2771 // Scale temperature every N steps, N=
2772 new_field -> equi_opts[3] = 1.0;
2773 // Cap forces fmax=
2774 new_field -> equi_opts[5] = 1000.0;
2775 // Resample the momenta distribution during N steps N=
2776 new_field -> equi_opts[7] = 1.0;
2777 // 8 (1/0) = Minimize system configuration
2778 // 9 Minimize target (0=Force, 1=Energy, 2= Distance)
2779 // Minimization threshold:
2780 new_field -> equi_opts[10] = 50.0;
2781 // Minimization to be done every N steps, N=
2782 new_field -> equi_opts[11] = 1.0;
2783 // 12 (1/0) = Optimize system configuration from start
2784 // 13 Minimize target (0=Force, 1=Energy, 2= Distance)
2785 // Minimization threshold:
2786 new_field -> equi_opts[14] = 50.0;
2787 // 15 (1/0) = Perform zero temperature (10K) optimization
2788 // 16 (1/0) = Include equilibration data in statistics
2789
2790 /*
2791 new_field -> vdw_opts[]
2792 */
2793 // Vdw cutoff
2794 // 0 (1/0) Use VdW interactions
2795 new_field -> vdw_opts[0] = 1.0;
2796 // Vdw Cutoff radius
2797 new_field -> vdw_opts[1] = 12.0;
2798 // 2 (1/0) Enforce direct VdW calculation
2799 // 3 (1/0) Apply force shifting
2800 // 4 (1/0) Apply mixing rule
2801 // 5 Mixing rule (0= LB, 1= FH, 2= H, 3= HHG, 4= TT, 5= F)
2802
2803 /*
2804 new_field -> elec_opts[]
2805 */
2806 // 0 (1/0) Evaluate electrostatics
2807 new_field -> elec_opts[0] = 1.0;
2808 // Long range cutoff:
2809 new_field -> elec_opts[1] = 12.0;
2810 // 2 (1/0) Add padding to rcut
2811 // 3 Padding in angstrom
2812 // 4 (1/0) Use extended coulombic exclusion
2813 // 5 Electrostatic calculation method (0= Direct, ...)
2814 new_field -> elec_opts[5] = 2.0;
2815 // Ewald precision:
2816 new_field -> elec_opts[6] = 1e-20;
2817 // 7 to 9 method options
2818 // 10 Evaluation every N steps, N=
2819 tmp_field -> elec_opts[10] = 1.0;
2820
2821 /*
2822 new_field -> met_opts[]
2823 */
2824 // 0 (1/0) Enable direct calculation of metal interactions
2825 // 1 (1/0) Switch the default embeding fuctions
2826
2827 // Ensemble (0= NVE, 1= NVT, 2= NPT, 4= NST)
2828 // new_field -> ensemble
2829 // Thermostat = f(ensemble):
2830 // new_field -> thermostat
2831 /*
2832 new_field -> thermo_opts[]
2833 */
2834 // 0-5 ensemble + thermostat related options
2835 // 6 (1/0) Use pseudo thermal bath
2836 // 7 Thermostat type (0= L+D, 1= L, 2= G, 3= D)
2837 // Thickness:
2838 new_field -> thermo_opts[8] = 2.0;
2839 // 9 = Target temperature
2840
2841 /*
2842 new_field -> md_opts[]
2843 */
2844 // MD target temperature (K)
2845 new_field -> md_opts[0] = 300.0;
2846 // 1 Verlet integrator (0 = velocity, 1 = leapfrog)
2847 // Number of MD steps:
2848 new_field -> md_opts[2] = 100.0;
2849 // 3 Time step (0 = fixed, 1 = variable)
2850 // Time step delta(t)= (ps)
2851 new_field -> md_opts[4] = 0.001;
2852 // Target pressure (katms):
2853 new_field -> md_opts[5] = 0.001;
2854 // Variable time step, maximum time step:
2855 new_field -> md_opts[6] = 0.001;
2856 // Variable time step, maximum distance allowed:
2857 new_field -> md_opts[7] = 0.5;
2858 // Variable time step, minimun distance allowed:
2859 new_field -> md_opts[8] = 0.01;
2860 // Shake iterations limit, in cycles:
2861 new_field -> md_opts[9] = 250.0;
2862 // Shake tolerance:
2863 new_field -> md_opts[10] = 1e-6;
2864 // FIQA iterations limit, in cycles:
2865 new_field -> md_opts[11] = 100;
2866 // FIQA quaternions tolerance:
2867 new_field -> md_opts[12] = 1e-8;
2868 // 13 (1/0) Initiate impact on particles
2869 // Index if the particle to impact:
2870 new_field -> md_opts[14] = 1.0;
2871 // Time step of the impact:
2872 new_field -> md_opts[15] = 1.0;
2873 // Energy of the imapct (keV):
2874 new_field -> md_opts[16] = 1.0;
2875 // Direction of the impact (from center of mass) ...
2876 // ... on x:
2877 new_field -> md_opts[17] = 1.0;
2878 // ... on y:
2879 new_field -> md_opts[18] = 1.0;
2880 // ... on z:
2881 new_field -> md_opts[19] = 1.0;
2882
2883 /*
2884 new_field -> out_opts[]
2885 */
2886 // 0 (1/0) Print Defect file ...
2887 // 1 ... starting from step
2888 // ... every N steps, N=
2889 new_field -> out_opts[2] = 1.0;
2890 // ... with a site inteerstitial cutoff:
2891 new_field -> out_opts[3] = 0.75;
2892 // 4 (1/0) Print Displacement file ..
2893 // 5 ... starting from step
2894 // ... every N steps, N=
2895 new_field -> out_opts[6] = 1.0;
2896 // ... with a qualifying cutoff:
2897 new_field -> out_opts[7] = 0.15;
2898 // 8 (1/0) Print HISTORY file ...
2899 // 9 ... starting from step
2900 // ... every N steps, N=
2901 new_field -> out_opts[10] = 1.0;
2902 // 11 Data level (0, 1, 2)
2903 // 12 (1/0) Print Atom's MSD ...
2904 // 13 ... starting from step
2905 // ... every N steps, N=
2906 new_field -> out_opts[14] = 1.0;
2907 // 15 (1/0) Print system data ...
2908 // ... every N steps, N=
2909 new_field -> out_opts[16] = 1.0;
2910 // 17 (1/0) Accumulate statistics ...
2911 // ... every N steps, N=
2912 new_field -> out_opts[18] = 1.0;
2913 // 19 (1/0) Rolling average stack ...
2914 // ... to N steps, N=
2915 new_field -> out_opts[20] = 1.0;
2916 // 21 (1/0) RDF ...
2917 // ... every N steps, N=
2918 new_field -> out_opts[22] = 1.0;
2919 // ... with a bin size=
2920 new_field -> out_opts[23] = 0.05;
2921 // 24 (1/0) Overall VAF analysis ...
2922 // ... every N steps, N=
2923 new_field -> out_opts[25] = 1.0;
2924 // ... with a bin size=
2925 new_field -> out_opts[26] = 50.0;
2926 // 27 (1/0) Z density profile ...
2927 // ... every N steps, N=
2928 new_field -> out_opts[28] = 1.0;
2929 // 29 (1/0) ignoring time-av. VAFs
2930 // Dump restart information every N steps, N=
2931 new_field -> out_opts[30] = 1000.0;
2932
2933 /*
2934 new_field -> ana_opts[]
2935 */
2936 // 0 (1/0) Activate all intra-molecular PDF analysis ...
2937 // ... every N steps, N=
2938 new_field -> ana_opts[1] = 1.0;
2939 // ... number of delta r in:
2940 new_field -> ana_opts[2] = 100.0;
2941 // ... Cutoff:
2942 new_field -> ana_opts[3] = 2.0;
2943 // 4 (1/0) Activate bonds PDF analysis ...
2944 // ... every N steps, N=
2945 new_field -> ana_opts[5] = 1.0;
2946 // ... number of delta r in:
2947 new_field -> ana_opts[6] = 100.0;
2948 // ... Cutoff:
2949 new_field -> ana_opts[7] = 2.0;
2950 // 8 (1/0) Activate angles analysis ...
2951 // ... every N steps, N=
2952 new_field -> ana_opts[9] = 1.0;
2953 // ... number of delta angle in:
2954 new_field -> ana_opts[10] = 90.0;
2955 // 11 (1/0) Activate dihedrals analysis ...
2956 // ... every N steps, N=
2957 new_field -> ana_opts[12] = 1.0;
2958 // ... number of delta angle in:
2959 new_field -> ana_opts[13] = 90.0;
2960 // 14 (1/0) Activate inversions analysis ...
2961 // ... every N steps, N=
2962 new_field -> ana_opts[15] = 1.0;
2963 // ... number of delta angle in:
2964 new_field -> ana_opts[16] = 90.0;
2965}
2966
2974classical_field * create_force_field_data_structure (int ai)
2975{
2976 int i;
2977 classical_field * new_field;
2978 new_field = g_malloc0 (sizeof*new_field);
2979 tmp_field = new_field;
2980 new_field -> energy_unit = 1;
2981 new_field -> atom_init = ai;
2982 switch (activef)
2983 {
2984 case 0:
2985 init_dlpoly_field (new_field);
2986 break;
2987 case 1:
2988 init_lammps_field (new_field);
2989 break;
2990 }
2991 new_field -> molecules = tmp_proj -> modelfc -> mol_by_step[0];
2992 new_field -> first_molecule = g_malloc0 (sizeof*new_field -> first_molecule);
2993 tmp_fmol = new_field -> first_molecule;
2994#ifdef DEBUG
2995 g_debug ("Total num of molecules: %d", tmp_proj -> modelfc -> mol_by_step[0]);
2996#endif // DEBUG
2997 for (i=0; i < tmp_proj -> modelfc -> mol_by_step[0]; i++)
2998 {
2999 tmp_mol = & tmp_proj -> modelfc -> mols[0][i];
3000#ifdef DEBUG
3001 g_debug ("Init FORCE FIELD for: ");
3002 g_debug (" tmp_mol -> id= %d", tmp_mol -> id);
3003 g_debug (" tmp_mol -> multi= %d", tmp_mol -> multiplicity);
3004 g_debug (" tmp_mol -> natoms= %d", tmp_mol -> natoms);
3005#endif // DEBUG
3006 init_field_molecule (i);
3007 if (i < tmp_proj -> modelfc -> mol_by_step[0]-1)
3008 {
3009 tmp_fmol -> next = g_malloc0 (sizeof*tmp_fmol -> next);
3010 tmp_fmol = tmp_fmol -> next;
3011 }
3012 }
3013
3014 for (i=0; i<5; i++)
3015 {
3016 tmp_field -> nbody[i] = 0;
3017 tmp_field -> first_body[i] = NULL;
3018 }
3019 // Tersoff cross terms
3020 tmp_field -> cross = NULL;
3021
3022
3023 // vdW start
3024 tmp_field -> nbody[0] = init_vdw (TRUE);
3025 // vdW end
3026
3027 // External field(s) start
3028 new_field -> extern_fields = 0;
3029 new_field -> first_external = NULL;
3030 // External field(s) end
3031
3032 // Non-bonded end
3033 return new_field;
3034}
limit
float limit[2]
Definition atom_edit.c:59
mol
insertion_menu mol[]
Definition w_library.c:193
search_type
int search_type
Definition calc_menu.c:92
chemistry
integer(kind=c_int) function chemistry()
Definition chemistry.F90:22
label
void label(cairo_t *cr, double val, int axe, int p, project *this_proj)
draw axis label
Definition labels.c:56
atoms
int atoms[NUM_STYLES][2]
Definition d_selection.c:361
get_active_atom
field_atom * get_active_atom(int a, int b)
retrieve field atom
Definition dlp_active.c:145
comp_fbody
field_nth_body * comp_fbody
Definition dlp_comp.c:44
sys_opts
gchar * sys_opts[10]
Definition dlp_control.c:2367
duplicate_field_struct_list
field_struct * duplicate_field_struct_list(field_struct *list_str, gboolean init)
create copy of list of field structural element(s)
Definition dlp_copy.c:325
atoms_id
int * atoms_id
Definition force_fields.c:238
cross
double *** cross
Definition dlp_edit.c:418
fvalues
int fvalues[2][15][21]
Definition dlp_field.c:255
tmp_fbody
field_nth_body * tmp_fbody
Definition dlp_field.c:965
tmp_proj
project * tmp_proj
Definition dlp_field.c:953
body_at
int body_at(int b)
find the number of atom(s) in a non bonded interaction
Definition dlp_field.c:1022
tmp_fat
field_atom * tmp_fat
Definition dlp_field.c:957
struct_id
int struct_id(int f)
number of atoms in a structural element
Definition dlp_field.c:999
tmp_fstr
field_struct * tmp_fstr
Definition dlp_field.c:964
tmp_fmol
field_molecule * tmp_fmol
Definition dlp_field.c:955
atom_init
int atom_init
Definition dlp_field.c:950
tmp_fdt
field_atom * tmp_fdt
Definition dlp_field.c:957
tmp_fct
field_atom * tmp_fct
Definition dlp_field.c:957
tmp_fbt
field_atom * tmp_fbt
Definition dlp_field.c:957
tmp_coord
coord_info * tmp_coord
Definition dlp_field.c:954
tmp_mol
molecule * tmp_mol
Definition dlp_field.c:956
tmp_field
classical_field * tmp_field
Definition dlp_field.c:951
dlp_field.h
Variable declarations for the creation of the DL_POLY input file(s)
dihedral_3d
angle dihedral_3d(cell_info *cell, int mdstep, atom *at, atom *bt, atom *ct, atom *dt)
dihedral between atom a, b, c and d in 3D
Definition ogl_utils.c:204
inversion_3d
angle inversion_3d(cell_info *cell, int mdstep, atom *at, atom *bt, atom *ct, atom *dt)
inversion angle between atom a, b, c and d in 3D
Definition ogl_utils.c:240
angle_3d
angle angle_3d(cell_info *cell, int mdstep, atom *at, atom *bt, atom *ct)
angle between atom a, b and c in 3D
Definition ogl_utils.c:179
TOTAL_COORD_AND_SPEC
@ TOTAL_COORD_AND_SPEC
Definition dlp_field.h:146
SPEC_ONLY
@ SPEC_ONLY
Definition dlp_field.h:144
PARTIAL_COORD_AND_SPEC
@ PARTIAL_COORD_AND_SPEC
Definition dlp_field.h:148
OTHER
@ OTHER
Definition dlp_field.h:149
distance_3d
distance distance_3d(cell_info *cell, int mdstep, atom *at, atom *bt)
distance between atom a and b in 3D
Definition ogl_utils.c:81
angles_from_bonds
int angles_from_bonds(int n, field_atom *at, field_atom *bt, field_atom *ct)
find, and initialize, angles using bonds
Definition dlp_init.c:1124
init_field_shell
field_shell * init_field_shell(int id, int ia, int ib)
initialize field core shell interaction
Definition dlp_init.c:394
init_field_tethered
field_tethered * init_field_tethered(int id, int num)
initialize new field tethered potential
Definition dlp_init.c:510
id_atom
int id_atom
Definition dlp_init.c:120
clean_field_struct_list
void clean_field_struct_list(field_struct *stru)
clean force field structural property(ies)
Definition dlp_init.c:2134
find_neighbor_loop
int find_neighbor_loop(int frag, int aid, int fai, field_atom *fat, field_neighbor *ngmb, int sid, gboolean save_it)
find / set up field atom data for a fragment
Definition dlp_init.c:1579
init_ngb
field_neighbor * init_ngb
Definition dlp_init.c:1462
init_dlpoly_field
void init_dlpoly_field(classical_field *new_field)
initialize DL-POLY classical force field
Definition dlp_init.c:2395
init_field_struct
field_struct * init_field_struct(int st, int ai, int an, int *aid)
initialize field molecule new structural element
Definition dlp_init.c:704
duplicate_int
int * duplicate_int(int num, int *old_val)
copy a list of int
Definition global.c:560
prepare_field_atom
int prepare_field_atom(int i, int j, int k, int l, int m)
initialiaze new type of field atom
Definition dlp_init.c:814
test_for_angles
int test_for_angles(field_atom *at, field_atom *bt, field_atom *ct)
search for angle(s) between these field atoms
Definition dlp_init.c:1079
assigned
int assigned
Definition dlp_init.c:1459
val
float val
Definition dlp_init.c:117
get_struct_id_from_atom_id
int get_struct_id_from_atom_id(int ids, int *aid)
retrieve field structural element id from a list of atom id
Definition dlp_init.c:559
get_atom_id_in_fragment_from_model_id
int get_atom_id_in_fragment_from_model_id(int frag, int at)
retrieve field atom id in fragment from model id
Definition dlp_init.c:197
id_n_fold_atoms_in_fragment
gboolean id_n_fold_atoms_in_fragment(int frag, int limit, int num_ngb, int search_type, int init, field_neighbor *ngma_init)
find atom id in fragment using the neighbor list
Definition dlp_init.c:1626
get_init_neighbor
field_neighbor * get_init_neighbor(int a)
retrieve neighbor data structure for atom a
Definition dlp_init.c:1471
print_all_field_struct
void print_all_field_struct(field_molecule *mol, int str)
print all field structural element(s)
Definition dlp_print.c:195
set_atom_id
int set_atom_id(field_atom *at, int c, int p, int id)
set new fragment/molecule atom parameters
Definition dlp_init.c:224
coord_sphere_multiplicity
int coord_sphere_multiplicity(atom *at, int id, gboolean set_atom_id)
get coord sphere multiplicity in the model
Definition dlp_init.c:1408
get_weight
void get_weight(int seq)
compute field weigth for the atom(s) in the model
Definition dlp_init.c:1820
init_all_impropers_inversions
void init_all_impropers_inversions(int stru)
initialize all impropers and inversions
Definition dlp_init.c:1331
set_field_atom_name
gchar * set_field_atom_name(field_atom *ato, field_molecule *mol)
get name string for field atom
Definition dlp_init.c:247
init_field_prop
field_prop * init_field_prop(int ti, int key, gboolean show, gboolean use)
initialize new field molecule structural property
Definition dlp_init.c:535
get_field_atom_id_from_model_id
int get_field_atom_id_from_model_id(field_molecule *fmol, int at)
retrieve field atom id using atom id in model
Definition dlp_init.c:151
init_field_molecule
void init_field_molecule(int i)
initialiaze of field molecule data structurre
Definition dlp_init.c:2227
multi
int multi
Definition dlp_init.c:121
are_in_bond
gboolean are_in_bond(atom ato, int at)
is at in ato neighbors ?
Definition dlp_init.c:930
was_not_created_struct
gboolean was_not_created_struct(int ids, int num, int *aid)
was this structural element already created ?
Definition dlp_init.c:603
astr
int * astr
Definition dlp_init.c:116
a_multi
int a_multi
Definition dlp_init.c:122
test_for_bonds
int test_for_bonds(field_atom *at, field_atom *bt)
search for bond(s) between 2 field atoms
Definition dlp_init.c:948
get_fragment_atom_id_from_model_id
int get_fragment_atom_id_from_model_id(field_molecule *fmol, int at)
retrieve fragment atom id using atom id in model
Definition dlp_init.c:174
i_weight
int * i_weight
Definition dlp_init.c:1460
bonds_between_atoms
int bonds_between_atoms(int n, field_atom *at, field_atom *bt, int a, int b)
search for bond(s) between two field atoms
Definition dlp_init.c:1015
setup_atomic_weight
int setup_atomic_weight(int seq)
setup atom(s) field weight
Definition dlp_init.c:1865
init_all_dihedrals
void init_all_dihedrals()
find, and intialiaze, all dihedral(s) using bonds
Definition dlp_init.c:1247
dihedrals_from_angles
int dihedrals_from_angles(int n, field_atom *at, field_atom *bt, field_atom *ct, field_atom *dt)
find, and initialize, dihedrals using angles
Definition dlp_init.c:1228
impropers_inversion
int impropers_inversion(int n, int stru, int at, int bt, int ct, int dt, int a, int b, int c, int d)
prepare new improper / inversion
Definition dlp_init.c:1289
init_field_pmf
field_pmf * init_field_pmf(int id, int num[2], int *list[2], float *w[2])
initialize new field mean force potential
Definition dlp_init.c:445
duplicate_float
float * duplicate_float(int num, float *old_val)
copy a list of float
Definition global.c:592
test_for_dihedrals
int test_for_dihedrals(field_atom *at, field_atom *bt, field_atom *ct, field_atom *dt)
search for dihedral(s) between these field atoms
Definition dlp_init.c:1173
create_force_field_data_structure
classical_field * create_force_field_data_structure(int ai)
initialize classical force field
Definition dlp_init.c:2974
n_weight
int * n_weight
Definition dlp_init.c:1461
get_position_in_field_atom_from_model_id
int get_position_in_field_atom_from_model_id(int fat, int at)
retrieve atom position id in field atom id using atom id in model
Definition dlp_init.c:132
prepare_field_struct
int prepare_field_struct(int ids, int sid, int yes_no_num, int *aid)
prepare the creation of a field structural element
Definition dlp_init.c:981
init_field_external
field_external * init_field_external(int bi)
initialize new field external potential
Definition dlp_init.c:790
find_atom_id_in_field_molecule
void find_atom_id_in_field_molecule()
find field atom id in field molecule (fragment multiplicity > 1)
Definition dlp_init.c:1960
are_neighbors
gboolean are_neighbors(field_neighbor *ngb, int at)
test if 2 atoms are neighbors
Definition dlp_init.c:912
init_all_atoms
void init_all_atoms(int i)
initialize all field atom(s)
Definition dlp_init.c:838
id_atoms_in_fragment
gboolean id_atoms_in_fragment(int frag, int seed)
find the atom field id in the fragment
Definition dlp_init.c:1777
in_bond
gboolean in_bond(int at, int bd[2])
is atom at in bond bd
Definition dlp_init.c:892
init_field_constraint
field_constraint * init_field_constraint(int id, int ia, int ib)
initialize field bond constraint
Definition dlp_init.c:419
init_all_angles
void init_all_angles()
find, and intialiaze, all angle(s) using bonds
Definition dlp_init.c:1140
init_field_nth_body
field_nth_body * init_field_nth_body(int bi, int bd, int *na, int **ma, int **ba)
initialize new field non bonded interaction
Definition dlp_init.c:740
init_field_atom
field_atom * init_field_atom(int id, int type, int at, int nat, int coo, int *list)
initialize a new type of field atom
Definition dlp_init.c:304
init_all_bonds
void init_all_bonds()
find, and initialize all bond(s)
Definition dlp_init.c:1031
init_field_rigid
field_rigid * init_field_rigid(int id, int num, int *list)
initialize new field rigid constraint
Definition dlp_init.c:487
init_lammps_field
void init_lammps_field(classical_field *new_field)
initialize a LAMMPS classical force field
Definition dlp_init.c:2685
init_all_field_struct
void init_all_field_struct(gboolean init)
initialize the force field structural properties for a field molecule
Definition dlp_init.c:2158
is_this_numbering_possible_for_this_atom
gboolean is_this_numbering_possible_for_this_atom(int frag, field_neighbor *ngma, int atom)
check if the id in fragment for this atom is possible or not
Definition dlp_init.c:1538
is_numbering_possible
gboolean is_numbering_possible(int frag)
compare if the atom numbering of two fragment(s) are identical, and it should
Definition dlp_init.c:1490
init_vdw
int init_vdw(gboolean init)
initialize the VdW section of a classical force field
Definition dlp_init.c:2287
allocdint
int ** allocdint(int xal, int yal)
allocate an int ** pointer
Definition global.c:330
allocbool
gboolean * allocbool(int val)
allocate a gboolean * pointer
Definition global.c:254
activef
int activef
Definition global.c:161
allocdouble
double * allocdouble(int val)
allocate a double * pointer
Definition global.c:459
allocint
int * allocint(int val)
allocate an int * pointer
Definition global.c:314
allocfloat
float * allocfloat(int val)
allocate a float * pointer
Definition global.c:398
CHEM_M
#define CHEM_M
Definition global.h:300
exact_name
gchar * exact_name(gchar *name)
short cut to print string without spaces
Definition interface.c:370
interface.h
Messaging function declarations.
molecules
integer(kind=c_int) function molecules(frag_and_mol, allbonds)
Definition molecules.F90:155
angle::angle
double angle
Definition glwin.h:116
atom
Definition global.h:891
atom::numv
int numv
Definition global.h:898
atom::vois
int * vois
Definition global.h:899
classical_field
Definition global.h:707
distance::length
double length
Definition glwin.h:124
field_atom
Definition dlp_field.h:168
field_constraint
Definition dlp_field.h:210
field_external
Definition dlp_field.h:299
field_molecule
Definition dlp_field.h:321
field_neighbor
Definition dlp_field.h:311
field_nth_body
Definition dlp_field.h:281
field_pmf
Definition dlp_field.h:223
field_prop
Definition dlp_field.h:154
field_rigid
Definition dlp_field.h:238
field_shell
Definition dlp_field.h:192
field_struct
Definition dlp_field.h:264
field_tethered
Definition dlp_field.h:250
project::coord
coord_info * coord
Definition global.h:961
b
int b
Definition tab-1.c:95
c
int c
Definition tab-1.c:95
d
int d
Definition tab-1.c:95
a
int a
Definition tab-1.c:95
res
GtkWidget * res[2]
Definition w_encode.c:212