Actual source code: aodata.c
1: /*
2: Defines the abstract operations on AOData
3: */
4: #include src/dm/ao/aoimpl.h
9: /*@C
10: AODataGetInfo - Gets the number of keys and their names in a database.
12: Not collective
14: Input Parameter:
15: . ao - the AOData database
17: Output Parameters:
18: + nkeys - the number of keys
19: - keys - the names of the keys (or PETSC_NULL)
21: Level: advanced
23: .keywords: application ordering
25: .seealso: AODataSegmentGetInfo()
26: @*/
27: PetscErrorCode AODataGetInfo(AOData ao,PetscInt *nkeys,char ***keys)
28: {
30: PetscInt n,i;
31: AODataKey *key = ao->keys;
36: *nkeys = n = ao->nkeys;
37: if (keys) {
38: PetscMalloc((n+1)*sizeof(char *),&keys);
39: for (i=0; i<n; i++) {
40: if (!key) SETERRQ(PETSC_ERR_COR,"Less keys in database then indicated");
41: (*keys)[i] = key->name;
42: key = key->next;
43: }
44: }
45: return(0);
46: }
50: /*
51: AODataKeyFind_Private - Given a keyname finds the key. Generates a flag if not found.
53: Not collective
55: Input Parameters:
56: . keyname - string name of key
58: Output Parameter:
59: + flag - PETSC_TRUE if found, PETSC_FALSE if not found
60: - key - the associated key
62: Level: advanced
64: */
65: PetscErrorCode AODataKeyFind_Private(AOData aodata,const char keyname[],PetscTruth *flag,AODataKey **key)
66: {
67: PetscTruth match;
69: AODataAlias *t = aodata->aliases;
70: const char *name = keyname;
71: AODataKey *nkey;
74: *key = PETSC_NULL;
75: *flag = PETSC_FALSE;
76: while (name) {
77: nkey = aodata->keys;
78: while (nkey) {
79: PetscStrcmp(nkey->name,name,&match);
80: if (match) {
81: /* found the key */
82: *key = nkey;
83: *flag = PETSC_TRUE;
84: return(0);
85: }
86: *key = nkey;
87: nkey = nkey->next;
88: }
89: name = 0;
90: while (t) {
91: PetscStrcmp(keyname,t->alias,&match);
92: if (match) {
93: name = t->name;
94: t = t->next;
95: break;
96: }
97: t = t->next;
98: }
99: }
100: return(0);
101: }
105: /*@C
106: AODataKeyExists - Determines if a key exists in the database.
108: Not collective
110: Input Parameters:
111: . keyname - string name of key
113: Output Parameter:
114: . flag - PETSC_TRUE if found, otherwise PETSC_FALSE
116: Level: advanced
118: @*/
119: PetscErrorCode AODataKeyExists(AOData aodata,const char keyname[],PetscTruth *flag)
120: {
122: PetscTruth iflag;
123: AODataKey *ikey;
127: AODataKeyFind_Private(aodata,keyname,&iflag,&ikey);
128: if (iflag) *flag = PETSC_TRUE;
129: else *flag = PETSC_FALSE;
130: return(0);
131: }
136: /*
137: AODataSegmentFind_Private - Given a key and segment finds the int key, segment
138: coordinates. Generates a flag if not found.
140: Not collective
142: Input Parameters:
143: + keyname - string name of key
144: - segname - string name of segment
146: Output Parameter:
147: + flag - PETSC_TRUE if found, PETSC_FALSE if not
148: . key - integer of keyname
149: - segment - integer of segment
151: If it doesn't find it, it returns the last seg in the key (if the key exists)
153: Level: advanced
155: */
156: PetscErrorCode AODataSegmentFind_Private(AOData aodata,const char keyname[],const char segname[],PetscTruth *flag,AODataKey **key,AODataSegment **seg)
157: {
159: PetscTruth keyflag,match;
160: AODataAlias *t = aodata->aliases;
161: const char *name;
162: AODataSegment *nseg;
165: *seg = PETSC_NULL;
166: *flag = PETSC_FALSE;
167: AODataKeyFind_Private(aodata,keyname,&keyflag,key);
168: if (keyflag) { /* found key now look for segment */
169: name = segname;
170: while (name) {
171: nseg = (*key)->segments;
172: while (nseg) {
173: PetscStrcmp(nseg->name,name,&match);
174: if (match) {
175: /* found the segment */
176: *seg = nseg;
177: *flag = PETSC_TRUE;
178: return(0);
179: }
180: *seg = nseg;
181: nseg = nseg->next;
182: }
183: name = 0;
184: while (t) {
185: PetscStrcmp(segname,t->alias,&match);
186: if (match) {
187: name = t->name;
188: t = t->next;
189: break;
190: }
191: t = t->next;
192: }
193: }
194: }
195: return(0);
196: }
200: /*@C
201: AODataSegmentExists - Determines if a key and segment exists in the database.
203: Not collective
205: Input Parameters:
206: + keyname - string name of key
207: - segname - string name of segment
209: Output Parameter:
210: . flag - PETSC_TRUE if found, else PETSC_FALSE
212: Level: advanced
214: @*/
215: PetscErrorCode AODataSegmentExists(AOData aodata,const char keyname[],const char segname[],PetscTruth *flag)
216: {
218: PetscTruth iflag;
219: AODataKey *ikey;
220: AODataSegment *iseg;
224: AODataSegmentFind_Private(aodata,keyname,segname,&iflag,&ikey,&iseg);
225: if (iflag) *flag = PETSC_TRUE;
226: else *flag = PETSC_FALSE;
227: return(0);
228: }
230: /* ------------------------------------------------------------------------------------*/
234: /*@C
235: AODataKeyGetActive - Get a sublist of key indices that have a logical flag on.
237: Collective on AOData
239: Input Parameters:
240: + aodata - the database
241: . name - the name of the key
242: . segment - the name of the segment
243: . n - the number of key indices provided by this processor
244: . keys - the keys provided by this processor
245: - wl - which logical key in the block (for block size 1 this is always 0)
247: Output Parameters:
248: . is - the list of key indices
250: Level: advanced
252: .keywords: database transactions
254: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
255: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
256: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
257: @*/
258: PetscErrorCode AODataKeyGetActive(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,PetscInt wl,IS *is)
259: {
264: (*aodata->ops->keygetactive)(aodata,name,segment,n,keys,wl,is);
265: return(0);
266: }
270: /*@C
271: AODataKeyGetActiveIS - Get a sublist of key indices that have a logical flag on.
273: Collective on AOData
275: Input Parameters:
276: + aodata - the database
277: . name - the name of the key
278: . segment - the name of the segment
279: . in - the key indices we are checking
280: - wl - which logical key in the block (for block size 1 this is always 0)
282: Output Parameters:
283: . IS - the list of key indices
285: Level: advanced
287: .keywords: database transactions
289: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
290: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
291: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
292: @*/
293: PetscErrorCode AODataKeyGetActiveIS(AOData aodata,const char name[],const char segname[],IS in,PetscInt wl,IS *is)
294: {
296: PetscInt n,*keys;
299: ISGetLocalSize(in,&n);
300: ISGetIndices(in,&keys);
301: AODataKeyGetActive(aodata,name,segname,n,keys,wl,is);
302: ISRestoreIndices(in,&keys);
303: return(0);
304: }
308: /*@C
309: AODataKeyGetActiveLocal - Get a sublist of key indices that have a logical flag on.
311: Collective on AOData
313: Input Parameters:
314: + aodata - the database
315: . name - the name of the key
316: . segment - the name of the segment
317: . n - the number of key indices provided by this processor
318: . keys - the keys provided by this processor
319: - wl - which logical key in the block (for block size 1 this is always 0)
321: Output Parameters:
322: . IS - the list of key indices
324: Level: advanced
326: .keywords: database transactions
328: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
329: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
330: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
331: @*/
332: PetscErrorCode AODataKeyGetActiveLocal(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,PetscInt wl,IS *is)
333: {
338: (*aodata->ops->keygetactivelocal)(aodata,name,segment,n,keys,wl,is);
339: return(0);
340: }
344: /*@C
345: AODataKeyGetActiveLocalIS - Get a sublist of key indices that have a logical flag on.
347: Collective on AOData
349: Input Parameters:
350: + aodata - the database
351: . name - the name of the key
352: . segment - the name of the segment
353: . in - the key indices we are checking
354: - wl - which logical key in the block (for block size 1 this is always 0)
356: Output Parameters:
357: . IS - the list of key indices
359: Level: advanced
361: .keywords: database transactions
363: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
364: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
365: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
366: @*/
367: PetscErrorCode AODataKeyGetActiveLocalIS(AOData aodata,const char name[],const char segname[],IS in,PetscInt wl,IS *is)
368: {
370: PetscInt n,*keys;
373: ISGetLocalSize(in,&n);
374: ISGetIndices(in,&keys);
375: AODataKeyGetActiveLocal(aodata,name,segname,n,keys,wl,is);
376: ISRestoreIndices(in,&keys);
377: return(0);
378: }
380: /* ------------------------------------------------------------------------------------*/
384: /*@C
385: AODataSegmentGet - Get data from a particular segment of a database.
387: Collective on AOData
389: Input Parameters:
390: + aodata - the database
391: . name - the name of the key
392: . segment - the name of the segment
393: . n - the number of data items needed by this processor
394: - keys - the keys provided by this processor
396: Output Parameters:
397: . data - the actual data
399: Level: advanced
401: .keywords: database transactions
403: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
404: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
405: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
406: @*/
407: PetscErrorCode AODataSegmentGet(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
408: {
413: (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
414: return(0);
415: }
419: /*@C
420: AODataSegmentRestore - Restores data from a particular segment of a database.
422: Collective on AOData
424: Input Parameters:
425: + aodata - the database
426: . name - the name of the key
427: . segment - the name of the segment
428: . n - the number of data items needed by this processor
429: - keys - the keys provided by this processor
431: Output Parameters:
432: . data - the actual data
434: Level: advanced
436: .keywords: database transactions
438: .seealso: AODataSegmentRestoreIS()
439: @*/
440: PetscErrorCode AODataSegmentRestore(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
441: {
446: (*aodata->ops->segmentrestore)(aodata,name,segment,n,keys,data);
447: return(0);
448: }
452: /*@C
453: AODataSegmentGetIS - Get data from a particular segment of a database.
455: Collective on AOData and IS
457: Input Parameters:
458: + aodata - the database
459: . name - the name of the key
460: . segment - the name of the segment
461: - is - the keys for data requested on this processor
463: Output Parameters:
464: . data - the actual data
466: Level: advanced
468: .keywords: database transactions
470: @*/
471: PetscErrorCode AODataSegmentGetIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
472: {
474: PetscInt n,*keys;
480: ISGetLocalSize(is,&n);
481: ISGetIndices(is,&keys);
482: (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
483: ISRestoreIndices(is,&keys);
484: return(0);
485: }
489: /*@C
490: AODataSegmentRestoreIS - Restores data from a particular segment of a database.
492: Collective on AOData and IS
494: Input Parameters:
495: + aodata - the database
496: . name - the name of the data key
497: . segment - the name of the segment
498: - is - the keys provided by this processor
500: Output Parameters:
501: . data - the actual data
503: Level: advanced
505: .keywords: database transactions
507: .seealso: AODataSegmentRestore()
508: @*/
509: PetscErrorCode AODataSegmentRestoreIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
510: {
516: (*aodata->ops->segmentrestore)(aodata,name,segment,0,0,data);
518: return(0);
519: }
521: /* ------------------------------------------------------------------------------------*/
524: /*@C
525: AODataSegmentGetLocal - Get data from a particular segment of a database. Returns the
526: values in the local numbering; valid only for integer segments.
528: Collective on AOData
530: Input Parameters:
531: + aodata - the database
532: . name - the name of the key
533: . segment - the name of the segment
534: . n - the number of data items needed by this processor
535: - keys - the keys provided by this processor in local numbering
537: Output Parameters:
538: . data - the actual data
540: Level: advanced
542: .keywords: database transactions
544: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
545: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
546: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
547: @*/
548: PetscErrorCode AODataSegmentGetLocal(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
549: {
554: (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
555: return(0);
556: }
560: /*@C
561: AODataSegmentRestoreLocal - Restores data from a particular segment of a database.
563: Collective on AOData
565: Input Parameters:
566: + aodata - the database
567: . name - the name of the key
568: . segment - the name of the segment
569: . n - the number of data items needed by this processor
570: - keys - the keys provided by this processor
572: Output Parameters:
573: . data - the actual data
575: Level: advanced
577: .keywords: database transactions
579: @*/
580: PetscErrorCode AODataSegmentRestoreLocal(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,void **data)
581: {
586: (*aodata->ops->segmentrestorelocal)(aodata,name,segment,n,keys,data);
587: return(0);
588: }
592: /*@C
593: AODataSegmentGetLocalIS - Get data from a particular segment of a database. Returns the
594: values in the local numbering; valid only for integer segments.
596: Collective on AOData and IS
598: Input Parameters:
599: + aodata - the database
600: . name - the name of the key
601: . segment - the name of the segment
602: - is - the keys for data requested on this processor
604: Output Parameters:
605: . data - the actual data
607: Level: advanced
609: .keywords: database transactions
611: .seealso: AODataSegmentRestoreLocalIS()
612: @*/
613: PetscErrorCode AODataSegmentGetLocalIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
614: {
616: PetscInt n,*keys;
622: ISGetLocalSize(is,&n);
623: ISGetIndices(is,&keys);
624: (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
625: ISRestoreIndices(is,&keys);
626: return(0);
627: }
631: /*@C
632: AODataSegmentRestoreLocalIS - Restores data from a particular segment of a database.
634: Collective on AOData and IS
636: Input Parameters:
637: + aodata - the database
638: . name - the name of the data key
639: . segment - the name of the segment
640: - is - the keys provided by this processor
642: Output Parameters:
643: . data - the actual data
645: Level: advanced
647: .keywords: database transactions
649: .seealso: AODataSegmentGetLocalIS()
650: @*/
651: PetscErrorCode AODataSegmentRestoreLocalIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
652: {
658: (*aodata->ops->segmentrestorelocal)(aodata,name,segment,0,0,data);
659: return(0);
660: }
662: /* ------------------------------------------------------------------------------------*/
666: /*@C
667: AODataKeyGetNeighbors - Given a list of keys generates a new list containing
668: those keys plus neighbors found in a neighbors list.
670: Collective on AOData
672: Input Parameters:
673: + aodata - the database
674: . name - the name of the key
675: . n - the number of data items needed by this processor
676: - keys - the keys provided by this processor
678: Output Parameters:
679: . is - the indices retrieved
681: Level: advanced
683: .keywords: database transactions
685: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
686: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
687: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(),
688: AODataKeyGetNeighborsIS()
689: @*/
690: PetscErrorCode AODataKeyGetNeighbors(AOData aodata,const char name[],PetscInt n,PetscInt *keys,IS *is)
691: {
693: IS reduced,input;
697:
698: /* get the list of neighbors */
699: AODataSegmentGetReduced(aodata,name,name,n,keys,&reduced);
701: ISCreateGeneral(aodata->comm,n,keys,&input);
702: ISExpand(input,reduced,is);
703: ISDestroy(input);
704: ISDestroy(reduced);
706: return(0);
707: }
711: /*@C
712: AODataKeyGetNeighborsIS - Given a list of keys generates a new list containing
713: those keys plus neighbors found in a neighbors list.
715: Collective on AOData and IS
717: Input Parameters:
718: + aodata - the database
719: . name - the name of the key
720: . n - the number of data items needed by this processor
721: - keys - the keys provided by this processor
723: Output Parameters:
724: . is - the indices retrieved
726: Level: advanced
728: .keywords: database transactions
730: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
731: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
732: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(),
733: AODataKeyGetNeighbors()
734: @*/
735: PetscErrorCode AODataKeyGetNeighborsIS(AOData aodata,const char name[],IS keys,IS *is)
736: {
738: IS reduced;
742:
743: /* get the list of neighbors */
744: AODataSegmentGetReducedIS(aodata,name,name,keys,&reduced);
745: /* combine keys and reduced is */
746: ISExpand(keys,reduced,is);
747: ISDestroy(reduced);
748: return(0);
749: }
753: /*@C
754: AODataSegmentGetReduced - Gets the unique list of segment values, by removing
755: duplicates.
757: Collective on AOData and IS
759: Input Parameters:
760: + aodata - the database
761: . name - the name of the key
762: . segment - the name of the segment
763: . n - the number of data items needed by this processor
764: - keys - the keys provided by this processor
766: Output Parameters:
767: . is - the indices retrieved
769: Level: advanced
771: Example:
772: .vb
773: keys -> 0 1 2 3 4 5 6 7
774: if the segment contains -> 1 2 1 3 1 4 2 0
775: and you request keys 0 1 2 5 7 it will return 1 2 4 0
776: .ve
778: .keywords: database transactions
780: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
781: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
782: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
783: @*/
784: PetscErrorCode AODataSegmentGetReduced(AOData aodata,const char name[],const char segment[],PetscInt n,PetscInt *keys,IS *is)
785: {
790: (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,is);
791: return(0);
792: }
796: /*@C
797: AODataSegmentGetExtrema - Gets the largest and smallest values for each entry in the block
799: Collective on AOData
801: Input Parameters:
802: + aodata - the database
803: . name - the name of the key
804: - segment - the name of the segment
806: Output Parameters:
807: + vmax - the maximum values (user must provide enough space)
808: - vmin - the minimum values (user must provide enough space)
810: Level: advanced
812: .keywords: database transactions
814: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
815: AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(),
816: AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
817: @*/
818: PetscErrorCode AODataSegmentGetExtrema(AOData aodata,const char name[],const char segment[],void *vmax,void *vmin)
819: {
824: (*aodata->ops->segmentgetextrema)(aodata,name,segment,vmax,vmin);
825: return(0);
826: }
830: /*@C
831: AODataSegmentGetReducedIS - Gets the unique list of segment values, by removing
832: duplicates.
834: Collective on AOData and IS
836: Input Parameters:
837: + aodata - the database
838: . name - the name of the key
839: . segment - the name of the segment
840: - is - the keys for data requested on this processor
842: Output Parameters:
843: . isout - the indices retreived
845: Level: advanced
847: Example:
848: .vb
849: keys -> 0 1 2 3 4 5 6 7
850: if the segment contains -> 1 2 1 3 1 4 2 0
852: and you request keys 0 1 2 5 7, AODataSegmentGetReducedIS() will return 1 2 4 0
853: .ve
855: .keywords: database transactions
857: .seealso:
858: @*/
859: PetscErrorCode AODataSegmentGetReducedIS(AOData aodata,const char name[],const char segment[],IS is,IS *isout)
860: {
862: PetscInt n,*keys;
868: ISGetLocalSize(is,&n);
869: ISGetIndices(is,&keys);
870: (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,isout);
871: ISRestoreIndices(is,&keys);
872: return(0);
873: }
875: /* ------------------------------------------------------------------------------------*/
879: /*@C
880: AODataKeySetLocalToGlobalMapping - Add a local to global mapping for a key in the
881: in the database
883: Not collective
885: Input Parameters:
886: + aodata - the database
887: . name - the name of the key
888: - map - local to global mapping
890: Level: advanced
892: .keywords: database additions
894: .seealso: AODataKeyGetLocalToGlobalMapping()
895: @*/
896: PetscErrorCode AODataKeySetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping map)
897: {
899: PetscTruth flag;
900: AODataKey *ikey;
905: AODataKeyFind_Private(aodata,name,&flag,&ikey);
906: if (!flag) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Key does not exist");
908: if (ikey->ltog) {
909: SETERRQ1(PETSC_ERR_ARG_WRONGSTATE,"Database key %s already has local to global mapping",name);
910: }
912: ikey->ltog = map;
913: PetscObjectReference((PetscObject)map);
915: return(0);
917: }
921: /*@C
922: AODataKeyGetLocalToGlobalMapping - gets a local to global mapping from a database
924: Not collective
926: Input Parameters:
927: + aodata - the database
928: - name - the name of the key
930: Output Parameters:
931: . map - local to global mapping
933: Level: advanced
935: .keywords: database additions
937: .seealso: AODataKeySetLocalToGlobalMapping()
938: @*/
939: PetscErrorCode AODataKeyGetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping *map)
940: {
942: PetscTruth flag;
943: AODataKey *ikey;
948: AODataKeyFind_Private(aodata,name,&flag,&ikey);
949: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key does not exist: %s",name);
951: *map = ikey->ltog;
952: return(0);
953: }
958: /*@C
959: AODataKeyGetOwnershipRange - Gets the ownership range to this key type.
961: Not collective
963: Input Parameters:
964: + aodata - the database
965: - name - the name of the key
967: Output Parameters:
968: + rstart - first key owned locally (or PETSC_NULL if not needed)
969: - rend - last key owned locally + 1 (or PETSC_NULL if not needed)
971: Level: advanced
973: .keywords: database accessing
975: .seealso: AODataKeyGetInfo()
976: @*/
977: PetscErrorCode AODataKeyGetOwnershipRange(AOData aodata,const char name[],PetscInt *rstart,PetscInt *rend)
978: {
980: PetscTruth flag;
981: AODataKey *key;
986: AODataKeyFind_Private(aodata,name,&flag,&key);
987: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key never created: %s",name);
989: if (rstart) *rstart = key->rstart;
990: if (rend) *rend = key->rend;
992: return(0);
993: }
997: /*@C
998: AODataKeyGetInfo - Gets the global size, local size and number of segments in a key.
1000: Not collective
1002: Input Parameters:
1003: + aodata - the database
1004: - name - the name of the key
1006: Output Parameters:
1007: + nglobal - global number of keys
1008: . nlocal - local number of keys
1009: . nsegments - number of segments associated with key
1010: - segnames - names of the segments or PETSC_NULL
1012: Level: advanced
1014: .keywords: database accessing
1016: .seealso: AODataKeyGetOwnershipRange()
1017: @*/
1018: PetscErrorCode AODataKeyGetInfo(AOData aodata,const char name[],PetscInt *nglobal,PetscInt *nlocal,PetscInt *nsegments,char ***segnames)
1019: {
1021: PetscInt i,n=0;
1022: AODataKey *key;
1023: AODataSegment *seg;
1024: PetscTruth flag;
1029: AODataKeyFind_Private(aodata,name,&flag,&key);
1030: if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key never created: %s",name);
1032: if (nglobal) *nglobal = key->N;
1033: if (nlocal) *nlocal = key->nlocal;
1034: if (nsegments) *nsegments = n = key->nsegments;
1035: if (nsegments && segnames) {
1036: PetscMalloc((n+1)*sizeof(char *),&segnames);
1037: seg = key->segments;
1038: for (i=0; i<n; i++) {
1039: if (!seg) SETERRQ(PETSC_ERR_COR,"Less segments in database then indicated");
1040: (*segnames)[i] = seg->name;
1041: seg = seg->next;
1042: }
1043: }
1045: return(0);
1046: }
1050: /*@C
1051: AODataSegmentGetInfo - Gets the blocksize and type of a data segment
1053: Not collective
1055: Input Parameters:
1056: + aodata - the database
1057: . keyname - the name of the key
1058: - segname - the name of the segment
1060: Output Parameters:
1061: + bs - the blocksize
1062: - dtype - the datatype
1064: Level: advanced
1066: .keywords: database accessing
1068: .seealso: AODataGetInfo()
1069: @*/
1070: PetscErrorCode AODataSegmentGetInfo(AOData aodata,const char keyname[],const char segname[],PetscInt *bs,PetscDataType *dtype)
1071: {
1073: PetscTruth flag;
1074: AODataKey *key;
1075: AODataSegment *seg;
1080: AODataSegmentFind_Private(aodata,keyname,segname,&flag,&key,&seg);
1081: if (flag == PETSC_FALSE) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Segment never created: %s",segname);
1082: if (bs) *bs = seg->bs;
1083: if (dtype) *dtype = seg->datatype;
1085: return(0);
1086: }
1090: /*@C
1091: AODataView - Displays an application ordering.
1093: Collective on AOData and PetscViewer
1095: Input Parameters:
1096: + aodata - the database
1097: - viewer - viewer used for display
1099: Level: intermediate
1101: The available visualization contexts include
1102: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
1103: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
1104: output where only the first processor opens
1105: the file. All other processors send their
1106: data to the first processor to print.
1108: The user can open an alternative visualization context with
1109: PetscViewerASCIIOpen() - output to a specified file.
1112: .keywords: database viewing
1114: .seealso: PetscViewerASCIIOpen()
1115: @*/
1116: PetscErrorCode AODataView(AOData aodata,PetscViewer viewer)
1117: {
1122: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(aodata->comm);
1124: (*aodata->ops->view)(aodata,viewer);
1125: return(0);
1126: }
1130: static PetscErrorCode AODataAliasDestroy_Private(AODataAlias *aliases)
1131: {
1132: AODataAlias *t = aliases;
1136: if (t) {
1137: t = aliases->next;
1138: PetscFree(aliases->name);
1139: PetscFree(aliases->alias);
1140: PetscFree(aliases);
1141: while (t) {
1142: aliases = t;
1143: t = t->next;
1144: PetscFree(aliases->name);
1145: PetscFree(aliases->alias);
1146: PetscFree(aliases);
1147: }
1148: }
1149: return(0);
1150: }
1154: PetscErrorCode AODataAliasAdd(AOData aodata,const char alias[],const char name[])
1155: {
1156: AODataAlias *t = aodata->aliases;
1160: if (t) {
1161: while (t->next) t = t->next;
1162: PetscNew(AODataAlias,&t->next);
1163: t = t->next;
1164: } else {
1165: PetscNew(AODataAlias,&t);
1166: aodata->aliases = t;
1167: }
1168: PetscStrallocpy(alias,&t->alias);
1169: PetscStrallocpy(name,&t->name);
1170: t->next = 0;
1171: return(0);
1172: }
1176: /*@C
1177: AODataDestroy - Destroys an application ordering set.
1179: Collective on AOData
1181: Input Parameters:
1182: . aodata - the database
1184: Level: intermediate
1186: .keywords: destroy, database
1188: .seealso: AODataCreateBasic()
1189: @*/
1190: PetscErrorCode AODataDestroy(AOData aodata)
1191: {
1196: if (!aodata) return(0);
1198: if (--aodata->refct > 0) return(0);
1199:
1200: AODataAliasDestroy_Private(aodata->aliases);
1201: (*aodata->ops->destroy)(aodata);
1203: return(0);
1204: }
1208: /*@C
1209: AODataKeyRemap - Remaps a key and all references to a key to a new numbering
1210: scheme where each processor indicates its new nodes by listing them in the
1211: previous numbering scheme.
1213: Collective on AOData and AO
1215: Input Parameters:
1216: + aodata - the database
1217: . key - the key to remap
1218: - ao - the old to new ordering
1220: Level: advanced
1222: .keywords: database remapping
1224: .seealso: AODataKeyGetAdjacency()
1225: @*/
1226: PetscErrorCode AODataKeyRemap(AOData aodata,const char key[],AO ao)
1227: {
1233: (*aodata->ops->keyremap)(aodata,key,ao);
1234: return(0);
1235: }
1239: /*@C
1240: AODataKeyGetAdjacency - Gets the adjacency graph for a key.
1242: Collective on AOData
1244: Input Parameters:
1245: + aodata - the database
1246: - key - the key
1248: Output Parameter:
1249: . adj - the adjacency graph
1251: Level: advanced
1253: .keywords: database, adjacency graph
1255: .seealso: AODataKeyRemap()
1256: @*/
1257: PetscErrorCode AODataKeyGetAdjacency(AOData aodata,const char key[],Mat *adj)
1258: {
1263: (*aodata->ops->keygetadjacency)(aodata,key,adj);
1264: return(0);
1265: }
1269: /*@C
1270: AODataSegmentPartition - Partitions a segment type across processors
1271: relative to a key that is partitioned. This will try to keep as
1272: many elements of the segment on the same processor as corresponding
1273: neighboring key elements are.
1275: Collective on AOData
1277: Input Parameters:
1278: + aodata - the database
1279: - key - the key to be partitioned and renumbered
1281: Level: advanced
1283: .seealso: AODataKeyPartition(), AODataPartitionAndSetupLocal()
1285: @*/
1286: PetscErrorCode AODataSegmentPartition(AOData aodata,const char key[],const char seg[])
1287: {
1292: (*aodata->ops->segmentpartition)(aodata,key,seg);
1293: return(0);
1294: }
1298: PetscErrorCode AODataPublish_Petsc(PetscObject obj)
1299: {
1301: return(0);
1302: }
1306: /*@C
1307: AODataKeyRemove - Remove a data key from a AOData database.
1309: Collective on AOData
1311: Input Parameters:
1312: + aodata - the database
1313: - name - the name of the key
1315: Level: advanced
1317: .keywords: database removal
1319: .seealso:
1320: @*/
1321: PetscErrorCode AODataKeyRemove(AOData aodata,const char name[])
1322: {
1327: (*aodata->ops->keyremove)(aodata,name);
1328: return(0);
1329: }
1333: /*@C
1334: AODataSegmentRemove - Remove a data segment from a AOData database.
1336: Collective on AOData
1338: Input Parameters:
1339: + aodata - the database
1340: . name - the name of the key
1341: - segname - name of the segment
1343: Level: advanced
1345: .keywords: database removal
1347: .seealso:
1348: @*/
1349: PetscErrorCode AODataSegmentRemove(AOData aodata,const char name[],const char segname[])
1350: {
1355: (*aodata->ops->segmentremove)(aodata,name,segname);
1356: return(0);
1357: }
1361: /*@C
1362: AODataKeyAdd - Add another data key to a AOData database.
1364: Collective on AOData
1366: Input Parameters:
1367: + aodata - the database
1368: . name - the name of the key
1369: . nlocal - number of indices to be associated with this processor
1370: - N - the number of indices in the key
1372: Level: advanced
1374: .keywords: database additions
1376: .seealso:
1377: @*/
1378: PetscErrorCode AODataKeyAdd(AOData aodata,const char name[],PetscInt nlocal,PetscInt N)
1379: {
1381: PetscMPIInt size,rank;
1382: PetscInt i;
1383: size_t len;
1384: AODataKey *key,*oldkey;
1385: MPI_Comm comm = aodata->comm;
1386: PetscTruth flag;
1391: AODataKeyFind_Private(aodata,name,&flag,&oldkey);
1392: if (flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key already exists with given name: %s",name);
1394: PetscNew(AODataKey,&key);
1395: if (oldkey) { oldkey->next = key;}
1396: else { aodata->keys = key;}
1397: PetscStrlen(name,&len);
1398: PetscMalloc((len+1)*sizeof(char),&key->name);
1399: PetscStrcpy(key->name,name);
1400: key->N = N;
1401: key->nsegments = 0;
1402: key->segments = 0;
1403: key->ltog = 0;
1404: key->next = 0;
1406: MPI_Comm_rank(comm,&rank);
1407: MPI_Comm_size(comm,&size);
1409: /* Set nlocal and ownership ranges */
1410: PetscSplitOwnership(comm,&nlocal,&N);
1411: PetscMalloc((size+1)*sizeof(PetscInt),&key->rowners);
1412: MPI_Allgather(&nlocal,1,MPIU_INT,key->rowners+1,1,MPIU_INT,comm);
1413: key->rowners[0] = 0;
1414: for (i=2; i<=size; i++) {
1415: key->rowners[i] += key->rowners[i-1];
1416: }
1417: key->rstart = key->rowners[rank];
1418: key->rend = key->rowners[rank+1];
1420: key->nlocal = nlocal;
1422: aodata->nkeys++;
1423: return(0);
1424: }
1428: /*@C
1429: AODataSegmentAdd - Adds another data segment to a AOData database.
1431: Collective on AOData
1433: Input Parameters:
1434: + aodata - the database
1435: . name - the name of the key
1436: . segment - the name of the data segment
1437: . bs - the fundamental blocksize of the data
1438: . n - the number of data items contributed by this processor
1439: . keys - the keys provided by this processor
1440: . data - the actual data
1441: - dtype - the data type (one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.)
1443: Level: advanced
1445: .keywords: database additions
1447: .seealso: AODataSegmentAddIS()
1448: @*/
1449: PetscErrorCode AODataSegmentAdd(AOData aodata,const char name[],const char segment[],PetscInt bs,PetscInt n,PetscInt *keys,void *data,PetscDataType dtype)
1450: {
1456: (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1458: /*
1459: PetscOptionsHasName(PETSC_NULL,"-ao_data_view",&flg1);
1460: if (flg1) {
1461: AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1462: }
1463: PetscOptionsHasName(PETSC_NULL,"-ao_data_view_info",&flg1);
1464: if (flg1) {
1465: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(comm),PETSC_VIEWER_ASCII_INFO);
1466: AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1467: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(comm));
1468: }
1469: */
1470: return(0);
1471: }
1475: /*@C
1476: AODataSegmentAddIS - Add another data segment to a AOData database.
1478: Collective on AOData and IS
1480: Input Parameters:
1481: + aodata - the database
1482: . name - the name of the key
1483: . segment - name of segment
1484: . bs - the fundamental blocksize of the data
1485: . is - the keys provided by this processor
1486: . data - the actual data
1487: - dtype - the data type, one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.
1489: Level: advanced
1491: .keywords: database additions
1493: .seealso: AODataSegmentAdd()
1494: @*/
1495: PetscErrorCode AODataSegmentAddIS(AOData aodata,const char name[],const char segment[],PetscInt bs,IS is,void *data,PetscDataType dtype)
1496: {
1498: PetscInt n,*keys;
1504: ISGetLocalSize(is,&n);
1505: ISGetIndices(is,&keys);
1506: (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1507: ISRestoreIndices(is,&keys);
1508: return(0);
1509: }