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: #undef __FUNCT__  
 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: }

 48: #undef __FUNCT__  
 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,char *keyname,PetscTruth *flag,AODataKey **key)
 66: {
 67:   PetscTruth  match;
 68:   int         ierr;
 69:   AODataAlias *t = aodata->aliases;
 70:   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: }

103: #undef __FUNCT__  
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,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: }


134: #undef __FUNCT__  
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,char *keyname,char *segname,PetscTruth *flag,AODataKey **key,AODataSegment **seg)
157: {
158:   int           ierr;
159:   PetscTruth    keyflag,match;
160:   AODataAlias   *t = aodata->aliases;
161:   char          *name;
162:   AODataSegment *nseg;

165:   *seg  = PETSC_NULL;
166:   *flag = PETSC_FALSE;
167:   ierr  = 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: }

198: #undef __FUNCT__  
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,char *keyname,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: /* ------------------------------------------------------------------------------------*/

232: #undef __FUNCT__  
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,char *name,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: }

268: #undef __FUNCT__  
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,char *name,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: }

305: #undef __FUNCT__  
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,char *name,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: }

341: #undef __FUNCT__  
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,char *name,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: /* ------------------------------------------------------------------------------------*/

380: #undef __FUNCT__  
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,char *name,char *segment,int n,int *keys,void **data)
406: {

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

415: #undef __FUNCT__  
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,char *name,char *segment,int n,int *keys,void **data)
439: {

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

448: #undef __FUNCT__  
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,char *name,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: }

484: #undef __FUNCT__  
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,char *name,char *segment,IS is,void **data)
507: {


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

515:   return(0);
516: }

518: /* ------------------------------------------------------------------------------------*/
519: #undef __FUNCT__  
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,char *name,char *segment,int n,int *keys,void **data)
546: {

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

555: #undef __FUNCT__  
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,char *name,char *segment,int n,int *keys,void **data)
578: {

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

587: #undef __FUNCT__  
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,char *name,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: }

625: #undef __FUNCT__  
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,char *name,char *segment,IS is,void **data)
648: {

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

657: /* ------------------------------------------------------------------------------------*/

659: #undef __FUNCT__  
661: /*@C
662:    AODataKeyGetNeighbors - Given a list of keys generates a new list containing
663:    those keys plus neighbors found in a neighbors list.

665:    Collective on AOData

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

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

676:    Level: advanced

678: .keywords: database transactions

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

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

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

701:   return(0);
702: }

704: #undef __FUNCT__  
706: /*@C
707:    AODataKeyGetNeighborsIS - Given a list of keys generates a new list containing
708:    those keys plus neighbors found in a neighbors list.

710:    Collective on AOData and IS

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

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

721:    Level: advanced

723: .keywords: database transactions

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

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

743:   return(0);
744: }

746: #undef __FUNCT__  
748: /*@C
749:    AODataSegmentGetReduced - Gets the unique list of segment values, by removing 
750:    duplicates.

752:    Collective on AOData and IS

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

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

764:    Level: advanced

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

773: .keywords: database transactions

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

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

789: #undef __FUNCT__  
791: /*@C
792:    AODataSegmentGetExtrema - Gets the largest and smallest values for each entry in the block

794:    Collective on AOData

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

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

805:    Level: advanced

807: .keywords: database transactions

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

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

823: #undef __FUNCT__  
825: /*@C
826:    AODataSegmentGetReducedIS -  Gets the unique list of segment values, by removing 
827:    duplicates.

829:    Collective on AOData and IS

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

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

840:    Level: advanced

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

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

850: .keywords: database transactions

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


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

869: /* ------------------------------------------------------------------------------------*/

871: #undef __FUNCT__  
873: /*@C
874:    AODataKeySetLocalToGlobalMapping - Add a local to global mapping for a key in the 
875:      in the database

877:    Not collective

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

884:    Level: advanced

886: .keywords: database additions

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


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

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

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

909:   return(0);

911: }

913: #undef __FUNCT__  
915: /*@C
916:    AODataKeyGetLocalToGlobalMapping - gets a local to global mapping from a database

918:    Not collective

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

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

927:    Level: advanced

929: .keywords: database additions

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


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

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


950: #undef __FUNCT__  
952: /*@C
953:    AODataKeyGetOwnershipRange - Gets the ownership range to this key type.

955:    Not collective

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

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

965:    Level: advanced

967: .keywords: database accessing

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


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

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

986:   return(0);
987: }

989: #undef __FUNCT__  
991: /*@C
992:    AODataKeyGetInfo - Gets the global size, local size and number of segments in a key.

994:    Not collective

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

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

1006:    Level: advanced

1008: .keywords: database accessing

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


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

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

1038:   return(0);
1039: }

1041: #undef __FUNCT__  
1043: /*@C
1044:    AODataSegmentGetInfo - Gets the blocksize and type of a data segment

1046:    Not collective

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

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

1057:    Level: advanced

1059: .keywords: database accessing

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


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

1078:   return(0);
1079: }

1081: #undef __FUNCT__  
1083: /*@C
1084:    AODataView - Displays an application ordering.

1086:    Collective on AOData and PetscViewer

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

1092:    Level: intermediate

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

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


1105: .keywords: database viewing

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

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

1121: #undef __FUNCT__  
1123: static int AODataAliasDestroy_Private(AODataAlias *aliases)
1124: {
1125:   AODataAlias *t = aliases;
1126:   int         ierr;

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

1145: #undef __FUNCT__  
1147: int AODataAliasAdd(AOData aodata,char *alias,char *name)
1148: {
1149:   AODataAlias *t = aodata->aliases;
1150:   int         ierr;

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

1167: #undef __FUNCT__  
1169: /*@C
1170:    AODataDestroy - Destroys an application ordering set.

1172:    Collective on AOData

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

1177:    Level: intermediate

1179: .keywords: destroy, database

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


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

1196:   return(0);
1197: }

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

1206:    Collective on AOData and AO

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

1213:    Level: advanced

1215: .keywords: database remapping

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

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

1230: #undef __FUNCT__  
1232: /*@C
1233:    AODataKeyGetAdjacency - Gets the adjacency graph for a key.

1235:    Collective on AOData

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

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

1244:    Level: advanced

1246: .keywords: database, adjacency graph

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

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

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

1268:     Collective on AOData

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

1274:    Level: advanced

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

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

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

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


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

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

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

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

1338:   PetscObjectPublishBaseEnd(obj);
1339: #endif

1341:   return(0);
1342: }

1344: #undef __FUNCT__  
1346: /*@C
1347:    AODataKeyRemove - Remove a data key from a AOData database.

1349:    Collective on AOData

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

1355:    Level: advanced

1357: .keywords: database removal

1359: .seealso:
1360: @*/
1361: int AODataKeyRemove(AOData aodata,char *name)
1362: {

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

1371: #undef __FUNCT__  
1373: /*@C
1374:    AODataSegmentRemove - Remove a data segment from a AOData database.

1376:    Collective on AOData

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

1383:    Level: advanced

1385: .keywords: database removal

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

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

1399: #undef __FUNCT__  
1401: /*@C
1402:    AODataKeyAdd - Add another data key to a AOData database.

1404:    Collective on AOData

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

1412:    Level: advanced

1414: .keywords: database additions

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


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

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

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

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

1457:   key->nlocal        = nlocal;

1459:   aodata->nkeys++;

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

1471:   return(0);
1472: }

1474: #undef __FUNCT__  
1476: /*@C
1477:    AODataSegmentAdd - Adds another data segment to a AOData database.

1479:    Collective on AOData

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

1491:    Level: advanced

1493: .keywords: database additions

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


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

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

1521: #undef __FUNCT__  
1523: /*@C
1524:    AODataSegmentAddIS - Add another data segment to a AOData database.

1526:    Collective on AOData and IS

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

1537:    Level: advanced

1539: .keywords: database additions

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


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