Branch data Line data Source code
1 : : /*
2 : : * International Chemical Identifier (InChI)
3 : : * Version 1
4 : : * Software version 1.07
5 : : * April 30, 2024
6 : : *
7 : : * MIT License
8 : : *
9 : : * Copyright (c) 2024 IUPAC and InChI Trust
10 : : *
11 : : * Permission is hereby granted, free of charge, to any person obtaining a copy
12 : : * of this software and associated documentation files (the "Software"), to deal
13 : : * in the Software without restriction, including without limitation the rights
14 : : * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
15 : : * copies of the Software, and to permit persons to whom the Software is
16 : : * furnished to do so, subject to the following conditions:
17 : : *
18 : : * The above copyright notice and this permission notice shall be included in all
19 : : * copies or substantial portions of the Software.
20 : : *
21 : : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
22 : : * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
24 : : * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
26 : : * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
27 : : * SOFTWARE.
28 : : *
29 : : * The InChI library and programs are free software developed under the
30 : : * auspices of the International Union of Pure and Applied Chemistry (IUPAC).
31 : : * Originally developed at NIST.
32 : : * Modifications and additions by IUPAC and the InChI Trust.
33 : : * Some portions of code were developed/changed by external contributors
34 : : * (either contractor or volunteer) which are listed in the file
35 : : * 'External-contributors' included in this distribution.
36 : : *
37 : : * info@inchi-trust.org
38 : : *
39 : : */
40 : :
41 : : #include <string.h>
42 : :
43 : : #include "mode.h"
44 : : #include "ichicomn.h"
45 : : #include "ichicant.h"
46 : :
47 : : #include "bcf_s.h"
48 : :
49 : : #if 0
50 : : #define RET_MAX 32767
51 : : #define CUTOFF 8 /* testing shows that this is good value */
52 : : #endif
53 : :
54 : : /* Note: the theoretical number of stack entries required is
55 : : no more than 1 + log2(num). But we switch to insertion
56 : : sort for CUTOFF elements or less, so we really only need
57 : : 1 + log2(num) - log2(CUTOFF) stack entries. For a CUTOFF
58 : : of 8, that means we need no more than 30 stack entries for
59 : : 32 bit platforms, and 62 for 64-bit platforms. */
60 : : #define STKSIZ (8*sizeof(void*) - 2)
61 : :
62 : :
63 : : /****************************************************************************
64 : : inchi's qsort
65 : : ****************************************************************************/
66 : 338 : void inchi_qsort( void *pParam,
67 : : void *base,
68 : : size_t num,
69 : : size_t width,
70 : : int( *comp )( const void *, const void *, void * ) )
71 : : {
72 : : char *lo, *hi; /* ends of sub-array currently sorting */
73 : : char *mid; /* points to middle of subarray */
74 : : char *loguy, *higuy; /* traveling pointers for partition step */
75 : : size_t size; /* size of the sub-array */
76 : : char *lostk[STKSIZ], *histk[STKSIZ];
77 : : int stkptr; /* stack for saving sub-array to be processed */
78 : :
79 [ + + ]: 338 : if (num < 2)
80 : 26 : return; /* nothing to do */
81 : :
82 : 312 : stkptr = 0; /* initialize stack */
83 : :
84 : 312 : lo = (char *) base;
85 : 312 : hi = (char *) base + width * ( num - 1 ); /* initialize limits */
86 : :
87 : : /* This entry point is for pseudo-recursion calling: setting
88 : : lo and hi and jumping to here is like recursion, but stkptr is
89 : : preserved, locals aren't, so we preserve stuff on the stack */
90 : 1347 : recurse:
91 : :
92 : 1659 : size = ( hi - lo ) / width + 1; /* number of el's to sort */
93 : :
94 : : /* First we pick a partitioning element. The efficiency of the
95 : : algorithm demands that we find one that is approximately the median
96 : : of the values, but also that we select one fast. We choose the
97 : : median of the first, middle, and last elements, to avoid bad
98 : : performance in the face of already sorted data, or data that is made
99 : : up of multiple sorted runs appended together. Testing shows that a
100 : : median-of-three algorithm provides better performance than simply
101 : : picking the middle element for the latter case. */
102 : :
103 : 1659 : mid = lo + ( size / 2 ) * width; /* find middle element */
104 : :
105 : : /* Sort the first, middle, last elements into order */
106 [ + + ]: 1659 : if (comp( lo, mid, pParam ) > 0)
107 : : {
108 : 137 : inchi_swap( lo, mid, width );
109 : : }
110 [ + + ]: 1659 : if (comp( lo, hi, pParam ) > 0)
111 : : {
112 : 77 : inchi_swap( lo, hi, width );
113 : : }
114 [ + + ]: 1659 : if (comp( mid, hi, pParam ) > 0)
115 : : {
116 : 204 : inchi_swap( mid, hi, width );
117 : : }
118 : :
119 : : /* We now wish to partition the array into three pieces, one consisting
120 : : of elements <= partition element, one of elements equal to the
121 : : partition element, and one of elements > than it. This is done
122 : : below; comments indicate conditions established at every step. */
123 : :
124 : 1659 : loguy = lo;
125 : 1659 : higuy = hi;
126 : :
127 : : /* Note that higuy decreases and loguy increases on every iteration,
128 : : so loop must terminate. */
129 : : for (;;)
130 : : {
131 : : /* lo <= loguy < hi, lo < higuy <= hi,
132 : : A[i] <= A[mid] for lo <= i <= loguy,
133 : : A[i] > A[mid] for higuy <= i < hi,
134 : : A[hi] >= A[mid] */
135 : :
136 : : /* The doubled loop is to avoid calling comp(mid,mid), since some
137 : : existing comparison funcs don't work when passed the same
138 : : value for both pointers. */
139 : :
140 [ + + ]: 1840 : if (mid > loguy)
141 : : {
142 : : do
143 : : {
144 : 3493 : loguy += width;
145 : : }
146 [ + + + + ]: 3493 : while (loguy < mid && comp( loguy, mid, pParam ) <= 0);
147 : : }
148 [ + + ]: 1840 : if (mid <= loguy)
149 : : {
150 : : do
151 : : {
152 : 1816 : loguy += width;
153 : : }
154 [ + + + + ]: 1816 : while (loguy <= hi && comp( loguy, mid, pParam ) <= 0);
155 : : }
156 : :
157 : : /* lo < loguy <= hi+1, A[i] <= A[mid] for lo <= i < loguy,
158 : : either loguy > hi or A[loguy] > A[mid] */
159 : :
160 : : do
161 : : {
162 : 3157 : higuy -= width;
163 : : }
164 [ + + + + ]: 3157 : while (higuy > mid && comp( higuy, mid, pParam ) > 0);
165 : :
166 : : /* lo <= higuy < hi, A[i] > A[mid] for higuy < i < hi,
167 : : either higuy == lo or A[higuy] <= A[mid] */
168 : :
169 [ + + ]: 1840 : if (higuy < loguy)
170 : : {
171 : 1659 : break;
172 : : }
173 : :
174 : : /* if loguy > hi or higuy == lo, then we would have exited, so
175 : : A[loguy] > A[mid], A[higuy] <= A[mid],
176 : : loguy <= hi, higuy > lo */
177 : :
178 : 181 : inchi_swap( loguy, higuy, width );
179 : :
180 : : /* If the partition element was moved, follow it. Only need
181 : : to check for mid == higuy, since before the swap,
182 : : A[loguy] > A[mid] implies loguy != mid. */
183 : :
184 [ + + ]: 181 : if (mid == higuy)
185 : : {
186 : 93 : mid = loguy;
187 : : }
188 : :
189 : : /* A[loguy] <= A[mid], A[higuy] > A[mid]; so condition at top
190 : : of loop is re-established */
191 : : }
192 : :
193 : : /* A[i] <= A[mid] for lo <= i < loguy,
194 : : A[i] > A[mid] for higuy < i < hi,
195 : : A[hi] >= A[mid]
196 : : higuy < loguy
197 : : implying:
198 : : higuy == loguy-1
199 : : or higuy == hi - 1, loguy == hi + 1, A[hi] == A[mid] */
200 : :
201 : : /* Find adjacent elements equal to the partition element. The
202 : : doubled loop is to avoid calling comp(mid,mid), since some
203 : : existing comparison funcs don't work when passed the same value
204 : : for both pointers. */
205 : :
206 : 1659 : higuy += width;
207 [ + + ]: 1659 : if (mid < higuy)
208 : : {
209 : : do
210 : : {
211 : 1140 : higuy -= width;
212 : : }
213 [ + + - + ]: 1140 : while (higuy > mid && comp( higuy, mid, pParam ) == 0);
214 : : }
215 [ + + ]: 1659 : if (mid >= higuy)
216 : : {
217 : : do
218 : : {
219 : 1568 : higuy -= width;
220 : : }
221 [ + + - + ]: 1568 : while (higuy > lo && comp( higuy, mid, pParam ) == 0);
222 : : }
223 : :
224 : : /* OK, now we have the following:
225 : : higuy < loguy
226 : : lo <= higuy <= hi
227 : : A[i] <= A[mid] for lo <= i <= higuy
228 : : A[i] == A[mid] for higuy < i < loguy
229 : : A[i] > A[mid] for loguy <= i < hi
230 : : A[hi] >= A[mid] */
231 : :
232 : : /* We've finished the partition, now we want to sort the subarrays
233 : : [lo, higuy] and [loguy, hi].
234 : : We do the smaller one first to minimize stack usage.
235 : : We only sort arrays of length 2 or more.*/
236 : :
237 [ + + ]: 1659 : if (higuy - lo >= hi - loguy)
238 : : {
239 [ + + ]: 1569 : if (lo < higuy)
240 : : {
241 : 804 : lostk[stkptr] = lo;
242 : 804 : histk[stkptr] = higuy;
243 : 804 : ++stkptr;
244 : : } /* save big recursion for later */
245 : :
246 [ + + ]: 1569 : if (loguy < hi)
247 : : {
248 : 424 : lo = loguy;
249 : 424 : goto recurse; /* do small recursion */
250 : : }
251 : : }
252 : : else
253 : : {
254 [ + - ]: 90 : if (loguy < hi)
255 : : {
256 : 90 : lostk[stkptr] = loguy;
257 : 90 : histk[stkptr] = hi;
258 : 90 : ++stkptr; /* save big recursion for later */
259 : : }
260 : :
261 [ + + ]: 90 : if (lo < higuy)
262 : : {
263 : 29 : hi = higuy;
264 : 29 : goto recurse; /* do small recursion */
265 : : }
266 : : }
267 : :
268 : : /* We have sorted the array, except for any pending sorts on the stack.
269 : : Check if there are any, and do them. */
270 : :
271 : 1206 : --stkptr;
272 [ + + ]: 1206 : if (stkptr >= 0)
273 : : {
274 : 894 : lo = lostk[stkptr];
275 : 894 : hi = histk[stkptr];
276 : 894 : goto recurse; /* pop subarray from stack */
277 : : }
278 : : else
279 : : {
280 : 312 : return; /* all subarrays done */
281 : : }
282 : : }
283 : :
284 : :
285 : : /****************************************************************************/
286 : 3422 : void inchi_swap( char *a, char *b, size_t width )
287 : : {
288 : : char tmp;
289 [ + - ]: 3422 : if (a != b)
290 : : {
291 [ + + ]: 11628 : while (width--)
292 : : {
293 : 8206 : tmp = *a;
294 : 8206 : *a++ = *b;
295 : 8206 : *b++ = tmp;
296 : : }
297 : : }
298 : 3422 : }
299 : :
300 : :
301 : : /****************************************************************************
302 : : Sort by insertions
303 : : ****************************************************************************/
304 : 3363 : int insertions_sort( void *pCG,
305 : : void *base,
306 : : size_t num, size_t width,
307 : : int( *compare )( const void *, const void *, void * ) ) /* djb-rwth: types of variables are sufficient */
308 : : {
309 : 3363 : char *i, *j, *pk = (char*) base;
310 : 3363 : int num_trans = 0;
311 : : size_t k;
312 [ + + ]: 10617 : for (k = 1; k < num; k++, pk += width)
313 : : {
314 : : /*for( i = pk, j = pk + width; j > (char*)base && (*compare)(i,j) > 0; j=i, i -= width )*/
315 : 7254 : for (i = j = pk + width;
316 [ + + + + ]: 9839 : j > ( char* )base && ( i -= width, ( *compare )( i, j, pCG ) ) > 0;
317 : 2585 : j = i) /* changed to keep BoundsChecker happy 2007-09-24 DT */
318 : : {
319 : 2585 : inchi_swap( i, j, width );
320 : 2585 : num_trans++;
321 : : }
322 : : }
323 : :
324 : 3363 : return num_trans;
325 : : }
326 : :
327 : :
328 : : /****************************************************************************
329 : : Sort by insertions
330 : : ****************************************************************************/
331 : 108 : int insertions_sort_AT_NUMBERS( void *pCG,
332 : : AT_NUMB *base,
333 : : int num,
334 : : int( *compare )( const void *e1, const void *e2, void * ) )
335 : : {
336 : : AT_NUMB *i, *j, *pk, tmp;
337 : 108 : int k, num_trans = 0;
338 [ + + ]: 244 : for (k = 1, pk = base; k < num; k++, pk++)
339 : : {
340 [ + + + + ]: 163 : for (j = ( i = pk ) + 1, tmp = *j; j > base && ( *compare )( i, &tmp, pCG ) > 0; j = i, i--)
341 : : {
342 : 27 : *j = *i;
343 : 27 : num_trans++;
344 : : }
345 : 136 : *j = tmp;
346 : : }
347 : :
348 : 108 : return num_trans;
349 : : }
350 : :
351 : :
352 : : /****************************************************************************
353 : : Sort neighbors according to ranks in ascending order
354 : : ****************************************************************************/
355 : 1512 : void insertions_sort_NeighList_AT_NUMBERS( NEIGH_LIST base, AT_RANK *nRank )
356 : : {
357 : : AT_NUMB *i, *j, *pk, tmp;
358 : : AT_RANK rj; /* optimization */
359 : 1512 : int k, num = (int) *base++;
360 [ + + ]: 4135 : for (k = 1, pk = base; k < num; k++, pk++)
361 : : {
362 [ + + + + ]: 2966 : for (j = ( i = pk ) + 1, rj = nRank[(int) *j]; j > base && nRank[(int) *i] > rj; j = i, i--)
363 : : {
364 : 343 : tmp = *i;
365 : 343 : *i = *j;
366 : 343 : *j = tmp;
367 : : }
368 : : }
369 : 1512 : }
370 : :
371 : :
372 : : /****************************************************************************
373 : : Sort neighbors according to ranks in ascending order
374 : : ****************************************************************************/
375 : 0 : int insertions_sort_AT_RANK( AT_RANK *base, int num )
376 : : {
377 : : AT_RANK *i, *j, *pk, tmp;
378 : 0 : int k, num_trans = 0;
379 [ # # ]: 0 : for (k = 1, pk = base; k < num; k++, pk++)
380 : : {
381 [ # # # # ]: 0 : for (j = ( i = pk ) + 1, tmp = *j; j > base && *i > tmp; j = i, i--)
382 : : {
383 : 0 : *j = *i;
384 : 0 : num_trans++;
385 : : }
386 : 0 : *j = tmp;
387 : : }
388 : :
389 : 0 : return num_trans;
390 : : }
391 : :
392 : :
393 : : /****************************************************************************
394 : : Sort neighbors according to ranks in ascending order
395 : : ****************************************************************************/
396 : 120 : int insertions_sort_NeighList_AT_NUMBERS3( NEIGH_LIST base, AT_RANK *nRank )
397 : : {
398 : : AT_NUMB *i, *j, *pk, tmp;
399 : : AT_RANK rj;
400 : 120 : int k, n, num = (int) *base++;
401 [ + + ]: 292 : for (k = 1, pk = base, n = 0; k < num; k++, pk++)
402 : : {
403 : 172 : for (j = ( i = pk ) + 1, rj = nRank[(int) ( tmp = *j )];
404 [ + + + + ]: 179 : j > base && nRank[(int) *i] > rj;
405 : 7 : j = i, i--)
406 : : {
407 : 7 : *j = *i;
408 : 7 : n++;
409 : : }
410 : 172 : *j = tmp;
411 : : }
412 : :
413 : 120 : return n;
414 : : }
415 : :
416 : :
417 : : /****************************************************************************
418 : : Sort neighbors according to symm. ranks (primary key) and canon.
419 : : ranks (secondary key), in descending order
420 : : ****************************************************************************/
421 : 3258 : void insertions_sort_NeighListBySymmAndCanonRank( NEIGH_LIST base,
422 : : const AT_RANK *nSymmRank,
423 : : const AT_RANK *nCanonRank )
424 : : {
425 : : AT_NUMB *i, *j, *pk, tmp;
426 : : int diff;
427 : 3258 : int k, num = (int) *base++;
428 [ + + ]: 6012 : for (k = 1, pk = base; k < num; k++, pk++)
429 : : {
430 : 2754 : for (j = ( i = pk ) + 1;
431 [ + + ]: 4857 : j > base && /* always j > i */
432 [ + + + + ]: 3639 : ( 0 > ( diff = (int) nSymmRank[(int) *i] - (int) nSymmRank[(int) *j] ) ||
433 [ + - ]: 60 : (!diff && nCanonRank[(int) *i] < nCanonRank[(int) *j]) ); /* djb-rwth: addressing LLVM warning */
434 : 2103 : j = i, i--)
435 : : {
436 : 2103 : tmp = *i;
437 : 2103 : *i = *j;
438 : 2103 : *j = tmp;
439 : : }
440 : : }
441 : :
442 : 3258 : }
443 : :
444 : :
445 : : /****************************************************************************
446 : : *
447 : : * Comparison functions
448 : : *
449 : : ****************************************************************************/
450 : :
451 : :
452 : : /****************************************************************************/
453 : 1952 : int CompNeighborsAT_NUMBER( const void* a1, const void* a2, void *p )
454 : : {
455 : 1952 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
456 : : #ifdef CT_NEIGH_INCREASE
457 : 3904 : return (int) pCG->m_pn_RankForSort[pCG->m_pNeighborsForSort[( int )*(const AT_NUMB*) a1]] -
458 : 1952 : (int) pCG->m_pn_RankForSort[pCG->m_pNeighborsForSort[( int )*(const AT_NUMB*) a2]];
459 : : #else
460 : : return (int) ( (CANON_GLOBALS *) pCG )->m_pn_RankForSort[pNeighborsForSort[( int )*(const AT_NUMB*) a2]] -
461 : : (int) ( (CANON_GLOBALS *) pCG )->m_pn_RankForSort[pNeighborsForSort[( int )*(const AT_NUMB*) a1]];
462 : : #endif
463 : : }
464 : :
465 : :
466 : : /**********************************************************************************/
467 : 741 : int comp_AT_RANK( const void* a1, const void* a2, void *p )
468 : : {
469 : 741 : return ( int )*(const AT_RANK*) a1 - ( int )*(const AT_RANK*) a2;
470 : : }
471 : :
472 : :
473 : : /**********************************************************************************/
474 : : /* Compare for sorting Ranks only */
475 : 1394 : int CompRank( const void* a1, const void* a2, void *p )
476 : : {
477 : 1394 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
478 : 1394 : int ret = (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a1] -
479 : 1394 : (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a2];
480 : :
481 : 1394 : return ret;
482 : : }
483 : :
484 : :
485 : : /**********************************************************************************/
486 : 8028 : int CompRanksOrd( const void* a1, const void* a2, void *p )
487 : : {
488 : : int ret;
489 : 8028 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
490 : 8028 : ret = (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a1] -
491 : 8028 : (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a2];
492 [ + + ]: 8028 : if (!ret)
493 : : {
494 : 1649 : ret = ( int )*(const AT_RANK*) a1 - ( int )*(const AT_RANK*) a2;
495 : : }
496 : :
497 : 8028 : return ret;
498 : : }
499 : :
500 : :
501 : : /**********************************************************************************/
502 : 3116 : int CompAtomInvariants2Only( const void* a1, const void* a2, void *p )
503 : : {
504 : 3116 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
505 : 3116 : const ATOM_INVARIANT2 *pAI1 = pCG->m_pAtomInvariant2ForSort + ( int )*(const AT_RANK*) a1;
506 : 3116 : const ATOM_INVARIANT2 *pAI2 = pCG->m_pAtomInvariant2ForSort + ( int )*(const AT_RANK*) a2;
507 : : int i;
508 [ + + ]: 15047 : for (i = 0; i < AT_INV_BREAK1; i++)
509 : : {
510 [ + + ]: 13431 : if (pAI1->val[i] == pAI2->val[i])
511 : 11931 : continue;
512 : 1500 : return (int) pAI1->val[i] - (int) pAI2->val[i];
513 : : }
514 [ - + ]: 1616 : if (pAI1->iso_sort_key != pAI2->iso_sort_key)
515 : : {
516 [ # # ]: 0 : return ( pAI1->iso_sort_key > pAI2->iso_sort_key ) ? 1 : -1;
517 : : }
518 [ + - ]: 1616 : for (; i < AT_INV_LENGTH; i++)
519 : : {
520 [ - + ]: 1616 : if (pAI1->val[i] != pAI2->val[i])
521 : : {
522 : 0 : continue;
523 : : }
524 : 1616 : return (int) pAI1->val[i] - (int) pAI2->val[i];
525 : : }
526 [ # # ]: 0 : if (pAI1->iso_aux_key != pAI2->iso_aux_key)
527 : : {
528 [ # # ]: 0 : return ( pAI1->iso_aux_key > pAI2->iso_aux_key ) ? 1 : -1;
529 : : }
530 : :
531 : 0 : return 0;
532 : : }
533 : :
534 : :
535 : : /**********************************************************************************/
536 : 2566 : int CompAtomInvariants2( const void* a1, const void* a2, void *p )
537 : : {
538 : : /* Warning: the following line may be compiler implementation dependent */
539 : 2566 : int ret = CompAtomInvariants2Only( a1, a2, p );
540 [ + + ]: 2566 : if (!ret)
541 : : {
542 : 1249 : ret = ( int )*(const AT_RANK*) a1 - ( int )*(const AT_RANK*) a2;
543 : : }
544 : :
545 : 2566 : return ret;
546 : : }
547 : :
548 : :
549 : : /**********************************************************************************/
550 : : /* Compare two elements lexicographically */
551 : 46 : int CompChemElemLex( const void *a1, const void *a2 )
552 : : {
553 : 46 : return memcmp( a1, a2, 2 );
554 : : }
555 : :
556 : :
557 : : /****************************************************************************
558 : : Lexicographic compare
559 : : ****************************************************************************/
560 : 2402 : int CompareNeighListLex( NEIGH_LIST pp1, NEIGH_LIST pp2, const AT_RANK *nRank )
561 : : {
562 : 2402 : int len1 = (int) *pp1++;
563 : 2402 : int len2 = (int) *pp2++;
564 : 2402 : int len = inchi_min( len1, len2 );
565 : 2402 : int diff = 0;
566 : : int ret;
567 : : /* djb-rwth: fixing oss-fuzz issue #25642 */
568 [ + + + + ]: 5718 : while ((len > 0) && !diff)
569 : : {
570 : 3316 : len--;
571 : 3316 : diff = (int)nRank[*pp1++] - (int)nRank[*pp2++];
572 : : };
573 : :
574 [ + + ]: 2402 : ret = diff ? diff : (len1 - len2);
575 : 2402 : return ret;
576 : : }
577 : :
578 : :
579 : : /****************************************************************************
580 : : Lexicographic compare
581 : : ****************************************************************************/
582 : 0 : int CompareNeighListLexUpToMaxRank( NEIGH_LIST pp1, NEIGH_LIST pp2, const AT_RANK *nRank, AT_RANK nMaxAtNeighRank )
583 : : {
584 : 0 : int len1 = (int) *pp1++;
585 : 0 : int len2 = (int) *pp2++;
586 : 0 : int diff = 0;
587 : : int len;
588 [ # # # # ]: 0 : while (0 < len1 && nRank[pp1[len1 - 1]] > nMaxAtNeighRank)
589 : : {
590 : 0 : len1--;
591 : : }
592 [ # # # # ]: 0 : while (0 < len2 && nRank[pp2[len2 - 1]] > nMaxAtNeighRank)
593 : : {
594 : 0 : len2--;
595 : : }
596 : 0 : len = inchi_min( len1, len2 );
597 [ # # # # ]: 0 : while (len-- > 0 && !( diff = (int) nRank[*pp1++] - (int) nRank[*pp2++] ))
598 : : {
599 : : ;
600 : : }
601 : :
602 [ # # ]: 0 : return diff ? diff : ( len1 - len2 );
603 : : }
604 : :
605 : :
606 : : /****************************************************************************/
607 : 2266 : int compare_NeighLists( const NEIGH_LIST *op1, const NEIGH_LIST *op2, void *p )
608 : : {
609 : 2266 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
610 : 2266 : return CompareNeighListLex( *op1, *op2, pCG->m_pn_RankForSort );
611 : : }
612 : :
613 : :
614 : : /****************************************************************************/
615 : 6057 : int CompNeighListRanks( const void* a1, const void* a2, void *p )
616 : : {
617 : 6057 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
618 : : int ret;
619 : 6057 : ret = (int) pCG->m_pn_RankForSort[*( (const AT_RANK*) a1 )] -
620 : 6057 : (int) pCG->m_pn_RankForSort[*( (const AT_RANK*) a2 )];
621 [ + + ]: 6057 : if (!ret)
622 : : {
623 : 2124 : ret = compare_NeighLists( pCG->m_pNeighList_RankForSort + *( (const AT_RANK*) a1 ),
624 : 2124 : pCG->m_pNeighList_RankForSort + *( (const AT_RANK*) a2 ), p );
625 : : }
626 : :
627 : 6057 : return ret;
628 : : }
629 : :
630 : :
631 : : /****************************************************************************/
632 : 142 : int CompNeighLists( const void* a1, const void* a2, void *p )
633 : : {
634 : 142 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
635 : : int ret;
636 : 142 : ret = compare_NeighLists( pCG->m_pNeighList_RankForSort + *( (const AT_RANK*) a1 ),
637 : 142 : pCG->m_pNeighList_RankForSort + *( (const AT_RANK*) a2 ), p );
638 : :
639 : 142 : return ret;
640 : : }
641 : :
642 : :
643 : : /****************************************************************************/
644 : 0 : int CompNeighListsUpToMaxRank( const void* a1, const void* a2, void *p )
645 : : {
646 : 0 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
647 : : int ret;
648 : 0 : ret = CompareNeighListLexUpToMaxRank( pCG->m_pNeighList_RankForSort[*( (const AT_RANK*) a1 )],
649 : 0 : pCG->m_pNeighList_RankForSort[*( (const AT_RANK*) a2 )],
650 : 0 : pCG->m_pn_RankForSort, pCG->m_nMaxAtNeighRankForSort );
651 : :
652 : 0 : return ret;
653 : : }
654 : :
655 : :
656 : : /****************************************************************************/
657 : 3275 : int CompNeighListRanksOrd( const void* a1, const void* a2, void *p )
658 : : {
659 : 3275 : int ret = CompNeighListRanks( a1, a2, p );
660 [ + + ]: 3275 : if (!ret)
661 : : {
662 : 532 : ret = ( int )*( (const AT_RANK*) a1 ) - ( int )*( (const AT_RANK*) a2 ); /* keep original order if identical */
663 : : }
664 : :
665 : 3275 : return ret;
666 : : }
667 : :
668 : :
669 : : /****************************************************************************/
670 : 0 : int CompRanksInvOrd( const void* a1, const void* a2, void *p )
671 : : {
672 : 0 : return ( int )*(const AT_RANK*) a2 - ( int )*(const AT_RANK*) a1;
673 : : }
674 : :
675 : :
676 : : /****************************************************************************/
677 : 0 : int CompNeighborsRanksCountEql( const void* a1, const void* a2, void *p )
678 : : {
679 : 0 : CANON_GLOBALS *pCG = (CANON_GLOBALS *) p;
680 : : #ifdef CT_NEIGH_INCREASE
681 : 0 : int ret = (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a1] -
682 : 0 : (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a2];
683 : : #else
684 : : int ret = (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a2] -
685 : : (int) pCG->m_pn_RankForSort[( int )*(const AT_RANK*) a1];
686 : : #endif
687 : 0 : pCG->m_nNumCompNeighborsRanksCountEql += !ret;
688 : :
689 : 0 : return ret;
690 : : }
691 : :
692 : :
693 : : /****************************************************************************
694 : : *
695 : : * In this neighbor list the (vertex number) = (canonical number) - 1
696 : : * Since LinearCT is sorted so that parents are in ascending order
697 : : * and all neighbors of a parent are smaller than the parent and are
698 : : * in ascending order, the neighbors in the NEIGH_LIST are automatically
699 : : * sorted in ascending order
700 : : ****************************************************************************/
701 : 58 : NEIGH_LIST *CreateNeighListFromLinearCT( AT_NUMB *LinearCT, int nLenCT, int num_atoms )
702 : : {
703 : : /* Atom numbers in LinearCT are canonical numbers
704 : : * order: parent[i] > neigh[i][0] < neigh[i][1]...<neigh[i][n] < parent[i+1] > neigh[i+1][0] < ...
705 : : * parent[i] < parent[i+1]
706 : : */
707 : : int i, j;
708 : 58 : S_CHAR *valence = NULL;
709 : 58 : NEIGH_LIST *pp = NULL;
710 : 58 : AT_NUMB *pAtList = NULL;
711 : : AT_RANK n_vertex, n_neigh;
712 : 58 : int err = 1, num_bonds;
713 : : int length, start;
714 : :
715 [ - + ]: 58 : if ((int) LinearCT[0] > num_atoms)
716 : : {
717 : 0 : goto exit_function;
718 : : }
719 : 58 : valence = (S_CHAR*)inchi_calloc((long long)num_atoms + 1, sizeof(valence[0]));
720 [ - + ]: 58 : if (!valence) /* djb-rwth: cast operator added */
721 : : {
722 : 0 : goto exit_function;
723 : : }
724 : :
725 [ + + ]: 1018 : for (i = 1, num_bonds = 0, n_vertex = LinearCT[0]; i < nLenCT; i++)
726 : : {
727 [ + + ]: 960 : if (( n_neigh = LinearCT[i] ) < n_vertex)
728 : : {
729 : 488 : valence[n_neigh] ++;
730 : 488 : valence[n_vertex] ++;
731 : 488 : num_bonds += 2;
732 : : }
733 : : else
734 : : {
735 [ - + ]: 472 : if ((int) ( n_vertex = n_neigh ) > num_atoms)
736 : : {
737 : 0 : goto exit_function;
738 : : }
739 : : }
740 : : }
741 [ - + ]: 58 : if ((int) n_vertex != num_atoms)
742 : : {
743 : 0 : goto exit_function;
744 : : }
745 : 58 : length = num_bonds + num_atoms + 1;
746 : 58 : pp = (NEIGH_LIST*)inchi_calloc(((long long)num_atoms + 1), sizeof(NEIGH_LIST));
747 : 58 : pAtList = (AT_NUMB*)inchi_malloc(length * sizeof(AT_NUMB));
748 [ - + - + ]: 58 : if (pp && pAtList) /* djb-rwth: cast operator added; addressing LLVM warning */
749 : : {
750 : : /* Create empty connection table */
751 [ + + ]: 588 : for (i = 1, length = 0; i <= num_atoms; i++)
752 : : {
753 : 530 : start = length;
754 : 530 : length += ( valence[i] + 1 );
755 : 530 : pp[i - 1] = pAtList + start;
756 : 530 : pp[i - 1][0] = 0;
757 : : }
758 : : /* Fill out the CT */
759 [ + + ]: 1018 : for (i = 1, n_vertex = LinearCT[0] - 1; i < nLenCT; i++)
760 : : {
761 [ + + ]: 960 : if (( n_neigh = LinearCT[i] - 1 ) < n_vertex)
762 : : {
763 : : /* Vertex - neighbor connection */
764 : 488 : j = (int) ( ++pp[(int) n_vertex][0] );
765 : 488 : pp[(int) n_vertex][j] = n_neigh;
766 : : /* neighbor - vertex connection */
767 : 488 : j = (int) ( ++pp[(int) n_neigh][0] );
768 : 488 : pp[(int) n_neigh][j] = n_vertex;
769 : : }
770 : : else
771 : : {
772 [ - + ]: 472 : if ((int) ( n_vertex = n_neigh ) >= num_atoms)
773 : : {
774 : 0 : goto exit_function;
775 : : }
776 : : }
777 : : }
778 : 58 : err = 0;
779 : : }
780 : :
781 : 0 : exit_function:
782 [ + - ]: 58 : if (valence)
783 : : {
784 [ + - ]: 58 : inchi_free( valence );
785 : : }
786 [ - + ]: 58 : if (err) /* djb-rwth: ignoring LLVM warning */
787 : : {
788 [ # # ]: 0 : if (pAtList)
789 : : {
790 [ # # ]: 0 : inchi_free( pAtList );
791 : : }
792 [ # # ]: 0 : if (pp)
793 : : {
794 [ # # ]: 0 : inchi_free( pp );
795 : 0 : pp = NULL;
796 : : }
797 : : }
798 : :
799 : 58 : return pp; /* djb-rwth: ignoring LLVM warning: since a pointer is returned, memory should be freed in a function which calls *CreateNeighListFromLinearCT */
800 : : }
801 : :
802 : :
803 : : /****************************************************************************
804 : : NEIGH_LIST pp[] is an array of pointers to the lists of neighboring
805 : : atoms numbers. The first number in each list is a number of neighbors.
806 : : In case of bDoubleBondSquare != 0 neighbors connected by the double bond
807 : : appear 2 times. The first element pp[0] is a pointer to be deallocated
808 : : to free all the lists.
809 : : ****************************************************************************/
810 : 393 : NEIGH_LIST *CreateNeighList( int num_atoms,
811 : : int num_at_tg,
812 : : sp_ATOM* at,
813 : : int bDoubleBondSquare,
814 : : T_GROUP_INFO *t_group_info )
815 : : {
816 : : /* +1 to add NULL termination */
817 : 393 : NEIGH_LIST *pp = (NEIGH_LIST *) inchi_calloc( ( (long long)num_at_tg + 1 ), sizeof( NEIGH_LIST ) ); /* djb-rwth: cast operator added */
818 : 393 : T_GROUP *t_group = NULL;
819 : 393 : AT_NUMB *nEndpointAtomNumber = NULL;
820 : 393 : int num_t_groups = 0;
821 : : int nFirstEndpointAtNoPos;
822 : :
823 : 393 : AT_NUMB *pAtList = NULL;
824 : : int length, start, val, i, j;
825 [ + - ]: 393 : if (pp)
826 : : {
827 [ - + ]: 393 : if (num_at_tg > num_atoms)
828 : : {
829 : 0 : t_group = t_group_info->t_group;
830 : 0 : num_t_groups = t_group_info->num_t_groups;
831 : 0 : nEndpointAtomNumber = t_group_info->nEndpointAtomNumber;
832 : : }
833 : :
834 [ + - ]: 393 : if (!bDoubleBondSquare)
835 : : {
836 [ + + ]: 4273 : for (i = 0, length = 0; i < num_atoms; i++)
837 : : {
838 [ - + - - ]: 3880 : length += (int) at[i].valence + ( num_t_groups && at[i].endpoint );
839 : : }
840 : 393 : length += num_atoms;
841 [ - + ]: 393 : for (i = 0; i < num_t_groups; i++)
842 : : {
843 : 0 : length += (int) t_group[i].nNumEndpoints;
844 : : }
845 : 393 : length += num_t_groups;
846 : : }
847 : : else
848 : : {
849 [ # # ]: 0 : for (i = 0, length = 0; i < num_atoms; i++)
850 : : {
851 : 0 : val = (int) at[i].valence;
852 [ # # ]: 0 : for (j = 0; j < val; j++)
853 : : {
854 [ # # # # ]: 0 : length += 1 + ( bDoubleBondSquare && BOND_DOUBLE == at[i].bond_type[j] );
855 : : }
856 [ # # # # ]: 0 : length += ( num_t_groups && at[i].endpoint );
857 : : }
858 : 0 : length += num_atoms;
859 [ # # ]: 0 : for (i = 0; i < num_t_groups; i++)
860 : : {
861 : 0 : length += (int) t_group[i].nNumEndpoints;
862 : : }
863 : 0 : length += num_t_groups;
864 : : }
865 : 393 : length++; /* +1 to save number of neighbors */
866 : 393 : pAtList = (AT_NUMB*)inchi_malloc(length * sizeof(AT_NUMB));
867 [ + - ]: 393 : if (pAtList) /* djb-rwth: addressing LLVM warning */
868 : : {
869 [ + - ]: 393 : if (!bDoubleBondSquare)
870 : : {
871 [ + + ]: 4273 : for (i = 0, length = 0; i < num_atoms; i++)
872 : : {
873 : 3880 : val = at[i].valence;
874 : 3880 : start = length++;
875 [ + + ]: 11026 : for (j = 0; j < val; j++)
876 : : {
877 : 7146 : pAtList[length++] = at[i].neighbor[j];
878 : : }
879 : : /* add endpoint */
880 [ - + - - ]: 3880 : if (num_t_groups && at[i].endpoint)
881 : : {
882 : 0 : pAtList[length++] = num_atoms + (int) at[i].endpoint - 1;
883 : : }
884 : 3880 : pAtList[start] = length - start - 1; /* number of neighbors before the list of neighbors */
885 : 3880 : pp[i] = pAtList + start; /* pointer to the <num.neigh.><list of neigh> */
886 : : }
887 : : }
888 : : else
889 : : {
890 [ # # ]: 0 : for (i = 0, length = 0; i < num_atoms; i++)
891 : : {
892 : 0 : val = at[i].valence;
893 : 0 : start = length++;
894 [ # # ]: 0 : for (j = 0; j < val; j++)
895 : : {
896 : 0 : pAtList[length++] = at[i].neighbor[j]; /* djb-rwth: buffer overrun avoided implicitly */
897 [ # # # # ]: 0 : if (bDoubleBondSquare && BOND_DOUBLE == at[i].bond_type[j])
898 : : {
899 : 0 : pAtList[length++] = at[i].neighbor[j]; /* a list of neighbor orig. numbers */
900 : : }
901 : : }
902 : : /* Add endpoint */
903 [ # # # # ]: 0 : if (num_t_groups && at[i].endpoint)
904 : : {
905 : 0 : pAtList[length++] = num_atoms + (int) at[i].endpoint - 1;
906 : : }
907 : 0 : pAtList[start] = length - start - 1; /* number of neighbors before the list of neighbors */
908 : 0 : pp[i] = pAtList + start; /* pointer to the <num.neigh.><list of neigh> */
909 : : }
910 : : }
911 : :
912 : : /* Add t-groups */
913 [ - + ]: 393 : for (i = 0; i < num_t_groups; i++)
914 : : {
915 : 0 : val = (int) t_group[i].nNumEndpoints;
916 : 0 : start = length++;
917 : 0 : nFirstEndpointAtNoPos = (int) t_group[i].nFirstEndpointAtNoPos;
918 [ # # ]: 0 : for (j = 0; j < val; j++)
919 : : {
920 : 0 : pAtList[length++] = nEndpointAtomNumber[nFirstEndpointAtNoPos + j];
921 : : }
922 : 0 : pAtList[start] = length - start - 1; /* number of neighbors before the list of neighbors */
923 : 0 : pp[num_atoms + i] = pAtList + start; /* pointer to the <num.neigh.><list of neigh> */
924 : : }
925 : : }
926 : : else
927 : : {
928 [ # # ]: 0 : inchi_free(pAtList); /* djb-rwth: fixing coverity ID #499598 */
929 [ # # ]: 0 : inchi_free( pp );
930 : 0 : return NULL;
931 : : }
932 : : } /* djb-rwth: ignoring LLVM warning */
933 : :
934 : : /* djb-rwth: fixing coverity ID #499598 -- pp uses pAtList values */
935 : :
936 : 393 : return pp; /* djb-rwth: ignoring LLVM warning: since a pointer is returned, memory should be freed in a function which calls *CreateNeighList */
937 : : }
938 : :
939 : :
940 : : /****************************************************************************/
941 : 793 : void FreeNeighList( NEIGH_LIST *pp )
942 : : {
943 [ + + ]: 793 : if (pp)
944 : : {
945 [ + - ]: 451 : if (pp[0])
946 : : {
947 [ + - ]: 451 : inchi_free( pp[0] );
948 : : }
949 [ + - ]: 451 : inchi_free( pp );
950 : : }
951 : :
952 : 793 : }
953 : :
954 : :
955 : : /****************************************************************************/
956 : 214 : int BreakAllTies( CANON_GLOBALS *pCG,
957 : : int num_atoms,
958 : : int num_max,
959 : : AT_RANK **pRankStack,
960 : : NEIGH_LIST *NeighList,
961 : : AT_RANK *nTempRank,
962 : : CANON_STAT *pCS )
963 : : {
964 : 214 : int i, nRet = -1, nNumRanks = 1 /* value does not matter*/;
965 : :
966 : 214 : AT_RANK *nPrevRank = *pRankStack++;
967 : 214 : AT_RANK *nPrevAtomNumber = *pRankStack++;
968 : :
969 : 214 : AT_RANK *nNewRank = NULL;
970 : 214 : AT_RANK *nNewAtomNumber = NULL;
971 : :
972 [ + + ]: 214 : if (!pRankStack[0])
973 : : {
974 : 96 : pRankStack[0] = (AT_RANK *) inchi_malloc( num_max * sizeof( *nNewRank ) );
975 : : }
976 [ + + ]: 214 : if (!pRankStack[1])
977 : : {
978 : 96 : pRankStack[1] = (AT_RANK *) inchi_malloc( num_max * sizeof( *nNewAtomNumber ) );
979 : : }
980 [ + - - + ]: 214 : if (!pRankStack[0] || !pRankStack[1])
981 : : {
982 : 0 : return CT_OUT_OF_RAM; /* <BRKPT> */
983 : : }
984 : 214 : nNewRank = pRankStack[0];
985 : 214 : nNewAtomNumber = pRankStack[1];
986 : :
987 [ + - + - ]: 214 : if (nNewRank && nNewAtomNumber)
988 : : {
989 : 214 : memcpy(nNewAtomNumber, nPrevAtomNumber, num_atoms * sizeof(nNewAtomNumber[0]));
990 : 214 : memcpy(nNewRank, nPrevRank, num_atoms * sizeof(nNewRank[0]));
991 [ + + ]: 2172 : for (i = 1, nRet = 0; i < num_atoms; i++)
992 : : {
993 : : /* 12-12-2001: replaced Prev... with New... */
994 [ + + ]: 1958 : if (nNewRank[(int) nNewAtomNumber[i - 1]] == nNewRank[(int) nNewAtomNumber[i]])
995 : : {
996 : 38 : nNewRank[nNewAtomNumber[i - 1]] = (AT_RANK) i;
997 : 38 : nNumRanks = DifferentiateRanks2( pCG, num_atoms, NeighList,
998 : : nNumRanks, nNewRank, nTempRank,
999 : : nNewAtomNumber,
1000 : : &pCS->lNumNeighListIter, 1 );
1001 : 38 : pCS->lNumBreakTies++;
1002 : 38 : nRet++;
1003 : : }
1004 : : }
1005 : : }
1006 : :
1007 : 214 : return nRet;
1008 : : }
1009 : :
1010 : :
1011 : : /****************************************************************************
1012 : : Int insertions sort
1013 : : ****************************************************************************/
1014 : 0 : int * iisort( int *list, int num )
1015 : : {
1016 : : int i;
1017 [ # # ]: 0 : for (i = 1; i < num; i++)
1018 : : {
1019 : 0 : int tmp = list[i];
1020 : 0 : int j = i - 1;
1021 [ # # # # ]: 0 : while (j >= 0 && list[j] > tmp)
1022 : : {
1023 : 0 : list[j + 1] = list[j];
1024 : 0 : j--;
1025 : : }
1026 : 0 : list[j + 1] = tmp;
1027 : : }
1028 : :
1029 : 0 : return list;
1030 : : }
|
