Actual source code: dupl.c
1: #define PETSC_DLL
3: #include private/viewerimpl.h
7: /*@
8: PetscViewerGetSingleton - Creates a new PetscViewer (same type as the old)
9: that lives on a single processor (with MPI_comm PETSC_COMM_SELF)
11: Collective on PetscViewer
13: Input Parameter:
14: . viewer - the PetscViewer to be duplicated
16: Output Parameter:
17: . outviewer - new PetscViewer
19: Level: advanced
21: Notes: Call PetscViewerRestoreSingleton() to return this PetscViewer, NOT PetscViewerDestroy()
23: This is most commonly used to view a sequential object that is part of a
24: parallel object. For example block Jacobi PC view could use this to obtain a
25: PetscViewer that is used with the sequential KSP on one block of the preconditioner.
27: Concepts: PetscViewer^sequential version
29: .seealso: PetscViewerSocketOpen(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerRestoreSingleton()
30: @*/
31: PetscErrorCode PetscViewerGetSingleton(PetscViewer viewer,PetscViewer *outviewer)
32: {
34: PetscMPIInt size;
39: MPI_Comm_size(((PetscObject)viewer)->comm,&size);
40: if (size == 1) {
41: PetscObjectReference((PetscObject)viewer);
42: *outviewer = viewer;
43: } else if (viewer->ops->getsingleton) {
44: (*viewer->ops->getsingleton)(viewer,outviewer);
45: } else {
46: SETERRQ1(PETSC_ERR_SUP,"Cannot get singleton PetscViewer for type %s",((PetscObject)viewer)->type_name);
47: }
48: return(0);
49: }
53: /*@
54: PetscViewerRestoreSingleton - Restores a new PetscViewer obtained with PetscViewerGetSingleton().
56: Collective on PetscViewer
58: Input Parameters:
59: + viewer - the PetscViewer to be duplicated
60: - outviewer - new PetscViewer
62: Level: advanced
64: Notes: Call PetscViewerGetSingleton() to get this PetscViewer, NOT PetscViewerCreate()
66: .seealso: PetscViewerSocketOpen(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerGetSingleton()
67: @*/
68: PetscErrorCode PetscViewerRestoreSingleton(PetscViewer viewer,PetscViewer *outviewer)
69: {
71: PetscMPIInt size;
76: MPI_Comm_size(((PetscObject)viewer)->comm,&size);
77: if (size == 1) {
78: PetscObjectDereference((PetscObject)viewer);
79: if (outviewer) *outviewer = 0;
80: } else if (viewer->ops->restoresingleton) {
81: (*viewer->ops->restoresingleton)(viewer,outviewer);
82: }
83: return(0);
84: }
88: /*@
89: PetscViewerGetSubcomm - Creates a new PetscViewer (same type as the old)
90: that lives on a subgroup of processors
92: Collective on PetscViewer
94: Input Parameter:
95: + viewer - the PetscViewer to be duplicated
96: - subcomm - MPI communicator
98: Output Parameter:
99: . outviewer - new PetscViewer
101: Level: advanced
103: Notes: Call PetscViewerRestoreSubcomm() to return this PetscViewer, NOT PetscViewerDestroy()
105: This is used to view a sequential or a parallel object that is part of a larger
106: parallel object. For example redundant PC view could use this to obtain a
107: PetscViewer that is used within a subcommunicator on one duplicated preconditioner.
109: Concepts: PetscViewer^sequential version
111: .seealso: PetscViewerSocketOpen(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerRestoreSubcomm()
112: @*/
113: PetscErrorCode PetscViewerGetSubcomm(PetscViewer viewer,MPI_Comm subcomm,PetscViewer *outviewer)
114: {
116: PetscMPIInt size;
121: MPI_Comm_size(((PetscObject)viewer)->comm,&size);
122: if (size == 1) {
123: PetscObjectReference((PetscObject)viewer);
124: *outviewer = viewer;
125: } else if (viewer->ops->getsubcomm) {
126: (*viewer->ops->getsubcomm)(viewer,subcomm,outviewer);
127: } else {
128: SETERRQ1(PETSC_ERR_SUP,"Cannot get subcommunicator PetscViewer for type %s",((PetscObject)viewer)->type_name);
129: }
130: return(0);
131: }
135: /*@
136: PetscViewerRestoreSubcomm - Restores a new PetscViewer obtained with PetscViewerGetSubcomm().
138: Collective on PetscViewer
140: Input Parameters:
141: + viewer - the PetscViewer to be duplicated
142: . subcomm - MPI communicator
143: - outviewer - new PetscViewer
145: Level: advanced
147: Notes: Call PetscViewerGetSubcomm() to get this PetscViewer, NOT PetscViewerCreate()
149: .seealso: PetscViewerSocketOpen(), PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscViewerGetSubcomm()
150: @*/
151: PetscErrorCode PetscViewerRestoreSubcomm(PetscViewer viewer,MPI_Comm subcomm,PetscViewer *outviewer)
152: {
154: PetscMPIInt size;
159: MPI_Comm_size(((PetscObject)viewer)->comm,&size);
160: if (size == 1) {
161: PetscObjectDereference((PetscObject)viewer);
162: if (outviewer) *outviewer = 0;
163: } else if (viewer->ops->restoresubcomm) {
164: (*viewer->ops->restoresubcomm)(viewer,subcomm,outviewer);
165: }
166: return(0);
167: }