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:   }
913: 
914:   PetscObjectReference((PetscObject)map);
915:   if (ikey->ltog) { ISLocalToGlobalMappingDestroy(ikey->ltog); }
916:   ikey->ltog = map;

918:   return(0);

920: }

924: /*@C
925:    AODataKeyGetLocalToGlobalMapping - gets a local to global mapping from a database

927:    Not collective

929:    Input Parameters:
930: +  aodata - the database
931: -  name - the name of the key

933:    Output Parameters:
934: .  map - local to global mapping

936:    Level: deprecated

938: .keywords: database additions

940: .seealso: AODataKeySetLocalToGlobalMapping()
941: @*/
942: PetscErrorCode  AODataKeyGetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping *map)
943: {
945:   PetscTruth     flag;
946:   AODataKey      *ikey;


951:   AODataKeyFind_Private(aodata,name,&flag,&ikey);
952:   if (!flag)  SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key does not exist: %s",name);

954:   *map = ikey->ltog;
955:   return(0);
956: }


961: /*@C
962:    AODataKeyGetOwnershipRange - Gets the ownership range to this key type.

964:    Not collective

966:    Input Parameters:
967: +  aodata - the database
968: -  name - the name of the key

970:    Output Parameters:
971: +  rstart - first key owned locally (or PETSC_NULL if not needed) 
972: -  rend - last key owned locally + 1 (or PETSC_NULL if not needed)

974:    Level: deprecated

976: .keywords: database accessing

978: .seealso: AODataKeyGetInfo()
979: @*/
980: PetscErrorCode  AODataKeyGetOwnershipRange(AOData aodata,const char name[],PetscInt *rstart,PetscInt *rend)
981: {
983:   PetscTruth     flag;
984:   AODataKey      *key;


989:   AODataKeyFind_Private(aodata,name,&flag,&key);
990:   if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key never created: %s",name);

992:   if (rstart) *rstart = key->rstart;
993:   if (rend)   *rend   = key->rend;

995:   return(0);
996: }

1000: /*@C
1001:    AODataKeyGetInfo - Gets the global size, local size and number of segments in a key.

1003:    Not collective

1005:    Input Parameters:
1006: +  aodata - the database
1007: -  name - the name of the key

1009:    Output Parameters:
1010: +  nglobal - global number of keys
1011: .  nlocal - local number of keys
1012: .  nsegments - number of segments associated with key
1013: -  segnames - names of the segments or PETSC_NULL

1015:    Level: deprecated

1017: .keywords: database accessing

1019: .seealso: AODataKeyGetOwnershipRange()
1020: @*/
1021: PetscErrorCode  AODataKeyGetInfo(AOData aodata,const char name[],PetscInt *nglobal,PetscInt *nlocal,PetscInt *nsegments,char ***segnames)
1022: {
1024:   PetscInt       i,n=0;
1025:   AODataKey     *key;
1026:   AODataSegment *seg;
1027:   PetscTruth    flag;


1032:   AODataKeyFind_Private(aodata,name,&flag,&key);
1033:   if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key never created: %s",name);

1035:   if (nglobal)   *nglobal   = key->N;
1036:   if (nlocal)    *nlocal    = key->nlocal;
1037:   if (nsegments) *nsegments = n = key->nsegments;
1038:   if (nsegments && segnames) {
1039:     PetscMalloc((n+1)*sizeof(char *),&segnames);
1040:     seg  = key->segments;
1041:     for (i=0; i<n; i++) {
1042:       if (!seg) SETERRQ(PETSC_ERR_COR,"Less segments in database then indicated");
1043:       (*segnames)[i] = seg->name;
1044:       seg            = seg->next;
1045:     }
1046:   }

1048:   return(0);
1049: }

1053: /*@C
1054:    AODataSegmentGetInfo - Gets the blocksize and type of a data segment

1056:    Not collective

1058:    Input Parameters:
1059: +  aodata - the database
1060: .  keyname - the name of the key
1061: -  segname - the name of the segment

1063:    Output Parameters:
1064: +  bs - the blocksize
1065: -  dtype - the datatype

1067:    Level: deprecated

1069: .keywords: database accessing

1071: .seealso:  AODataGetInfo()
1072: @*/
1073: PetscErrorCode  AODataSegmentGetInfo(AOData aodata,const char keyname[],const char segname[],PetscInt *bs,PetscDataType *dtype)
1074: {
1076:   PetscTruth    flag;
1077:   AODataKey     *key;
1078:   AODataSegment *seg;


1083:   AODataSegmentFind_Private(aodata,keyname,segname,&flag,&key,&seg);
1084:   if (!flag) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Segment never created: %s",segname);
1085:   if (bs)        *bs        = seg->bs;
1086:   if (dtype)     *dtype     = seg->datatype;

1088:   return(0);
1089: }

1093: /*@C
1094:    AODataView - Displays an application ordering.

1096:    Collective on AOData and PetscViewer

1098:    Input Parameters:
1099: +  aodata - the database
1100: -  viewer - viewer used for display

1102:    Level: intermediate

1104:    The available visualization contexts include
1105: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
1106: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
1107:          output where only the first processor opens
1108:          the file.  All other processors send their 
1109:          data to the first processor to print. 

1111:    The user can open an alternative visualization context with
1112:    PetscViewerASCIIOpen() - output to a specified file.


1115: .keywords: database viewing

1117: .seealso: PetscViewerASCIIOpen()
1118: @*/
1119: PetscErrorCode  AODataView(AOData aodata,PetscViewer viewer)
1120: {

1125:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(aodata->comm);
1127:   (*aodata->ops->view)(aodata,viewer);
1128:   return(0);
1129: }

1133: static PetscErrorCode AODataAliasDestroy_Private(AODataAlias *aliases)
1134: {
1135:   AODataAlias    *t = aliases;

1139:   if (t) {
1140:     t = aliases->next;
1141:     PetscFree(aliases->name);
1142:     PetscFree(aliases->alias);
1143:     PetscFree(aliases);
1144:     while (t) {
1145:       aliases = t;
1146:       t       = t->next;
1147:       PetscFree(aliases->name);
1148:       PetscFree(aliases->alias);
1149:       PetscFree(aliases);
1150:     }
1151:   }
1152:   return(0);
1153: }

1157: /*@C
1158:   AODataAliasAdd - Man page needed.

1160:   Level: deprecated

1162: @*/
1163: PetscErrorCode  AODataAliasAdd(AOData aodata,const char alias[],const char name[])
1164: {
1165:   AODataAlias *t = aodata->aliases;

1169:   if (t) {
1170:     while (t->next) t = t->next;
1171:     PetscNew(AODataAlias,&t->next);
1172:     t    = t->next;
1173:   } else {
1174:     PetscNew(AODataAlias,&t);
1175:     aodata->aliases = t;
1176:   }
1177:   PetscStrallocpy(alias,&t->alias);
1178:   PetscStrallocpy(name,&t->name);
1179:   t->next = 0;
1180:   return(0);
1181: }

1185: /*@C
1186:    AODataDestroy - Destroys an application ordering set.

1188:    Collective on AOData

1190:    Input Parameters:
1191: .  aodata - the database

1193:    Level: deprecated

1195: .keywords: destroy, database

1197: .seealso: AODataCreateBasic()
1198: @*/
1199: PetscErrorCode  AODataDestroy(AOData aodata)
1200: {


1205:   if (!aodata) return(0);
1207:   if (--aodata->refct > 0) return(0);
1208: 
1209:   AODataAliasDestroy_Private(aodata->aliases);
1210:   (*aodata->ops->destroy)(aodata);

1212:   return(0);
1213: }

1217: /*@C
1218:    AODataKeyRemap - Remaps a key and all references to a key to a new numbering 
1219:    scheme where each processor indicates its new nodes by listing them in the
1220:    previous numbering scheme.

1222:    Collective on AOData and AO

1224:    Input Parameters:
1225: +  aodata - the database
1226: .  key  - the key to remap
1227: -  ao - the old to new ordering

1229:    Level: deprecated

1231: .keywords: database remapping

1233: .seealso: AODataKeyGetAdjacency()
1234: @*/
1235: PetscErrorCode  AODataKeyRemap(AOData aodata,const char key[],AO ao)
1236: {

1242:   (*aodata->ops->keyremap)(aodata,key,ao);
1243:   return(0);
1244: }

1248: /*@C
1249:    AODataKeyGetAdjacency - Gets the adjacency graph for a key.

1251:    Collective on AOData

1253:    Input Parameters:
1254: +  aodata - the database
1255: -  key  - the key

1257:    Output Parameter:
1258: .  adj - the adjacency graph

1260:    Level: deprecated

1262: .keywords: database, adjacency graph

1264: .seealso: AODataKeyRemap()
1265: @*/
1266: PetscErrorCode  AODataKeyGetAdjacency(AOData aodata,const char key[],Mat *adj)
1267: {

1272:   (*aodata->ops->keygetadjacency)(aodata,key,adj);
1273:   return(0);
1274: }

1278: /*@C
1279:     AODataSegmentPartition - Partitions a segment type across processors 
1280:     relative to a key that is partitioned. This will try to keep as
1281:     many elements of the segment on the same processor as corresponding
1282:     neighboring key elements are.

1284:     Collective on AOData

1286:     Input Parameters:
1287: +   aodata - the database
1288: -   key - the key to be partitioned and renumbered

1290:    Level: deprecated

1292: .seealso: AODataKeyPartition(), AODataPartitionAndSetupLocal()

1294: @*/
1295: PetscErrorCode  AODataSegmentPartition(AOData aodata,const char key[],const char seg[])
1296: {

1301:   (*aodata->ops->segmentpartition)(aodata,key,seg);
1302:   return(0);
1303: }

1307: PetscErrorCode AODataPublish_Petsc(PetscObject obj)
1308: {
1310:   return(0);
1311: }

1315: /*@C
1316:    AODataKeyRemove - Remove a data key from a AOData database.

1318:    Collective on AOData

1320:    Input Parameters:
1321: +  aodata - the database
1322: -  name - the name of the key

1324:    Level: deprecated

1326: .keywords: database removal

1328: .seealso:
1329: @*/
1330: PetscErrorCode  AODataKeyRemove(AOData aodata,const char name[])
1331: {

1336:   (*aodata->ops->keyremove)(aodata,name);
1337:   return(0);
1338: }

1342: /*@C
1343:    AODataSegmentRemove - Remove a data segment from a AOData database.

1345:    Collective on AOData

1347:    Input Parameters:
1348: +  aodata - the database
1349: .  name - the name of the key
1350: -  segname - name of the segment

1352:    Level: deprecated

1354: .keywords: database removal

1356: .seealso:
1357: @*/
1358: PetscErrorCode  AODataSegmentRemove(AOData aodata,const char name[],const char segname[])
1359: {

1364:   (*aodata->ops->segmentremove)(aodata,name,segname);
1365:   return(0);
1366: }

1370: /*@C
1371:    AODataKeyAdd - Add another data key to a AOData database.

1373:    Collective on AOData

1375:    Input Parameters:
1376: +  aodata - the database
1377: .  name - the name of the key
1378: .  nlocal - number of indices to be associated with this processor
1379: -  N - the number of indices in the key

1381:    Level: deprecated

1383: .keywords: database additions

1385: .seealso:
1386: @*/
1387: PetscErrorCode  AODataKeyAdd(AOData aodata,const char name[],PetscInt nlocal,PetscInt N)
1388: {
1390:   PetscMPIInt    size,rank;
1391:   PetscInt       i;
1392:   size_t         len;
1393:   AODataKey      *key,*oldkey;
1394:   MPI_Comm       comm = aodata->comm;
1395:   PetscTruth     flag;


1400:   AODataKeyFind_Private(aodata,name,&flag,&oldkey);
1401:   if (flag)  SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Key already exists with given name: %s",name);

1403:   PetscNew(AODataKey,&key);
1404:   if (oldkey) { oldkey->next = key;}
1405:   else        { aodata->keys = key;}
1406:   PetscStrlen(name,&len);
1407:   PetscMalloc((len+1)*sizeof(char),&key->name);
1408:   PetscStrcpy(key->name,name);
1409:   key->N         = N;
1410:   key->nsegments = 0;
1411:   key->segments  = 0;
1412:   key->ltog      = 0;
1413:   key->next      = 0;

1415:   MPI_Comm_rank(comm,&rank);
1416:   MPI_Comm_size(comm,&size);

1418:   /*  Set nlocal and ownership ranges */
1419:   PetscSplitOwnership(comm,&nlocal,&N);
1420:   PetscMalloc((size+1)*sizeof(PetscInt),&key->rowners);
1421:   MPI_Allgather(&nlocal,1,MPIU_INT,key->rowners+1,1,MPIU_INT,comm);
1422:   key->rowners[0] = 0;
1423:   for (i=2; i<=size; i++) {
1424:     key->rowners[i] += key->rowners[i-1];
1425:   }
1426:   key->rstart        = key->rowners[rank];
1427:   key->rend          = key->rowners[rank+1];

1429:   key->nlocal        = nlocal;

1431:   aodata->nkeys++;
1432:   return(0);
1433: }

1437: /*@C
1438:    AODataSegmentAdd - Adds another data segment to a AOData database.

1440:    Collective on AOData

1442:    Input Parameters:
1443: +  aodata  - the database
1444: .  name    - the name of the key
1445: .  segment - the name of the data segment
1446: .  bs      - the fundamental blocksize of the data
1447: .  n       - the number of data items contributed by this processor
1448: .  keys    - the keys provided by this processor
1449: .  data    - the actual data
1450: -  dtype   - the data type (one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.)

1452:    Level: deprecated

1454: .keywords: database additions

1456: .seealso: AODataSegmentAddIS()
1457: @*/
1458: PetscErrorCode  AODataSegmentAdd(AOData aodata,const char name[],const char segment[],PetscInt bs,PetscInt n,PetscInt *keys,void *data,PetscDataType dtype)
1459: {


1465:   (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);

1467:   /*
1468:   PetscOptionsHasName(PETSC_NULL,"-ao_data_view",&flg1);
1469:   if (flg1) {
1470:     AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1471:   }
1472:   PetscOptionsHasName(PETSC_NULL,"-ao_data_view_info",&flg1);
1473:   if (flg1) {
1474:     PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(comm),PETSC_VIEWER_ASCII_INFO);
1475:     AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1476:     PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(comm));
1477:   }
1478:   */
1479:   return(0);
1480: }

1484: /*@C
1485:    AODataSegmentAddIS - Add another data segment to a AOData database.

1487:    Collective on AOData and IS

1489:    Input Parameters:
1490: +  aodata - the database
1491: .  name - the name of the key
1492: .  segment - name of segment
1493: .  bs - the fundamental blocksize of the data
1494: .  is - the keys provided by this processor
1495: .  data - the actual data
1496: -  dtype - the data type, one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.

1498:    Level: deprecated

1500: .keywords: database additions

1502: .seealso: AODataSegmentAdd()
1503: @*/
1504: PetscErrorCode  AODataSegmentAddIS(AOData aodata,const char name[],const char segment[],PetscInt bs,IS is,void *data,PetscDataType dtype)
1505: {
1507:   PetscInt       n,*keys;


1513:   ISGetLocalSize(is,&n);
1514:   ISGetIndices(is,&keys);
1515:   (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1516:   ISRestoreIndices(is,&keys);
1517:   return(0);
1518: }