Actual source code: pmap.c
petsc-3.7.3 2016-08-01
2: /*
3: This file contains routines for basic map object implementation.
4: */
6: #include <petscis.h> /*I "petscis.h" I*/
7: #include <petscsf.h>
11: /*@
12: PetscLayoutCreate - Allocates PetscLayout space and sets the map contents to the default.
14: Collective on MPI_Comm
16: Input Parameters:
17: + comm - the MPI communicator
18: - map - pointer to the map
20: Level: advanced
22: Notes:
23: Typical calling sequence
24: .vb
25: PetscLayoutCreate(MPI_Comm,PetscLayout *);
26: PetscLayoutSetBlockSize(PetscLayout,1);
27: PetscLayoutSetSize(PetscLayout,N) // or PetscLayoutSetLocalSize(PetscLayout,n);
28: PetscLayoutSetUp(PetscLayout);
29: .ve
30: Optionally use any of the following:
32: + PetscLayoutGetSize(PetscLayout,PetscInt *);
33: . PetscLayoutGetLocalSize(PetscLayout,PetscInt *);
34: . PetscLayoutGetRange(PetscLayout,PetscInt *rstart,PetscInt *rend);
35: . PetscLayoutGetRanges(PetscLayout,const PetscInt *range[]);
36: - PetscLayoutDestroy(PetscLayout*);
38: The PetscLayout object and methods are intended to be used in the PETSc Vec and Mat implementions; it is often not needed in
39: user codes unless you really gain something in their use.
41: .seealso: PetscLayoutSetLocalSize(), PetscLayoutSetSize(), PetscLayoutGetSize(), PetscLayoutGetLocalSize(), PetscLayout, PetscLayoutDestroy(),
42: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize(), PetscLayoutSetUp()
44: @*/
45: PetscErrorCode PetscLayoutCreate(MPI_Comm comm,PetscLayout *map)
46: {
50: PetscNew(map);
52: (*map)->comm = comm;
53: (*map)->bs = -1;
54: (*map)->n = -1;
55: (*map)->N = -1;
56: (*map)->range = 0;
57: (*map)->rstart = 0;
58: (*map)->rend = 0;
59: return(0);
60: }
64: /*@
65: PetscLayoutDestroy - Frees a map object and frees its range if that exists.
67: Collective on MPI_Comm
69: Input Parameters:
70: . map - the PetscLayout
72: Level: developer
74: Note:
75: The PetscLayout object and methods are intended to be used in the PETSc Vec and Mat implementions; it is
76: recommended they not be used in user codes unless you really gain something in their use.
78: .seealso: PetscLayoutSetLocalSize(), PetscLayoutSetSize(), PetscLayoutGetSize(), PetscLayoutGetLocalSize(), PetscLayout, PetscLayoutCreate(),
79: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize(), PetscLayoutSetUp()
81: @*/
82: PetscErrorCode PetscLayoutDestroy(PetscLayout *map)
83: {
87: if (!*map) return(0);
88: if (!(*map)->refcnt--) {
89: PetscFree((*map)->range);
90: ISLocalToGlobalMappingDestroy(&(*map)->mapping);
91: PetscFree((*map));
92: }
93: *map = NULL;
94: return(0);
95: }
99: /*@
100: PetscLayoutSetUp - given a map where you have set either the global or local
101: size sets up the map so that it may be used.
103: Collective on MPI_Comm
105: Input Parameters:
106: . map - pointer to the map
108: Level: developer
110: Notes: Typical calling sequence
111: $ PetscLayoutCreate(MPI_Comm,PetscLayout *);
112: $ PetscLayoutSetBlockSize(PetscLayout,1);
113: $ PetscLayoutSetSize(PetscLayout,n) or PetscLayoutSetLocalSize(PetscLayout,N); or both
114: $ PetscLayoutSetUp(PetscLayout);
115: $ PetscLayoutGetSize(PetscLayout,PetscInt *);
118: If the local size, global size are already set and range exists then this does nothing.
120: .seealso: PetscLayoutSetLocalSize(), PetscLayoutSetSize(), PetscLayoutGetSize(), PetscLayoutGetLocalSize(), PetscLayout, PetscLayoutDestroy(),
121: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize(), PetscLayoutCreate()
122: @*/
123: PetscErrorCode PetscLayoutSetUp(PetscLayout map)
124: {
125: PetscMPIInt rank,size;
126: PetscInt p;
130: if ((map->n >= 0) && (map->N >= 0) && (map->range)) return(0);
132: if (map->n > 0 && map->bs > 1) {
133: if (map->n % map->bs) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Local matrix size %D must be divisible by blocksize %D",map->n,map->bs);
134: }
135: if (map->N > 0 && map->bs > 1) {
136: if (map->N % map->bs) SETERRQ2(map->comm,PETSC_ERR_PLIB,"Global matrix size %D must be divisible by blocksize %D",map->N,map->bs);
137: }
139: MPI_Comm_size(map->comm, &size);
140: MPI_Comm_rank(map->comm, &rank);
141: if (map->n > 0) map->n = map->n/PetscAbs(map->bs);
142: if (map->N > 0) map->N = map->N/PetscAbs(map->bs);
143: PetscSplitOwnership(map->comm,&map->n,&map->N);
144: map->n = map->n*PetscAbs(map->bs);
145: map->N = map->N*PetscAbs(map->bs);
146: if (!map->range) {
147: PetscMalloc1(size+1, &map->range);
148: }
149: MPI_Allgather(&map->n, 1, MPIU_INT, map->range+1, 1, MPIU_INT, map->comm);
151: map->range[0] = 0;
152: for (p = 2; p <= size; p++) map->range[p] += map->range[p-1];
154: map->rstart = map->range[rank];
155: map->rend = map->range[rank+1];
156: return(0);
157: }
161: /*@
162: PetscLayoutDuplicate - creates a new PetscLayout with the same information as a given one. If the PetscLayout already exists it is destroyed first.
164: Collective on PetscLayout
166: Input Parameter:
167: . in - input PetscLayout to be duplicated
169: Output Parameter:
170: . out - the copy
172: Level: developer
174: Notes: PetscLayoutSetUp() does not need to be called on the resulting PetscLayout
176: .seealso: PetscLayoutCreate(), PetscLayoutDestroy(), PetscLayoutSetUp(), PetscLayoutReference()
177: @*/
178: PetscErrorCode PetscLayoutDuplicate(PetscLayout in,PetscLayout *out)
179: {
180: PetscMPIInt size;
182: MPI_Comm comm = in->comm;
185: PetscLayoutDestroy(out);
186: PetscLayoutCreate(comm,out);
187: MPI_Comm_size(comm,&size);
188: PetscMemcpy(*out,in,sizeof(struct _n_PetscLayout));
189: PetscMalloc1(size+1,&(*out)->range);
190: PetscMemcpy((*out)->range,in->range,(size+1)*sizeof(PetscInt));
192: (*out)->refcnt = 0;
193: return(0);
194: }
198: /*@
199: PetscLayoutReference - Causes a PETSc Vec or Mat to share a PetscLayout with one that already exists. Used by Vec/MatDuplicate_XXX()
201: Collective on PetscLayout
203: Input Parameter:
204: . in - input PetscLayout to be copied
206: Output Parameter:
207: . out - the reference location
209: Level: developer
211: Notes: PetscLayoutSetUp() does not need to be called on the resulting PetscLayout
213: If the out location already contains a PetscLayout it is destroyed
215: .seealso: PetscLayoutCreate(), PetscLayoutDestroy(), PetscLayoutSetUp(), PetscLayoutDuplicate()
216: @*/
217: PetscErrorCode PetscLayoutReference(PetscLayout in,PetscLayout *out)
218: {
222: in->refcnt++;
223: PetscLayoutDestroy(out);
224: *out = in;
225: return(0);
226: }
230: /*@
231: PetscLayoutSetISLocalToGlobalMapping - sets a ISLocalGlobalMapping into a PetscLayout
233: Collective on PetscLayout
235: Input Parameter:
236: + in - input PetscLayout
237: - ltog - the local to global mapping
240: Level: developer
242: Notes: PetscLayoutSetUp() does not need to be called on the resulting PetscLayout
244: If the ltog location already contains a PetscLayout it is destroyed
246: .seealso: PetscLayoutCreate(), PetscLayoutDestroy(), PetscLayoutSetUp(), PetscLayoutDuplicate()
247: @*/
248: PetscErrorCode PetscLayoutSetISLocalToGlobalMapping(PetscLayout in,ISLocalToGlobalMapping ltog)
249: {
251: PetscInt bs;
254: ISLocalToGlobalMappingGetBlockSize(ltog,&bs);
255: if (in->bs > 0 && in->bs != bs) SETERRQ2(in->comm,PETSC_ERR_PLIB,"Blocksize of layout %D must match that of mapping %D",in->bs,bs);
256: PetscObjectReference((PetscObject)ltog);
257: ISLocalToGlobalMappingDestroy(&in->mapping);
258: in->mapping = ltog;
259: return(0);
260: }
264: /*@
265: PetscLayoutSetLocalSize - Sets the local size for a PetscLayout object.
267: Collective on PetscLayout
269: Input Parameters:
270: + map - pointer to the map
271: - n - the local size
273: Level: developer
275: Notes:
276: Call this after the call to PetscLayoutCreate()
278: .seealso: PetscLayoutCreate(), PetscLayoutSetSize(), PetscLayoutGetSize(), PetscLayoutGetLocalSize(), PetscLayoutSetUp()
279: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize()
280: @*/
281: PetscErrorCode PetscLayoutSetLocalSize(PetscLayout map,PetscInt n)
282: {
284: if (map->bs > 1 && n % map->bs) SETERRQ2(map->comm,PETSC_ERR_ARG_INCOMP,"Local size %D not compatible with block size %D",n,map->bs);
285: map->n = n;
286: return(0);
287: }
289: /*@C
290: PetscLayoutGetLocalSize - Gets the local size for a PetscLayout object.
292: Not Collective
294: Input Parameters:
295: . map - pointer to the map
297: Output Parameters:
298: . n - the local size
300: Level: developer
302: Notes:
303: Call this after the call to PetscLayoutSetUp()
305: Fortran Notes:
306: Not available from Fortran
308: .seealso: PetscLayoutCreate(), PetscLayoutSetSize(), PetscLayoutGetSize(), PetscLayoutGetLocalSize(), PetscLayoutSetUp()
309: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize()
311: @*/
314: PetscErrorCode PetscLayoutGetLocalSize(PetscLayout map,PetscInt *n)
315: {
317: *n = map->n;
318: return(0);
319: }
323: /*@
324: PetscLayoutSetSize - Sets the global size for a PetscLayout object.
326: Logically Collective on PetscLayout
328: Input Parameters:
329: + map - pointer to the map
330: - n - the global size
332: Level: developer
334: Notes:
335: Call this after the call to PetscLayoutCreate()
337: .seealso: PetscLayoutCreate(), PetscLayoutSetLocalSize(), PetscLayoutGetLocalSize(), PetscLayoutGetSize(), PetscLayoutSetUp()
338: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize()
339: @*/
340: PetscErrorCode PetscLayoutSetSize(PetscLayout map,PetscInt n)
341: {
343: map->N = n;
344: return(0);
345: }
349: /*@
350: PetscLayoutGetSize - Gets the global size for a PetscLayout object.
352: Not Collective
354: Input Parameters:
355: . map - pointer to the map
357: Output Parameters:
358: . n - the global size
360: Level: developer
362: Notes:
363: Call this after the call to PetscLayoutSetUp()
365: .seealso: PetscLayoutCreate(), PetscLayoutSetLocalSize(), PetscLayoutGetLocalSize(), PetscLayoutSetSize(), PetscLayoutSetUp()
366: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetBlockSize()
367: @*/
368: PetscErrorCode PetscLayoutGetSize(PetscLayout map,PetscInt *n)
369: {
371: *n = map->N;
372: return(0);
373: }
377: /*@
378: PetscLayoutSetBlockSize - Sets the block size for a PetscLayout object.
380: Logically Collective on PetscLayout
382: Input Parameters:
383: + map - pointer to the map
384: - bs - the size
386: Level: developer
388: Notes:
389: Call this after the call to PetscLayoutCreate()
391: .seealso: PetscLayoutCreate(), PetscLayoutSetLocalSize(), PetscLayoutGetLocalSize(), PetscLayoutGetBlockSize(),
392: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetSize(), PetscLayoutGetSize(), PetscLayoutSetUp()
393: @*/
394: PetscErrorCode PetscLayoutSetBlockSize(PetscLayout map,PetscInt bs)
395: {
397: if (bs < 0) return(0);
398: if (map->n > 0 && map->n % bs) SETERRQ2(map->comm,PETSC_ERR_ARG_INCOMP,"Local size %D not compatible with block size %D",map->n,bs);
399: if (map->range && map->bs > 0 && map->bs != bs) SETERRQ2(map->comm,PETSC_ERR_ARG_INCOMP,"Cannot change block size %D to %D",map->bs,bs);
400: if (map->mapping) {
401: PetscInt lbs;
404: ISLocalToGlobalMappingGetBlockSize(map->mapping,&lbs);
405: if (lbs != bs) SETERRQ2(map->comm,PETSC_ERR_PLIB,"Blocksize of localtoglobalmapping %D must match that of layout %D",lbs,bs);
406: }
407: map->bs = bs;
408: return(0);
409: }
413: /*@
414: PetscLayoutGetBlockSize - Gets the block size for a PetscLayout object.
416: Not Collective
418: Input Parameters:
419: . map - pointer to the map
421: Output Parameters:
422: . bs - the size
424: Level: developer
426: Notes:
427: Call this after the call to PetscLayoutSetUp()
429: .seealso: PetscLayoutCreate(), PetscLayoutSetLocalSize(), PetscLayoutGetLocalSize(), PetscLayoutSetSize(), PetscLayoutSetUp()
430: PetscLayoutGetRange(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetSize()
431: @*/
432: PetscErrorCode PetscLayoutGetBlockSize(PetscLayout map,PetscInt *bs)
433: {
435: *bs = PetscAbs(map->bs);
436: return(0);
437: }
441: /*@
442: PetscLayoutGetRange - gets the range of values owned by this process
444: Not Collective
446: Input Parameters:
447: . map - pointer to the map
449: Output Parameters:
450: + rstart - first index owned by this process
451: - rend - one more than the last index owned by this process
453: Level: developer
455: Notes:
456: Call this after the call to PetscLayoutSetUp()
458: .seealso: PetscLayoutCreate(), PetscLayoutSetLocalSize(), PetscLayoutGetLocalSize(), PetscLayoutSetSize(),
459: PetscLayoutGetSize(), PetscLayoutGetRanges(), PetscLayoutSetBlockSize(), PetscLayoutGetSize(), PetscLayoutSetUp()
460: @*/
461: PetscErrorCode PetscLayoutGetRange(PetscLayout map,PetscInt *rstart,PetscInt *rend)
462: {
464: if (rstart) *rstart = map->rstart;
465: if (rend) *rend = map->rend;
466: return(0);
467: }
469: /*@C
470: PetscLayoutGetRanges - gets the range of values owned by all processes
472: Not Collective
474: Input Parameters:
475: . map - pointer to the map
477: Output Parameters:
478: . range - start of each processors range of indices (the final entry is one more then the
479: last index on the last process)
481: Level: developer
483: Notes:
484: Call this after the call to PetscLayoutSetUp()
486: Fortran Notes:
487: Not available from Fortran
489: .seealso: PetscLayoutCreate(), PetscLayoutSetLocalSize(), PetscLayoutGetLocalSize(), PetscLayoutSetSize(),
490: PetscLayoutGetSize(), PetscLayoutGetRange(), PetscLayoutSetBlockSize(), PetscLayoutGetSize(), PetscLayoutSetUp()
492: @*/
495: PetscErrorCode PetscLayoutGetRanges(PetscLayout map,const PetscInt *range[])
496: {
498: *range = map->range;
499: return(0);
500: }
504: /*@C
505: PetscSFSetGraphLayout - Set a parallel star forest via global indices and a PetscLayout
507: Collective
509: Input Arguments:
510: + sf - star forest
511: . layout - PetscLayout defining the global space
512: . nleaves - number of leaf vertices on the current process, each of these references a root on any process
513: . ilocal - locations of leaves in leafdata buffers, pass NULL for contiguous storage
514: - iremote - remote locations of root vertices for each leaf on the current process
516: Level: intermediate
518: .seealso: PetscSFCreate(), PetscSFView(), PetscSFSetGraph(), PetscSFGetGraph()
519: @*/
520: PetscErrorCode PetscSFSetGraphLayout(PetscSF sf,PetscLayout layout,PetscInt nleaves,const PetscInt *ilocal,PetscCopyMode localmode,const PetscInt *iremote)
521: {
523: PetscInt i,nroots;
524: PetscSFNode *remote;
527: PetscLayoutGetLocalSize(layout,&nroots);
528: PetscMalloc1(nleaves,&remote);
529: for (i=0; i<nleaves; i++) {
530: PetscInt owner = -1;
531: PetscLayoutFindOwner(layout,iremote[i],&owner);
532: remote[i].rank = owner;
533: remote[i].index = iremote[i] - layout->range[owner];
534: }
535: PetscSFSetGraph(sf,nroots,nleaves,ilocal,localmode,remote,PETSC_OWN_POINTER);
536: return(0);
537: }