Actual source code: aodata.c

  1: /*$Id: aodata.c,v 1.55 2001/05/17 15:14:27 bsmith Exp $*/
  2: /*  
  3:    Defines the abstract operations on AOData
  4: */
 5:  #include src/dm/ao/aoimpl.h


  8: /*@C
  9:     AODataGetInfo - Gets the number of keys and their names in a database.

 11:     Not collective

 13:     Input Parameter:
 14: .   ao - the AOData database

 16:     Output Parameters:
 17: +   nkeys - the number of keys
 18: -   keys - the names of the keys (or PETSC_NULL)

 20:    Level: advanced

 22: .keywords: application ordering

 24: .seealso:  AODataSegmentGetInfo()
 25: @*/
 26: int AODataGetInfo(AOData ao,int *nkeys,char ***keys)
 27: {
 28:   int       n,i,ierr;
 29:   AODataKey *key = ao->keys;


 34:   *nkeys = n = ao->nkeys;
 35:   if (keys) {
 36:     PetscMalloc((n+1)*sizeof(char *),&keys);
 37:     for (i=0; i<n; i++) {
 38:       if (!key) SETERRQ(PETSC_ERR_COR,"Less keys in database then indicated");
 39:       (*keys)[i] = key->name;
 40:       key        = key->next;
 41:     }
 42:   }
 43:   return(0);
 44: }

 46: /*
 47:    AODataKeyFind_Private - Given a keyname  finds the key. Generates a flag if not found.

 49:    Not collective

 51:    Input Parameters:
 52: .  keyname - string name of key

 54:    Output Parameter:
 55: +  flag - PETSC_TRUE if found, PETSC_FALSE if not found
 56: -  key - the associated key

 58:    Level: advanced

 60: */
 61: int AODataKeyFind_Private(AOData aodata,char *keyname,PetscTruth *flag,AODataKey **key)
 62: {
 63:   PetscTruth  match;
 64:   int         ierr;
 65:   AODataAlias *t = aodata->aliases;
 66:   char        *name = keyname;
 67:   AODataKey   *nkey;

 70:   *key   = PETSC_NULL;
 71:   *flag  = PETSC_FALSE;
 72:   while (name) {
 73:     nkey  = aodata->keys;
 74:     while (nkey) {
 75:       PetscStrcmp(nkey->name,name,&match);
 76:       if (match) {
 77:         /* found the key */
 78:         *key   = nkey;
 79:         *flag  = PETSC_TRUE;
 80:         return(0);
 81:       }
 82:       *key = nkey;
 83:       nkey = nkey->next;
 84:     }
 85:     name = 0;
 86:     while (t) {
 87:       PetscStrcmp(keyname,t->alias,&match);
 88:       if (match) {
 89:         name = t->name;
 90:         t    = t->next;
 91:         break;
 92:       }
 93:       t = t->next;
 94:     }
 95:   }
 96:   return(0);
 97: }

 99: /*@C
100:    AODataKeyExists - Determines if a key exists in the database.

102:    Not collective

104:    Input Parameters:
105: .  keyname - string name of key

107:    Output Parameter:
108: .  flag - PETSC_TRUE if found, otherwise PETSC_FALSE

110:    Level: advanced

112: @*/
113: int AODataKeyExists(AOData aodata,char *keyname,PetscTruth *flag)
114: {
115:   int        ierr;
116:   PetscTruth iflag;
117:   AODataKey  *ikey;

121:   AODataKeyFind_Private(aodata,keyname,&iflag,&ikey);
122:   if (iflag) *flag = PETSC_TRUE;
123:   else       *flag = PETSC_FALSE;
124:   return(0);
125: }


128: /*
129:    AODataSegmentFind_Private - Given a key and segment finds the int key, segment
130:    coordinates. Generates a flag if not found.

132:    Not collective

134:    Input Parameters:
135: +  keyname - string name of key
136: -  segname - string name of segment

138:    Output Parameter:
139: +  flag - PETSC_TRUE if found, PETSC_FALSE if not
140: .  key - integer of keyname
141: -  segment - integer of segment

143:    If it doesn't find it, it returns the last seg in the key (if the key exists)

145:    Level: advanced

147: */
148: int AODataSegmentFind_Private(AOData aodata,char *keyname,char *segname,PetscTruth *flag,AODataKey **key,AODataSegment **seg)
149: {
150:   int           ierr;
151:   PetscTruth    keyflag,match;
152:   AODataAlias   *t = aodata->aliases;
153:   char          *name;
154:   AODataSegment *nseg;

157:   *seg  = PETSC_NULL;
158:   *flag = PETSC_FALSE;
159:   ierr  = AODataKeyFind_Private(aodata,keyname,&keyflag,key);
160:   if (keyflag) { /* found key now look for segment */
161:     name = segname;
162:     while (name) {
163:       nseg = (*key)->segments;
164:       while (nseg) {
165:         PetscStrcmp(nseg->name,name,&match);
166:         if (match) {
167:           /* found the segment */
168:           *seg   = nseg;
169:           *flag  = PETSC_TRUE;
170:           return(0);
171:         }
172:         *seg = nseg;
173:         nseg = nseg->next;
174:       }
175:       name = 0;
176:       while (t) {
177:         PetscStrcmp(segname,t->alias,&match);
178:         if (match) {
179:           name = t->name;
180:           t    = t->next;
181:           break;
182:         }
183:         t = t->next;
184:       }
185:     }
186:   }
187:   return(0);
188: }

190: /*@C
191:    AODataSegmentExists - Determines if a key  and segment exists in the database.

193:    Not collective

195:    Input Parameters:
196: +  keyname - string name of key
197: -  segname - string name of segment

199:    Output Parameter:
200: .  flag - PETSC_TRUE if found, else PETSC_FALSE

202:    Level: advanced

204: @*/
205: int AODataSegmentExists(AOData aodata,char *keyname,char *segname,PetscTruth *flag)
206: {
207:   int           ierr;
208:   PetscTruth    iflag;
209:   AODataKey     *ikey;
210:   AODataSegment *iseg;

214:   AODataSegmentFind_Private(aodata,keyname,segname,&iflag,&ikey,&iseg);
215:   if (iflag) *flag = PETSC_TRUE;
216:   else       *flag = PETSC_FALSE;
217:   return(0);
218: }

220: /* ------------------------------------------------------------------------------------*/

222: /*@C
223:    AODataKeyGetActive - Get a sublist of key indices that have a logical flag on.

225:    Collective on AOData

227:    Input Parameters:
228: +  aodata - the database
229: .  name - the name of the key
230: .  segment - the name of the segment
231: .  n - the number of key indices provided by this processor
232: .  keys - the keys provided by this processor
233: -  wl - which logical key in the block (for block size 1 this is always 0)

235:    Output Parameters:
236: .  IS - the list of key indices

238:    Level: advanced

240: .keywords: database transactions

242: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
243:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
244:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
245: @*/
246: int AODataKeyGetActive(AOData aodata,char *name,char *segment,int n,int *keys,int wl,IS *is)
247: {

252:   (*aodata->ops->keygetactive)(aodata,name,segment,n,keys,wl,is);
253:   return(0);
254: }

256: /*@C
257:    AODataKeyGetActiveIS - Get a sublist of key indices that have a logical flag on.

259:    Collective on AOData

261:    Input Parameters:
262: +  aodata - the database
263: .  name - the name of the key
264: .  segment - the name of the segment
265: .  in - the key indices we are checking
266: -  wl - which logical key in the block (for block size 1 this is always 0)

268:    Output Parameters:
269: .  IS - the list of key indices

271:    Level: advanced

273: .keywords: database transactions

275: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
276:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
277:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
278: @*/
279: int AODataKeyGetActiveIS(AOData aodata,char *name,char *segname,IS in,int wl,IS *is)
280: {
281:   int ierr,n,*keys;

284:   ISGetLocalSize(in,&n);
285:   ISGetIndices(in,&keys);
286:   AODataKeyGetActive(aodata,name,segname,n,keys,wl,is);
287:   ISRestoreIndices(in,&keys);
288:   return(0);
289: }

291: /*@C
292:    AODataKeyGetActiveLocal - Get a sublist of key indices that have a logical flag on.

294:    Collective on AOData

296:    Input Parameters:
297: +  aodata - the database
298: .  name - the name of the key
299: .  segment - the name of the segment
300: .  n - the number of key indices provided by this processor
301: .  keys - the keys provided by this processor
302: -  wl - which logical key in the block (for block size 1 this is always 0)

304:    Output Parameters:
305: .  IS - the list of key indices

307:    Level: advanced

309: .keywords: database transactions

311: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
312:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
313:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
314: @*/
315: int AODataKeyGetActiveLocal(AOData aodata,char *name,char *segment,int n,int *keys,int wl,IS *is)
316: {

321:   (*aodata->ops->keygetactivelocal)(aodata,name,segment,n,keys,wl,is);
322:   return(0);
323: }

325: /*@C
326:    AODataKeyGetActiveLocalIS - Get a sublist of key indices that have a logical flag on.

328:    Collective on AOData

330:    Input Parameters:
331: +  aodata - the database
332: .  name - the name of the key
333: .  segment - the name of the segment
334: .  in - the key indices we are checking
335: -  wl - which logical key in the block (for block size 1 this is always 0)

337:    Output Parameters:
338: .  IS - the list of key indices

340:    Level: advanced

342: .keywords: database transactions

344: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
345:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
346:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
347: @*/
348: int AODataKeyGetActiveLocalIS(AOData aodata,char *name,char *segname,IS in,int wl,IS *is)
349: {
350:   int ierr,n,*keys;

353:   ISGetLocalSize(in,&n);
354:   ISGetIndices(in,&keys);
355:   AODataKeyGetActiveLocal(aodata,name,segname,n,keys,wl,is);
356:   ISRestoreIndices(in,&keys);
357:   return(0);
358: }

360: /* ------------------------------------------------------------------------------------*/

362: /*@C
363:    AODataSegmentGet - Get data from a particular segment of a database.

365:    Collective on AOData

367:    Input Parameters:
368: +  aodata - the database
369: .  name - the name of the key
370: .  segment - the name of the segment
371: .  n - the number of data items needed by this processor
372: -  keys - the keys provided by this processor

374:    Output Parameters:
375: .  data - the actual data

377:    Level: advanced

379: .keywords: database transactions

381: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
382:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
383:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
384: @*/
385: int AODataSegmentGet(AOData aodata,char *name,char *segment,int n,int *keys,void **data)
386: {

391:   (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
392:   return(0);
393: }

395: /*@C
396:    AODataSegmentRestore - Restores data from a particular segment of a database.

398:    Collective on AOData

400:    Input Parameters:
401: +  aodata - the database
402: .  name - the name of the key
403: .  segment - the name of the segment
404: .  n - the number of data items needed by this processor
405: -  keys - the keys provided by this processor

407:    Output Parameters:
408: .  data - the actual data

410:    Level: advanced

412: .keywords: database transactions

414: .seealso: AODataSegmentRestoreIS()
415: @*/
416: int AODataSegmentRestore(AOData aodata,char *name,char *segment,int n,int *keys,void **data)
417: {

422:   (*aodata->ops->segmentrestore)(aodata,name,segment,n,keys,data);
423:   return(0);
424: }

426: /*@C
427:    AODataSegmentGetIS - Get data from a particular segment of a database.

429:    Collective on AOData and IS

431:    Input Parameters:
432: +  aodata - the database
433: .  name - the name of the key
434: .  segment - the name of the segment
435: -  is - the keys for data requested on this processor

437:    Output Parameters:
438: .  data - the actual data

440:    Level: advanced

442: .keywords: database transactions

444: @*/
445: int AODataSegmentGetIS(AOData aodata,char *name,char *segment,IS is,void **data)
446: {
447:   int ierr,n,*keys;


453:   ISGetLocalSize(is,&n);
454:   ISGetIndices(is,&keys);
455:   (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
456:   ISRestoreIndices(is,&keys);
457:   return(0);
458: }

460: /*@C
461:    AODataSegmentRestoreIS - Restores data from a particular segment of a database.

463:    Collective on AOData and IS

465:    Input Parameters:
466: +  aodata - the database
467: .  name - the name of the data key
468: .  segment - the name of the segment
469: -  is - the keys provided by this processor

471:    Output Parameters:
472: .  data - the actual data

474:    Level: advanced

476: .keywords: database transactions

478: .seealso: AODataSegmentRestore()
479: @*/
480: int AODataSegmentRestoreIS(AOData aodata,char *name,char *segment,IS is,void **data)
481: {


487:   (*aodata->ops->segmentrestore)(aodata,name,segment,0,0,data);

489:   return(0);
490: }

492: /* ------------------------------------------------------------------------------------*/
493: /*@C
494:    AODataSegmentGetLocal - Get data from a particular segment of a database. Returns the 
495:    values in the local numbering; valid only for integer segments.

497:    Collective on AOData

499:    Input Parameters:
500: +  aodata - the database
501: .  name - the name of the key
502: .  segment - the name of the segment
503: .  n - the number of data items needed by this processor
504: -  keys - the keys provided by this processor in local numbering

506:    Output Parameters:
507: .  data - the actual data

509:    Level: advanced

511: .keywords: database transactions

513: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
514:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
515:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
516: @*/
517: int AODataSegmentGetLocal(AOData aodata,char *name,char *segment,int n,int *keys,void **data)
518: {

523:   (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
524:   return(0);
525: }

527: /*@C
528:    AODataSegmentRestoreLocal - Restores data from a particular segment of a database.

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

539:    Output Parameters:
540: .  data - the actual data

542:    Level: advanced

544: .keywords: database transactions

546: @*/
547: int AODataSegmentRestoreLocal(AOData aodata,char *name,char *segment,int n,int *keys,void **data)
548: {

553:   (*aodata->ops->segmentrestorelocal)(aodata,name,segment,n,keys,data);
554:   return(0);
555: }

557: /*@C
558:    AODataSegmentGetLocalIS - Get data from a particular segment of a database. Returns the 
559:    values in the local numbering; valid only for integer segments.

561:    Collective on AOData and IS

563:    Input Parameters:
564: +  aodata - the database
565: .  name - the name of the key
566: .  segment - the name of the segment
567: -  is - the keys for data requested on this processor

569:    Output Parameters:
570: .  data - the actual data

572:    Level: advanced

574: .keywords: database transactions

576: .seealso: AODataSegmentRestoreLocalIS()
577: @*/
578: int AODataSegmentGetLocalIS(AOData aodata,char *name,char *segment,IS is,void **data)
579: {
580:   int ierr,n,*keys;


586:   ISGetLocalSize(is,&n);
587:   ISGetIndices(is,&keys);
588:   (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
589:   ISRestoreIndices(is,&keys);
590:   return(0);
591: }

593: /*@C
594:    AODataSegmentRestoreLocalIS - Restores data from a particular segment of a database.

596:    Collective on AOData and IS

598:    Input Parameters:
599: +  aodata - the database
600: .  name - the name of the data key
601: .  segment - the name of the segment
602: -  is - the keys provided by this processor

604:    Output Parameters:
605: .  data - the actual data

607:    Level: advanced

609: .keywords: database transactions

611: .seealso: AODataSegmentGetLocalIS()
612: @*/
613: int AODataSegmentRestoreLocalIS(AOData aodata,char *name,char *segment,IS is,void **data)
614: {

619:   (*aodata->ops->segmentrestorelocal)(aodata,name,segment,0,0,data);
620:   return(0);
621: }

623: /* ------------------------------------------------------------------------------------*/

625: /*@C
626:    AODataKeyGetNeighbors - Given a list of keys generates a new list containing
627:    those keys plus neighbors found in a neighbors list.

629:    Collective on AOData

631:    Input Parameters:
632: +  aodata - the database
633: .  name - the name of the key
634: .  n - the number of data items needed by this processor
635: -  keys - the keys provided by this processor

637:    Output Parameters:
638: .  is - the indices retrieved

640:    Level: advanced

642: .keywords: database transactions

644: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
645:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
646:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(), 
647:           AODataKeyGetNeighborsIS()
648: @*/
649: int AODataKeyGetNeighbors(AOData aodata,char *name,int n,int *keys,IS *is)
650: {
652:   IS  reduced,input;

656: 
657:   /* get the list of neighbors */
658:   AODataSegmentGetReduced(aodata,name,name,n,keys,&reduced);

660:   ISCreateGeneral(aodata->comm,n,keys,&input);
661:   ISSum(input,reduced,is);
662:   ISDestroy(input);
663:   ISDestroy(reduced);

665:   return(0);
666: }

668: /*@C
669:    AODataKeyGetNeighborsIS - Given a list of keys generates a new list containing
670:    those keys plus neighbors found in a neighbors list.

672:    Collective on AOData and IS

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: advanced

685: .keywords: database transactions

687: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
688:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
689:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(), 
690:           AODataKeyGetNeighbors()
691: @*/
692: int AODataKeyGetNeighborsIS(AOData aodata,char *name,IS keys,IS *is)
693: {
695:   IS  reduced;

699: 
700:   /* get the list of neighbors */
701:   AODataSegmentGetReducedIS(aodata,name,name,keys,&reduced);
702:   ISSum(keys,reduced,is);
703:   ISDestroy(reduced);

705:   return(0);
706: }

708: /*@C
709:    AODataSegmentGetReduced - Gets the unique list of segment values, by removing 
710:    duplicates.

712:    Collective on AOData and IS

714:    Input Parameters:
715: +  aodata - the database
716: .  name - the name of the key
717: .  segment - the name of the segment
718: .  n - the number of data items needed by this processor
719: -  keys - the keys provided by this processor

721:    Output Parameters:
722: .  is - the indices retrieved

724:    Level: advanced

726:    Example:
727: .vb
728:                       keys    ->      0  1  2  3  4   5  6  7
729:       if the segment contains ->      1  2  1  3  1   4  2  0
730:    and you request keys 0 1 2 5 7 it will return 1 2 4 0
731: .ve

733: .keywords: database transactions

735: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
736:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
737:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
738: @*/
739: int AODataSegmentGetReduced(AOData aodata,char *name,char *segment,int n,int *keys,IS *is)
740: {

745:   (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,is);
746:   return(0);
747: }

749: /*@C
750:    AODataSegmentGetExtrema - Gets the largest and smallest values for each entry in the block

752:    Collective on AOData

754:    Input Parameters:
755: +  aodata - the database
756: .  name - the name of the key
757: -  segment - the name of the segment

759:    Output Parameters:
760: +  vmax - the maximum values (user must provide enough space)
761: -  vmin - the minimum values (user must provide enough space)

763:    Level: advanced

765: .keywords: database transactions

767: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
768:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
769:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
770: @*/
771: int AODataSegmentGetExtrema(AOData aodata,char *name,char *segment,void *vmax,void *vmin)
772: {

777:   (*aodata->ops->segmentgetextrema)(aodata,name,segment,vmax,vmin);
778:   return(0);
779: }

781: /*@C
782:    AODataSegmentGetReducedIS -  Gets the unique list of segment values, by removing 
783:    duplicates.

785:    Collective on AOData and IS

787:    Input Parameters:
788: +  aodata - the database
789: .  name - the name of the key
790: .  segment - the name of the segment
791: -  is - the keys for data requested on this processor

793:    Output Parameters:
794: .  isout - the indices retreived

796:    Level: advanced

798:    Example:
799: .vb
800:                       keys    ->      0  1  2  3  4   5  6  7
801:       if the segment contains ->      1  2  1  3  1   4  2  0

803:   and you request keys 0 1 2 5 7, AODataSegmentGetReducedIS() will return 1 2 4 0
804: .ve

806: .keywords: database transactions

808: .seealso:
809: @*/
810: int AODataSegmentGetReducedIS(AOData aodata,char *name,char *segment,IS is,IS *isout)
811: {
812:   int ierr,n,*keys;


818:   ISGetLocalSize(is,&n);
819:   ISGetIndices(is,&keys);
820:   (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,isout);
821:   ISRestoreIndices(is,&keys);
822:   return(0);
823: }

825: /* ------------------------------------------------------------------------------------*/

827: /*@C
828:    AODataKeySetLocalToGlobalMapping - Add a local to global mapping for a key in the 
829:      in the database

831:    Not collective

833:    Input Parameters:
834: +  aodata - the database
835: .   name - the name of the key
836: -  map - local to global mapping

838:    Level: advanced

840: .keywords: database additions

842: .seealso: AODataKeyGetLocalToGlobalMapping()
843: @*/
844: int AODataKeySetLocalToGlobalMapping(AOData aodata,char *name,ISLocalToGlobalMapping map)
845: {
846:   int        ierr;
847:   PetscTruth flag;
848:   AODataKey  *ikey;


853:   AODataKeyFind_Private(aodata,name,&flag,&ikey);
854:   if (!flag)  SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Key does not exist");

856:   if (ikey->ltog) {
857:     SETERRQ1(1,"Database key %s already has local to global mapping",name);
858:   }

860:   ikey->ltog = map;
861:   PetscObjectReference((PetscObject)map);

863:   return(0);

865: }

867: /*@C
868:    AODataKeyGetLocalToGlobalMapping - gets a local to global mapping from a database

870:    Not collective

872:    Input Parameters:
873: +  aodata - the database
874: -  name - the name of the key

876:    Output Parameters:
877: .  map - local to global mapping

879:    Level: advanced

881: .keywords: database additions

883: .seealso: AODataKeySetLocalToGlobalMapping()
884: @*/
885: int AODataKeyGetLocalToGlobalMapping(AOData aodata,char *name,ISLocalToGlobalMapping *map)
886: {
887:   int        ierr;
888:   PetscTruth flag;
889:   AODataKey  *ikey;


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

897:   *map = ikey->ltog;
898:   return(0);
899: }


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

905:    Not collective

907:    Input Parameters:
908: +  aodata - the database
909: -  name - the name of the key

911:    Output Parameters:
912: +  rstart - first key owned locally (or PETSC_NULL if not needed) 
913: -  rend - last key owned locally + 1 (or PETSC_NULL if not needed)

915:    Level: advanced

917: .keywords: database accessing

919: .seealso: AODataKeyGetInfo()
920: @*/
921: int AODataKeyGetOwnershipRange(AOData aodata,char *name,int *rstart,int *rend)
922: {
923:   int        ierr;
924:   PetscTruth flag;
925:   AODataKey  *key;


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

933:   if (rstart) *rstart = key->rstart;
934:   if (rend)   *rend   = key->rend;

936:   return(0);
937: }

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

942:    Not collective

944:    Input Parameters:
945: +  aodata - the database
946: -  name - the name of the key

948:    Output Parameters:
949: +  nglobal - global number of keys
950: .  nlocal - local number of keys
951: .  nsegments - number of segments associated with key
952: -  segnames - names of the segments or PETSC_NULL

954:    Level: advanced

956: .keywords: database accessing

958: .seealso: AODataKeyGetOwnershipRange()
959: @*/
960: int AODataKeyGetInfo(AOData aodata,char *name,int *nglobal,int *nlocal,int *nsegments,char ***segnames)
961: {
962:   int           ierr,i,n=0;
963:   AODataKey     *key;
964:   AODataSegment *seg;
965:   PetscTruth    flag;


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

973:   if (nglobal)   *nglobal   = key->N;
974:   if (nlocal)    *nlocal    = key->nlocal;
975:   if (nsegments) *nsegments = n = key->nsegments;
976:   if (nsegments && segnames) {
977:     PetscMalloc((n+1)*sizeof(char *),&segnames);
978:     seg  = key->segments;
979:     for (i=0; i<n; i++) {
980:       if (!seg) SETERRQ(PETSC_ERR_COR,"Less segments in database then indicated");
981:       (*segnames)[i] = seg->name;
982:       seg            = seg->next;
983:     }
984:   }

986:   return(0);
987: }

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

992:    Not collective

994:    Input Parameters:
995: +  aodata - the database
996: .  keyname - the name of the key
997: -  segname - the name of the segment

999:    Output Parameters:
1000: +  bs - the blocksize
1001: -  dtype - the datatype

1003:    Level: advanced

1005: .keywords: database accessing

1007: .seealso:  AODataGetInfo()
1008: @*/
1009: int AODataSegmentGetInfo(AOData aodata,char *keyname,char *segname,int *bs,PetscDataType *dtype)
1010: {
1011:   int           ierr;
1012:   PetscTruth    flag;
1013:   AODataKey     *key;
1014:   AODataSegment *seg;


1019:   AODataSegmentFind_Private(aodata,keyname,segname,&flag,&key,&seg);
1020:   if (flag == PETSC_FALSE) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Segment never created: %s",segname);
1021:   if (bs)        *bs        = seg->bs;
1022:   if (dtype)     *dtype     = seg->datatype;

1024:   return(0);
1025: }

1027: /*@C
1028:    AODataView - Displays an application ordering.

1030:    Collective on AOData and PetscViewer

1032:    Input Parameters:
1033: +  aodata - the database
1034: -  viewer - viewer used for display

1036:    Level: intermediate

1038:    The available visualization contexts include
1039: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
1040: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
1041:          output where only the first processor opens
1042:          the file.  All other processors send their 
1043:          data to the first processor to print. 

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


1049: .keywords: database viewing

1051: .seealso: PetscViewerASCIIOpen()
1052: @*/
1053: int AODataView(AOData aodata,PetscViewer viewer)
1054: {

1059:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(aodata->comm);
1061:   (*aodata->ops->view)(aodata,viewer);
1062:   return(0);
1063: }

1065: static int AODataAliasDestroy_Private(AODataAlias *aliases)
1066: {
1067:   AODataAlias *t = aliases;
1068:   int         ierr;

1071:   if (t) {
1072:     t = aliases->next;
1073:     PetscFree(aliases->name);
1074:     PetscFree(aliases->alias);
1075:     PetscFree(aliases);
1076:     while (t) {
1077:       aliases = t;
1078:       t       = t->next;
1079:       PetscFree(aliases->name);
1080:       PetscFree(aliases->alias);
1081:       ierr    = PetscFree(aliases);
1082:     }
1083:   }
1084:   return(0);
1085: }

1087: int AODataAliasAdd(AOData aodata,char *alias,char *name)
1088: {
1089:   AODataAlias *t = aodata->aliases;
1090:   int         ierr;

1093:   if (t) {
1094:     while (t->next) t = t->next;
1095:     PetscNew(AODataAlias,&t->next);
1096:     t    = t->next;
1097:   } else {
1098:     ierr            = PetscNew(AODataAlias,&t);
1099:     aodata->aliases = t;
1100:   }
1101:   ierr    = PetscStrallocpy(alias,&t->alias);
1102:   ierr    = PetscStrallocpy(name,&t->name);
1103:   t->next = 0;
1104:   return(0);
1105: }

1107: /*@C
1108:    AODataDestroy - Destroys an application ordering set.

1110:    Collective on AOData

1112:    Input Parameters:
1113: .  aodata - the database

1115:    Level: intermediate

1117: .keywords: destroy, database

1119: .seealso: AODataCreateBasic()
1120: @*/
1121: int AODataDestroy(AOData aodata)
1122: {


1127:   if (!aodata) return(0);
1129:   if (--aodata->refct > 0) return(0);
1130: 
1131:   AODataAliasDestroy_Private(aodata->aliases);
1132:   (*aodata->ops->destroy)(aodata);

1134:   return(0);
1135: }

1137: /*@C
1138:    AODataKeyRemap - Remaps a key and all references to a key to a new numbering 
1139:    scheme where each processor indicates its new nodes by listing them in the
1140:    previous numbering scheme.

1142:    Collective on AOData and AO

1144:    Input Parameters:
1145: +  aodata - the database
1146: .  key  - the key to remap
1147: -  ao - the old to new ordering

1149:    Level: advanced

1151: .keywords: database remapping

1153: .seealso: AODataKeyGetAdjacency()
1154: @*/
1155: int AODataKeyRemap(AOData aodata,char *key,AO ao)
1156: {

1162:   (*aodata->ops->keyremap)(aodata,key,ao);
1163:   return(0);
1164: }

1166: /*@C
1167:    AODataKeyGetAdjacency - Gets the adjacency graph for a key.

1169:    Collective on AOData

1171:    Input Parameters:
1172: +  aodata - the database
1173: -  key  - the key

1175:    Output Parameter:
1176: .  adj - the adjacency graph

1178:    Level: advanced

1180: .keywords: database, adjacency graph

1182: .seealso: AODataKeyRemap()
1183: @*/
1184: int AODataKeyGetAdjacency(AOData aodata,char *key,Mat *adj)
1185: {

1190:   (*aodata->ops->keygetadjacency)(aodata,key,adj);
1191:   return(0);
1192: }

1194: /*@C
1195:     AODataSegmentPartition - Partitions a segment type across processors 
1196:     relative to a key that is partitioned. This will try to keep as
1197:     many elements of the segment on the same processor as corresponding
1198:     neighboring key elements are.

1200:     Collective on AOData

1202:     Input Parameters:
1203: +   aodata - the database
1204: -   key - the key to be partitioned and renumbered

1206:    Level: advanced

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

1210: @*/
1211: int AODataSegmentPartition(AOData aodata,char *key,char *seg)
1212: {
1213:   int             ierr;

1217:   (*aodata->ops->segmentpartition)(aodata,key,seg);
1218:   return(0);
1219: }

1221: int AODataPublish_Petsc(PetscObject obj)
1222: {
1223: #if defined(PETSC_HAVE_AMS)
1224:   AOData        ao = (AOData) obj;
1225:   AODataKey     *key = 0;
1226:   AODataSegment *segment = 0;
1227:   int           ierr,keys,segments;
1228:   char          tmp[1024];
1229: #endif


1233: #if defined(PETSC_HAVE_AMS)
1234:   /* if it is already published then return */
1235:   if (ao->amem >=0) return(0);

1237:   PetscObjectPublishBaseBegin(obj);
1238:   AMS_Memory_add_field((AMS_Memory)ao->amem,"Number_of_Keys",&ao->nkeys,1,AMS_INT,AMS_READ,
1239:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1240:   /* Loop over keys publishing info on each */
1241:   for (keys=0; keys<ao->nkeys; keys++) {
1242:     if (!keys) key = ao->keys;
1243:     else       key = key->next;

1245:     PetscStrcpy(tmp,key->name);
1246:     PetscStrcat(tmp,"_N");
1247:     AMS_Memory_add_field((AMS_Memory)ao->amem,tmp,&key->N,1,AMS_INT,AMS_READ,
1248:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1249: 
1250:     PetscStrcpy(tmp,key->name);
1251:     PetscStrcat(tmp,"_nsegments");
1252:     AMS_Memory_add_field((AMS_Memory)ao->amem,tmp,&key->nsegments,1,AMS_INT,AMS_READ,
1253:                                 AMS_COMMON,AMS_REDUCT_UNDEF);

1255:     for (segments=0; segments<key->nsegments; segments++) {
1256:       if (!segments) segment = key->segments;
1257:       else           segment = segment->next;
1258: 
1259:       PetscStrcpy(tmp,key->name);
1260:       PetscStrcat(tmp,"_");
1261:       PetscStrcat(tmp,segment->name);
1262:       PetscStrcat(tmp,"_bs");
1263:       AMS_Memory_add_field((AMS_Memory)ao->amem,tmp,&segment->bs,1,AMS_INT,AMS_READ,
1264:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1265:     }
1266:   }

1268:   PetscObjectPublishBaseEnd(obj);
1269: #endif

1271:   return(0);
1272: }

1274: /*@C
1275:    AODataKeyRemove - Remove a data key from a AOData database.

1277:    Collective on AOData

1279:    Input Parameters:
1280: +  aodata - the database
1281: -  name - the name of the key

1283:    Level: advanced

1285: .keywords: database removal

1287: .seealso:
1288: @*/
1289: int AODataKeyRemove(AOData aodata,char *name)
1290: {

1295:   (*aodata->ops->keyremove)(aodata,name);
1296:   return(0);
1297: }

1299: /*@C
1300:    AODataSegmentRemove - Remove a data segment from a AOData database.

1302:    Collective on AOData

1304:    Input Parameters:
1305: +  aodata - the database
1306: .  name - the name of the key
1307: -  segname - name of the segment

1309:    Level: advanced

1311: .keywords: database removal

1313: .seealso:
1314: @*/
1315: int AODataSegmentRemove(AOData aodata,char *name,char *segname)
1316: {

1321:   (*aodata->ops->segmentremove)(aodata,name,segname);
1322:   return(0);
1323: }

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

1328:    Collective on AOData

1330:    Input Parameters:
1331: +  aodata - the database
1332: .  name - the name of the key
1333: .  N - the number of indices in the key
1334: -  nlocal - number of indices to be associated with this processor

1336:    Level: advanced

1338: .keywords: database additions

1340: .seealso:
1341: @*/
1342: int AODataKeyAdd(AOData aodata,char *name,int nlocal,int N)
1343: {
1344:   int        ierr,size,rank,i,len;
1345:   AODataKey  *key,*oldkey;
1346:   MPI_Comm   comm = aodata->comm;
1347:   PetscTruth flag;


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

1355:   PetscNew(AODataKey,&key);
1356:   if (oldkey) { oldkey->next = key;}
1357:   else        { aodata->keys = key;}
1358:   ierr           = PetscStrlen(name,&len);
1359:   ierr           = PetscMalloc((len+1)*sizeof(char),&key->name);
1360:   ierr           = PetscStrcpy(key->name,name);
1361:   key->N         = N;
1362:   key->nsegments = 0;
1363:   key->segments  = 0;
1364:   key->ltog      = 0;
1365:   key->next      = 0;

1367:   MPI_Comm_rank(comm,&rank);
1368:   MPI_Comm_size(comm,&size);

1370:   /*  Set nlocal and ownership ranges */
1371:   ierr            = PetscSplitOwnership(comm,&nlocal,&N);
1372:   ierr            = PetscMalloc((size+1)*sizeof(int),&key->rowners);
1373:   ierr            = MPI_Allgather(&nlocal,1,MPI_INT,key->rowners+1,1,MPI_INT,comm);
1374:   key->rowners[0] = 0;
1375:   for (i=2; i<=size; i++) {
1376:     key->rowners[i] += key->rowners[i-1];
1377:   }
1378:   key->rstart        = key->rowners[rank];
1379:   key->rend          = key->rowners[rank+1];

1381:   key->nlocal        = nlocal;

1383:   aodata->nkeys++;

1385: #if defined(PETSC_HAVE_AMS)
1386:   if (aodata->amem >=0) {
1387:     char namesize[1024];
1388:     PetscStrcpy(namesize,name);
1389:     PetscStrcat(namesize,"_N");
1390:     AMS_Memory_add_field((AMS_Memory)aodata->amem,namesize,&key->N,1,AMS_INT,AMS_READ,
1391:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1392:   }
1393: #endif

1395:   return(0);
1396: }

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

1401:    Collective on AOData

1403:    Input Parameters:
1404: +  aodata  - the database
1405: .  name    - the name of the key
1406: .  segment - the name of the data segment
1407: .  bs      - the fundamental blocksize of the data
1408: .  n       - the number of data items contributed by this processor
1409: .  keys    - the keys provided by this processor
1410: .  data    - the actual data
1411: -  dtype   - the data type (one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.)

1413:    Level: advanced

1415: .keywords: database additions

1417: .seealso: AODataSegmentAddIS()
1418: @*/
1419: int AODataSegmentAdd(AOData aodata,char *name,char *segment,int bs,int n,int *keys,void *data,PetscDataType dtype)
1420: {
1421:   int      ierr;


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

1428:   /*
1429:   PetscOptionsHasName(PETSC_NULL,"-ao_data_view",&flg1);
1430:   if (flg1) {
1431:     AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1432:   }
1433:   PetscOptionsHasName(PETSC_NULL,"-ao_data_view_info",&flg1);
1434:   if (flg1) {
1435:     PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(comm),PETSC_VIEWER_ASCII_INFO);
1436:     AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1437:     PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(comm));
1438:   }
1439:   */
1440:   return(0);
1441: }

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

1446:    Collective on AOData and IS

1448:    Input Parameters:
1449: +  aodata - the database
1450: .  name - the name of the key
1451: .  segment - name of segment
1452: .  bs - the fundamental blocksize of the data
1453: .  is - the keys provided by this processor
1454: .  data - the actual data
1455: -  dtype - the data type, one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.

1457:    Level: advanced

1459: .keywords: database additions

1461: .seealso: AODataSegmentAdd()
1462: @*/
1463: int AODataSegmentAddIS(AOData aodata,char *name,char *segment,int bs,IS is,void *data,PetscDataType dtype)
1464: {
1465:   int n,*keys,ierr;


1471:   ISGetLocalSize(is,&n);
1472:   ISGetIndices(is,&keys);
1473:   (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1474:   ISRestoreIndices(is,&keys);
1475:   return(0);
1476: }