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


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

 13:     Not collective

 15:     Input Parameter:
 16: .   ao - the AOData database

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

 22:    Level: advanced

 24: .keywords: application ordering

 26: .seealso:  AODataSegmentGetInfo()
 27: @*/
 28: int AODataGetInfo(AOData ao,int *nkeys,char ***keys)
 29: {
 30:   int       n,i,ierr;
 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: int AODataKeyFind_Private(AOData aodata,const char keyname[],PetscTruth *flag,AODataKey **key)
 66: {
 67:   PetscTruth  match;
 68:   int         ierr;
 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: int AODataKeyExists(AOData aodata,const char keyname[],PetscTruth *flag)
120: {
121:   int        ierr;
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: int AODataSegmentFind_Private(AOData aodata,const char keyname[],const char segname[],PetscTruth *flag,AODataKey **key,AODataSegment **seg)
157: {
158:   int           ierr;
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: int AODataSegmentExists(AOData aodata,const char keyname[],const char segname[],PetscTruth *flag)
216: {
217:   int           ierr;
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: int AODataKeyGetActive(AOData aodata,const char name[],const char segment[],int n,int *keys,int 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: int AODataKeyGetActiveIS(AOData aodata,const char name[],const char segname[],IS in,int wl,IS *is)
294: {
295:   int ierr,n,*keys;

298:   ISGetLocalSize(in,&n);
299:   ISGetIndices(in,&keys);
300:   AODataKeyGetActive(aodata,name,segname,n,keys,wl,is);
301:   ISRestoreIndices(in,&keys);
302:   return(0);
303: }

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

310:    Collective on AOData

312:    Input Parameters:
313: +  aodata - the database
314: .  name - the name of the key
315: .  segment - the name of the segment
316: .  n - the number of key indices provided by this processor
317: .  keys - the keys provided by this processor
318: -  wl - which logical key in the block (for block size 1 this is always 0)

320:    Output Parameters:
321: .  IS - the list of key indices

323:    Level: advanced

325: .keywords: database transactions

327: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
328:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
329:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
330: @*/
331: int AODataKeyGetActiveLocal(AOData aodata,const char name[],const char segment[],int n,int *keys,int wl,IS *is)
332: {

337:   (*aodata->ops->keygetactivelocal)(aodata,name,segment,n,keys,wl,is);
338:   return(0);
339: }

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

346:    Collective on AOData

348:    Input Parameters:
349: +  aodata - the database
350: .  name - the name of the key
351: .  segment - the name of the segment
352: .  in - the key indices we are checking
353: -  wl - which logical key in the block (for block size 1 this is always 0)

355:    Output Parameters:
356: .  IS - the list of key indices

358:    Level: advanced

360: .keywords: database transactions

362: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
363:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
364:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
365: @*/
366: int AODataKeyGetActiveLocalIS(AOData aodata,const char name[],const char segname[],IS in,int wl,IS *is)
367: {
368:   int ierr,n,*keys;

371:   ISGetLocalSize(in,&n);
372:   ISGetIndices(in,&keys);
373:   AODataKeyGetActiveLocal(aodata,name,segname,n,keys,wl,is);
374:   ISRestoreIndices(in,&keys);
375:   return(0);
376: }

378: /* ------------------------------------------------------------------------------------*/

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

385:    Collective on AOData

387:    Input Parameters:
388: +  aodata - the database
389: .  name - the name of the key
390: .  segment - the name of the segment
391: .  n - the number of data items needed by this processor
392: -  keys - the keys provided by this processor

394:    Output Parameters:
395: .  data - the actual data

397:    Level: advanced

399: .keywords: database transactions

401: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
402:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
403:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
404: @*/
405: int AODataSegmentGet(AOData aodata,const char name[],const char segment[],int n,int *keys,void **data)
406: {

411:   (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
412:   return(0);
413: }

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

420:    Collective on AOData

422:    Input Parameters:
423: +  aodata - the database
424: .  name - the name of the key
425: .  segment - the name of the segment
426: .  n - the number of data items needed by this processor
427: -  keys - the keys provided by this processor

429:    Output Parameters:
430: .  data - the actual data

432:    Level: advanced

434: .keywords: database transactions

436: .seealso: AODataSegmentRestoreIS()
437: @*/
438: int AODataSegmentRestore(AOData aodata,const char name[],const char segment[],int n,int *keys,void **data)
439: {

444:   (*aodata->ops->segmentrestore)(aodata,name,segment,n,keys,data);
445:   return(0);
446: }

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

453:    Collective on AOData and IS

455:    Input Parameters:
456: +  aodata - the database
457: .  name - the name of the key
458: .  segment - the name of the segment
459: -  is - the keys for data requested on this processor

461:    Output Parameters:
462: .  data - the actual data

464:    Level: advanced

466: .keywords: database transactions

468: @*/
469: int AODataSegmentGetIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
470: {
471:   int ierr,n,*keys;


477:   ISGetLocalSize(is,&n);
478:   ISGetIndices(is,&keys);
479:   (*aodata->ops->segmentget)(aodata,name,segment,n,keys,data);
480:   ISRestoreIndices(is,&keys);
481:   return(0);
482: }

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

489:    Collective on AOData and IS

491:    Input Parameters:
492: +  aodata - the database
493: .  name - the name of the data key
494: .  segment - the name of the segment
495: -  is - the keys provided by this processor

497:    Output Parameters:
498: .  data - the actual data

500:    Level: advanced

502: .keywords: database transactions

504: .seealso: AODataSegmentRestore()
505: @*/
506: int AODataSegmentRestoreIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
507: {


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

515:   return(0);
516: }

518: /* ------------------------------------------------------------------------------------*/
521: /*@C
522:    AODataSegmentGetLocal - Get data from a particular segment of a database. Returns the 
523:    values in the local numbering; valid only for integer segments.

525:    Collective on AOData

527:    Input Parameters:
528: +  aodata - the database
529: .  name - the name of the key
530: .  segment - the name of the segment
531: .  n - the number of data items needed by this processor
532: -  keys - the keys provided by this processor in local numbering

534:    Output Parameters:
535: .  data - the actual data

537:    Level: advanced

539: .keywords: database transactions

541: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
542:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
543:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
544: @*/
545: int AODataSegmentGetLocal(AOData aodata,const char name[],const char segment[],int n,int *keys,void **data)
546: {

551:   (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
552:   return(0);
553: }

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

560:    Collective on AOData

562:    Input Parameters:
563: +  aodata - the database
564: .  name - the name of the key
565: .  segment - the name of the segment
566: .  n - the number of data items needed by this processor
567: -  keys - the keys provided by this processor

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

572:    Level: advanced

574: .keywords: database transactions

576: @*/
577: int AODataSegmentRestoreLocal(AOData aodata,const char name[],const char segment[],int n,int *keys,void **data)
578: {

583:   (*aodata->ops->segmentrestorelocal)(aodata,name,segment,n,keys,data);
584:   return(0);
585: }

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

593:    Collective on AOData and IS

595:    Input Parameters:
596: +  aodata - the database
597: .  name - the name of the key
598: .  segment - the name of the segment
599: -  is - the keys for data requested on this processor

601:    Output Parameters:
602: .  data - the actual data

604:    Level: advanced

606: .keywords: database transactions

608: .seealso: AODataSegmentRestoreLocalIS()
609: @*/
610: int AODataSegmentGetLocalIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
611: {
612:   int ierr,n,*keys;


618:   ISGetLocalSize(is,&n);
619:   ISGetIndices(is,&keys);
620:   (*aodata->ops->segmentgetlocal)(aodata,name,segment,n,keys,data);
621:   ISRestoreIndices(is,&keys);
622:   return(0);
623: }

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

630:    Collective on AOData and IS

632:    Input Parameters:
633: +  aodata - the database
634: .  name - the name of the data key
635: .  segment - the name of the segment
636: -  is - the keys provided by this processor

638:    Output Parameters:
639: .  data - the actual data

641:    Level: advanced

643: .keywords: database transactions

645: .seealso: AODataSegmentGetLocalIS()
646: @*/
647: int AODataSegmentRestoreLocalIS(AOData aodata,const char name[],const char segment[],IS is,void **data)
648: {

654:   (*aodata->ops->segmentrestorelocal)(aodata,name,segment,0,0,data);
655:   return(0);
656: }

658: /* ------------------------------------------------------------------------------------*/

662: /*@C
663:    AODataKeyGetNeighbors - Given a list of keys generates a new list containing
664:    those keys plus neighbors found in a neighbors list.

666:    Collective on AOData

668:    Input Parameters:
669: +  aodata - the database
670: .  name - the name of the key
671: .  n - the number of data items needed by this processor
672: -  keys - the keys provided by this processor

674:    Output Parameters:
675: .  is - the indices retrieved

677:    Level: advanced

679: .keywords: database transactions

681: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
682:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
683:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(), 
684:           AODataKeyGetNeighborsIS()
685: @*/
686: int AODataKeyGetNeighbors(AOData aodata,const char name[],int n,int *keys,IS *is)
687: {
689:   IS  reduced,input;

693: 
694:   /* get the list of neighbors */
695:   AODataSegmentGetReduced(aodata,name,name,n,keys,&reduced);

697:   ISCreateGeneral(aodata->comm,n,keys,&input);
698:   ISSum(input,reduced,is);
699:   ISDestroy(input);
700:   ISDestroy(reduced);

702:   return(0);
703: }

707: /*@C
708:    AODataKeyGetNeighborsIS - Given a list of keys generates a new list containing
709:    those keys plus neighbors found in a neighbors list.

711:    Collective on AOData and IS

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

719:    Output Parameters:
720: .  is - the indices retrieved

722:    Level: advanced

724: .keywords: database transactions

726: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
727:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
728:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd(), 
729:           AODataKeyGetNeighbors()
730: @*/
731: int AODataKeyGetNeighborsIS(AOData aodata,const char name[],IS keys,IS *is)
732: {
734:   IS  reduced;

738: 
739:   /* get the list of neighbors */
740:   AODataSegmentGetReducedIS(aodata,name,name,keys,&reduced);
741:   ISSum(keys,reduced,is);
742:   ISDestroy(reduced);

744:   return(0);
745: }

749: /*@C
750:    AODataSegmentGetReduced - Gets the unique list of segment values, by removing 
751:    duplicates.

753:    Collective on AOData and IS

755:    Input Parameters:
756: +  aodata - the database
757: .  name - the name of the key
758: .  segment - the name of the segment
759: .  n - the number of data items needed by this processor
760: -  keys - the keys provided by this processor

762:    Output Parameters:
763: .  is - the indices retrieved

765:    Level: advanced

767:    Example:
768: .vb
769:                       keys    ->      0  1  2  3  4   5  6  7
770:       if the segment contains ->      1  2  1  3  1   4  2  0
771:    and you request keys 0 1 2 5 7 it will return 1 2 4 0
772: .ve

774: .keywords: database transactions

776: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
777:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
778:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
779: @*/
780: int AODataSegmentGetReduced(AOData aodata,const char name[],const char segment[],int n,int *keys,IS *is)
781: {

786:   (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,is);
787:   return(0);
788: }

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

795:    Collective on AOData

797:    Input Parameters:
798: +  aodata - the database
799: .  name - the name of the key
800: -  segment - the name of the segment

802:    Output Parameters:
803: +  vmax - the maximum values (user must provide enough space)
804: -  vmin - the minimum values (user must provide enough space)

806:    Level: advanced

808: .keywords: database transactions

810: .seealso: AODataCreateBasic(), AODataDestroy(), AODataKeyAdd(), AODataSegmentRestore(),
811:           AODataSegmentGetIS(), AODataSegmentRestoreIS(), AODataSegmentAdd(), 
812:           AODataKeyGetInfo(), AODataSegmentGetInfo(), AODataSegmentAdd()
813: @*/
814: int AODataSegmentGetExtrema(AOData aodata,const char name[],const char segment[],void *vmax,void *vmin)
815: {

820:   (*aodata->ops->segmentgetextrema)(aodata,name,segment,vmax,vmin);
821:   return(0);
822: }

826: /*@C
827:    AODataSegmentGetReducedIS -  Gets the unique list of segment values, by removing 
828:    duplicates.

830:    Collective on AOData and IS

832:    Input Parameters:
833: +  aodata - the database
834: .  name - the name of the key
835: .  segment - the name of the segment
836: -  is - the keys for data requested on this processor

838:    Output Parameters:
839: .  isout - the indices retreived

841:    Level: advanced

843:    Example:
844: .vb
845:                       keys    ->      0  1  2  3  4   5  6  7
846:       if the segment contains ->      1  2  1  3  1   4  2  0

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

851: .keywords: database transactions

853: .seealso:
854: @*/
855: int AODataSegmentGetReducedIS(AOData aodata,const char name[],const char segment[],IS is,IS *isout)
856: {
857:   int ierr,n,*keys;


863:   ISGetLocalSize(is,&n);
864:   ISGetIndices(is,&keys);
865:   (*aodata->ops->segmentgetreduced)(aodata,name,segment,n,keys,isout);
866:   ISRestoreIndices(is,&keys);
867:   return(0);
868: }

870: /* ------------------------------------------------------------------------------------*/

874: /*@C
875:    AODataKeySetLocalToGlobalMapping - Add a local to global mapping for a key in the 
876:      in the database

878:    Not collective

880:    Input Parameters:
881: +  aodata - the database
882: .   name - the name of the key
883: -  map - local to global mapping

885:    Level: advanced

887: .keywords: database additions

889: .seealso: AODataKeyGetLocalToGlobalMapping()
890: @*/
891: int AODataKeySetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping map)
892: {
893:   int        ierr;
894:   PetscTruth flag;
895:   AODataKey  *ikey;


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

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

907:   ikey->ltog = map;
908:   PetscObjectReference((PetscObject)map);

910:   return(0);

912: }

916: /*@C
917:    AODataKeyGetLocalToGlobalMapping - gets a local to global mapping from a database

919:    Not collective

921:    Input Parameters:
922: +  aodata - the database
923: -  name - the name of the key

925:    Output Parameters:
926: .  map - local to global mapping

928:    Level: advanced

930: .keywords: database additions

932: .seealso: AODataKeySetLocalToGlobalMapping()
933: @*/
934: int AODataKeyGetLocalToGlobalMapping(AOData aodata,const char name[],ISLocalToGlobalMapping *map)
935: {
936:   int        ierr;
937:   PetscTruth flag;
938:   AODataKey  *ikey;


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

946:   *map = ikey->ltog;
947:   return(0);
948: }


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

956:    Not collective

958:    Input Parameters:
959: +  aodata - the database
960: -  name - the name of the key

962:    Output Parameters:
963: +  rstart - first key owned locally (or PETSC_NULL if not needed) 
964: -  rend - last key owned locally + 1 (or PETSC_NULL if not needed)

966:    Level: advanced

968: .keywords: database accessing

970: .seealso: AODataKeyGetInfo()
971: @*/
972: int AODataKeyGetOwnershipRange(AOData aodata,const char name[],int *rstart,int *rend)
973: {
974:   int        ierr;
975:   PetscTruth flag;
976:   AODataKey  *key;


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

984:   if (rstart) *rstart = key->rstart;
985:   if (rend)   *rend   = key->rend;

987:   return(0);
988: }

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

995:    Not collective

997:    Input Parameters:
998: +  aodata - the database
999: -  name - the name of the key

1001:    Output Parameters:
1002: +  nglobal - global number of keys
1003: .  nlocal - local number of keys
1004: .  nsegments - number of segments associated with key
1005: -  segnames - names of the segments or PETSC_NULL

1007:    Level: advanced

1009: .keywords: database accessing

1011: .seealso: AODataKeyGetOwnershipRange()
1012: @*/
1013: int AODataKeyGetInfo(AOData aodata,const char name[],int *nglobal,int *nlocal,int *nsegments,char ***segnames)
1014: {
1015:   int           ierr,i,n=0;
1016:   AODataKey     *key;
1017:   AODataSegment *seg;
1018:   PetscTruth    flag;


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

1026:   if (nglobal)   *nglobal   = key->N;
1027:   if (nlocal)    *nlocal    = key->nlocal;
1028:   if (nsegments) *nsegments = n = key->nsegments;
1029:   if (nsegments && segnames) {
1030:     PetscMalloc((n+1)*sizeof(char *),&segnames);
1031:     seg  = key->segments;
1032:     for (i=0; i<n; i++) {
1033:       if (!seg) SETERRQ(PETSC_ERR_COR,"Less segments in database then indicated");
1034:       (*segnames)[i] = seg->name;
1035:       seg            = seg->next;
1036:     }
1037:   }

1039:   return(0);
1040: }

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

1047:    Not collective

1049:    Input Parameters:
1050: +  aodata - the database
1051: .  keyname - the name of the key
1052: -  segname - the name of the segment

1054:    Output Parameters:
1055: +  bs - the blocksize
1056: -  dtype - the datatype

1058:    Level: advanced

1060: .keywords: database accessing

1062: .seealso:  AODataGetInfo()
1063: @*/
1064: int AODataSegmentGetInfo(AOData aodata,const char keyname[],const char segname[],int *bs,PetscDataType *dtype)
1065: {
1066:   int           ierr;
1067:   PetscTruth    flag;
1068:   AODataKey     *key;
1069:   AODataSegment *seg;


1074:   AODataSegmentFind_Private(aodata,keyname,segname,&flag,&key,&seg);
1075:   if (flag == PETSC_FALSE) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Segment never created: %s",segname);
1076:   if (bs)        *bs        = seg->bs;
1077:   if (dtype)     *dtype     = seg->datatype;

1079:   return(0);
1080: }

1084: /*@C
1085:    AODataView - Displays an application ordering.

1087:    Collective on AOData and PetscViewer

1089:    Input Parameters:
1090: +  aodata - the database
1091: -  viewer - viewer used for display

1093:    Level: intermediate

1095:    The available visualization contexts include
1096: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
1097: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
1098:          output where only the first processor opens
1099:          the file.  All other processors send their 
1100:          data to the first processor to print. 

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


1106: .keywords: database viewing

1108: .seealso: PetscViewerASCIIOpen()
1109: @*/
1110: int AODataView(AOData aodata,PetscViewer viewer)
1111: {

1116:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(aodata->comm);
1118:   (*aodata->ops->view)(aodata,viewer);
1119:   return(0);
1120: }

1124: static int AODataAliasDestroy_Private(AODataAlias *aliases)
1125: {
1126:   AODataAlias *t = aliases;
1127:   int         ierr;

1130:   if (t) {
1131:     t = aliases->next;
1132:     PetscFree(aliases->name);
1133:     PetscFree(aliases->alias);
1134:     PetscFree(aliases);
1135:     while (t) {
1136:       aliases = t;
1137:       t       = t->next;
1138:       PetscFree(aliases->name);
1139:       PetscFree(aliases->alias);
1140:       PetscFree(aliases);
1141:     }
1142:   }
1143:   return(0);
1144: }

1148: int AODataAliasAdd(AOData aodata,const char alias[],const char name[])
1149: {
1150:   AODataAlias *t = aodata->aliases;
1151:   int         ierr;

1154:   if (t) {
1155:     while (t->next) t = t->next;
1156:     PetscNew(AODataAlias,&t->next);
1157:     t    = t->next;
1158:   } else {
1159:     PetscNew(AODataAlias,&t);
1160:     aodata->aliases = t;
1161:   }
1162:   PetscStrallocpy(alias,&t->alias);
1163:   PetscStrallocpy(name,&t->name);
1164:   t->next = 0;
1165:   return(0);
1166: }

1170: /*@C
1171:    AODataDestroy - Destroys an application ordering set.

1173:    Collective on AOData

1175:    Input Parameters:
1176: .  aodata - the database

1178:    Level: intermediate

1180: .keywords: destroy, database

1182: .seealso: AODataCreateBasic()
1183: @*/
1184: int AODataDestroy(AOData aodata)
1185: {


1190:   if (!aodata) return(0);
1192:   if (--aodata->refct > 0) return(0);
1193: 
1194:   AODataAliasDestroy_Private(aodata->aliases);
1195:   (*aodata->ops->destroy)(aodata);

1197:   return(0);
1198: }

1202: /*@C
1203:    AODataKeyRemap - Remaps a key and all references to a key to a new numbering 
1204:    scheme where each processor indicates its new nodes by listing them in the
1205:    previous numbering scheme.

1207:    Collective on AOData and AO

1209:    Input Parameters:
1210: +  aodata - the database
1211: .  key  - the key to remap
1212: -  ao - the old to new ordering

1214:    Level: advanced

1216: .keywords: database remapping

1218: .seealso: AODataKeyGetAdjacency()
1219: @*/
1220: int AODataKeyRemap(AOData aodata,const char key[],AO ao)
1221: {

1227:   (*aodata->ops->keyremap)(aodata,key,ao);
1228:   return(0);
1229: }

1233: /*@C
1234:    AODataKeyGetAdjacency - Gets the adjacency graph for a key.

1236:    Collective on AOData

1238:    Input Parameters:
1239: +  aodata - the database
1240: -  key  - the key

1242:    Output Parameter:
1243: .  adj - the adjacency graph

1245:    Level: advanced

1247: .keywords: database, adjacency graph

1249: .seealso: AODataKeyRemap()
1250: @*/
1251: int AODataKeyGetAdjacency(AOData aodata,const char key[],Mat *adj)
1252: {

1257:   (*aodata->ops->keygetadjacency)(aodata,key,adj);
1258:   return(0);
1259: }

1263: /*@C
1264:     AODataSegmentPartition - Partitions a segment type across processors 
1265:     relative to a key that is partitioned. This will try to keep as
1266:     many elements of the segment on the same processor as corresponding
1267:     neighboring key elements are.

1269:     Collective on AOData

1271:     Input Parameters:
1272: +   aodata - the database
1273: -   key - the key to be partitioned and renumbered

1275:    Level: advanced

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

1279: @*/
1280: int AODataSegmentPartition(AOData aodata,const char key[],const char seg[])
1281: {
1282:   int             ierr;

1286:   (*aodata->ops->segmentpartition)(aodata,key,seg);
1287:   return(0);
1288: }

1292: int AODataPublish_Petsc(PetscObject obj)
1293: {
1294: #if defined(PETSC_HAVE_AMS)
1295:   AOData        ao = (AOData) obj;
1296:   AODataKey     *key = 0;
1297:   AODataSegment *segment = 0;
1298:   int           ierr,keys,segments;
1299:   char          tmp[1024];
1300: #endif


1304: #if defined(PETSC_HAVE_AMS)
1305:   /* if it is already published then return */
1306:   if (ao->amem >=0) return(0);

1308:   PetscObjectPublishBaseBegin(obj);
1309:   AMS_Memory_add_field((AMS_Memory)ao->amem,"Number_of_Keys",&ao->nkeys,1,AMS_INT,AMS_READ,
1310:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1311:   /* Loop over keys publishing info on each */
1312:   for (keys=0; keys<ao->nkeys; keys++) {
1313:     if (!keys) key = ao->keys;
1314:     else       key = key->next;

1316:     PetscStrcpy(tmp,key->name);
1317:     PetscStrcat(tmp,"_N");
1318:     AMS_Memory_add_field((AMS_Memory)ao->amem,tmp,&key->N,1,AMS_INT,AMS_READ,
1319:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1320: 
1321:     PetscStrcpy(tmp,key->name);
1322:     PetscStrcat(tmp,"_nsegments");
1323:     AMS_Memory_add_field((AMS_Memory)ao->amem,tmp,&key->nsegments,1,AMS_INT,AMS_READ,
1324:                                 AMS_COMMON,AMS_REDUCT_UNDEF);

1326:     for (segments=0; segments<key->nsegments; segments++) {
1327:       if (!segments) segment = key->segments;
1328:       else           segment = segment->next;
1329: 
1330:       PetscStrcpy(tmp,key->name);
1331:       PetscStrcat(tmp,"_");
1332:       PetscStrcat(tmp,segment->name);
1333:       PetscStrcat(tmp,"_bs");
1334:       AMS_Memory_add_field((AMS_Memory)ao->amem,tmp,&segment->bs,1,AMS_INT,AMS_READ,
1335:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1336:     }
1337:   }

1339:   PetscObjectPublishBaseEnd(obj);
1340: #endif

1342:   return(0);
1343: }

1347: /*@C
1348:    AODataKeyRemove - Remove a data key from a AOData database.

1350:    Collective on AOData

1352:    Input Parameters:
1353: +  aodata - the database
1354: -  name - the name of the key

1356:    Level: advanced

1358: .keywords: database removal

1360: .seealso:
1361: @*/
1362: int AODataKeyRemove(AOData aodata,const char name[])
1363: {

1368:   (*aodata->ops->keyremove)(aodata,name);
1369:   return(0);
1370: }

1374: /*@C
1375:    AODataSegmentRemove - Remove a data segment from a AOData database.

1377:    Collective on AOData

1379:    Input Parameters:
1380: +  aodata - the database
1381: .  name - the name of the key
1382: -  segname - name of the segment

1384:    Level: advanced

1386: .keywords: database removal

1388: .seealso:
1389: @*/
1390: int AODataSegmentRemove(AOData aodata,const char name[],const char segname[])
1391: {

1396:   (*aodata->ops->segmentremove)(aodata,name,segname);
1397:   return(0);
1398: }

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

1405:    Collective on AOData

1407:    Input Parameters:
1408: +  aodata - the database
1409: .  name - the name of the key
1410: .  nlocal - number of indices to be associated with this processor
1411: -  N - the number of indices in the key

1413:    Level: advanced

1415: .keywords: database additions

1417: .seealso:
1418: @*/
1419: int AODataKeyAdd(AOData aodata,const char name[],int nlocal,int N)
1420: {
1421:   int        ierr,size,rank,i,len;
1422:   AODataKey  *key,*oldkey;
1423:   MPI_Comm   comm = aodata->comm;
1424:   PetscTruth flag;


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

1432:   PetscNew(AODataKey,&key);
1433:   if (oldkey) { oldkey->next = key;}
1434:   else        { aodata->keys = key;}
1435:   PetscStrlen(name,&len);
1436:   PetscMalloc((len+1)*sizeof(char),&key->name);
1437:   PetscStrcpy(key->name,name);
1438:   key->N         = N;
1439:   key->nsegments = 0;
1440:   key->segments  = 0;
1441:   key->ltog      = 0;
1442:   key->next      = 0;

1444:   MPI_Comm_rank(comm,&rank);
1445:   MPI_Comm_size(comm,&size);

1447:   /*  Set nlocal and ownership ranges */
1448:   PetscSplitOwnership(comm,&nlocal,&N);
1449:   PetscMalloc((size+1)*sizeof(int),&key->rowners);
1450:   MPI_Allgather(&nlocal,1,MPI_INT,key->rowners+1,1,MPI_INT,comm);
1451:   key->rowners[0] = 0;
1452:   for (i=2; i<=size; i++) {
1453:     key->rowners[i] += key->rowners[i-1];
1454:   }
1455:   key->rstart        = key->rowners[rank];
1456:   key->rend          = key->rowners[rank+1];

1458:   key->nlocal        = nlocal;

1460:   aodata->nkeys++;

1462: #if defined(PETSC_HAVE_AMS)
1463:   if (aodata->amem >=0) {
1464:     char namesize[1024];
1465:     PetscStrcpy(namesize,name);
1466:     PetscStrcat(namesize,"_N");
1467:     AMS_Memory_add_field((AMS_Memory)aodata->amem,namesize,&key->N,1,AMS_INT,AMS_READ,
1468:                                 AMS_COMMON,AMS_REDUCT_UNDEF);
1469:   }
1470: #endif

1472:   return(0);
1473: }

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

1480:    Collective on AOData

1482:    Input Parameters:
1483: +  aodata  - the database
1484: .  name    - the name of the key
1485: .  segment - the name of the data segment
1486: .  bs      - the fundamental blocksize of the data
1487: .  n       - the number of data items contributed by this processor
1488: .  keys    - the keys provided by this processor
1489: .  data    - the actual data
1490: -  dtype   - the data type (one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.)

1492:    Level: advanced

1494: .keywords: database additions

1496: .seealso: AODataSegmentAddIS()
1497: @*/
1498: int AODataSegmentAdd(AOData aodata,const char name[],const char segment[],int bs,int n,int *keys,void *data,PetscDataType dtype)
1499: {
1500:   int      ierr;


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

1507:   /*
1508:   PetscOptionsHasName(PETSC_NULL,"-ao_data_view",&flg1);
1509:   if (flg1) {
1510:     AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1511:   }
1512:   PetscOptionsHasName(PETSC_NULL,"-ao_data_view_info",&flg1);
1513:   if (flg1) {
1514:     PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(comm),PETSC_VIEWER_ASCII_INFO);
1515:     AODataView(aodata,PETSC_VIEWER_STDOUT_(comm));
1516:     PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(comm));
1517:   }
1518:   */
1519:   return(0);
1520: }

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

1527:    Collective on AOData and IS

1529:    Input Parameters:
1530: +  aodata - the database
1531: .  name - the name of the key
1532: .  segment - name of segment
1533: .  bs - the fundamental blocksize of the data
1534: .  is - the keys provided by this processor
1535: .  data - the actual data
1536: -  dtype - the data type, one of PETSC_INT, PETSC_DOUBLE, PETSC_SCALAR, etc.

1538:    Level: advanced

1540: .keywords: database additions

1542: .seealso: AODataSegmentAdd()
1543: @*/
1544: int AODataSegmentAddIS(AOData aodata,const char name[],const char segment[],int bs,IS is,void *data,PetscDataType dtype)
1545: {
1546:   int n,*keys,ierr;


1552:   ISGetLocalSize(is,&n);
1553:   ISGetIndices(is,&keys);
1554:   (*aodata->ops->segmentadd)(aodata,name,segment,bs,n,keys,data,dtype);
1555:   ISRestoreIndices(is,&keys);
1556:   return(0);
1557: }