Actual source code: aodata.c
1: #define PETSCDM_DLL
3: /*
4: Defines the abstract operations on AOData
5: */
6: #include src/dm/ao/aoimpl.h
11: /*@C
12: AODataGetInfo - Gets the number of keys and their names in a database.
14: Not collective
16: Input Parameter:
17: . ao - the AOData database
19: Output Parameters:
20: + nkeys - the number of keys
21: - keys - the names of the keys (or PETSC_NULL)
23: Level: deprecated
25: .keywords: application ordering
27: .seealso: AODataSegmentGetInfo()
28: @*/
29: PetscErrorCode AODataGetInfo(AOData ao,PetscInt *nkeys,char ***keys)
30: {
32: PetscInt n,i;
33: AODataKey *key = ao->keys;
38: *nkeys = n = ao->nkeys;
39: if (keys) {
40: PetscMalloc((n+1)*sizeof(char *),&keys);
41: for (i=0; i<n; i++) {
42: if (!key) SETERRQ(PETSC_ERR_COR,"Less keys in database then indicated");
43: (*keys)[i] = key->name;
44: key = key->next;
45: }
46: }
47: return(0);
48: }
52: /*
53: AODataKeyFind_Private - Given a keyname finds the key. Generates a flag if not found.
55: Not collective
57: Input Parameters:
58: . keyname - string name of key
60: Output Parameter:
61: + flag - PETSC_TRUE if found, PETSC_FALSE if not found
62: - key - the associated key
64: Level: deprecated
66: */
67: PetscErrorCode AODataKeyFind_Private(AOData aodata,const char keyname[],PetscTruth *flag,AODataKey **key)
68: {
69: PetscTruth match;
71: AODataAlias *t = aodata->aliases;
72: const char *name = keyname;
73: AODataKey *nkey;
76: *key = PETSC_NULL;
77: *flag = PETSC_FALSE;
78: while (name) {
79: nkey = aodata->keys;
80: while (nkey) {
81: PetscStrcmp(nkey->name,name,&match);
82: if (match) {
83: /* found the key */
84: *key = nkey;
85: *flag = PETSC_TRUE;
86: return(0);
87: }
88: *key = nkey;
89: nkey = nkey->next;
90: }
91: name = 0;
92: while (t) {
93: PetscStrcmp(keyname,t->alias,&match);
94: if (match) {
95: name = t->name;
96: t = t->next;
97: break;
98: }
99: t = t->next;
100: }
101: }
102: return(0);
103: }
107: /*@C
108: AODataKeyExists - Determines if a key exists in the database.
110: Not collective
112: Input Parameters:
113: . keyname - string name of key
115: Output Parameter:
116: . flag - PETSC_TRUE if found, otherwise PETSC_FALSE
118: Level: deprecated
120: @*/
121: PetscErrorCode AODataKeyExists(AOData aodata,const char keyname[],PetscTruth *flag)
122: {
124: PetscTruth iflag;
125: AODataKey *ikey;
129: AODataKeyFind_Private(aodata,keyname,&iflag,&ikey);
130: if (iflag) *flag = PETSC_TRUE;
131: else *flag = PETSC_FALSE;
132: return(0);
133: }
138: /*
139: AODataSegmentFind_Private - Given a key and segment finds the int key, segment
140: coordinates. Generates a flag if not found.
142: Not collective
144: Input Parameters:
145: + keyname - string name of key
146: - segname - string name of segment
148: Output Parameter:
149: + flag - PETSC_TRUE if found, PETSC_FALSE if not
150: . key - integer of keyname
151: - segment - integer of segment
153: If it doesn't find it, it returns the last seg in the key (if the key exists)
155: Level: deprecated
157: */
158: PetscErrorCode AODataSegmentFind_Private(AOData aodata,const char keyname[],const char segname[],PetscTruth *flag,AODataKey **key,AODataSegment **seg)
159: {
161: PetscTruth keyflag,match;
162: AODataAlias *t = aodata->aliases;
163: const char *name;
164: AODataSegment *nseg;
167: *seg = PETSC_NULL;
168: *flag = PETSC_FALSE;
169: AODataKeyFind_Private(aodata,keyname,&keyflag,key);
170: if (keyflag) { /* found key now look for segment */
171: name = segname;
172: while (name) {
173: nseg = (*key)->segments;
174: while (nseg) {
175: PetscStrcmp(nseg->name,name,&match);
176: if (match) {
177: /* found the segment */
178: *seg = nseg;
179: *flag = PETSC_TRUE;
180: return(0);
181: }
182: *seg = nseg;
183: nseg = nseg->next;
184: }
185: name = 0;
186: while (t) {
187: PetscStrcmp(segname,t->alias,&match);
188: if (match) {
189: name = t->name;
190: t = t->next;
191: break;
192: }
193: t = t->next;
194: }
195: }
196: }
197: return(0);
198: }
202: /*@C
203: AODataSegmentExists - Determines if a key and segment exists in the database.
205: Not collective
207: Input Parameters:
208: + keyname - string name of key
209: - segname - string name of segment
211: Output Parameter:
212: . flag - PETSC_TRUE if found, else PETSC_FALSE
214: Level: deprecated
216: @*/
217: PetscErrorCode AODataSegmentExists(AOData aodata,const char keyname[],const char segname[],PetscTruth *flag)
218: {
220: PetscTruth iflag;
221: AODataKey *ikey;
222: AODataSegment *iseg;
226: AODataSegmentFind_Private(aodata,keyname,segname,&iflag,&ikey,&iseg);
227: if (iflag) *flag = PETSC_TRUE;
228: else *flag = PETSC_FALSE;
229: return(0);
230: }
232: /* ------------------------------------------------------------------------------------*/
236: /*@C
237: AODataKeyGetActive - Get a sublist of key indices that have a logical flag on.
239: Collective on AOData
241: Input Parameters:
242: + aodata - the database
243: . name - the name of the key
244: . segment - the name of the segment
245: . n - the number of key indices provided by this processor
246: . keys - the keys provided by this processor
247: - wl - which logical key in the block (for block size 1 this is always 0)
249: Output Parameters:
250: . is - the list of key indices
252: Level: deprecated
254: .keywords: database transactions
256: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
257: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
258: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
259: @*/
260: PetscErrorCode AODataKeyGetActive(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,PetscInt wl,IS *is)
261: {
266: (*aodata->ops->keygetactive)(aodata,name,segment,n,keys,wl,is);
267: return(0);
268: }
272: /*@C
273: AODataKeyGetActiveIS - Get a sublist of key indices that have a logical flag on.
275: Collective on AOData
277: Input Parameters:
278: + aodata - the database
279: . name - the name of the key
280: . segment - the name of the segment
281: . in - the key indices we are checking
282: - wl - which logical key in the block (for block size 1 this is always 0)
284: Output Parameters:
285: . IS - the list of key indices
287: Level: deprecated
289: .keywords: database transactions
291: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
292: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
293: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
294: @*/
295: PetscErrorCode AODataKeyGetActiveIS(AOData aodata,const char name[],const char segname[],IS in,PetscInt wl,IS *is)
296: {
298: PetscInt n,*keys;
301: ISGetLocalSize(in,&n);
302: ISGetIndices(in,&keys);
303: AODataKeyGetActive(aodata,name,segname,n,keys,wl,is);
304: ISRestoreIndices(in,&keys);
305: return(0);
306: }
310: /*@C
311: AODataKeyGetActiveLocal - Get a sublist of key indices that have a logical flag on.
313: Collective on AOData
315: Input Parameters:
316: + aodata - the database
317: . name - the name of the key
318: . segment - the name of the segment
319: . n - the number of key indices provided by this processor
320: . keys - the keys provided by this processor
321: - wl - which logical key in the block (for block size 1 this is always 0)
323: Output Parameters:
324: . IS - the list of key indices
326: Level: deprecated
328: .keywords: database transactions
330: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
331: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
332: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
333: @*/
334: PetscErrorCode AODataKeyGetActiveLocal(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,PetscInt wl,IS *is)
335: {
340: (*aodata->ops->keygetactivelocal)(aodata,name,segment,n,keys,wl,is);
341: return(0);
342: }
346: /*@C
347: AODataKeyGetActiveLocalIS - Get a sublist of key indices that have a logical flag on.
349: Collective on AOData
351: Input Parameters:
352: + aodata - the database
353: . name - the name of the key
354: . segment - the name of the segment
355: . in - the key indices we are checking
356: - wl - which logical key in the block (for block size 1 this is always 0)
358: Output Parameters:
359: . IS - the list of key indices
361: Level: advanced
363: .keywords: database transactions
365: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
366: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
367: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
368: @*/
369: PetscErrorCode AODataKeyGetActiveLocalIS(AOData aodata,const char name[],const char segname[],IS in,PetscInt wl,IS *is)
370: {
372: PetscInt n,*keys;
375: ISGetLocalSize(in,&n);
376: ISGetIndices(in,&keys);
377: AODataKeyGetActiveLocal(aodata,name,segname,n,keys,wl,is);
378: ISRestoreIndices(in,&keys);
379: return(0);
380: }
382: /* ------------------------------------------------------------------------------------*/
386: /*@C
387: AODataSegmentGet - Get data from a particular segment of a database.
389: Collective on AOData
391: Input Parameters:
392: + aodata - the database
393: . name - the name of the key
394: . segment - the name of the segment
395: . n - the number of data items needed by this processor
396: - keys - the keys provided by this processor
398: Output Parameters:
399: . data - the actual data
401: Level: deprecated
403: .keywords: database transactions
405: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
406: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
407: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
408: @*/
409: PetscErrorCode AODataSegmentGet(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
410: {
415: (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
416: return(0);
417: }
421: /*@C
422: AODataSegmentRestore - Restores data from a particular segment of a database.
424: Collective on AOData
426: Input Parameters:
427: + aodata - the database
428: . name - the name of the key
429: . segment - the name of the segment
430: . n - the number of data items needed by this processor
431: - keys - the keys provided by this processor
433: Output Parameters:
434: . data - the actual data
436: Level: deprecated
438: .keywords: database transactions
440: .seealso: AODataSegmentRestoreIS()
441: @*/
442: PetscErrorCode AODataSegmentRestore(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
443: {
448: (*aodata->ops->segmentrestore)(aodata,name,segment,n,keys,data);
449: return(0);
450: }
454: /*@C
455: AODataSegmentGetIS - Get data from a particular segment of a database.
457: Collective on AOData and IS
459: Input Parameters:
460: + aodata - the database
461: . name - the name of the key
462: . segment - the name of the segment
463: - is - the keys for data requested on this processor
465: Output Parameters:
466: . data - the actual data
468: Level: deprecated
470: .keywords: database transactions
472: @*/
473: PetscErrorCode AODataSegmentGetIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
474: {
476: PetscInt n,*keys;
482: ISGetLocalSize(is,&n);
483: ISGetIndices(is,&keys);
484: (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
485: ISRestoreIndices(is,&keys);
486: return(0);
487: }
491: /*@C
492: AODataSegmentRestoreIS - Restores data from a particular segment of a database.
494: Collective on AOData and IS
496: Input Parameters:
497: + aodata - the database
498: . name - the name of the data key
499: . segment - the name of the segment
500: - is - the keys provided by this processor
502: Output Parameters:
503: . data - the actual data
505: Level: deprecated
507: .keywords: database transactions
509: .seealso: AODataSegmentRestore()
510: @*/
511: PetscErrorCode AODataSegmentRestoreIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
512: {
518: (*aodata->ops->segmentrestore)(aodata,name,segment,0,0,data);
520: return(0);
521: }
523: /* ------------------------------------------------------------------------------------*/
526: /*@C
527: AODataSegmentGetLocal - Get data from a particular segment of a database. Returns the
528: values in the local numbering; valid only for integer segments.
530: Collective on AOData
532: Input Parameters:
533: + aodata - the database
534: . name - the name of the key
535: . segment - the name of the segment
536: . n - the number of data items needed by this processor
537: - keys - the keys provided by this processor in local numbering
539: Output Parameters:
540: . data - the actual data
542: Level: deprecated
544: .keywords: database transactions
546: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
547: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
548: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
549: @*/
550: PetscErrorCode AODataSegmentGetLocal(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
551: {
556: (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
557: return(0);
558: }
562: /*@C
563: AODataSegmentRestoreLocal - Restores data from a particular segment of a database.
565: Collective on AOData
567: Input Parameters:
568: + aodata - the database
569: . name - the name of the key
570: . segment - the name of the segment
571: . n - the number of data items needed by this processor
572: - keys - the keys provided by this processor
574: Output Parameters:
575: . data - the actual data
577: Level: deprecated
579: .keywords: database transactions
581: @*/
582: PetscErrorCode AODataSegmentRestoreLocal(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
583: {
588: (*aodata->ops->segmentrestorelocal)(aodata,name,segment,n,keys,data);
589: return(0);
590: }
594: /*@C
595: AODataSegmentGetLocalIS - Get data from a particular segment of a database. Returns the
596: values in the local numbering; valid only for integer segments.
598: Collective on AOData and IS
600: Input Parameters:
601: + aodata - the database
602: . name - the name of the key
603: . segment - the name of the segment
604: - is - the keys for data requested on this processor
606: Output Parameters:
607: . data - the actual data
609: Level: deprecated
611: .keywords: database transactions
613: .seealso: AODataSegmentRestoreLocalIS()
614: @*/
615: PetscErrorCode AODataSegmentGetLocalIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
616: {
618: PetscInt n,*keys;
624: ISGetLocalSize(is,&n);
625: ISGetIndices(is,&keys);
626: (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
627: ISRestoreIndices(is,&keys);
628: return(0);
629: }
633: /*@C
634: AODataSegmentRestoreLocalIS - Restores data from a particular segment of a database.
636: Collective on AOData and IS
638: Input Parameters:
639: + aodata - the database
640: . name - the name of the data key
641: . segment - the name of the segment
642: - is - the keys provided by this processor
644: Output Parameters:
645: . data - the actual data
647: Level: deprecated
649: .keywords: database transactions
651: .seealso: AODataSegmentGetLocalIS()
652: @*/
653: PetscErrorCode AODataSegmentRestoreLocalIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
654: {
660: (*aodata->ops->segmentrestorelocal)(aodata,name,segment,0,0,data);
661: return(0);
662: }
664: /* ------------------------------------------------------------------------------------*/
668: /*@C
669: AODataKeyGetNeighbors - Given a list of keys generates a new list containing
670: those keys plus neighbors found in a neighbors list.
672: Collective on AOData
674: Input Parameters:
675: + aodata - the database
676: . name - the name of the key
677: . n - the number of data items needed by this processor
678: - keys - the keys provided by this processor
680: Output Parameters:
681: . is - the indices retrieved
683: Level: deprecated
685: .keywords: database transactions
687: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
688: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
689: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(),
690: AODataKeyGetNeighborsIS()
691: @*/
692: PetscErrorCode AODataKeyGetNeighbors(AOData aodata,const char name[],PetscInt n,PetscInt *keys,IS *is)
693: {
695: IS reduced,input;
699:
700: /* get the list of neighbors */
701: AODataSegmentGetReduced(aodata,name,name,n,keys,&reduced);
703: ISCreateGeneral(aodata->comm,n,keys,&input);
704: ISExpand(input,reduced,is);
705: ISDestroy(input);
706: ISDestroy(reduced);
708: return(0);
709: }
713: /*@C
714: AODataKeyGetNeighborsIS - Given a list of keys generates a new list containing
715: those keys plus neighbors found in a neighbors list.
717: Collective on AOData and IS
719: Input Parameters:
720: + aodata - the database
721: . name - the name of the key
722: . n - the number of data items needed by this processor
723: - keys - the keys provided by this processor
725: Output Parameters:
726: . is - the indices retrieved
728: Level: deprecated
730: .keywords: database transactions
732: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
733: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
734: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(),
735: AODataKeyGetNeighbors()
736: @*/
737: PetscErrorCode AODataKeyGetNeighborsIS(AOData aodata,const char name[],IS keys,IS *is)
738: {
740: IS reduced;
744:
745: /* get the list of neighbors */
746: AODataSegmentGetReducedIS(aodata,name,name,keys,&reduced);
747: /* combine keys and reduced is */
748: ISExpand(keys,reduced,is);
749: ISDestroy(reduced);
750: return(0);
751: }
755: /*@C
756: AODataSegmentGetReduced - Gets the unique list of segment values, by removing
757: duplicates.
759: Collective on AOData and IS
761: Input Parameters:
762: + aodata - the database
763: . name - the name of the key
764: . segment - the name of the segment
765: . n - the number of data items needed by this processor
766: - keys - the keys provided by this processor
768: Output Parameters:
769: . is - the indices retrieved
771: Level: deprecated
773: Example:
774: .vb
775: keys -> 0 1 2 3 4 5 6 7
776: if the segment contains -> 1 2 1 3 1 4 2 0
777: and you request keys 0 1 2 5 7 it will return 1 2 4 0
778: .ve
780: .keywords: database transactions
782: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
783: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
784: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
785: @*/
786: PetscErrorCode AODataSegmentGetReduced(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,IS *is)
787: {
792: (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,is);
793: return(0);
794: }
798: /*@C
799: AODataSegmentGetExtrema - Gets the largest and smallest values for each entry in the block
801: Collective on AOData
803: Input Parameters:
804: + aodata - the database
805: . name - the name of the key
806: - segment - the name of the segment
808: Output Parameters:
809: + vmax - the maximum values (user must provide enough space)
810: - vmin - the minimum values (user must provide enough space)
812: Level: deprecated
814: .keywords: database transactions
816: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
817: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
818: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
819: @*/
820: PetscErrorCode AODataSegmentGetExtrema(AOData aodata,const char name[],const char segment[],void *vmax,void *vmin)
821: {
826: (*aodata->ops->segmentgetextrema)(aodata,name,segment,vmax,vmin);
827: return(0);
828: }
832: /*@C
833: AODataSegmentGetReducedIS - Gets the unique list of segment values, by removing
834: duplicates.
836: Collective on AOData and IS
838: Input Parameters:
839: + aodata - the database
840: . name - the name of the key
841: . segment - the name of the segment
842: - is - the keys for data requested on this processor
844: Output Parameters:
845: . isout - the indices retreived
847: Level: deprecated
849: Example:
850: .vb
851: keys -> 0 1 2 3 4 5 6 7
852: if the segment contains -> 1 2 1 3 1 4 2 0
854: and you request keys 0 1 2 5 7, AODataSegmentGetReducedIS() will return 1 2 4 0
855: .ve
857: .keywords: database transactions
859: .seealso:
860: @*/
861: PetscErrorCode AODataSegmentGetReducedIS(AOData aodata,const char name[],const char segment[],IS is,IS *isout)
862: {
864: PetscInt n,*keys;
870: ISGetLocalSize(is,&n);
871: ISGetIndices(is,&keys);
872: (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,isout);
873: ISRestoreIndices(is,&keys);
874: return(0);
875: }
877: /* ------------------------------------------------------------------------------------*/
881: /*@C
882: AODataKeySetLocalToGlobalMapping - Add a local to global mapping for a key in the
883: in the database
885: Not collective
887: Input Parameters:
888: + aodata - the database
889: . name - the name of the key
890: - map - local to global mapping
892: Level: deprecated
894: .keywords: database additions
896: .seealso: AODataKeyGetLocalToGlobalMapping()
897: @*/
898: PetscErrorCode AODataKeySetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping map)
899: {
901: PetscTruth flag;
902: AODataKey *ikey;
907: AODataKeyFind_Private(aodata,name,&flag,&ikey);
908: if (!flag) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Key does not exist");
910: if (ikey->ltog) {
911: SETERRQ1(PETSC_ERR_ARG_WRONGSTATE,"Database key %s already has local to global mapping",name);
912: }
914: ikey->ltog = map;
915: PetscObjectReference((PetscObject)map);
917: return(0);
919: }
923: /*@C
924: AODataKeyGetLocalToGlobalMapping - gets a local to global mapping from a database
926: Not collective
928: Input Parameters:
929: + aodata - the database
930: - name - the name of the key
932: Output Parameters:
933: . map - local to global mapping
935: Level: deprecated
937: .keywords: database additions
939: .seealso: AODataKeySetLocalToGlobalMapping()
940: @*/
941: PetscErrorCode AODataKeyGetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping *map)
942: {
944: PetscTruth flag;
945: AODataKey *ikey;
950: AODataKeyFind_Private(aodata,name,&flag,&ikey);
951: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key does not exist: %s",name);
953: *map = ikey->ltog;
954: return(0);
955: }
960: /*@C
961: AODataKeyGetOwnershipRange - Gets the ownership range to this key type.
963: Not collective
965: Input Parameters:
966: + aodata - the database
967: - name - the name of the key
969: Output Parameters:
970: + rstart - first key owned locally (or PETSC_NULL if not needed)
971: - rend - last key owned locally + 1 (or PETSC_NULL if not needed)
973: Level: deprecated
975: .keywords: database accessing
977: .seealso: AODataKeyGetInfo()
978: @*/
979: PetscErrorCode AODataKeyGetOwnershipRange(AOData aodata,const char name[],PetscInt *rstart,PetscInt *rend)
980: {
982: PetscTruth flag;
983: AODataKey *key;
988: AODataKeyFind_Private(aodata,name,&flag,&key);
989: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key never created: %s",name);
991: if (rstart) *rstart = key->rstart;
992: if (rend) *rend = key->rend;
994: return(0);
995: }
999: /*@C
1000: AODataKeyGetInfo - Gets the global size, local size and number of segments in a key.
1002: Not collective
1004: Input Parameters:
1005: + aodata - the database
1006: - name - the name of the key
1008: Output Parameters:
1009: + nglobal - global number of keys
1010: . nlocal - local number of keys
1011: . nsegments - number of segments associated with key
1012: - segnames - names of the segments or PETSC_NULL
1014: Level: deprecated
1016: .keywords: database accessing
1018: .seealso: AODataKeyGetOwnershipRange()
1019: @*/
1020: PetscErrorCode AODataKeyGetInfo(AOData aodata,const char name[],PetscInt *nglobal,PetscInt *nlocal,PetscInt *nsegments,char ***segnames)
1021: {
1023: PetscInt i,n=0;
1024: AODataKey *key;
1025: AODataSegment *seg;
1026: PetscTruth flag;
1031: AODataKeyFind_Private(aodata,name,&flag,&key);
1032: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key never created: %s",name);
1034: if (nglobal) *nglobal = key->N;
1035: if (nlocal) *nlocal = key->nlocal;
1036: if (nsegments) *nsegments = n = key->nsegments;
1037: if (nsegments && segnames) {
1038: PetscMalloc((n+1)*sizeof(char *),&segnames);
1039: seg = key->segments;
1040: for (i=0; i<n; i++) {
1041: if (!seg) SETERRQ(PETSC_ERR_COR,"Less segments in database then indicated");
1042: (*segnames)[i] = seg->name;
1043: seg = seg->next;
1044: }
1045: }
1047: return(0);
1048: }
1052: /*@C
1053: AODataSegmentGetInfo - Gets the blocksize and type of a data segment
1055: Not collective
1057: Input Parameters:
1058: + aodata - the database
1059: . keyname - the name of the key
1060: - segname - the name of the segment
1062: Output Parameters:
1063: + bs - the blocksize
1064: - dtype - the datatype
1066: Level: deprecated
1068: .keywords: database accessing
1070: .seealso: AODataGetInfo()
1071: @*/
1072: PetscErrorCode AODataSegmentGetInfo(AOData aodata,const char keyname[],const char segname[],PetscInt *bs,PetscDataType *dtype)
1073: {
1075: PetscTruth flag;
1076: AODataKey *key;
1077: AODataSegment *seg;
1082: AODataSegmentFind_Private(aodata,keyname,segname,&flag,&key,&seg);
1083: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Segment never created: %s",segname);
1084: if (bs) *bs = seg->bs;
1085: if (dtype) *dtype = seg->datatype;
1087: return(0);
1088: }
1092: /*@C
1093: AODataView - Displays an application ordering.
1095: Collective on AOData and PetscViewer
1097: Input Parameters:
1098: + aodata - the database
1099: - viewer - viewer used for display
1101: Level: intermediate
1103: The available visualization contexts include
1104: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
1105: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
1106: output where only the first processor opens
1107: the file. All other processors send their
1108: data to the first processor to print.
1110: The user can open an alternative visualization context with
1111: PetscViewerASCIIOpen() - output to a specified file.
1114: .keywords: database viewing
1116: .seealso: PetscViewerASCIIOpen()
1117: @*/
1118: PetscErrorCode AODataView(AOData aodata,PetscViewer viewer)
1119: {
1124: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(aodata->comm);
1126: (*aodata->ops->view)(aodata,viewer);
1127: return(0);
1128: }
1132: static PetscErrorCode AODataAliasDestroy_Private(AODataAlias *aliases)
1133: {
1134: AODataAlias *t = aliases;
1138: if (t) {
1139: t = aliases->next;
1140: PetscFree(aliases->name);
1141: PetscFree(aliases->alias);
1142: PetscFree(aliases);
1143: while (t) {
1144: aliases = t;
1145: t = t->next;
1146: PetscFree(aliases->name);
1147: PetscFree(aliases->alias);
1148: PetscFree(aliases);
1149: }
1150: }
1151: return(0);
1152: }
1156: /*@C
1157: AODataAliasAdd - Man page needed.
1159: Level: deprecated
1161: @*/
1162: PetscErrorCode AODataAliasAdd(AOData aodata,const char alias[],const char name[])
1163: {
1164: AODataAlias *t = aodata->aliases;
1168: if (t) {
1169: while (t->next) t = t->next;
1170: PetscNew(AODataAlias,&t->next);
1171: t = t->next;
1172: } else {
1173: PetscNew(AODataAlias,&t);
1174: aodata->aliases = t;
1175: }
1176: PetscStrallocpy(alias,&t->alias);
1177: PetscStrallocpy(name,&t->name);
1178: t->next = 0;
1179: return(0);
1180: }
1184: /*@C
1185: AODataDestroy - Destroys an application ordering set.
1187: Collective on AOData
1189: Input Parameters:
1190: . aodata - the database
1192: Level: deprecated
1194: .keywords: destroy, database
1196: .seealso: AODataCreateBasic()
1197: @*/
1198: PetscErrorCode AODataDestroy(AOData aodata)
1199: {
1204: if (!aodata) return(0);
1206: if (--aodata->refct > 0) return(0);
1207:
1208: AODataAliasDestroy_Private(aodata->aliases);
1209: (*aodata->ops->destroy)(aodata);
1211: return(0);
1212: }
1216: /*@C
1217: AODataKeyRemap - Remaps a key and all references to a key to a new numbering
1218: scheme where each processor indicates its new nodes by listing them in the
1219: previous numbering scheme.
1221: Collective on AOData and AO
1223: Input Parameters:
1224: + aodata - the database
1225: . key - the key to remap
1226: - ao - the old to new ordering
1228: Level: deprecated
1230: .keywords: database remapping
1232: .seealso: AODataKeyGetAdjacency()
1233: @*/
1234: PetscErrorCode AODataKeyRemap(AOData aodata,const char key[],AO ao)
1235: {
1241: (*aodata->ops->keyremap)(aodata,key,ao);
1242: return(0);
1243: }
1247: /*@C
1248: AODataKeyGetAdjacency - Gets the adjacency graph for a key.
1250: Collective on AOData
1252: Input Parameters:
1253: + aodata - the database
1254: - key - the key
1256: Output Parameter:
1257: . adj - the adjacency graph
1259: Level: deprecated
1261: .keywords: database, adjacency graph
1263: .seealso: AODataKeyRemap()
1264: @*/
1265: PetscErrorCode AODataKeyGetAdjacency(AOData aodata,const char key[],Mat *adj)
1266: {
1271: (*aodata->ops->keygetadjacency)(aodata,key,adj);
1272: return(0);
1273: }
1277: /*@C
1278: AODataSegmentPartition - Partitions a segment type across processors
1279: relative to a key that is partitioned. This will try to keep as
1280: many elements of the segment on the same processor as corresponding
1281: neighboring key elements are.
1283: Collective on AOData
1285: Input Parameters:
1286: + aodata - the database
1287: - key - the key to be partitioned and renumbered
1289: Level: deprecated
1291: .seealso: AODataKeyPartition(), AODataPartitionAndSetupLocal()
1293: @*/
1294: PetscErrorCode AODataSegmentPartition(AOData aodata,const char key[],const char seg[])
1295: {
1300: (*aodata->ops->segmentpartition)(aodata,key,seg);
1301: return(0);
1302: }
1306: PetscErrorCode AODataPublish_Petsc(PetscObject obj)
1307: {
1309: return(0);
1310: }
1314: /*@C
1315: AODataKeyRemove - Remove a data key from a AOData database.
1317: Collective on AOData
1319: Input Parameters:
1320: + aodata - the database
1321: - name - the name of the key
1323: Level: deprecated
1325: .keywords: database removal
1327: .seealso:
1328: @*/
1329: PetscErrorCode AODataKeyRemove(AOData aodata,const char name[])
1330: {
1335: (*aodata->ops->keyremove)(aodata,name);
1336: return(0);
1337: }
1341: /*@C
1342: AODataSegmentRemove - Remove a data segment from a AOData database.
1344: Collective on AOData
1346: Input Parameters:
1347: + aodata - the database
1348: . name - the name of the key
1349: - segname - name of the segment
1351: Level: deprecated
1353: .keywords: database removal
1355: .seealso:
1356: @*/
1357: PetscErrorCode AODataSegmentRemove(AOData aodata,const char name[],const char segname[])
1358: {
1363: (*aodata->ops->segmentremove)(aodata,name,segname);
1364: return(0);
1365: }
1369: /*@C
1370: AODataKeyAdd - Add another data key to a AOData database.
1372: Collective on AOData
1374: Input Parameters:
1375: + aodata - the database
1376: . name - the name of the key
1377: . nlocal - number of indices to be associated with this processor
1378: - N - the number of indices in the key
1380: Level: deprecated
1382: .keywords: database additions
1384: .seealso:
1385: @*/
1386: PetscErrorCode AODataKeyAdd(AOData aodata,const char name[],PetscInt nlocal,PetscInt N)
1387: {
1389: PetscMPIInt size,rank;
1390: PetscInt i;
1391: size_t len;
1392: AODataKey *key,*oldkey;
1393: MPI_Comm comm = aodata->comm;
1394: PetscTruth flag;
1399: AODataKeyFind_Private(aodata,name,&flag,&oldkey);
1400: if (flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key already exists with given name: %s",name);
1402: PetscNew(AODataKey,&key);
1403: if (oldkey) { oldkey->next = key;}
1404: else { aodata->keys = key;}
1405: PetscStrlen(name,&len);
1406: PetscMalloc((len+1)*sizeof(char),&key->name);
1407: PetscStrcpy(key->name,name);
1408: key->N = N;
1409: key->nsegments = 0;
1410: key->segments = 0;
1411: key->ltog = 0;
1412: key->next = 0;
1414: MPI_Comm_rank(comm,&rank);
1415: MPI_Comm_size(comm,&size);
1417: /* Set nlocal and ownership ranges */
1418: PetscSplitOwnership(comm,&nlocal,&N);
1419: PetscMalloc((size+1)*sizeof(PetscInt),&key->rowners);
1420: MPI_Allgather(&nlocal,1,MPIU_INT,key->rowners+1,1,MPIU_INT,comm);
1421: key->rowners[0] = 0;
1422: for (i=2; i<=size; i++) {
1423: key->rowners[i] += key->rowners[i-1];
1424: }
1425: key->rstart = key->rowners[rank];
1426: key->rend = key->rowners[rank+1];
1428: key->nlocal = nlocal;
1430: aodata->nkeys++;
1431: return(0);
1432: }
1436: /*@C
1437: AODataSegmentAdd - Adds another data segment to a AOData database.
1439: Collective on AOData
1441: Input Parameters:
1442: + aodata - the database
1443: . name - the name of the key
1444: . segment - the name of the data segment
1445: . bs - the fundamental blocksize of the data
1446: . n - the number of data items contributed by this processor
1447: . keys - the keys provided by this processor
1448: . data - the actual data
1449: - dtype - the data type (one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.)
1451: Level: deprecated
1453: .keywords: database additions
1455: .seealso: AODataSegmentAddIS()
1456: @*/
1457: PetscErrorCode AODataSegmentAdd(AOData aodata,const char name[],const char segment[],PetscInt bs,PetscInt n,PetscInt *keys,void *data,PetscDataType dtype)
1458: {
1464: (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1466: /*
1467: PetscOptionsHasName(PETSC_NULL,"-ao_data_view",&flg1);
1468: if (flg1) {
1469: AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1470: }
1471: PetscOptionsHasName(PETSC_NULL,"-ao_data_view_info",&flg1);
1472: if (flg1) {
1473: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(comm),PETSC_VIEWER_ASCII_INFO);
1474: AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1475: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(comm));
1476: }
1477: */
1478: return(0);
1479: }
1483: /*@C
1484: AODataSegmentAddIS - Add another data segment to a AOData database.
1486: Collective on AOData and IS
1488: Input Parameters:
1489: + aodata - the database
1490: . name - the name of the key
1491: . segment - name of segment
1492: . bs - the fundamental blocksize of the data
1493: . is - the keys provided by this processor
1494: . data - the actual data
1495: - dtype - the data type, one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.
1497: Level: deprecated
1499: .keywords: database additions
1501: .seealso: AODataSegmentAdd()
1502: @*/
1503: PetscErrorCode AODataSegmentAddIS(AOData aodata,const char name[],const char segment[],PetscInt bs,IS is,void *data,PetscDataType dtype)
1504: {
1506: PetscInt n,*keys;
1512: ISGetLocalSize(is,&n);
1513: ISGetIndices(is,&keys);
1514: (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1515: ISRestoreIndices(is,&keys);
1516: return(0);
1517: }