Actual source code: pmap.c
2: /*
3: This file contains routines for basic map object implementation.
4: */
6: #include <petsc/private/isimpl.h>
8: /*@
9: PetscLayoutCreate - Allocates PetscLayout space and sets the PetscLayout contents to the default.
11: Collective
13: Input Parameters:
14: . comm - the MPI communicator
16: Output Parameters:
17: . map - the new PetscLayout
19: Level: advanced
21: Notes:
22: Typical calling sequence
23: .vb
24: PetscLayoutCreate(MPI_Comm,PetscLayout *);
25: PetscLayoutSetBlockSize(PetscLayout,bs);
26: PetscLayoutSetSize(PetscLayout,N); // or PetscLayoutSetLocalSize(PetscLayout,n);
27: PetscLayoutSetUp(PetscLayout);
28: .ve
29: Alternatively,
30: .vb
31: PetscLayoutCreateFromSizes(comm,n,N,bs,&layout);
32: .ve
34: Optionally use any of the following
35: .vb
36: PetscLayoutGetSize(PetscLayout,PetscInt *);
37: PetscLayoutGetLocalSize(PetscLayout,PetscInt *);
38: PetscLayoutGetRange(PetscLayout,PetscInt *rstart,PetscInt *rend);
39: PetscLayoutGetRanges(PetscLayout,const PetscInt *range[]);
40: PetscLayoutDestroy(PetscLayout*);
41: .ve
43: The PetscLayout object and methods are intended to be used in the PETSc Vec and Mat implementations; it is often not needed in
44: user codes unless you really gain something in their use.
46: .seealso: `PetscLayoutSetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayout`, `PetscLayoutDestroy()`,
47: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`, `PetscLayoutSetUp()`,
48: `PetscLayoutCreateFromSizes()`
50: @*/
51: PetscErrorCode PetscLayoutCreate(MPI_Comm comm, PetscLayout *map)
52: {
53: PetscNew(map);
54: MPI_Comm_size(comm, &(*map)->size);
55: (*map)->comm = comm;
56: (*map)->bs = -1;
57: (*map)->n = -1;
58: (*map)->N = -1;
59: (*map)->range = NULL;
60: (*map)->range_alloc = PETSC_TRUE;
61: (*map)->rstart = 0;
62: (*map)->rend = 0;
63: (*map)->setupcalled = PETSC_FALSE;
64: (*map)->oldn = -1;
65: (*map)->oldN = -1;
66: (*map)->oldbs = -1;
67: return 0;
68: }
70: /*@
71: PetscLayoutCreateFromSizes - Allocates PetscLayout space, sets the layout sizes, and sets the layout up.
73: Collective
75: Input Parameters:
76: + comm - the MPI communicator
77: . n - the local size (or PETSC_DECIDE)
78: . N - the global size (or PETSC_DECIDE)
79: - bs - the block size (or PETSC_DECIDE)
81: Output Parameters:
82: . map - the new PetscLayout
84: Level: advanced
86: Notes:
87: $ PetscLayoutCreateFromSizes(comm,n,N,bs,&layout);
88: is a shorthand for
89: .vb
90: PetscLayoutCreate(comm,&layout);
91: PetscLayoutSetLocalSize(layout,n);
92: PetscLayoutSetSize(layout,N);
93: PetscLayoutSetBlockSize(layout,bs);
94: PetscLayoutSetUp(layout);
95: .ve
97: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayout`, `PetscLayoutDestroy()`,
98: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`, `PetscLayoutSetUp()`, `PetscLayoutCreateFromRanges()`
100: @*/
101: PetscErrorCode PetscLayoutCreateFromSizes(MPI_Comm comm, PetscInt n, PetscInt N, PetscInt bs, PetscLayout *map)
102: {
103: PetscLayoutCreate(comm, map);
104: PetscLayoutSetLocalSize(*map, n);
105: PetscLayoutSetSize(*map, N);
106: PetscLayoutSetBlockSize(*map, bs);
107: PetscLayoutSetUp(*map);
108: return 0;
109: }
111: /*@
112: PetscLayoutDestroy - Frees a map object and frees its range if that exists.
114: Collective
116: Input Parameters:
117: . map - the PetscLayout
119: Level: developer
121: Note:
122: The PetscLayout object and methods are intended to be used in the PETSc Vec and Mat implementations; it is
123: recommended they not be used in user codes unless you really gain something in their use.
125: .seealso: `PetscLayoutSetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayout`, `PetscLayoutCreate()`,
126: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`, `PetscLayoutSetUp()`
128: @*/
129: PetscErrorCode PetscLayoutDestroy(PetscLayout *map)
130: {
131: if (!*map) return 0;
132: if (!(*map)->refcnt--) {
133: if ((*map)->range_alloc) PetscFree((*map)->range);
134: ISLocalToGlobalMappingDestroy(&(*map)->mapping);
135: PetscFree((*map));
136: }
137: *map = NULL;
138: return 0;
139: }
141: /*@
142: PetscLayoutCreateFromRanges - Creates a new PetscLayout with the given ownership ranges and sets it up.
144: Collective
146: Input Parameters:
147: + comm - the MPI communicator
148: . range - the array of ownership ranges for each rank with length commsize+1
149: . mode - the copy mode for range
150: - bs - the block size (or PETSC_DECIDE)
152: Output Parameters:
153: . newmap - the new PetscLayout
155: Level: developer
157: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayout`, `PetscLayoutDestroy()`,
158: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`, `PetscLayoutSetUp()`, `PetscLayoutCreateFromSizes()`
160: @*/
161: PetscErrorCode PetscLayoutCreateFromRanges(MPI_Comm comm, const PetscInt range[], PetscCopyMode mode, PetscInt bs, PetscLayout *newmap)
162: {
163: PetscLayout map;
164: PetscMPIInt rank;
166: MPI_Comm_rank(comm, &rank);
167: PetscLayoutCreate(comm, &map);
168: PetscLayoutSetBlockSize(map, bs);
169: switch (mode) {
170: case PETSC_COPY_VALUES:
171: PetscMalloc1(map->size + 1, &map->range);
172: PetscArraycpy(map->range, range, map->size + 1);
173: break;
174: case PETSC_USE_POINTER:
175: map->range_alloc = PETSC_FALSE;
176: break;
177: default:
178: map->range = (PetscInt *)range;
179: break;
180: }
181: map->rstart = map->range[rank];
182: map->rend = map->range[rank + 1];
183: map->n = map->rend - map->rstart;
184: map->N = map->range[map->size];
185: if (PetscDefined(USE_DEBUG)) { /* just check that n, N and bs are consistent */
186: PetscInt tmp;
187: MPIU_Allreduce(&map->n, &tmp, 1, MPIU_INT, MPI_SUM, map->comm);
191: }
192: /* lock the layout */
193: map->setupcalled = PETSC_TRUE;
194: map->oldn = map->n;
195: map->oldN = map->N;
196: map->oldbs = map->bs;
197: *newmap = map;
198: return 0;
199: }
201: /*@
202: PetscLayoutSetUp - given a map where you have set either the global or local
203: size sets up the map so that it may be used.
205: Collective
207: Input Parameters:
208: . map - pointer to the map
210: Level: developer
212: Notes:
213: Typical calling sequence
214: $ PetscLayoutCreate(MPI_Comm,PetscLayout *);
215: $ PetscLayoutSetBlockSize(PetscLayout,1);
216: $ PetscLayoutSetSize(PetscLayout,n) or PetscLayoutSetLocalSize(PetscLayout,N); or both
217: $ PetscLayoutSetUp(PetscLayout);
218: $ PetscLayoutGetSize(PetscLayout,PetscInt *);
220: If range exists, and local size is not set, everything gets computed from the range.
222: If the local size, global size are already set and range exists then this does nothing.
224: .seealso: `PetscLayoutSetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayout`, `PetscLayoutDestroy()`,
225: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`, `PetscLayoutCreate()`, `PetscSplitOwnership()`
226: @*/
227: PetscErrorCode PetscLayoutSetUp(PetscLayout map)
228: {
229: PetscMPIInt rank;
230: PetscInt p;
233: map->oldn, map->oldN, map->n, map->N);
234: if (map->setupcalled) return 0;
239: MPI_Comm_rank(map->comm, &rank);
240: if (map->n > 0) map->n = map->n / PetscAbs(map->bs);
241: if (map->N > 0) map->N = map->N / PetscAbs(map->bs);
242: PetscSplitOwnership(map->comm, &map->n, &map->N);
243: map->n = map->n * PetscAbs(map->bs);
244: map->N = map->N * PetscAbs(map->bs);
245: if (!map->range) PetscMalloc1(map->size + 1, &map->range);
246: MPI_Allgather(&map->n, 1, MPIU_INT, map->range + 1, 1, MPIU_INT, map->comm);
248: map->range[0] = 0;
249: for (p = 2; p <= map->size; p++) map->range[p] += map->range[p - 1];
251: map->rstart = map->range[rank];
252: map->rend = map->range[rank + 1];
254: /* lock the layout */
255: map->setupcalled = PETSC_TRUE;
256: map->oldn = map->n;
257: map->oldN = map->N;
258: map->oldbs = map->bs;
259: return 0;
260: }
262: /*@
263: PetscLayoutDuplicate - creates a new PetscLayout with the same information as a given one. If the PetscLayout already exists it is destroyed first.
265: Collective on PetscLayout
267: Input Parameter:
268: . in - input PetscLayout to be duplicated
270: Output Parameter:
271: . out - the copy
273: Level: developer
275: Notes:
276: PetscLayoutSetUp() does not need to be called on the resulting PetscLayout
278: .seealso: `PetscLayoutCreate()`, `PetscLayoutDestroy()`, `PetscLayoutSetUp()`, `PetscLayoutReference()`
279: @*/
280: PetscErrorCode PetscLayoutDuplicate(PetscLayout in, PetscLayout *out)
281: {
282: MPI_Comm comm = in->comm;
284: PetscLayoutDestroy(out);
285: PetscLayoutCreate(comm, out);
286: PetscMemcpy(*out, in, sizeof(struct _n_PetscLayout));
287: if (in->range) {
288: PetscMalloc1((*out)->size + 1, &(*out)->range);
289: PetscArraycpy((*out)->range, in->range, (*out)->size + 1);
290: }
291: (*out)->refcnt = 0;
292: return 0;
293: }
295: /*@
296: PetscLayoutReference - Causes a PETSc Vec or Mat to share a PetscLayout with one that already exists. Used by Vec/MatDuplicate_XXX()
298: Collective on PetscLayout
300: Input Parameter:
301: . in - input PetscLayout to be copied
303: Output Parameter:
304: . out - the reference location
306: Level: developer
308: Notes:
309: PetscLayoutSetUp() does not need to be called on the resulting PetscLayout
311: If the out location already contains a PetscLayout it is destroyed
313: .seealso: `PetscLayoutCreate()`, `PetscLayoutDestroy()`, `PetscLayoutSetUp()`, `PetscLayoutDuplicate()`
314: @*/
315: PetscErrorCode PetscLayoutReference(PetscLayout in, PetscLayout *out)
316: {
317: in->refcnt++;
318: PetscLayoutDestroy(out);
319: *out = in;
320: return 0;
321: }
323: /*@
324: PetscLayoutSetISLocalToGlobalMapping - sets a ISLocalGlobalMapping into a PetscLayout
326: Collective on PetscLayout
328: Input Parameters:
329: + in - input PetscLayout
330: - ltog - the local to global mapping
332: Level: developer
334: Notes:
335: PetscLayoutSetUp() does not need to be called on the resulting PetscLayout
337: If the ltog location already contains a PetscLayout it is destroyed
339: .seealso: `PetscLayoutCreate()`, `PetscLayoutDestroy()`, `PetscLayoutSetUp()`, `PetscLayoutDuplicate()`
340: @*/
341: PetscErrorCode PetscLayoutSetISLocalToGlobalMapping(PetscLayout in, ISLocalToGlobalMapping ltog)
342: {
343: if (ltog) {
344: PetscInt bs;
346: ISLocalToGlobalMappingGetBlockSize(ltog, &bs);
348: PetscObjectReference((PetscObject)ltog);
349: }
350: ISLocalToGlobalMappingDestroy(&in->mapping);
351: in->mapping = ltog;
352: return 0;
353: }
355: /*@
356: PetscLayoutSetLocalSize - Sets the local size for a PetscLayout object.
358: Collective on PetscLayout
360: Input Parameters:
361: + map - pointer to the map
362: - n - the local size
364: Level: developer
366: Notes:
367: Call this after the call to PetscLayoutCreate()
369: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutSetUp()`
370: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`
371: @*/
372: PetscErrorCode PetscLayoutSetLocalSize(PetscLayout map, PetscInt n)
373: {
375: map->n = n;
376: return 0;
377: }
379: /*@
380: PetscLayoutGetLocalSize - Gets the local size for a PetscLayout object.
382: Not Collective
384: Input Parameters:
385: . map - pointer to the map
387: Output Parameters:
388: . n - the local size
390: Level: developer
392: Notes:
393: Call this after the call to PetscLayoutSetUp()
395: Fortran Notes:
396: Not available from Fortran
398: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutSetUp()`
399: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`
401: @*/
402: PetscErrorCode PetscLayoutGetLocalSize(PetscLayout map, PetscInt *n)
403: {
404: *n = map->n;
405: return 0;
406: }
408: /*@
409: PetscLayoutSetSize - Sets the global size for a PetscLayout object.
411: Logically Collective on PetscLayout
413: Input Parameters:
414: + map - pointer to the map
415: - n - the global size
417: Level: developer
419: Notes:
420: Call this after the call to PetscLayoutCreate()
422: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutGetSize()`, `PetscLayoutSetUp()`
423: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`
424: @*/
425: PetscErrorCode PetscLayoutSetSize(PetscLayout map, PetscInt n)
426: {
427: map->N = n;
428: return 0;
429: }
431: /*@
432: PetscLayoutGetSize - Gets the global size for a PetscLayout object.
434: Not Collective
436: Input Parameters:
437: . map - pointer to the map
439: Output Parameters:
440: . n - the global size
442: Level: developer
444: Notes:
445: Call this after the call to PetscLayoutSetUp()
447: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutSetUp()`
448: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetBlockSize()`
449: @*/
450: PetscErrorCode PetscLayoutGetSize(PetscLayout map, PetscInt *n)
451: {
452: *n = map->N;
453: return 0;
454: }
456: /*@
457: PetscLayoutSetBlockSize - Sets the block size for a PetscLayout object.
459: Logically Collective on PetscLayout
461: Input Parameters:
462: + map - pointer to the map
463: - bs - the size
465: Level: developer
467: Notes:
468: Call this after the call to PetscLayoutCreate()
470: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutGetBlockSize()`,
471: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutSetUp()`
472: @*/
473: PetscErrorCode PetscLayoutSetBlockSize(PetscLayout map, PetscInt bs)
474: {
475: if (bs < 0) return 0;
477: if (map->mapping) {
478: PetscInt obs;
480: ISLocalToGlobalMappingGetBlockSize(map->mapping, &obs);
481: if (obs > 1) ISLocalToGlobalMappingSetBlockSize(map->mapping, bs);
482: }
483: map->bs = bs;
484: return 0;
485: }
487: /*@
488: PetscLayoutGetBlockSize - Gets the block size for a PetscLayout object.
490: Not Collective
492: Input Parameters:
493: . map - pointer to the map
495: Output Parameters:
496: . bs - the size
498: Level: developer
500: Notes:
501: Call this after the call to PetscLayoutSetUp()
503: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutSetSize()`, `PetscLayoutSetUp()`
504: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetSize()`
505: @*/
506: PetscErrorCode PetscLayoutGetBlockSize(PetscLayout map, PetscInt *bs)
507: {
508: *bs = PetscAbs(map->bs);
509: return 0;
510: }
512: /*@
513: PetscLayoutGetRange - gets the range of values owned by this process
515: Not Collective
517: Input Parameter:
518: . map - pointer to the map
520: Output Parameters:
521: + rstart - first index owned by this process
522: - rend - one more than the last index owned by this process
524: Level: developer
526: Notes:
527: Call this after the call to PetscLayoutSetUp()
529: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutSetSize()`,
530: `PetscLayoutGetSize()`, `PetscLayoutGetRanges()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetSize()`, `PetscLayoutSetUp()`
531: @*/
532: PetscErrorCode PetscLayoutGetRange(PetscLayout map, PetscInt *rstart, PetscInt *rend)
533: {
534: if (rstart) *rstart = map->rstart;
535: if (rend) *rend = map->rend;
536: return 0;
537: }
539: /*@C
540: PetscLayoutGetRanges - gets the range of values owned by all processes
542: Not Collective
544: Input Parameters:
545: . map - pointer to the map
547: Output Parameters:
548: . range - start of each processors range of indices (the final entry is one more than the
549: last index on the last process)
551: Level: developer
553: Notes:
554: Call this after the call to PetscLayoutSetUp()
556: Fortran Notes:
557: Not available from Fortran
559: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutSetSize()`,
560: `PetscLayoutGetSize()`, `PetscLayoutGetRange()`, `PetscLayoutSetBlockSize()`, `PetscLayoutGetSize()`, `PetscLayoutSetUp()`
562: @*/
563: PetscErrorCode PetscLayoutGetRanges(PetscLayout map, const PetscInt *range[])
564: {
565: *range = map->range;
566: return 0;
567: }
569: /*@
570: PetscLayoutCompare - Compares two layouts
572: Not Collective
574: Input Parameters:
575: + mapa - pointer to the first map
576: - mapb - pointer to the second map
578: Output Parameters:
579: . congruent - PETSC_TRUE if the two layouts are congruent, PETSC_FALSE otherwise
581: Level: beginner
583: Notes:
585: .seealso: `PetscLayoutCreate()`, `PetscLayoutSetLocalSize()`, `PetscLayoutGetLocalSize()`, `PetscLayoutGetBlockSize()`,
586: `PetscLayoutGetRange()`, `PetscLayoutGetRanges()`, `PetscLayoutSetSize()`, `PetscLayoutGetSize()`, `PetscLayoutSetUp()`
587: @*/
588: PetscErrorCode PetscLayoutCompare(PetscLayout mapa, PetscLayout mapb, PetscBool *congruent)
589: {
590: *congruent = PETSC_FALSE;
591: if (mapa->N == mapb->N && mapa->range && mapb->range && mapa->size == mapb->size) PetscArraycmp(mapa->range, mapb->range, mapa->size + 1, congruent);
592: return 0;
593: }