Actual source code: section.c

  1: #include "private/meshimpl.h"   /*I      "petscmesh.h"   I*/
  2: #include <petscmesh_viewers.hh>

  4: /* Logging support */
  5: PetscCookie  SECTIONREAL_COOKIE;
  6: PetscLogEvent  SectionReal_View;
  7: PetscCookie  SECTIONINT_COOKIE;
  8: PetscLogEvent  SectionInt_View;
  9: PetscCookie  SECTIONPAIR_COOKIE;
 10: PetscLogEvent  SectionPair_View;

 14: PetscErrorCode SectionRealView_Sieve(SectionReal section, PetscViewer viewer)
 15: {
 16:   PetscTruth     iascii, isbinary, isdraw;

 20:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &iascii);
 21:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &isbinary);
 22:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_DRAW, &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_ERR_SUP, "Binary viewer not implemented for Section");
 35:   } else if (isdraw){
 36:     SETERRQ(PETSC_ERR_SUP, "Draw viewer not implemented for Section");
 37:   } else {
 38:     SETERRQ1(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_COOKIE,0,"SectionReal",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: /*@C
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: {

298:   if (--((PetscObject)section)->refct > 0) return(0);
299:   section->s = PETSC_NULL;
300:   PetscHeaderDestroy(section);
301:   return(0);
302: }

306: /*@C
307:   SectionRealDistribute - Distributes the sections.

309:   Not Collective

311:   Input Parameters:
312: + serialSection - The original Section object
313: - parallelMesh - The parallel Mesh

315:   Output Parameter:
316: . parallelSection - The distributed Section object

318:   Level: intermediate

320: .keywords: mesh, section, distribute
321: .seealso: MeshCreate()
322: @*/
323: PetscErrorCode SectionRealDistribute(SectionReal serialSection, Mesh 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:   MeshGetMesh(parallelMesh, m);
332:   SectionRealCreate(oldSection->comm(), parallelSection);
333: #ifdef PETSC_OPT_SIEVE
334:   SETERRQ(PETSC_ERR_SUP, "I am being lazy, bug me.");
335: #else
336:   ALE::Obj<PETSC_MESH_TYPE::real_section_type> newSection = ALE::Distribution<PETSC_MESH_TYPE>::distributeSection(oldSection, m, m->getDistSendOverlap(), m->getDistRecvOverlap());
337:   SectionRealSetSection(*parallelSection, newSection);
338: #endif
339:   return(0);
340: }

344: /*@C
345:   SectionRealRestrict - Restricts the SectionReal to a subset of the topology, returning an array of values.

347:   Not collective

349:   Input Parameters:
350: + section - the section object
351: - point - the Sieve point

353:   Output Parameter:
354: . values - The values associated with the submesh

356:   Level: advanced

358: .seealso SectionUpdate(), SectionCreate(), SectionView()
359: @*/
360: PetscErrorCode  SectionRealRestrict(SectionReal section, PetscInt point, PetscScalar *values[])
361: {
365:   *values = (PetscScalar *) section->b->restrictClosure(section->s, point);
366:   return(0);
367: }

371: /*@C
372:   SectionRealUpdate - Updates the array of values associated to a subset of the topology in this Section.

374:   Not collective

376:   Input Parameters:
377: + section - the section object
378: . point - the Sieve point
379: - values - The values associated with the submesh

381:   Level: advanced

383: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
384: @*/
385: PetscErrorCode  SectionRealUpdate(SectionReal section, PetscInt point, const PetscScalar values[])
386: {
390:   section->b->update(section->s, point, values);
391:   return(0);
392: }

396: /*@C
397:   SectionRealUpdateAdd - Updates the array of values associated to a subset of the topology in this Section.

399:   Not collective

401:   Input Parameters:
402: + section - the section object
403: . point - the Sieve point
404: - values - The values associated with the submesh

406:   Level: advanced

408: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
409: @*/
410: PetscErrorCode  SectionRealUpdateAdd(SectionReal section, PetscInt point, const PetscScalar values[])
411: {
415:   section->b->updateAdd(section->s, point, values);
416:   return(0);
417: }

421: /*@C
422:   SectionRealComplete - Exchanges data across the mesh overlap.

424:   Not collective

426:   Input Parameter:
427: . section - the section object

429:   Level: advanced

431: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
432: @*/
433: PetscErrorCode SectionRealComplete(SectionReal section)
434: {
435:   Obj<PETSC_MESH_TYPE::real_section_type> s;
436:   Obj<PETSC_MESH_TYPE>                    b;

440:   SectionRealGetSection(section, s);
441:   SectionRealGetBundle(section, b);
442:   ALE::Distribution<PETSC_MESH_TYPE>::completeSection(b, s);
443:   return(0);
444: }

448: /*@C
449:   SectionRealZero - Zero out the entries

451:   Not collective

453:   Input Parameter:
454: . section - the section object

456:   Level: advanced

458: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
459: @*/
460: PetscErrorCode SectionRealZero(SectionReal section)
461: {
462:   Obj<PETSC_MESH_TYPE::real_section_type> s;

466:   SectionRealGetSection(section, s);
467:   s->zero();
468:   return(0);
469: }

473: /*@C
474:   SectionRealSetFiberDimension - Set the size of the vector space attached to the point

476:   Not collective

478:   Input Parameters:
479: + section - the section object
480: . point - the Sieve point
481: - size - The fiber dimension

483:   Level: advanced

485: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
486: @*/
487: PetscErrorCode  SectionRealSetFiberDimension(SectionReal section, PetscInt point, const PetscInt size)
488: {
491:   section->s->setFiberDimension(point, size);
492:   return(0);
493: }

497: /*@C
498:   SectionRealAllocate - Allocate storage for this section

500:   Not collective

502:   Input Parameter:
503: . section - the section object

505:   Level: advanced

507: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
508: @*/
509: PetscErrorCode  SectionRealAllocate(SectionReal section)
510: {
513:   section->b->allocate(section->s);
514:   return(0);
515: }

519: /*@C
520:   SectionRealCreateLocalVector - Creates a vector with the local piece of the Section

522:   Collective on Mesh

524:   Input Parameter:
525: . section - the Section  

527:   Output Parameter:
528: . localVec - the local vector

530:   Level: advanced

532: .seealso MeshDestroy(), MeshCreate()
533: @*/
534: PetscErrorCode  SectionRealCreateLocalVector(SectionReal section, Vec *localVec)
535: {
536:   ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;

540:   SectionRealGetSection(section, s);
541:   VecCreateSeqWithArray(PETSC_COMM_SELF, s->getStorageSize(), s->restrictSpace(), localVec);
542:   return(0);
543: }

547: /*@C
548:   SectionRealToVec - Maps the given section to a Vec

550:   Collective on Section

552:   Input Parameters:
553: + section - the real Section
554: - mesh - The Mesh

556:   Output Parameter:
557: . vec - the Vec 

559:   Level: intermediate

561: .seealso VecCreate(), SectionRealCreate()
562: @*/
563: PetscErrorCode  SectionRealToVec(SectionReal section, Mesh mesh, ScatterMode mode, Vec vec)
564: {
565:   Vec            localVec;
566:   VecScatter     scatter;

571:   SectionRealCreateLocalVector(section, &localVec);
572:   MeshGetGlobalScatter(mesh, &scatter);
573:   if (mode == SCATTER_FORWARD) {
574:     VecScatterBegin(scatter, localVec, vec, INSERT_VALUES, mode);
575:     VecScatterEnd(scatter, localVec, vec, INSERT_VALUES, mode);
576:   } else {
577:     VecScatterBegin(scatter, vec, localVec, INSERT_VALUES, mode);
578:     VecScatterEnd(scatter, vec, localVec, INSERT_VALUES, mode);
579:   }
580:   VecDestroy(localVec);
581:   return(0);
582: }

584: PetscErrorCode  SectionRealToVec(SectionReal section, Mesh mesh, VecScatter scatter, ScatterMode mode, Vec vec)
585: {
586:   Vec            localVec;

591:   SectionRealCreateLocalVector(section, &localVec);
592:   if (mode == SCATTER_FORWARD) {
593:     VecScatterBegin(scatter, localVec, vec, INSERT_VALUES, mode);
594:     VecScatterEnd(scatter, localVec, vec, INSERT_VALUES, mode);
595:   } else {
596:     VecScatterBegin(scatter, vec, localVec, INSERT_VALUES, mode);
597:     VecScatterEnd(scatter, vec, localVec, INSERT_VALUES, mode);
598:   }
599:   VecDestroy(localVec);
600:   return(0);
601: }

605: /*@C
606:   SectionRealClear - Dellocate storage for this section

608:   Not collective

610:   Input Parameter:
611: . section - the section object

613:   Level: advanced

615: .seealso SectionRealRestrict(), SectionRealCreate(), SectionRealView()
616: @*/
617: PetscErrorCode  SectionRealClear(SectionReal section)
618: {
621:   section->s->clear();
622:   return(0);
623: }

627: /*@C
628:   SectionRealNorm - Computes the vector norm.

630:   Collective on Section

632:   Input Parameters:
633: +  section - the real Section
634: -  type - one of NORM_1, NORM_2, NORM_INFINITY.  Also available
635:           NORM_1_AND_2, which computes both norms and stores them
636:           in a two element array.

638:   Output Parameter:
639: . val - the norm 

641:   Notes:
642: $     NORM_1 denotes sum_i |x_i|
643: $     NORM_2 denotes sqrt(sum_i (x_i)^2)
644: $     NORM_INFINITY denotes max_i |x_i|

646:   Level: intermediate

648: .seealso VecNorm(), SectionRealCreate()
649: @*/
650: PetscErrorCode  SectionRealNorm(SectionReal section, Mesh mesh, NormType type, PetscReal *val)
651: {
652:   Obj<PETSC_MESH_TYPE> m;
653:   Obj<PETSC_MESH_TYPE::real_section_type> s;
654:   Vec            v;

659:   MeshGetMesh(mesh, m);
660:   SectionRealGetSection(section, s);
661:   const ALE::Obj<PETSC_MESH_TYPE::order_type>& order = m->getFactory()->getGlobalOrder(m, s->getName(), s);
662:   VecCreate(m->comm(), &v);
663:   VecSetSizes(v, order->getLocalSize(), order->getGlobalSize());
664:   VecSetFromOptions(v);
665:   SectionRealToVec(section, mesh, SCATTER_FORWARD, v);
666:   VecNorm(v, type, val);
667:   VecDestroy(v);
668:   return(0);
669: }

673: /*@C
674:   SectionRealAXPY - 

676:   Collective on Section

678:   Input Parameters:
679: +  section - the real Section
680: .  alpha - a scalar
681: -  X - the other real Section

683:   Output Parameter:
684: . section - the difference 

686:   Level: intermediate

688: .seealso VecNorm(), SectionRealCreate()
689: @*/
690: PetscErrorCode  SectionRealAXPY(SectionReal section, Mesh mesh, PetscScalar alpha, SectionReal X)
691: {
692:   Obj<PETSC_MESH_TYPE> m;
693:   Obj<PETSC_MESH_TYPE::real_section_type> s;
694:   Obj<PETSC_MESH_TYPE::real_section_type> sX;
695:   Vec            v, x;

700:   MeshGetMesh(mesh, m);
701:   SectionRealGetSection(section, s);
702:   SectionRealGetSection(X, sX);
703:   const ALE::Obj<PETSC_MESH_TYPE::order_type>& order = m->getFactory()->getGlobalOrder(m, s->getName(), s);
704:   VecCreate(m->comm(), &v);
705:   VecSetSizes(v, order->getLocalSize(), order->getGlobalSize());
706:   VecSetFromOptions(v);
707:   VecDuplicate(v, &x);
708:   SectionRealToVec(section, mesh, SCATTER_FORWARD, v);
709:   SectionRealToVec(X,       mesh, SCATTER_FORWARD, x);
710:   VecAXPY(v, alpha, x);
711:   SectionRealToVec(section, mesh, SCATTER_REVERSE, v);
712:   VecDestroy(v);
713:   VecDestroy(x);
714:   return(0);
715: }

719: /*@C
720:   MeshGetVertexSectionReal - Create a Section over the vertices with the specified fiber dimension

722:   Collective on Mesh

724:   Input Parameters:
725: + mesh - The Mesh object
726: - fiberDim - The number of degrees of freedom per vertex

728:   Output Parameter:
729: . section - The section

731:   Level: intermediate

733: .keywords: mesh, section, vertex
734: .seealso: MeshCreate(), SectionRealCreate()
735: @*/
736: PetscErrorCode MeshGetVertexSectionReal(Mesh mesh, PetscInt fiberDim, SectionReal *section)
737: {
738:   ALE::Obj<PETSC_MESH_TYPE> m;
739:   ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
740:   PetscErrorCode      ierr;

743:   MeshGetMesh(mesh, m);
744:   SectionRealCreate(m->comm(), section);
745:   SectionRealSetBundle(*section, m);
746:   SectionRealGetSection(*section, s);
747: #ifdef PETSC_OPT_SIEVE
748:   s->setChart(m->getSieve()->getChart());
749: #endif
750:   s->setFiberDimension(m->depthStratum(0), fiberDim);
751:   m->allocate(s);
752:   return(0);
753: }

757: /*@C
758:   MeshGetCellSectionReal - Create a Section over the cells with the specified fiber dimension

760:   Collective on Mesh

762:   Input Parameters:
763: + mesh - The Mesh object
764: - fiberDim - The section name

766:   Output Parameter:
767: . section - The section

769:   Level: intermediate

771: .keywords: mesh, section, cell
772: .seealso: MeshCreate(), SectionRealCreate()
773: @*/
774: PetscErrorCode MeshGetCellSectionReal(Mesh mesh, PetscInt fiberDim, SectionReal *section)
775: {
776:   ALE::Obj<PETSC_MESH_TYPE> m;
777:   ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
778:   PetscErrorCode      ierr;

781:   MeshGetMesh(mesh, m);
782:   SectionRealCreate(m->comm(), section);
783:   SectionRealSetBundle(*section, m);
784:   SectionRealGetSection(*section, s);
785:   s->setFiberDimension(m->heightStratum(0), fiberDim);
786:   m->allocate(s);
787:   return(0);
788: }

792: /*@C
793:     MeshCreateGlobalRealVector - Creates a vector of the correct size to be gathered into 
794:         by the mesh.

796:     Collective on Mesh

798:     Input Parameters:
799: +    mesh - the mesh object
800: -    section - The SectionReal

802:     Output Parameters:
803: .   gvec - the global vector

805:     Level: advanced

807: .seealso MeshDestroy(), MeshCreate(), MeshCreateGlobalVector()

809: @*/
810: PetscErrorCode  MeshCreateGlobalRealVector(Mesh mesh, SectionReal section, Vec *gvec)
811: {
812:   ALE::Obj<PETSC_MESH_TYPE> m;
813:   ALE::Obj<PETSC_MESH_TYPE::real_section_type> s;
814:   const char    *name;

818:   MeshGetMesh(mesh, m);
819:   SectionRealGetSection(section, s);
820:   PetscObjectGetName((PetscObject) section, &name);
821:   const ALE::Obj<PETSC_MESH_TYPE::order_type>& order = m->getFactory()->getGlobalOrder(m, name, s);

823:   VecCreate(m->comm(), gvec);
824:   VecSetSizes(*gvec, order->getLocalSize(), order->getGlobalSize());
825:   VecSetFromOptions(*gvec);
826:   return(0);
827: }

831: PetscErrorCode SectionIntView_Sieve(SectionInt section, PetscViewer viewer)
832: {
833:   PetscTruth     iascii, isbinary, isdraw;

837:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &iascii);
838:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &isbinary);
839:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_DRAW, &isdraw);

841:   if (iascii){
842:     ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
843:     ALE::Obj<PETSC_MESH_TYPE>                   b;
844:     const char                                  *name;

846:     SectionIntGetSection(section, s);
847:     SectionIntGetBundle(section, b);
848:     PetscObjectGetName((PetscObject) section, &name);
849:     SectionView_Sieve_Ascii(b, s, name, viewer);
850:   } else if (isbinary) {
851:     SETERRQ(PETSC_ERR_SUP, "Binary viewer not implemented for Section");
852:   } else if (isdraw){
853:     SETERRQ(PETSC_ERR_SUP, "Draw viewer not implemented for Section");
854:   } else {
855:     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported by this section object", ((PetscObject)viewer)->type_name);
856:   }
857:   return(0);
858: }

862: /*@C
863:    SectionIntView - Views a Section object. 

865:    Collective on Section

867:    Input Parameters:
868: +  section - the Section
869: -  viewer - an optional visualization context

871:    Notes:
872:    The available visualization contexts include
873: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
874: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
875:          output where only the first processor opens
876:          the file.  All other processors send their 
877:          data to the first processor to print. 

879:    You can change the format the section is printed using the 
880:    option PetscViewerSetFormat().

882:    The user can open alternative visualization contexts with
883: +    PetscViewerASCIIOpen() - Outputs section to a specified file
884: .    PetscViewerBinaryOpen() - Outputs section in binary to a
885:          specified file; corresponding input uses SectionLoad()
886: .    PetscViewerDrawOpen() - Outputs section to an X window display

888:    The user can call PetscViewerSetFormat() to specify the output
889:    format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
890:    PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen).  Available formats include
891: +    PETSC_VIEWER_DEFAULT - default, prints section information
892: -    PETSC_VIEWER_ASCII_VTK - outputs a VTK file describing the section

894:    Level: beginner

896:    Concepts: section^printing
897:    Concepts: section^saving to disk

899: .seealso: VecView(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerBinaryOpen(), PetscViewerCreate()
900: @*/
901: PetscErrorCode SectionIntView(SectionInt section, PetscViewer viewer)
902: {

908:   if (!viewer) {
909:     PetscViewerASCIIGetStdout(((PetscObject)section)->comm,&viewer);
910:   }

914:   PetscLogEventBegin(SectionInt_View,0,0,0,0);
915:   (*section->ops->view)(section, viewer);
916:   PetscLogEventEnd(SectionInt_View,0,0,0,0);
917:   return(0);
918: }

922: /*@C
923:   SectionIntGetSection - Gets the internal section object

925:   Not collective

927:   Input Parameter:
928: . section - the section object

930:   Output Parameter:
931: . s - the internal section object
932:  
933:   Level: advanced

935: .seealso SectionIntCreate(), SectionIntSetSection()
936: @*/
937: PetscErrorCode  SectionIntGetSection(SectionInt section, ALE::Obj<PETSC_MESH_TYPE::int_section_type>& s)
938: {
941:   s = section->s;
942:   return(0);
943: }

947: /*@C
948:   SectionIntSetSection - Sets the internal section object

950:   Not collective

952:   Input Parameters:
953: + section - the section object
954: - s - the internal section object
955:  
956:   Level: advanced

958: .seealso SectionIntCreate(), SectionIntGetSection()
959: @*/
960: PetscErrorCode  SectionIntSetSection(SectionInt section, const ALE::Obj<PETSC_MESH_TYPE::int_section_type>& s)
961: {

966:   if (!s.isNull()) {PetscObjectSetName((PetscObject) section, s->getName().c_str());}
967:   section->s = s;
968:   return(0);
969: }

973: /*@C
974:   SectionIntGetBundle - Gets the section bundle

976:   Not collective

978:   Input Parameter:
979: . section - the section object

981:   Output Parameter:
982: . b - the section bundle
983:  
984:   Level: advanced

986: .seealso SectionIntCreate(), SectionIntGetSection(), SectionIntSetSection()
987: @*/
988: PetscErrorCode  SectionIntGetBundle(SectionInt section, ALE::Obj<PETSC_MESH_TYPE>& b)
989: {
992:   b = section->b;
993:   return(0);
994: }

998: /*@C
999:   SectionIntSetBundle - Sets the section bundle

1001:   Not collective

1003:   Input Parameters:
1004: + section - the section object
1005: - b - the section bundle
1006:  
1007:   Level: advanced

1009: .seealso SectionIntCreate(), SectionIntGetSection(), SectionIntSetSection()
1010: @*/
1011: PetscErrorCode  SectionIntSetBundle(SectionInt section, const ALE::Obj<PETSC_MESH_TYPE>& b)
1012: {
1015:   section->b = b;
1016:   return(0);
1017: }

1021: /*@C
1022:   SectionIntCreate - Creates a Section object, used to manage data for an unstructured problem
1023:   described by a Sieve.

1025:   Collective on MPI_Comm

1027:   Input Parameter:
1028: . comm - the processors that will share the global section

1030:   Output Parameters:
1031: . section - the section object

1033:   Level: advanced

1035: .seealso SectionIntDestroy(), SectionIntView()
1036: @*/
1037: PetscErrorCode  SectionIntCreate(MPI_Comm comm, SectionInt *section)
1038: {
1040:   SectionInt    s;

1044:   *section = PETSC_NULL;

1046:   PetscHeaderCreate(s,_p_SectionInt,struct _SectionIntOps,SECTIONINT_COOKIE,0,"SectionInt",comm,SectionIntDestroy,0);
1047:   s->ops->view     = SectionIntView_Sieve;
1048:   s->ops->restrictClosure = SectionIntRestrict;
1049:   s->ops->update   = SectionIntUpdate;

1051:   PetscObjectChangeTypeName((PetscObject) s, "sieve");

1053:   new(&s->s) ALE::Obj<PETSC_MESH_TYPE::int_section_type>(PETSC_MESH_TYPE::int_section_type(comm));
1054:   new(&s->b) ALE::Obj<PETSC_MESH_TYPE>(PETSC_NULL);
1055:   *section = s;
1056:   return(0);
1057: }

1061: /*@C
1062:   SectionIntDestroy - Destroys a section.

1064:   Collective on Section

1066:   Input Parameter:
1067: . section - the section object

1069:   Level: advanced

1071: .seealso SectionIntCreate(), SectionIntView()
1072: @*/
1073: PetscErrorCode  SectionIntDestroy(SectionInt section)
1074: {

1079:   if (--((PetscObject)section)->refct > 0) return(0);
1080:   section->s = PETSC_NULL;
1081:   PetscHeaderDestroy(section);
1082:   return(0);
1083: }

1087: /*@C
1088:   SectionIntDistribute - Distributes the sections.

1090:   Not Collective

1092:   Input Parameters:
1093: + serialSection - The original Section object
1094: - parallelMesh - The parallel Mesh

1096:   Output Parameter:
1097: . parallelSection - The distributed Section object

1099:   Level: intermediate

1101: .keywords: mesh, section, distribute
1102: .seealso: MeshCreate()
1103: @*/
1104: PetscErrorCode SectionIntDistribute(SectionInt serialSection, Mesh parallelMesh, SectionInt *parallelSection)
1105: {
1106:   ALE::Obj<PETSC_MESH_TYPE::int_section_type> oldSection;
1107:   ALE::Obj<PETSC_MESH_TYPE> m;
1108:   PetscErrorCode      ierr;

1111:   SectionIntGetSection(serialSection, oldSection);
1112:   MeshGetMesh(parallelMesh, m);
1113:   SectionIntCreate(oldSection->comm(), parallelSection);
1114: #ifdef PETSC_OPT_SIEVE
1115:   SETERRQ(PETSC_ERR_SUP, "I am being lazy, bug me.");
1116: #else
1117:   ALE::Obj<PETSC_MESH_TYPE::int_section_type> newSection = ALE::Distribution<PETSC_MESH_TYPE>::distributeSection(oldSection, m, m->getDistSendOverlap(), m->getDistRecvOverlap());
1118:   SectionIntSetSection(*parallelSection, newSection);
1119: #endif
1120:   return(0);
1121: }

1125: /*@C
1126:   SectionIntRestrict - Restricts the SectionInt to a subset of the topology, returning an array of values.

1128:   Not collective

1130:   Input Parameters:
1131: + section - the section object
1132: - point - the Sieve point

1134:   Output Parameter:
1135: . values - The values associated with the submesh

1137:   Level: advanced

1139: .seealso SectionIntUpdate(), SectionIntCreate(), SectionIntView()
1140: @*/
1141: PetscErrorCode  SectionIntRestrict(SectionInt section, PetscInt point, PetscInt *values[])
1142: {
1146:   *values = (PetscInt *) section->b->restrictClosure(section->s, point);
1147:   return(0);
1148: }

1152: /*@C
1153:   SectionIntUpdate - Updates the array of values associated to a subset of the topology in this Section.

1155:   Not collective

1157:   Input Parameters:
1158: + section - the section object
1159: . point - the Sieve point
1160: - values - The values associated with the submesh

1162:   Level: advanced

1164: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1165: @*/
1166: PetscErrorCode  SectionIntUpdate(SectionInt section, PetscInt point, const PetscInt values[])
1167: {
1171:   section->b->update(section->s, point, values);
1172:   return(0);
1173: }

1177: /*@C
1178:   SectionIntUpdateAdd - Updates the array of values associated to a subset of the topology in this Section.

1180:   Not collective

1182:   Input Parameters:
1183: + section - the section object
1184: . point - the Sieve point
1185: - values - The values associated with the submesh

1187:   Level: advanced

1189: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1190: @*/
1191: PetscErrorCode  SectionIntUpdateAdd(SectionInt section, PetscInt point, const PetscInt values[])
1192: {
1196:   section->b->updateAdd(section->s, point, values);
1197:   return(0);
1198: }

1202: /*@C
1203:   SectionIntComplete - Exchanges data across the mesh overlap.

1205:   Not collective

1207:   Input Parameter:
1208: . section - the section object

1210:   Level: advanced

1212: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1213: @*/
1214: PetscErrorCode SectionIntComplete(SectionInt section)
1215: {
1216:   Obj<PETSC_MESH_TYPE::int_section_type> s;
1217:   Obj<PETSC_MESH_TYPE>                   b;

1221:   SectionIntGetSection(section, s);
1222:   SectionIntGetBundle(section, b);
1223:   ALE::Distribution<PETSC_MESH_TYPE>::completeSection(b, s);
1224:   return(0);
1225: }

1229: /*@C
1230:   SectionIntSetFiberDimension - Set the size of the vector space attached to the point

1232:   Not collective

1234:   Input Parameters:
1235: + section - the section object
1236: . point - the Sieve point
1237: - size - The fiber dimension

1239:   Level: advanced

1241: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1242: @*/
1243: PetscErrorCode  SectionIntSetFiberDimension(SectionInt section, PetscInt point, const PetscInt size)
1244: {
1247:   section->s->setFiberDimension(point, size);
1248:   return(0);
1249: }

1253: /*@C
1254:   SectionIntAllocate - Allocate storage for this section

1256:   Not collective

1258:   Input Parameter:
1259: . section - the section object

1261:   Level: advanced

1263: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1264: @*/
1265: PetscErrorCode  SectionIntAllocate(SectionInt section)
1266: {
1269:   section->b->allocate(section->s);
1270:   return(0);
1271: }

1275: /*@C
1276:   SectionIntClear - Dellocate storage for this section

1278:   Not collective

1280:   Input Parameter:
1281: . section - the section object

1283:   Level: advanced

1285: .seealso SectionIntRestrict(), SectionIntCreate(), SectionIntView()
1286: @*/
1287: PetscErrorCode  SectionIntClear(SectionInt section)
1288: {
1291:   section->s->clear();
1292:   return(0);
1293: }

1297: /*@C
1298:   MeshGetVertexSectionInt - Create a Section over the vertices with the specified fiber dimension

1300:   Collective on Mesh

1302:   Input Parameters:
1303: + mesh - The Mesh object
1304: - fiberDim - The section name

1306:   Output Parameter:
1307: . section - The section

1309:   Level: intermediate

1311: .keywords: mesh, section, vertex
1312: .seealso: MeshCreate(), SectionIntCreate()
1313: @*/
1314: PetscErrorCode MeshGetVertexSectionInt(Mesh mesh, PetscInt fiberDim, SectionInt *section)
1315: {
1316:   ALE::Obj<PETSC_MESH_TYPE> m;
1317:   ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1318:   PetscErrorCode      ierr;

1321:   MeshGetMesh(mesh, m);
1322:   SectionIntCreate(m->comm(), section);
1323:   SectionIntSetBundle(*section, m);
1324:   SectionIntGetSection(*section, s);
1325:   s->setFiberDimension(m->depthStratum(0), fiberDim);
1326:   m->allocate(s);
1327:   return(0);
1328: }

1332: /*@C
1333:   MeshGetCellSectionInt - Create a Section over the cells with the specified fiber dimension

1335:   Collective on Mesh

1337:   Input Parameters:
1338: + mesh - The Mesh object
1339: - fiberDim - The section name

1341:   Output Parameter:
1342: . section - The section

1344:   Level: intermediate

1346: .keywords: mesh, section, cell
1347: .seealso: MeshCreate(), SectionIntCreate()
1348: @*/
1349: PetscErrorCode MeshGetCellSectionInt(Mesh mesh, PetscInt fiberDim, SectionInt *section)
1350: {
1351:   ALE::Obj<PETSC_MESH_TYPE> m;
1352:   ALE::Obj<PETSC_MESH_TYPE::int_section_type> s;
1353:   PetscErrorCode      ierr;

1356:   MeshGetMesh(mesh, m);
1357:   SectionIntCreate(m->comm(), section);
1358:   SectionIntSetBundle(*section, m);
1359:   SectionIntGetSection(*section, s);
1360:   s->setFiberDimension(m->heightStratum(0), fiberDim);
1361:   m->allocate(s);
1362:   return(0);
1363: }