Actual source code: section.c
1: #include <private/meshimpl.h> /*I "petscdmmesh.h" I*/
2: #include <petscdmmesh_viewers.hh>
4: /* Logging support */
5: PetscClassId SECTIONREAL_CLASSID;
6: PetscLogEvent SectionReal_View;
7: PetscClassId SECTIONINT_CLASSID;
8: PetscLogEvent SectionInt_View;
9: PetscClassId SECTIONPAIR_CLASSID;
10: PetscLogEvent SectionPair_View;
14: PetscErrorCode SectionRealView_Sieve(SectionReal section, PetscViewer viewer)
15: {
16: PetscBool iascii, isbinary, isdraw;
20: PetscTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
21: PetscTypeCompare((PetscObject) viewer, PETSCVIEWERBINARY, &isbinary);
22: PetscTypeCompare((PetscObject) viewer, PETSCVIEWERDRAW, &isdraw);
24: if (iascii){
25: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
26: ALE::Obj<PETSC_MESH_TYPE> b;
27: const char *name;
29: SectionRealGetSection(section, s);
30: SectionRealGetBundle(section, b);
31: PetscObjectGetName((PetscObject) section, &name);
32: SectionView_Sieve_Ascii(b, s, name, viewer);
33: } else if (isbinary) {
34: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "Binary viewer not implemented for Section");
35: } else if (isdraw){
36: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "Draw viewer not implemented for Section");
37: } else {
38: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Viewer type %s not supported by this section object", ((PetscObject)viewer)->type_name);
39: }
40: return(0);
41: }
45: /*@C
46: SectionRealView - Views a Section object.
48: Collective on Section
50: Input Parameters:
51: + section - the Section
52: - viewer - an optional visualization context
54: Notes:
55: The available visualization contexts include
56: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
57: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
58: output where only the first processor opens
59: the file. All other processors send their
60: data to the first processor to print.
62: You can change the format the section is printed using the
63: option PetscViewerSetFormat().
65: The user can open alternative visualization contexts with
66: + PetscViewerASCIIOpen() - Outputs section to a specified file
67: . PetscViewerBinaryOpen() - Outputs section in binary to a
68: specified file; corresponding input uses SectionLoad()
69: . PetscViewerDrawOpen() - Outputs section to an X window display
71: The user can call PetscViewerSetFormat() to specify the output
72: format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
73: PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen). Available formats include
74: + PETSC_VIEWER_DEFAULT - default, prints section information
75: - PETSC_VIEWER_ASCII_VTK - outputs a VTK file describing the section
77: Level: beginner
79: Concepts: section^printing
80: Concepts: section^saving to disk
82: .seealso: VecView(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerBinaryOpen(), PetscViewerCreate()
83: @*/
84: PetscErrorCode SectionRealView(SectionReal section, PetscViewer viewer)
85: {
91: if (!viewer) {
92: PetscViewerASCIIGetStdout(((PetscObject)section)->comm,&viewer);
93: }
97: PetscLogEventBegin(SectionReal_View,0,0,0,0);
98: (*section->ops->view)(section, viewer);
99: PetscLogEventEnd(SectionReal_View,0,0,0,0);
100: return(0);
101: }
105: /*@C
106: SectionRealDuplicate - Create an equivalent Section object
108: Not collective
110: Input Parameter:
111: . section - the section object
113: Output Parameter:
114: . newSection - the duplicate
115:
116: Level: advanced
118: .seealso SectionRealCreate(), SectionRealSetSection()
119: @*/
120: PetscErrorCode SectionRealDuplicate(SectionReal section, SectionReal *newSection)
121: {
127: const ALE::Obj<PETSC_MESH_TYPE::real_section_type>& s = section->s;
128: ALE::Obj<PETSC_MESH_TYPE::real_section_type> t = new PETSC_MESH_TYPE::real_section_type(s->comm(), s->debug());
130: t->setAtlas(s->getAtlas());
131: t->allocateStorage();
132: t->copyBC(s);
133: SectionRealCreate(s->comm(), newSection);
134: SectionRealSetSection(*newSection, t);
135: SectionRealSetBundle(*newSection, section->b);
136: return(0);
137: }
141: /*@C
142: SectionRealGetSection - Gets the internal section object
144: Not collective
146: Input Parameter:
147: . section - the section object
149: Output Parameter:
150: . s - the internal section object
151:
152: Level: advanced
154: .seealso SectionRealCreate(), SectionRealSetSection()
155: @*/
156: PetscErrorCode SectionRealGetSection(SectionReal section, ALE::Obj<PETSC_MESH_TYPE::real_section_type>& s)
157: {
160: s = section->s;
161: return(0);
162: }
166: /*@C
167: SectionRealSetSection - Sets the internal section object
169: Not collective
171: Input Parameters:
172: + section - the section object
173: - s - the internal section object
174:
175: Level: advanced
177: .seealso SectionRealCreate(), SectionRealGetSection()
178: @*/
179: PetscErrorCode SectionRealSetSection(SectionReal section, const ALE::Obj<PETSC_MESH_TYPE::real_section_type>& s)
180: {
185: if (!s.isNull()) {PetscObjectSetName((PetscObject) section, s->getName().c_str());}
186: section->s = s;
187: return(0);
188: }
192: /*@C
193: SectionRealGetBundle - Gets the section bundle
195: Not collective
197: Input Parameter:
198: . section - the section object
200: Output Parameter:
201: . b - the section bundle
202:
203: Level: advanced
205: .seealso SectionRealCreate(), SectionRealGetSection(), SectionRealSetSection()
206: @*/
207: PetscErrorCode SectionRealGetBundle(SectionReal section, ALE::Obj<PETSC_MESH_TYPE>& b)
208: {
211: b = section->b;
212: return(0);
213: }
217: /*@C
218: SectionRealSetBundle - Sets the section bundle
220: Not collective
222: Input Parameters:
223: + section - the section object
224: - b - the section bundle
225:
226: Level: advanced
228: .seealso SectionRealCreate(), SectionRealGetSection(), SectionRealSetSection()
229: @*/
230: PetscErrorCode SectionRealSetBundle(SectionReal section, const ALE::Obj<PETSC_MESH_TYPE>& b)
231: {
234: section->b = b;
235: return(0);
236: }
240: /*@C
241: SectionRealCreate - Creates a Section object, used to manage data for an unstructured problem
242: described by a Sieve.
244: Collective on MPI_Comm
246: Input Parameter:
247: . comm - the processors that will share the global section
249: Output Parameters:
250: . section - the section object
252: Level: advanced
254: .seealso SectionRealDestroy(), SectionRealView()
255: @*/
256: PetscErrorCode SectionRealCreate(MPI_Comm comm, SectionReal *section)
257: {
259: SectionReal s;
263: *section = PETSC_NULL;
265: PetscHeaderCreate(s,_p_SectionReal,struct _SectionRealOps,SECTIONREAL_CLASSID,0,"SectionReal","Section","DM",comm,SectionRealDestroy,0);
266: s->ops->view = SectionRealView_Sieve;
267: s->ops->restrictClosure = SectionRealRestrict;
268: s->ops->update = SectionRealUpdate;
270: PetscObjectChangeTypeName((PetscObject) s, "sieve");
272: new(&s->s) ALE::Obj<PETSC_MESH_TYPE::real_section_type>(PETSC_MESH_TYPE::real_section_type(comm));
273: new(&s->b) ALE::Obj<PETSC_MESH_TYPE>(PETSC_NULL);
274: *section = s;
275: return(0);
276: }
280: /*@
281: SectionRealDestroy - Destroys a section.
283: Collective on Section
285: Input Parameter:
286: . section - the section object
288: Level: advanced
290: .seealso SectionRealCreate(), SectionRealView()
291: @*/
292: PetscErrorCode SectionRealDestroy(SectionReal *section)
293: {
297: if (!*section) return(0);
299: if (--((PetscObject)(*section))->refct > 0) return(0);
300: PetscHeaderDestroy(section);
301: return(0);
302: }
306: /*@
307: SectionRealDistribute - Distributes the sections.
309: Not Collective
311: Input Parameters:
312: + serialSection - The original Section object
313: - parallelMesh - The parallel DMMesh
315: Output Parameter:
316: . parallelSection - The distributed Section object
318: Level: intermediate
320: .keywords: mesh, section, distribute
321: .seealso: DMMeshCreate()
322: @*/
323: PetscErrorCode SectionRealDistribute(SectionReal serialSection, DM parallelMesh, SectionReal *parallelSection)
324: {
325: ALE::Obj<PETSC_MESH_TYPE::real_section_type> oldSection;
326: ALE::Obj<PETSC_MESH_TYPE> m;
330: SectionRealGetSection(serialSection, oldSection);
331: DMMeshGetMesh(parallelMesh, m);
332: SectionRealCreate(oldSection->comm(), parallelSection);
333: #ifdef PETSC_OPT_SIEVE
334: ALE::Obj<PETSC_MESH_TYPE::real_section_type> newSection;
336: // We assume all integer sections are complete sections
337: newSection->setName(oldSection->getName());
338: newSection->setChart(m->getSieve()->getChart());
339: //distributeSection(oldSection, partition, m->getRenumbering(), m->getDistSendOverlap(), m->getDistRecvOverlap(), newSection);
340: SectionRealSetSection(*parallelSection, newSection);
341: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "Not working because the partition is unavailable");
342: #else
343: ALE::Obj<PETSC_MESH_TYPE::real_section_type> newSection = ALE::Distribution<PETSC_MESH_TYPE>::distributeSection(oldSection, m, m->getDistSendOverlap(), m->getDistRecvOverlap());
344: SectionRealSetSection(*parallelSection, newSection);
345: #endif
346: return(0);
347: }
351: /*@C
352: SectionRealRestrict - Restricts the SectionReal to a subset of the topology, returning an array of values.
354: Not collective
356: Input Parameters:
357: + section - the section object
358: - point - the Sieve point
360: Output Parameter:
361: . values - The values associated with the submesh
363: Level: advanced
365: .seealso SectionUpdate(), SectionCreate(), SectionView()
366: @*/
367: PetscErrorCode SectionRealRestrict(SectionReal section, PetscInt point, PetscScalar *values[])
368: {
372: try {
373: *values = (PetscScalar *) section->s->restrictPoint(point);
374: } catch(ALE::Exception e) {
375: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid argument: %s", e.message());
376: }
377: return(0);
378: }
382: /*@C
383: SectionRealUpdate - Updates the array of values associated to a subset of the topology in this Section.
385: Not collective
387: Input Parameters:
388: + section - the section object
389: . point - the Sieve point
390: . values - The values associated with the submesh
391: - mode - The insertion mode
393: Level: advanced
395: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
396: @*/
397: PetscErrorCode SectionRealUpdate(SectionReal section, PetscInt point, const PetscScalar values[], InsertMode mode)
398: {
402: #ifdef PETSC_USE_COMPLEX
403: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "SectionReal does not support complex updates");
404: #else
405: try {
406: if (mode == INSERT_VALUES) {
407: section->b->update(section->s, point, values);
408: } else if (mode == ADD_VALUES) {
409: section->b->updateAdd(section->s, point, values);
410: } else {
411: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG, "Invalid insertion mode: %d", mode);
412: }
413: } catch(ALE::Exception e) {
414: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid argument: %s", e.message());
415: }
416: #endif
417: return(0);
418: }
422: /*@C
423: SectionRealRestrictClosure - Returns an array with the values in a given closure
425: Not Collective
427: Input Parameters:
428: + section - The section
429: . mesh - The DMMesh object
430: - point - The sieve point
432: Output Parameter:
433: . array - The array full of values in the closure
435: Level: intermediate
437: .keywords: mesh, elements
438: .seealso: DMMeshCreate()
439: @*/
440: PetscErrorCode SectionRealRestrictClosure(SectionReal section, DM dm, PetscInt point, const PetscScalar *values[])
441: {
442: ALE::Obj<PETSC_MESH_TYPE> m;
443: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
444: PetscErrorCode ierr;
447: DMMeshGetMesh(dm, m);
448: SectionRealGetSection(section, s);
449: #ifdef PETSC_USE_COMPLEX
450: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "SectionReal does not support complex restriction");
451: #else
452: try {
453: *values = m->restrictClosure(s, point);
454: } catch(ALE::Exception e) {
455: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid argument: %s", e.message());
456: }
457: #endif
458: return(0);
459: }
463: /*@C
464: SectionRealRestrictClosure - Returns an array with the values in a given closure
466: Not Collective
468: Input Parameters:
469: + section - The section
470: . mesh - The DMMesh object
471: . point - The sieve point
472: . n - The array size
473: - array - The array to fill up
475: Output Parameter:
476: . array - The array full of values in the closure
478: Level: intermediate
480: .keywords: mesh, elements
481: .seealso: DMMeshCreate()
482: @*/
483: PetscErrorCode SectionRealRestrictClosure(SectionReal section, DM dm, PetscInt point, PetscInt n, PetscScalar values[])
484: {
485: ALE::Obj<PETSC_MESH_TYPE> m;
486: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
487: PetscErrorCode ierr;
490: DMMeshGetMesh(dm, m);
491: SectionRealGetSection(section, s);
492: #ifdef PETSC_USE_COMPLEX
493: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "SectionReal does not support complex restriction");
494: #else
495: try {
496: m->restrictClosure(s, point, values, n);
497: } catch(ALE::Exception e) {
498: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid argument: %s", e.message());
499: }
500: #endif
501: return(0);
502: }
506: /*@C
507: SectionRealUpdateClosure - Updates the values in a given closure from the array
509: Not Collective
511: Input Parameters:
512: + section - The section
513: . mesh - The DMMesh object
514: . point - The sieve point
515: . array - The array to fill up
516: - mode - The insertion mode
518: Output Parameter:
519: . array - The array full of values in the closure
521: Level: intermediate
523: .keywords: mesh, elements
524: .seealso: DMMeshCreate()
525: @*/
526: PetscErrorCode SectionRealUpdateClosure(SectionReal section, DM dm, PetscInt point, PetscScalar values[], InsertMode mode)
527: {
528: ALE::Obj<PETSC_MESH_TYPE> m;
529: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
530: PetscErrorCode ierr;
533: DMMeshGetMesh(dm, m);
534: SectionRealGetSection(section, s);
535: #ifdef PETSC_USE_COMPLEX
536: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "SectionReal does not support complex update");
537: #else
538: try {
539: if (mode == INSERT_VALUES) {
540: m->update(s, point, values);
541: } else if (mode == ADD_VALUES) {
542: m->updateAdd(s, point, values);
543: } else {
544: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG, "Invalid insertion mode: %d", mode);
545: }
546: } catch(ALE::Exception e) {
547: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid argument: %s", e.message());
548: }
549: #endif
550: return(0);
551: }
555: /*@
556: SectionRealComplete - Exchanges data across the mesh overlap.
558: Not collective
560: Input Parameter:
561: . section - the section object
563: Level: advanced
565: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
566: @*/
567: PetscErrorCode SectionRealComplete(SectionReal section)
568: {
569: Obj<PETSC_MESH_TYPE::real_section_type> s;
570: Obj<PETSC_MESH_TYPE> b;
574: SectionRealGetSection(section, s);
575: SectionRealGetBundle(section, b);
576: #if 0
577: ALE::Distribution<PETSC_MESH_TYPE>::completeSection(b, s);
578: #else
579: ALE::Completion::completeSectionAdd(b->getSendOverlap(), b->getRecvOverlap(), s, s);
580: #endif
581: return(0);
582: }
586: /*@
587: SectionRealZero - Zero out the entries
589: Not collective
591: Input Parameter:
592: . section - the section object
594: Level: advanced
596: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
597: @*/
598: PetscErrorCode SectionRealZero(SectionReal section)
599: {
600: Obj<PETSC_MESH_TYPE::real_section_type> s;
604: SectionRealGetSection(section, s);
605: s->zero();
606: return(0);
607: }
611: /*@
612: SectionRealGetFiberDimension - Get the size of the vector space attached to the point
614: Not collective
616: Input Parameters:
617: + section - the section object
618: - point - the Sieve point
620: Output Parameters:
621: . size - The fiber dimension
623: Level: advanced
625: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
626: @*/
627: PetscErrorCode SectionRealGetFiberDimension(SectionReal section, PetscInt point, PetscInt *size)
628: {
631: *size = section->s->getFiberDimension(point);
632: return(0);
633: }
637: /*@
638: SectionRealSetFiberDimension - Set the size of the vector space attached to the point
640: Not collective
642: Input Parameters:
643: + section - the section object
644: . point - the Sieve point
645: - size - The fiber dimension
647: Level: advanced
649: .seealso SectionRealSetFiberDimensionField(), SectionRealRestrict(), SectionRealCreate(), SectionRealView()
650: @*/
651: PetscErrorCode SectionRealSetFiberDimension(SectionReal section, PetscInt point, const PetscInt size)
652: {
655: section->s->setFiberDimension(point, size);
656: return(0);
657: }
661: /*@
662: SectionRealSetFiberDimensionField - Set the size of the vector space attached to the point for a given field
664: Not collective
666: Input Parameters:
667: + section - the section object
668: . point - the Sieve point
669: . size - The fiber dimension
670: - field - The field number
672: Level: advanced
674: .seealso SectionRealSetFiberDimension(), SectionRealRestrict(), SectionRealCreate(), SectionRealView()
675: @*/
676: PetscErrorCode SectionRealSetFiberDimensionField(SectionReal section, PetscInt point, const PetscInt size, const PetscInt field)
677: {
680: section->s->setFiberDimension(point, size, field);
681: return(0);
682: }
686: /*@
687: SectionRealGetSize - Gets the number of local dofs in this Section
689: Not collective
691: Input Parameter:
692: . section - the section object
694: Output Parameter:
695: . size - the section size
697: Level: advanced
699: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
700: @*/
701: PetscErrorCode SectionRealGetSize(SectionReal section, PetscInt *size)
702: {
703: Obj<PETSC_MESH_TYPE::real_section_type> s;
708: SectionRealGetSection(section, s);
709: *size = s->size();
710: return(0);
711: }
715: /*@
716: SectionRealAllocate - Allocate storage for this section
718: Not collective
720: Input Parameter:
721: . section - the section object
723: Level: advanced
725: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
726: @*/
727: PetscErrorCode SectionRealAllocate(SectionReal section)
728: {
731: section->b->allocate(section->s);
732: return(0);
733: }
737: /*@
738: SectionRealCreateLocalVector - Creates a vector with the local piece of the Section
740: Collective on DMMesh
742: Input Parameter:
743: . section - the Section
745: Output Parameter:
746: . localVec - the local vector
748: Level: advanced
750: Notes: The vector can safely be destroyed using VecDestroy().
751: .seealso DMMeshDestroy(), DMMeshCreate()
752: @*/
753: PetscErrorCode SectionRealCreateLocalVector(SectionReal section, Vec *localVec)
754: {
755: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
759: #ifdef PETSC_USE_COMPLEX
760: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "SectionReal does not support complex Vec");
761: #else
762: SectionRealGetSection(section, s);
763: VecCreateSeqWithArray(PETSC_COMM_SELF, s->getStorageSize(), s->restrictSpace(), localVec);
764: #endif
765: return(0);
766: }
770: /*@
771: SectionRealAddSpace - Add another field to this section
773: Collective on DMMesh
775: Input Parameter:
776: . section - the Section
778: Level: advanced
780: .seealso SectionRealCreate(), SectionRealGetFibration()
781: @*/
782: PetscErrorCode SectionRealAddSpace(SectionReal section)
783: {
784: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
788: SectionRealGetSection(section, s);
789: s->addSpace();
790: return(0);
791: }
795: /*@
796: SectionRealGetFibration - Creates a section for only the data associated with the given field
798: Collective on DMMesh
800: Input Parameter:
801: + section - the Section
802: - field- The field number
804: Output Parameter:
805: . subsection - the section of the given field
807: Level: advanced
809: .seealso SectionRealCreate()
810: @*/
811: PetscErrorCode SectionRealGetFibration(SectionReal section, const PetscInt field, SectionReal *subsection)
812: {
813: ALE::Obj<PETSC_MESH_TYPE> b;
814: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
815: ALE::Obj<PETSC_MESH_TYPE::real_section_type> t;
816: MPI_Comm comm;
820: PetscObjectGetComm((PetscObject) section, &comm);
821: SectionRealGetBundle(section, b);
822: SectionRealGetSection(section, s);
823: SectionRealCreate(comm, subsection);
824: SectionRealSetBundle(*subsection, b);
825: t = s->getFibration(field);
826: SectionRealSetSection(*subsection, t);
827: return(0);
828: }
832: /*@C
833: SectionRealToVec - Maps the given section to a Vec
835: Collective on Section
837: Input Parameters:
838: + section - the real Section
839: - mesh - The DMMesh
841: Output Parameter:
842: . vec - the Vec
844: Level: intermediate
846: .seealso VecCreate(), SectionRealCreate()
847: @*/
848: PetscErrorCode SectionRealToVec(SectionReal section, DM dm, ScatterMode mode, Vec vec)
849: {
850: Vec localVec;
851: VecScatter scatter;
856: SectionRealCreateLocalVector(section, &localVec);
857: DMMeshGetGlobalScatter(dm, &scatter);
858: if (mode == SCATTER_FORWARD) {
859: VecScatterBegin(scatter, localVec, vec, INSERT_VALUES, mode);
860: VecScatterEnd(scatter, localVec, vec, INSERT_VALUES, mode);
861: } else {
862: VecScatterBegin(scatter, vec, localVec, INSERT_VALUES, mode);
863: VecScatterEnd(scatter, vec, localVec, INSERT_VALUES, mode);
864: }
865: VecDestroy(&localVec);
866: return(0);
867: }
871: /*@
872: SectionRealToVec - Map between unassembled local Section storage and a globally assembled Vec
874: Collective on VecScatter
876: Input Parameters:
877: + section - the Section
878: . scatter - the scatter from the Section to the Vec
879: . mode - the mode, SCATTER_FORWARD (Section to Vec) or SCATTER_REVERSE (Vec to Section)
880: - vec - the Vec
882: Level: advanced
884: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
885: @*/
886: PetscErrorCode SectionRealToVec(SectionReal section, VecScatter scatter, ScatterMode mode, Vec vec)
887: {
888: Vec localVec;
893: SectionRealCreateLocalVector(section, &localVec);
894: if (mode == SCATTER_FORWARD) {
895: VecScatterBegin(scatter, localVec, vec, INSERT_VALUES, mode);
896: VecScatterEnd(scatter, localVec, vec, INSERT_VALUES, mode);
897: } else {
898: VecScatterBegin(scatter, vec, localVec, INSERT_VALUES, mode);
899: VecScatterEnd(scatter, vec, localVec, INSERT_VALUES, mode);
900: }
901: VecDestroy(&localVec);
902: return(0);
903: }
907: /*@
908: SectionRealClear - Dellocate storage for this section
910: Not collective
912: Input Parameter:
913: . section - the section object
915: Level: advanced
917: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
918: @*/
919: PetscErrorCode SectionRealClear(SectionReal section)
920: {
923: section->s->clear();
924: return(0);
925: }
929: /*@
930: SectionRealSet - Sets all the values to the given value
932: Not collective
934: Input Parameters:
935: + section - the real Section
936: - val - the value
938: Level: intermediate
940: .seealso VecNorm(), SectionRealCreate()
941: @*/
942: PetscErrorCode SectionRealSet(SectionReal section, PetscReal val)
943: {
944: Obj<PETSC_MESH_TYPE::real_section_type> s;
948: SectionRealGetSection(section, s);
949: s->set(val);
950: return(0);
951: }
955: /*@C
956: SectionRealNorm - Computes the vector norm.
958: Collective on Section
960: Input Parameters:
961: + section - the real Section
962: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
963: NORM_1_AND_2, which computes both norms and stores them
964: in a two element array.
966: Output Parameter:
967: . val - the norm
969: Notes:
970: $ NORM_1 denotes sum_i |x_i|
971: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
972: $ NORM_INFINITY denotes max_i |x_i|
974: Level: intermediate
976: .seealso VecNorm(), SectionRealCreate()
977: @*/
978: PetscErrorCode SectionRealNorm(SectionReal section, DM dm, NormType type, PetscReal *val)
979: {
980: Obj<PETSC_MESH_TYPE> m;
981: Obj<PETSC_MESH_TYPE::real_section_type> s;
982: Vec v;
987: DMMeshGetMesh(dm, m);
988: SectionRealGetSection(section, s);
989: const ALE::Obj<PETSC_MESH_TYPE::order_type>& order = m->getFactory()->getGlobalOrder(m, s->getName(), s);
990: VecCreate(m->comm(), &v);
991: VecSetSizes(v, order->getLocalSize(), order->getGlobalSize());
992: VecSetFromOptions(v);
993: SectionRealToVec(section, dm, SCATTER_FORWARD, v);
994: VecNorm(v, type, val);
995: VecDestroy(&v);
996: return(0);
997: }
1001: /*@
1002: SectionRealAXPY -
1004: Collective on Section
1006: Input Parameters:
1007: + section - the real Section
1008: . alpha - a scalar
1009: - X - the other real Section
1011: Output Parameter:
1012: . section - the difference
1014: Level: intermediate
1016: .seealso VecNorm(), SectionRealCreate()
1017: @*/
1018: PetscErrorCode SectionRealAXPY(SectionReal section, DM dm, PetscScalar alpha, SectionReal X)
1019: {
1020: Obj<PETSC_MESH_TYPE> m;
1021: Obj<PETSC_MESH_TYPE::real_section_type> s;
1022: Obj<PETSC_MESH_TYPE::real_section_type> sX;
1023: Vec v, x;
1028: DMMeshGetMesh(dm, m);
1029: SectionRealGetSection(section, s);
1030: SectionRealGetSection(X, sX);
1031: const ALE::Obj<PETSC_MESH_TYPE::order_type>& order = m->getFactory()->getGlobalOrder(m, s->getName(), s);
1032: VecCreate(m->comm(), &v);
1033: VecSetSizes(v, order->getLocalSize(), order->getGlobalSize());
1034: VecSetFromOptions(v);
1035: VecDuplicate(v, &x);
1036: SectionRealToVec(section, dm, SCATTER_FORWARD, v);
1037: SectionRealToVec(X, dm, SCATTER_FORWARD, x);
1038: VecAXPY(v, alpha, x);
1039: SectionRealToVec(section, dm, SCATTER_REVERSE, v);
1040: VecDestroy(&v);
1041: VecDestroy(&x);
1042: return(0);
1043: }
1047: /*@C
1048: DMMeshGetVertexSectionReal - Create a Section over the vertices with the specified fiber dimension
1050: Collective on DMMesh
1052: Input Parameters:
1053: + mesh - The DMMesh object
1054: - fiberDim - The number of degrees of freedom per vertex
1056: Output Parameter:
1057: . section - The section
1059: Level: intermediate
1061: .keywords: mesh, section, vertex
1062: .seealso: DMMeshCreate(), SectionRealCreate()
1063: @*/
1064: PetscErrorCode DMMeshGetVertexSectionReal(DM dm, const char name[], PetscInt fiberDim, SectionReal *section)
1065: {
1066: ALE::Obj<PETSC_MESH_TYPE> m;
1067: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
1068: PetscErrorCode ierr;
1071: DMMeshGetMesh(dm, m);
1072: SectionRealCreate(m->comm(), section);
1073: PetscObjectSetName((PetscObject) *section, name);
1074: SectionRealSetBundle(*section, m);
1075: SectionRealGetSection(*section, s);
1076: s->setChart(m->getSieve()->getChart());
1077: s->setName(name);
1078: s->setFiberDimension(m->depthStratum(0), fiberDim);
1079: m->allocate(s);
1080: return(0);
1081: }
1085: /*@C
1086: DMMeshGetCellSectionReal - Create a Section over the cells with the specified fiber dimension
1088: Collective on DMMesh
1090: Input Parameters:
1091: + mesh - The DMMesh object
1092: - fiberDim - The section name
1094: Output Parameter:
1095: . section - The section
1097: Level: intermediate
1099: .keywords: mesh, section, cell
1100: .seealso: DMMeshCreate(), SectionRealCreate()
1101: @*/
1102: PetscErrorCode DMMeshGetCellSectionReal(DM dm, const char name[], PetscInt fiberDim, SectionReal *section)
1103: {
1104: ALE::Obj<PETSC_MESH_TYPE> m;
1105: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
1106: PetscErrorCode ierr;
1109: DMMeshGetMesh(dm, m);
1110: SectionRealCreate(m->comm(), section);
1111: PetscObjectSetName((PetscObject) *section, name);
1112: SectionRealSetBundle(*section, m);
1113: SectionRealGetSection(*section, s);
1114: s->setChart(m->getSieve()->getChart());
1115: s->setName(name);
1116: s->setFiberDimension(m->heightStratum(0), fiberDim);
1117: m->allocate(s);
1118: return(0);
1119: }
1123: /*@C
1124: DMMeshCreateGlobalRealVector - Creates a vector of the correct size to be gathered into by the mesh.
1126: Collective on DMMesh
1128: Input Parameters:
1129: + mesh - the mesh object
1130: - section - The SectionReal
1132: Output Parameters:
1133: . gvec - the global vector
1135: Level: advanced
1137: .seealso DMMeshDestroy(), DMMeshCreate(), DMMeshCreateGlobalRealVector()
1138: @*/
1139: PetscErrorCode DMMeshCreateGlobalRealVector(DM dm, SectionReal section, Vec *gvec)
1140: {
1141: ALE::Obj<PETSC_MESH_TYPE> m;
1142: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
1143: const char *name;
1147: DMMeshGetMesh(dm, m);
1148: SectionRealGetSection(section, s);
1149: PetscObjectGetName((PetscObject) section, &name);
1150: const ALE::Obj<PETSC_MESH_TYPE::order_type>& order = m->getFactory()->getGlobalOrder(m, name, s);
1152: VecCreate(m->comm(), gvec);
1153: VecSetSizes(*gvec, order->getLocalSize(), order->getGlobalSize());
1154: VecSetFromOptions(*gvec);
1155: return(0);
1156: }
1160: PetscErrorCode SectionIntView_Sieve(SectionInt section, PetscViewer viewer)
1161: {
1162: PetscBool iascii, isbinary, isdraw;
1166: PetscTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
1167: PetscTypeCompare((PetscObject) viewer, PETSCVIEWERBINARY, &isbinary);
1168: PetscTypeCompare((PetscObject) viewer, PETSCVIEWERDRAW, &isdraw);
1170: if (iascii){
1171: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1172: ALE::Obj<PETSC_MESH_TYPE> b;
1173: const char *name;
1175: SectionIntGetSection(section, s);
1176: SectionIntGetBundle(section, b);
1177: PetscObjectGetName((PetscObject) section, &name);
1178: SectionView_Sieve_Ascii(b, s, name, viewer);
1179: } else if (isbinary) {
1180: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "Binary viewer not implemented for Section");
1181: } else if (isdraw){
1182: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "Draw viewer not implemented for Section");
1183: } else {
1184: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Viewer type %s not supported by this section object", ((PetscObject)viewer)->type_name);
1185: }
1186: return(0);
1187: }
1191: /*@C
1192: SectionIntView - Views a Section object.
1194: Collective on Section
1196: Input Parameters:
1197: + section - the Section
1198: - viewer - an optional visualization context
1200: Notes:
1201: The available visualization contexts include
1202: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
1203: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
1204: output where only the first processor opens
1205: the file. All other processors send their
1206: data to the first processor to print.
1208: You can change the format the section is printed using the
1209: option PetscViewerSetFormat().
1211: The user can open alternative visualization contexts with
1212: + PetscViewerASCIIOpen() - Outputs section to a specified file
1213: . PetscViewerBinaryOpen() - Outputs section in binary to a
1214: specified file; corresponding input uses SectionLoad()
1215: . PetscViewerDrawOpen() - Outputs section to an X window display
1217: The user can call PetscViewerSetFormat() to specify the output
1218: format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
1219: PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen). Available formats include
1220: + PETSC_VIEWER_DEFAULT - default, prints section information
1221: - PETSC_VIEWER_ASCII_VTK - outputs a VTK file describing the section
1223: Level: beginner
1225: Concepts: section^printing
1226: Concepts: section^saving to disk
1228: .seealso: VecView(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerBinaryOpen(), PetscViewerCreate()
1229: @*/
1230: PetscErrorCode SectionIntView(SectionInt section, PetscViewer viewer)
1231: {
1237: if (!viewer) {
1238: PetscViewerASCIIGetStdout(((PetscObject)section)->comm,&viewer);
1239: }
1243: PetscLogEventBegin(SectionInt_View,0,0,0,0);
1244: (*section->ops->view)(section, viewer);
1245: PetscLogEventEnd(SectionInt_View,0,0,0,0);
1246: return(0);
1247: }
1251: /*@C
1252: SectionIntGetSection - Gets the internal section object
1254: Not collective
1256: Input Parameter:
1257: . section - the section object
1259: Output Parameter:
1260: . s - the internal section object
1261:
1262: Level: advanced
1264: .seealso SectionIntCreate(), SectionIntSetSection()
1265: @*/
1266: PetscErrorCode SectionIntGetSection(SectionInt section, ALE::Obj<PETSC_MESH_TYPE::int_section_type>& s)
1267: {
1270: s = section->s;
1271: return(0);
1272: }
1276: /*@C
1277: SectionIntSetSection - Sets the internal section object
1279: Not collective
1281: Input Parameters:
1282: + section - the section object
1283: - s - the internal section object
1284:
1285: Level: advanced
1287: .seealso SectionIntCreate(), SectionIntGetSection()
1288: @*/
1289: PetscErrorCode SectionIntSetSection(SectionInt section, const ALE::Obj<PETSC_MESH_TYPE::int_section_type>& s)
1290: {
1295: if (!s.isNull()) {PetscObjectSetName((PetscObject) section, s->getName().c_str());}
1296: section->s = s;
1297: return(0);
1298: }
1302: /*@C
1303: SectionIntGetBundle - Gets the section bundle
1305: Not collective
1307: Input Parameter:
1308: . section - the section object
1310: Output Parameter:
1311: . b - the section bundle
1312:
1313: Level: advanced
1315: .seealso SectionIntCreate(), SectionIntGetSection(), SectionIntSetSection()
1316: @*/
1317: PetscErrorCode SectionIntGetBundle(SectionInt section, ALE::Obj<PETSC_MESH_TYPE>& b)
1318: {
1321: b = section->b;
1322: return(0);
1323: }
1327: /*@C
1328: SectionIntSetBundle - Sets the section bundle
1330: Not collective
1332: Input Parameters:
1333: + section - the section object
1334: - b - the section bundle
1335:
1336: Level: advanced
1338: .seealso SectionIntCreate(), SectionIntGetSection(), SectionIntSetSection()
1339: @*/
1340: PetscErrorCode SectionIntSetBundle(SectionInt section, const ALE::Obj<PETSC_MESH_TYPE>& b)
1341: {
1344: section->b = b;
1345: return(0);
1346: }
1350: /*@C
1351: SectionIntCreate - Creates a Section object, used to manage data for an unstructured problem
1352: described by a Sieve.
1354: Collective on MPI_Comm
1356: Input Parameter:
1357: . comm - the processors that will share the global section
1359: Output Parameters:
1360: . section - the section object
1362: Level: advanced
1364: .seealso SectionIntDestroy(), SectionIntView()
1365: @*/
1366: PetscErrorCode SectionIntCreate(MPI_Comm comm, SectionInt *section)
1367: {
1369: SectionInt s;
1373: *section = PETSC_NULL;
1375: PetscHeaderCreate(s,_p_SectionInt,struct _SectionIntOps,SECTIONINT_CLASSID,0,"SectionInt","Section","DM",comm,SectionIntDestroy,0);
1376: s->ops->view = SectionIntView_Sieve;
1377: s->ops->restrictClosure = SectionIntRestrict;
1378: s->ops->update = SectionIntUpdate;
1380: PetscObjectChangeTypeName((PetscObject) s, "sieve");
1382: new(&s->s) ALE::Obj<PETSC_MESH_TYPE::int_section_type>(PETSC_MESH_TYPE::int_section_type(comm));
1383: new(&s->b) ALE::Obj<PETSC_MESH_TYPE>(PETSC_NULL);
1384: *section = s;
1385: return(0);
1386: }
1390: /*@
1391: SectionIntDestroy - Destroys a section.
1393: Collective on Section
1395: Input Parameter:
1396: . section - the section object
1398: Level: advanced
1400: .seealso SectionIntCreate(), SectionIntView()
1401: @*/
1402: PetscErrorCode SectionIntDestroy(SectionInt *section)
1403: {
1407: if (!*section) return(0);
1409: if (--((PetscObject)(*section))->refct > 0) return(0);
1410: PetscHeaderDestroy(section);
1411: return(0);
1412: }
1416: /*@
1417: SectionIntDistribute - Distributes the sections.
1419: Not Collective
1421: Input Parameters:
1422: + serialSection - The original Section object
1423: - parallelMesh - The parallel DMMesh
1425: Output Parameter:
1426: . parallelSection - The distributed Section object
1428: Level: intermediate
1430: .keywords: mesh, section, distribute
1431: .seealso: DMMeshCreate()
1432: @*/
1433: PetscErrorCode SectionIntDistribute(SectionInt serialSection, DM parallelMesh, SectionInt *parallelSection)
1434: {
1435: ALE::Obj<PETSC_MESH_TYPE::int_section_type> oldSection;
1436: ALE::Obj<PETSC_MESH_TYPE> m;
1437: PetscErrorCode ierr;
1440: SectionIntGetSection(serialSection, oldSection);
1441: DMMeshGetMesh(parallelMesh, m);
1442: SectionIntCreate(oldSection->comm(), parallelSection);
1443: #ifdef PETSC_OPT_SIEVE
1444: ALE::Obj<PETSC_MESH_TYPE::int_section_type> newSection;
1446: // We assume all integer sections are complete sections
1447: newSection->setName(oldSection->getName());
1448: newSection->setChart(m->getSieve()->getChart());
1449: //distributeSection(oldSection, partition, m->getRenumbering(), m->getDistSendOverlap(), m->getDistRecvOverlap(), newSection);
1450: SectionIntSetSection(*parallelSection, newSection);
1451: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP, "Not working because the partition is unavailable");
1452: #else
1453: ALE::Obj<PETSC_MESH_TYPE::int_section_type> newSection = ALE::Distribution<PETSC_MESH_TYPE>::distributeSection(oldSection, m, m->getDistSendOverlap(), m->getDistRecvOverlap());
1454: SectionIntSetSection(*parallelSection, newSection);
1455: #endif
1456: return(0);
1457: }
1461: /*@C
1462: SectionIntRestrict - Restricts the SectionInt to a subset of the topology, returning an array of values.
1464: Not collective
1466: Input Parameters:
1467: + section - the section object
1468: - point - the Sieve point
1470: Output Parameter:
1471: . values - The values associated with the submesh
1473: Level: advanced
1475: .seealso SectionIntUpdate(), SectionIntCreate(), SectionIntView()
1476: @*/
1477: PetscErrorCode SectionIntRestrict(SectionInt section, PetscInt point, PetscInt *values[])
1478: {
1482: *values = (PetscInt *) section->b->restrictClosure(section->s, point);
1483: return(0);
1484: }
1488: /*@C
1489: SectionIntUpdate - Updates the array of values associated to a subset of the topology in this Section.
1491: Not collective
1493: Input Parameters:
1494: + section - the section object
1495: . point - the Sieve point
1496: . values - The values associated with the submesh
1497: - mode - The insertion mode
1499: Level: advanced
1501: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1502: @*/
1503: PetscErrorCode SectionIntUpdate(SectionInt section, PetscInt point, const PetscInt values[], InsertMode mode)
1504: {
1508: if (mode == INSERT_VALUES) {
1509: section->b->update(section->s, point, values);
1510: } else if (mode == ADD_VALUES) {
1511: section->b->updateAdd(section->s, point, values);
1512: } else {
1513: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG, "Invalid insertion mode: %d", mode);
1514: }
1515: return(0);
1516: }
1520: /*@C
1521: SectionIntRestrictClosure - Returns an array with the values in a given closure
1523: Not Collective
1525: Input Parameters:
1526: + section - The section
1527: . mesh - The DMMesh object
1528: . point - The sieve point
1529: . n - The array size
1530: - array - The array to fill up
1532: Output Parameter:
1533: . array - The array full of values in the closure
1535: Level: intermediate
1537: .keywords: mesh, elements
1538: .seealso: DMMeshCreate()
1539: @*/
1540: PetscErrorCode SectionIntRestrictClosure(SectionInt section, DM dm, PetscInt point, PetscInt n, PetscInt values[])
1541: {
1542: ALE::Obj<PETSC_MESH_TYPE> m;
1543: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1544: PetscErrorCode ierr;
1547: DMMeshGetMesh(dm, m);
1548: SectionIntGetSection(section, s);
1549: m->restrictClosure(s, point, values, n);
1550: return(0);
1551: }
1555: /*@C
1556: SectionIntUpdateClosure - Updates the values in a given closure from the array
1558: Not Collective
1560: Input Parameters:
1561: + section - The section
1562: . mesh - The DMMesh object
1563: . point - The sieve point
1564: . array - The array to fill up
1565: - mode - The insertion mode
1567: Output Parameter:
1568: . array - The array full of values in the closure
1570: Level: intermediate
1572: .keywords: mesh, elements
1573: .seealso: DMMeshCreate()
1574: @*/
1575: PetscErrorCode SectionIntUpdateClosure(SectionInt section, DM dm, PetscInt point, PetscInt values[], InsertMode mode)
1576: {
1577: ALE::Obj<PETSC_MESH_TYPE> m;
1578: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1579: PetscErrorCode ierr;
1582: DMMeshGetMesh(dm, m);
1583: SectionIntGetSection(section, s);
1584: if (mode == INSERT_VALUES) {
1585: m->update(s, point, values);
1586: } else if (mode == ADD_VALUES) {
1587: m->updateAdd(s, point, values);
1588: } else {
1589: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG, "Invalid insertion mode: %d", mode);
1590: }
1591: return(0);
1592: }
1596: /*@
1597: SectionIntComplete - Exchanges data across the mesh overlap.
1599: Not collective
1601: Input Parameter:
1602: . section - the section object
1604: Level: advanced
1606: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1607: @*/
1608: PetscErrorCode SectionIntComplete(SectionInt section)
1609: {
1610: Obj<PETSC_MESH_TYPE::int_section_type> s;
1611: Obj<PETSC_MESH_TYPE> b;
1615: SectionIntGetSection(section, s);
1616: SectionIntGetBundle(section, b);
1617: #if 0
1618: ALE::Distribution<PETSC_MESH_TYPE>::completeSection(b, s);
1619: #else
1620: ALE::Completion::completeSectionAdd(b->getSendOverlap(), b->getRecvOverlap(), s, s);
1621: #endif
1622: return(0);
1623: }
1627: /*@
1628: SectionIntZero - Zero out the entries
1630: Not collective
1632: Input Parameter:
1633: . section - the section object
1635: Level: advanced
1637: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1638: @*/
1639: PetscErrorCode SectionIntZero(SectionInt section)
1640: {
1641: Obj<PETSC_MESH_TYPE::int_section_type> s;
1645: SectionIntGetSection(section, s);
1646: s->zero();
1647: return(0);
1648: }
1652: /*@
1653: SectionIntGetFiberDimension - Get the size of the vector space attached to the point
1655: Not collective
1657: Input Parameters:
1658: + section - the section object
1659: - point - the Sieve point
1661: Output Parameters:
1662: . size - The fiber dimension
1664: Level: advanced
1666: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
1667: @*/
1668: PetscErrorCode SectionIntGetFiberDimension(SectionInt section, PetscInt point, PetscInt *size)
1669: {
1672: *size = section->s->getFiberDimension(point);
1673: return(0);
1674: }
1678: /*@
1679: SectionIntSetFiberDimension - Set the size of the vector space attached to the point
1681: Not collective
1683: Input Parameters:
1684: + section - the section object
1685: . point - the Sieve point
1686: - size - The fiber dimension
1688: Level: advanced
1690: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1691: @*/
1692: PetscErrorCode SectionIntSetFiberDimension(SectionInt section, PetscInt point, const PetscInt size)
1693: {
1696: section->s->setFiberDimension(point, size);
1697: return(0);
1698: }
1702: /*@
1703: SectionIntSetFiberDimensionField - Set the size of the vector space attached to the point for a given field
1705: Not collective
1707: Input Parameters:
1708: + section - the section object
1709: . point - the Sieve point
1710: . size - The fiber dimension
1711: - field - The field number
1713: Level: advanced
1715: .seealso SectionIntSetFiberDimension(), SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1716: @*/
1717: PetscErrorCode SectionIntSetFiberDimensionField(SectionInt section, PetscInt point, const PetscInt size, const PetscInt field)
1718: {
1721: section->s->setFiberDimension(point, size, field);
1722: return(0);
1723: }
1727: /*@
1728: SectionIntGetSize - Gets the number of local dofs in this Section
1730: Not collective
1732: Input Parameter:
1733: . section - the section object
1735: Output Parameter:
1736: . size - the section size
1738: Level: advanced
1740: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1741: @*/
1742: PetscErrorCode SectionIntGetSize(SectionInt section, PetscInt *size)
1743: {
1744: Obj<PETSC_MESH_TYPE::int_section_type> s;
1749: SectionIntGetSection(section, s);
1750: *size = s->size();
1751: return(0);
1752: }
1756: /*@
1757: SectionIntAllocate - Allocate storage for this section
1759: Not collective
1761: Input Parameter:
1762: . section - the section object
1764: Level: advanced
1766: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1767: @*/
1768: PetscErrorCode SectionIntAllocate(SectionInt section)
1769: {
1772: section->b->allocate(section->s);
1773: return(0);
1774: }
1778: /*@C
1779: SectionIntClear - Dellocate storage for this section
1781: Not collective
1783: Input Parameter:
1784: . section - the section object
1786: Level: advanced
1788: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1789: @*/
1790: PetscErrorCode SectionIntClear(SectionInt section)
1791: {
1794: section->s->clear();
1795: return(0);
1796: }
1800: /*@
1801: SectionIntSet - Sets all the values to the given value
1803: Not collective
1805: Input Parameters:
1806: + section - the real Section
1807: - val - the value
1809: Level: intermediate
1811: .seealso VecNorm(), SectionIntCreate()
1812: @*/
1813: PetscErrorCode SectionIntSet(SectionInt section, PetscInt val)
1814: {
1815: Obj<PETSC_MESH_TYPE::int_section_type> s;
1819: SectionIntGetSection(section, s);
1820: s->set(val);
1821: return(0);
1822: }
1826: /*@
1827: SectionIntAddSpace - Add another field to this section
1829: Collective on DMMesh
1831: Input Parameter:
1832: . section - the Section
1834: Level: advanced
1836: .seealso SectionIntCreate(), SectionIntGetFibration()
1837: @*/
1838: PetscErrorCode SectionIntAddSpace(SectionInt section)
1839: {
1840: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1844: SectionIntGetSection(section, s);
1845: s->addSpace();
1846: return(0);
1847: }
1851: /*@
1852: SectionIntGetFibration - Creates a section for only the data associated with the given field
1854: Collective on DMMesh
1856: Input Parameter:
1857: + section - the Section
1858: - field- The field number
1860: Output Parameter:
1861: . subsection - the section of the given field
1863: Level: advanced
1865: .seealso SectionIntCreate(), SectionIntAddSpace()
1866: @*/
1867: PetscErrorCode SectionIntGetFibration(SectionInt section, const PetscInt field, SectionInt *subsection)
1868: {
1869: ALE::Obj<PETSC_MESH_TYPE> b;
1870: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1871: ALE::Obj<PETSC_MESH_TYPE::int_section_type> t;
1872: MPI_Comm comm;
1876: PetscObjectGetComm((PetscObject) section, &comm);
1877: SectionIntGetBundle(section, b);
1878: SectionIntGetSection(section, s);
1879: SectionIntCreate(comm, subsection);
1880: SectionIntSetBundle(*subsection, b);
1881: t = s->getFibration(field);
1882: SectionIntSetSection(*subsection, t);
1883: return(0);
1884: }
1888: /*@C
1889: DMMeshGetVertexSectionInt - Create a Section over the vertices with the specified fiber dimension
1891: Collective on DMMesh
1893: Input Parameters:
1894: + mesh - The DMMesh object
1895: - fiberDim - The section name
1897: Output Parameter:
1898: . section - The section
1900: Level: intermediate
1902: .keywords: mesh, section, vertex
1903: .seealso: DMMeshCreate(), SectionIntCreate()
1904: @*/
1905: PetscErrorCode DMMeshGetVertexSectionInt(DM dm, const char name[], PetscInt fiberDim, SectionInt *section)
1906: {
1907: ALE::Obj<PETSC_MESH_TYPE> m;
1908: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1909: PetscErrorCode ierr;
1912: DMMeshGetMesh(dm, m);
1913: SectionIntCreate(m->comm(), section);
1914: PetscObjectSetName((PetscObject) *section, name);
1915: SectionIntSetBundle(*section, m);
1916: SectionIntGetSection(*section, s);
1917: #ifdef PETSC_OPT_SIEVE
1918: s->setChart(m->getSieve()->getChart());
1919: #endif
1920: s->setName(name);
1921: s->setFiberDimension(m->depthStratum(0), fiberDim);
1922: m->allocate(s);
1923: return(0);
1924: }
1928: /*@C
1929: DMMeshGetCellSectionInt - Create a Section over the cells with the specified fiber dimension
1931: Collective on DMMesh
1933: Input Parameters:
1934: + mesh - The DMMesh object
1935: - fiberDim - The section name
1937: Output Parameter:
1938: . section - The section
1940: Level: intermediate
1942: .keywords: mesh, section, cell
1943: .seealso: DMMeshCreate(), SectionIntCreate()
1944: @*/
1945: PetscErrorCode DMMeshGetCellSectionInt(DM dm, const char name[], PetscInt fiberDim, SectionInt *section)
1946: {
1947: ALE::Obj<PETSC_MESH_TYPE> m;
1948: ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1949: PetscErrorCode ierr;
1952: DMMeshGetMesh(dm, m);
1953: SectionIntCreate(m->comm(), section);
1954: PetscObjectSetName((PetscObject) *section, name);
1955: SectionIntSetBundle(*section, m);
1956: SectionIntGetSection(*section, s);
1957: #ifdef PETSC_OPT_SIEVE
1958: s->setChart(m->getSieve()->getChart());
1959: #endif
1960: s->setName(name);
1961: s->setFiberDimension(m->heightStratum(0), fiberDim);
1962: m->allocate(s);
1963: return(0);
1964: }
1968: /*@C
1969: DMMeshSetupSection - Layout Section based upon discretization and boundary condition information in the Mesh
1971: Collective on DMMesh
1973: Input Parameters:
1974: . mesh - The DMMesh object
1976: Output Parameter:
1977: . section - The section
1979: Level: intermediate
1981: .keywords: mesh, section, cell
1982: .seealso: DMMeshCreate(), SectionRealCreate()
1983: @*/
1984: PetscErrorCode DMMeshSetupSection(DM dm, SectionReal section)
1985: {
1986: ALE::Obj<PETSC_MESH_TYPE> m;
1987: ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
1988: PetscErrorCode ierr;
1991: DMMeshGetMesh(dm, m);
1992: SectionRealGetSection(section, s);
1993: try {
1994: m->setupField(s);
1995: } catch(ALE::Exception e) {
1996: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "%s", e.message());
1997: }
1998: return(0);
1999: }