Actual source code: andor.c
1: #include <petsc/private/vecimpl.h>
2: #include "../src/vec/vec/utils/tagger/impls/andor.h"
4: static PetscErrorCode VecTaggerDestroy_AndOr(VecTagger tagger)
5: {
6: VecTagger_AndOr *andOr = (VecTagger_AndOr *)tagger->data;
7: PetscInt i;
9: for (i = 0; i < andOr->nsubs; i++) VecTaggerDestroy(&andOr->subs[i]);
10: if (andOr->mode == PETSC_OWN_POINTER) PetscFree(andOr->subs);
11: PetscFree(tagger->data);
12: return 0;
13: }
15: PetscErrorCode VecTaggerGetSubs_AndOr(VecTagger tagger, PetscInt *nsubs, VecTagger **subs)
16: {
17: VecTagger_AndOr *andOr = (VecTagger_AndOr *)tagger->data;
20: if (nsubs) {
22: *nsubs = andOr->nsubs;
23: }
24: if (subs) {
26: *subs = andOr->subs;
27: }
28: return 0;
29: }
31: PetscErrorCode VecTaggerSetSubs_AndOr(VecTagger tagger, PetscInt nsubs, VecTagger *subs, PetscCopyMode mode)
32: {
33: PetscInt i;
34: VecTagger_AndOr *andOr = (VecTagger_AndOr *)tagger->data;
38: if (nsubs == andOr->nsubs && subs == andOr->subs && mode != PETSC_COPY_VALUES) return 0;
39: if (subs) {
40: for (i = 0; i < nsubs; i++) PetscObjectReference((PetscObject)subs[i]);
41: }
42: for (i = 0; i < andOr->nsubs; i++) VecTaggerDestroy(&(andOr->subs[i]));
43: if (andOr->mode == PETSC_OWN_POINTER && andOr->subs != subs) PetscFree(andOr->subs);
44: andOr->nsubs = nsubs;
45: if (subs) {
46: if (mode == PETSC_COPY_VALUES) {
47: andOr->mode = PETSC_OWN_POINTER;
48: PetscMalloc1(nsubs, &(andOr->subs));
49: for (i = 0; i < nsubs; i++) andOr->subs[i] = subs[i];
50: } else {
51: andOr->subs = subs;
52: andOr->mode = mode;
53: for (i = 0; i < nsubs; i++) PetscObjectDereference((PetscObject)subs[i]);
54: }
55: } else {
56: MPI_Comm comm = PetscObjectComm((PetscObject)tagger);
57: PetscInt bs;
58: const char *prefix;
59: char tprefix[128];
61: VecTaggerGetBlockSize(tagger, &bs);
62: PetscObjectGetOptionsPrefix((PetscObject)tagger, &prefix);
63: andOr->mode = PETSC_OWN_POINTER;
64: PetscMalloc1(nsubs, &(andOr->subs));
65: for (i = 0; i < nsubs; i++) {
66: VecTagger sub;
68: PetscSNPrintf(tprefix, 128, "sub_%" PetscInt_FMT "_", i);
69: VecTaggerCreate(comm, &sub);
70: VecTaggerSetBlockSize(sub, bs);
71: PetscObjectSetOptionsPrefix((PetscObject)sub, prefix);
72: PetscObjectAppendOptionsPrefix((PetscObject)sub, tprefix);
73: andOr->subs[i] = sub;
74: }
75: }
76: return 0;
77: }
79: static PetscErrorCode VecTaggerSetFromOptions_AndOr(VecTagger tagger, PetscOptionItems *PetscOptionsObject)
80: {
81: PetscInt i, nsubs, nsubsOrig;
82: const char *name;
83: char headstring[BUFSIZ];
84: char funcstring[BUFSIZ];
85: char descstring[BUFSIZ];
86: VecTagger *subs;
88: PetscObjectGetType((PetscObject)tagger, &name);
89: VecTaggerGetSubs_AndOr(tagger, &nsubs, NULL);
90: nsubsOrig = nsubs;
91: PetscSNPrintf(headstring, sizeof(headstring), "VecTagger %s options", name);
92: PetscSNPrintf(funcstring, sizeof(funcstring), "VecTagger%sSetSubs()", name);
93: PetscSNPrintf(descstring, sizeof(descstring), "number of sub tags in %s tag", name);
94: PetscOptionsHeadBegin(PetscOptionsObject, headstring);
95: PetscOptionsInt("-vec_tagger_num_subs", descstring, funcstring, nsubs, &nsubs, NULL);
96: PetscOptionsHeadEnd();
97: if (nsubs != nsubsOrig) {
98: VecTaggerSetSubs_AndOr(tagger, nsubs, NULL, PETSC_OWN_POINTER);
99: VecTaggerGetSubs_AndOr(tagger, NULL, &subs);
100: for (i = 0; i < nsubs; i++) VecTaggerSetFromOptions(subs[i]);
101: }
102: return 0;
103: }
105: static PetscErrorCode VecTaggerSetUp_AndOr(VecTagger tagger)
106: {
107: PetscInt nsubs, i;
108: VecTagger *subs;
110: VecTaggerGetSubs_AndOr(tagger, &nsubs, &subs);
112: for (i = 0; i < nsubs; i++) VecTaggerSetUp(subs[i]);
113: return 0;
114: }
116: static PetscErrorCode VecTaggerView_AndOr(VecTagger tagger, PetscViewer viewer)
117: {
118: PetscBool iascii;
120: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
121: if (iascii) {
122: PetscInt i, nsubs;
123: VecTagger *subs;
124: const char *name;
126: VecTaggerGetSubs_AndOr(tagger, &nsubs, &subs);
127: PetscObjectGetType((PetscObject)tagger, &name);
128: PetscViewerASCIIPrintf(viewer, " %s of %" PetscInt_FMT " subtags:\n", name, nsubs);
129: PetscViewerASCIIPushTab(viewer);
130: for (i = 0; i < nsubs; i++) VecTaggerView(subs[i], viewer);
131: PetscViewerASCIIPopTab(viewer);
132: }
133: return 0;
134: }
136: PetscErrorCode VecTaggerCreate_AndOr(VecTagger tagger)
137: {
138: VecTagger_AndOr *andOr;
140: tagger->ops->destroy = VecTaggerDestroy_AndOr;
141: tagger->ops->setfromoptions = VecTaggerSetFromOptions_AndOr;
142: tagger->ops->setup = VecTaggerSetUp_AndOr;
143: tagger->ops->view = VecTaggerView_AndOr;
144: tagger->ops->computeis = VecTaggerComputeIS_FromBoxes;
145: PetscNew(&andOr);
146: tagger->data = andOr;
147: return 0;
148: }
150: PetscErrorCode VecTaggerAndOrIsSubBox_Private(PetscInt bs, const VecTaggerBox *superBox, const VecTaggerBox *subBox, PetscBool *isSub)
151: {
152: PetscInt i;
154: *isSub = PETSC_TRUE;
155: for (i = 0; i < bs; i++) {
156: #if !defined(PETSC_USE_COMPLEX)
157: if (superBox[i].min > subBox[i].min || superBox[i].max < subBox[i].max) {
158: *isSub = PETSC_FALSE;
159: break;
160: }
161: #else
162: if (PetscRealPart(superBox[i].min) > PetscRealPart(subBox[i].min) || PetscImaginaryPart(superBox[i].min) > PetscImaginaryPart(subBox[i].min) || PetscRealPart(superBox[i].max) < PetscRealPart(subBox[i].max) ||
163: PetscImaginaryPart(superBox[i].max) < PetscImaginaryPart(subBox[i].max)) {
164: *isSub = PETSC_FALSE;
165: break;
166: }
167: #endif
168: }
169: return 0;
170: }
172: PetscErrorCode VecTaggerAndOrIntersect_Private(PetscInt bs, const VecTaggerBox *a, const VecTaggerBox *b, VecTaggerBox *c, PetscBool *empty)
173: {
174: PetscInt i;
176: *empty = PETSC_FALSE;
177: for (i = 0; i < bs; i++) {
178: #if !defined(PETSC_USE_COMPLEX)
179: c[i].min = PetscMax(a[i].min, b[i].min);
180: c[i].max = PetscMin(a[i].max, b[i].max);
181: if (c[i].max < c[i].min) {
182: *empty = PETSC_TRUE;
183: break;
184: }
185: #else
186: {
187: PetscReal maxMinReal = PetscMax(PetscRealPart(a[i].min), PetscRealPart(b[i].min));
188: PetscReal maxMinImag = PetscMax(PetscImaginaryPart(a[i].min), PetscImaginaryPart(b[i].min));
189: PetscReal minMaxReal = PetscMin(PetscRealPart(a[i].max), PetscRealPart(b[i].max));
190: PetscReal minMaxImag = PetscMin(PetscImaginaryPart(a[i].max), PetscImaginaryPart(b[i].max));
192: c[i].min = PetscCMPLX(maxMinReal, maxMinImag);
193: c[i].max = PetscCMPLX(minMaxReal, minMaxImag);
194: if ((PetscRealPart(c[i].max - c[i].min) < 0.) || (PetscImaginaryPart(c[i].max - c[i].min) < 0.)) {
195: *empty = PETSC_TRUE;
196: break;
197: }
198: }
199: #endif
200: }
201: return 0;
202: }