Actual source code: ctable.c
petsc-3.7.3 2016-08-01
2: /* Contributed by - Mark Adams */
4: #include <petscsys.h>
5: #include <petscctable.h>
9: /*
10: PetscTableCreate Creates a PETSc look up table
12: Input Parameters:
13: + n - expected number of keys
14: - maxkey- largest possible key
16: Notes: keys are between 1 and maxkey inclusive
18: */
19: PetscErrorCode PetscTableCreate(const PetscInt n,PetscInt maxkey,PetscTable *rta)
20: {
21: PetscTable ta;
25: if (n < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"n < 0");
26: PetscNew(&ta);
27: ta->tablesize = 17 + PetscIntMultTruncate(3,n/2);
28: PetscCalloc1(ta->tablesize,&ta->keytable);
29: PetscMalloc1(ta->tablesize,&ta->table);
30: ta->head = 0;
31: ta->count = 0;
32: ta->maxkey = maxkey;
33: *rta = ta;
34: return(0);
35: }
39: /* PetscTableCreate() ********************************************
40: *
41: * hash table for non-zero data and keys
42: *
43: */
44: PetscErrorCode PetscTableCreateCopy(const PetscTable intable,PetscTable *rta)
45: {
47: PetscInt i;
48: PetscTable ta;
51: PetscNew(&ta);
52: ta->tablesize = intable->tablesize;
53: PetscMalloc1(ta->tablesize,&ta->keytable);
54: PetscMalloc1(ta->tablesize,&ta->table);
55: for (i = 0; i < ta->tablesize; i++) {
56: ta->keytable[i] = intable->keytable[i];
57: ta->table[i] = intable->table[i];
58: #if defined(PETSC_USE_DEBUG)
59: if (ta->keytable[i] < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_COR,"ta->keytable[i] < 0");
60: #endif
61: }
62: ta->head = 0;
63: ta->count = intable->count;
64: ta->maxkey = intable->maxkey;
65: *rta = ta;
66: return(0);
67: }
71: /* PetscTableDestroy() ********************************************
72: *
73: *
74: */
75: PetscErrorCode PetscTableDestroy(PetscTable *ta)
76: {
80: if (!*ta) return(0);
81: PetscFree((*ta)->keytable);
82: PetscFree((*ta)->table);
83: PetscFree(*ta);
84: return(0);
85: }
89: /* PetscTableGetCount() ********************************************
90: */
91: PetscErrorCode PetscTableGetCount(const PetscTable ta,PetscInt *count)
92: {
94: *count = ta->count;
95: return(0);
96: }
100: /* PetscTableIsEmpty() ********************************************
101: */
102: PetscErrorCode PetscTableIsEmpty(const PetscTable ta,PetscInt *flag)
103: {
105: *flag = !(ta->count);
106: return(0);
107: }
111: /*
112: PetscTableAddExpand - called by PetscTableAdd() if more space is needed
114: */
115: PetscErrorCode PetscTableAddExpand(PetscTable ta,PetscInt key,PetscInt data,InsertMode imode)
116: {
118: PetscInt ii = 0;
119: const PetscInt tsize = ta->tablesize,tcount = ta->count;
120: PetscInt *oldtab = ta->table,*oldkt = ta->keytable,newk,ndata;
123: ta->tablesize = PetscIntMultTruncate(2,ta->tablesize);
124: if (tsize == ta->tablesize) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Table is as large as possible; ./configure with the option --with-64-bit-integers to run this large case");
126: PetscMalloc1(ta->tablesize,&ta->table);
127: PetscCalloc1(ta->tablesize,&ta->keytable);
129: ta->count = 0;
130: ta->head = 0;
132: PetscTableAdd(ta,key,data,INSERT_VALUES);
133: /* rehash */
134: for (ii = 0; ii < tsize; ii++) {
135: newk = oldkt[ii];
136: if (newk) {
137: ndata = oldtab[ii];
138: PetscTableAdd(ta,newk,ndata,imode);
139: }
140: }
141: if (ta->count != tcount + 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_COR,"corrupted ta->count");
143: PetscFree(oldtab);
144: PetscFree(oldkt);
145: return(0);
146: }
151: /* PetscTableRemoveAll() ********************************************
152: *
153: *
154: */
155: PetscErrorCode PetscTableRemoveAll(PetscTable ta)
156: {
160: ta->head = 0;
161: if (ta->count) {
162: ta->count = 0;
164: PetscMemzero(ta->keytable,ta->tablesize*sizeof(PetscInt));
165: }
166: return(0);
167: }
173: /* PetscTableGetHeadPosition() ********************************************
174: *
175: */
176: PetscErrorCode PetscTableGetHeadPosition(PetscTable ta,PetscTablePosition *ppos)
177: {
178: PetscInt i = 0;
181: *ppos = NULL;
182: if (!ta->count) return(0);
184: /* find first valid place */
185: do {
186: if (ta->keytable[i]) {
187: *ppos = (PetscTablePosition)&ta->table[i];
188: break;
189: }
190: } while (i++ < ta->tablesize);
191: if (!*ppos) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_COR,"No head");
192: return(0);
193: }
197: /* PetscTableGetNext() ********************************************
198: *
199: * - iteration - PetscTablePosition is always valid (points to a data)
200: *
201: */
202: PetscErrorCode PetscTableGetNext(PetscTable ta,PetscTablePosition *rPosition,PetscInt *pkey,PetscInt *data)
203: {
204: PetscInt idex;
205: PetscTablePosition pos;
208: pos = *rPosition;
209: if (!pos) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Null position");
210: *data = *pos;
211: if (!*data) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Null data");
212: idex = pos - ta->table;
213: *pkey = ta->keytable[idex];
214: if (!*pkey) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Null key");
216: /* get next */
217: do {
218: pos++; idex++;
219: if (idex >= ta->tablesize) {
220: pos = 0; /* end of list */
221: break;
222: } else if (ta->keytable[idex]) {
223: pos = ta->table + idex;
224: break;
225: }
226: } while (idex < ta->tablesize);
227: *rPosition = pos;
228: return(0);
229: }
234: PetscErrorCode PetscTableAddCountExpand(PetscTable ta,PetscInt key)
235: {
237: PetscInt ii = 0,hash = PetscHash(ta,key);
238: const PetscInt tsize = ta->tablesize,tcount = ta->count;
239: PetscInt *oldtab = ta->table,*oldkt = ta->keytable,newk,ndata;
242: /* before making the table larger check if key is already in table */
243: while (ii++ < tsize) {
244: if (ta->keytable[hash] == key) return(0);
245: hash = (hash == (ta->tablesize-1)) ? 0 : hash+1;
246: }
248: ta->tablesize = PetscIntMultTruncate(2,ta->tablesize);
249: if (tsize == ta->tablesize) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Table is as large as possible; ./configure with the option --with-64-bit-integers to run this large case");
250: PetscMalloc1(ta->tablesize,&ta->table);
251: PetscCalloc1(ta->tablesize,&ta->keytable);
253: ta->count = 0;
254: ta->head = 0;
256: /* Build a new copy of the data */
257: for (ii = 0; ii < tsize; ii++) {
258: newk = oldkt[ii];
259: if (newk) {
260: ndata = oldtab[ii];
261: PetscTableAdd(ta,newk,ndata,INSERT_VALUES);
262: }
263: }
264: PetscTableAddCount(ta,key);
265: if (ta->count != tcount + 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_COR,"corrupted ta->count");
267: PetscFree(oldtab);
268: PetscFree(oldkt);
269: return(0);
270: }