File: | dm/field/impls/da/dmfieldda.c |
Warning: | line 361, column 21 Division by zero |
[?] Use j/k keys for keyboard navigation
1 | #include <petsc/private/dmfieldimpl.h> /*I "petscdmfield.h" I*/ | |||
2 | #include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/ | |||
3 | #include <petscdmda.h> | |||
4 | ||||
5 | typedef struct _n_DMField_DA | |||
6 | { | |||
7 | PetscScalar *cornerVals; | |||
8 | PetscScalar *cornerCoeffs; | |||
9 | PetscScalar *work; | |||
10 | PetscReal coordRange[3][2]; | |||
11 | } | |||
12 | DMField_DA; | |||
13 | ||||
14 | static PetscErrorCode DMFieldDestroy_DA(DMField field) | |||
15 | { | |||
16 | DMField_DA *dafield; | |||
17 | PetscErrorCode ierr; | |||
18 | ||||
19 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 19; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
20 | dafield = (DMField_DA *) field->data; | |||
21 | ierr = PetscFree3(dafield->cornerVals,dafield->cornerCoeffs,dafield->work)PetscFreeA(3,21,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,&(dafield->cornerVals),&(dafield->cornerCoeffs ),&(dafield->work));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),21,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
22 | ierr = PetscFree(dafield)((*PetscTrFree)((void*)(dafield),22,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ) || ((dafield) = 0,0));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),22,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
23 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
24 | } | |||
25 | ||||
26 | static PetscErrorCode DMFieldView_DA(DMField field,PetscViewer viewer) | |||
27 | { | |||
28 | DMField_DA *dafield = (DMField_DA *) field->data; | |||
29 | PetscBool iascii; | |||
30 | PetscErrorCode ierr; | |||
31 | ||||
32 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 32; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
33 | ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII"ascii",&iascii);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),33,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
34 | if (iascii) { | |||
35 | PetscInt i, c, dim; | |||
36 | PetscInt nc; | |||
37 | DM dm = field->dm; | |||
38 | ||||
39 | PetscViewerASCIIPrintf(viewer, "Field corner values:\n");CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),39,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
40 | ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),40,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
41 | ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),41,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
42 | nc = field->numComponents; | |||
43 | for (i = 0, c = 0; i < (1 << dim); i++) { | |||
44 | PetscInt j; | |||
45 | ||||
46 | for (j = 0; j < nc; j++, c++) { | |||
47 | PetscScalar val = dafield->cornerVals[nc * i + j]; | |||
48 | ||||
49 | #if !defined(PETSC_USE_COMPLEX) | |||
50 | ierr = PetscViewerASCIIPrintf(viewer,"%g ",(double) val);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),50,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
51 | #else | |||
52 | ierr = PetscViewerASCIIPrintf(viewer,"%g+i%g ",(double) PetscRealPart(val)(val),(double) PetscImaginaryPart(val)((PetscReal)0.));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),52,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
53 | #endif | |||
54 | } | |||
55 | ierr = PetscViewerASCIIPrintf(viewer,"\n");CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),55,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
56 | } | |||
57 | ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),57,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
58 | } | |||
59 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
60 | } | |||
61 | ||||
62 | #define MEdot(y,A,x,m,c,cast)do { PetscInt _k, _l; for (_k = 0; _k < (c); _k++) (y)[_k] = 0.; for (_l = 0; _l < (m); _l++) { for (_k = 0; _k < (c); _k++) { (y)[_k] += cast((A)[(c) * _l + _k]) * (x)[_l]; } } } while (0) \ | |||
63 | do { \ | |||
64 | PetscInt _k, _l; \ | |||
65 | for (_k = 0; _k < (c); _k++) (y)[_k] = 0.; \ | |||
66 | for (_l = 0; _l < (m); _l++) { \ | |||
67 | for (_k = 0; _k < (c); _k++) { \ | |||
68 | (y)[_k] += cast((A)[(c) * _l + _k]) * (x)[_l]; \ | |||
69 | } \ | |||
70 | } \ | |||
71 | } while (0) | |||
72 | ||||
73 | #define MEHess(out,cf,etaB,etaD,dim,nc,cast)do { PetscInt _m, _j, _k; for (_m = 0; _m < (nc) * (dim) * (dim); _m++) (out)[_m] = 0.; for (_j = 0; _j < (dim); _j++ ) { for (_k = _j + 1; _k < (dim); _k++) { PetscInt _ind = ( 1 << _j) + (1 << _k); for (_m = 0; _m < (nc); _m ++) { PetscScalar c = (cf)[_m] * (etaB)[_ind] * (etaD)[_ind]; (out)[(_m * (dim) + _k) * (dim) + _j] += cast(c); (out)[(_m * (dim) + _j) * (dim) + _k] += cast(c); } } } } while (0) \ | |||
74 | do { \ | |||
75 | PetscInt _m, _j, _k; \ | |||
76 | for (_m = 0; _m < (nc) * (dim) * (dim); _m++) (out)[_m] = 0.; \ | |||
77 | for (_j = 0; _j < (dim); _j++) { \ | |||
78 | for (_k = _j + 1; _k < (dim); _k++) { \ | |||
79 | PetscInt _ind = (1 << _j) + (1 << _k); \ | |||
80 | for (_m = 0; _m < (nc); _m++) { \ | |||
81 | PetscScalar c = (cf)[_m] * (etaB)[_ind] * (etaD)[_ind]; \ | |||
82 | (out)[(_m * (dim) + _k) * (dim) + _j] += cast(c); \ | |||
83 | (out)[(_m * (dim) + _j) * (dim) + _k] += cast(c); \ | |||
84 | } \ | |||
85 | } \ | |||
86 | } \ | |||
87 | } while (0) | |||
88 | ||||
89 | static void MultilinearEvaluate(PetscInt dim, PetscReal (*coordRange)[2], PetscInt nc, PetscScalar *cf, PetscScalar *cfWork, PetscInt nPoints, const PetscScalar *points, PetscDataType datatype, void *B, void *D, void *H) | |||
90 | { | |||
91 | PetscInt i, j, k, l, m; | |||
92 | PetscInt whol = 1 << dim; | |||
93 | PetscInt half = whol >> 1; | |||
94 | ||||
95 | PetscFunctionBeginHotdo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 95; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_TRUE || petscstack->hotdepth); } ; } while (0); ; } while (0); | |||
96 | if (!B && !D && !H) PetscFunctionReturnVoid()do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return;} while (0); | |||
97 | for (i = 0; i < nPoints; i++) { | |||
98 | const PetscScalar *point = &points[dim * i]; | |||
99 | PetscReal deta[3] = {0.}; | |||
100 | PetscReal etaB[8] = {1.,1.,1.,1.,1.,1.,1.,1.}; | |||
101 | PetscReal etaD[8] = {1.,1.,1.,1.,1.,1.,1.,1.}; | |||
102 | PetscReal work[8]; | |||
103 | ||||
104 | for (j = 0; j < dim; j++) { | |||
105 | PetscReal e, d; | |||
106 | ||||
107 | e = (PetscRealPart(point[j])(point[j]) - coordRange[j][0]) / coordRange[j][1]; | |||
108 | deta[j] = d = 1. / coordRange[j][1]; | |||
109 | for (k = 0; k < whol; k++) {work[k] = etaB[k];} | |||
110 | for (k = 0; k < half; k++) { | |||
111 | etaB[k] = work[2 * k] * e; | |||
112 | etaB[k + half] = work[2 * k + 1]; | |||
113 | } | |||
114 | if (H) { | |||
115 | for (k = 0; k < whol; k++) {work[k] = etaD[k];} | |||
116 | for (k = 0; k < half; k++) { | |||
117 | etaD[k + half] = work[2 * k]; | |||
118 | etaD[k ] = work[2 * k + 1] * d; | |||
119 | } | |||
120 | } | |||
121 | } | |||
122 | if (B) { | |||
123 | if (datatype == PETSC_SCALARPETSC_DOUBLE) { | |||
124 | PetscScalar *out = &((PetscScalar *)B)[nc * i]; | |||
125 | ||||
126 | MEdot(out,cf,etaB,(1 << dim),nc,(PetscScalar))do { PetscInt _k, _l; for (_k = 0; _k < (nc); _k++) (out)[ _k] = 0.; for (_l = 0; _l < ((1 << dim)); _l++) { for (_k = 0; _k < (nc); _k++) { (out)[_k] += (PetscScalar)((cf )[(nc) * _l + _k]) * (etaB)[_l]; } } } while (0); | |||
127 | } else { | |||
128 | PetscReal *out = &((PetscReal *)B)[nc * i]; | |||
129 | ||||
130 | MEdot(out,cf,etaB,(1 << dim),nc,PetscRealPart)do { PetscInt _k, _l; for (_k = 0; _k < (nc); _k++) (out)[ _k] = 0.; for (_l = 0; _l < ((1 << dim)); _l++) { for (_k = 0; _k < (nc); _k++) { (out)[_k] += ((cf)[(nc) * _l + _k]) * (etaB)[_l]; } } } while (0); | |||
131 | } | |||
132 | } | |||
133 | if (D) { | |||
134 | if (datatype == PETSC_SCALARPETSC_DOUBLE) { | |||
135 | PetscScalar *out = &((PetscScalar *)D)[nc * dim * i]; | |||
136 | ||||
137 | for (m = 0; m < nc * dim; m++) out[m] = 0.; | |||
138 | } else { | |||
139 | PetscReal *out = &((PetscReal *)D)[nc * dim * i]; | |||
140 | ||||
141 | for (m = 0; m < nc * dim; m++) out[m] = 0.; | |||
142 | } | |||
143 | for (j = 0; j < dim; j++) { | |||
144 | PetscReal d = deta[j]; | |||
145 | ||||
146 | for (k = 0; k < whol * nc; k++) {cfWork[k] = cf[k];} | |||
147 | for (k = 0; k < whol; k++) {work[k] = etaB[k];} | |||
148 | for (k = 0; k < half; k++) { | |||
149 | PetscReal e; | |||
150 | ||||
151 | etaB[k] = work[2 * k]; | |||
152 | etaB[k + half] = e = work[2 * k + 1]; | |||
153 | ||||
154 | for (l = 0; l < nc; l++) { | |||
155 | cf[ k * nc + l] = cfWork[ 2 * k * nc + l]; | |||
156 | cf[(k + half) * nc + l] = cfWork[(2 * k + 1) * nc + l]; | |||
157 | } | |||
158 | if (datatype == PETSC_SCALARPETSC_DOUBLE) { | |||
159 | PetscScalar *out = &((PetscScalar *)D)[nc * dim * i]; | |||
160 | ||||
161 | for (l = 0; l < nc; l++) { | |||
162 | out[l * dim + j] += d * e * cf[k * nc + l]; | |||
163 | } | |||
164 | } else { | |||
165 | PetscReal *out = &((PetscReal *)D)[nc * dim * i]; | |||
166 | ||||
167 | for (l = 0; l < nc; l++) { | |||
168 | out[l * dim + j] += d * e * PetscRealPart(cf[k * nc + l])(cf[k * nc + l]); | |||
169 | } | |||
170 | } | |||
171 | } | |||
172 | } | |||
173 | } | |||
174 | if (H) { | |||
175 | if (datatype == PETSC_SCALARPETSC_DOUBLE) { | |||
176 | PetscScalar *out = &((PetscScalar *)H)[nc * dim * dim * i]; | |||
177 | ||||
178 | MEHess(out,cf,etaB,etaD,dim,nc,(PetscScalar))do { PetscInt _m, _j, _k; for (_m = 0; _m < (nc) * (dim) * (dim); _m++) (out)[_m] = 0.; for (_j = 0; _j < (dim); _j++ ) { for (_k = _j + 1; _k < (dim); _k++) { PetscInt _ind = ( 1 << _j) + (1 << _k); for (_m = 0; _m < (nc); _m ++) { PetscScalar c = (cf)[_m] * (etaB)[_ind] * (etaD)[_ind]; (out)[(_m * (dim) + _k) * (dim) + _j] += (PetscScalar)(c); ( out)[(_m * (dim) + _j) * (dim) + _k] += (PetscScalar)(c); } } } } while (0); | |||
179 | } else { | |||
180 | PetscReal *out = &((PetscReal *)H)[nc * dim * dim * i]; | |||
181 | ||||
182 | MEHess(out,cf,etaB,etaD,dim,nc,PetscRealPart)do { PetscInt _m, _j, _k; for (_m = 0; _m < (nc) * (dim) * (dim); _m++) (out)[_m] = 0.; for (_j = 0; _j < (dim); _j++ ) { for (_k = _j + 1; _k < (dim); _k++) { PetscInt _ind = ( 1 << _j) + (1 << _k); for (_m = 0; _m < (nc); _m ++) { PetscScalar c = (cf)[_m] * (etaB)[_ind] * (etaD)[_ind]; (out)[(_m * (dim) + _k) * (dim) + _j] += (c); (out)[(_m * (dim ) + _j) * (dim) + _k] += (c); } } } } while (0); | |||
183 | } | |||
184 | } | |||
185 | } | |||
186 | PetscFunctionReturnVoid()do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return;} while (0); | |||
187 | } | |||
188 | ||||
189 | static PetscErrorCode DMFieldEvaluate_DA(DMField field, Vec points, PetscDataType datatype, void *B, void *D, void *H) | |||
190 | { | |||
191 | DM dm; | |||
192 | DMField_DA *dafield; | |||
193 | PetscInt dim; | |||
194 | PetscInt N, n, nc; | |||
195 | const PetscScalar *array; | |||
196 | PetscReal (*coordRange)[2]; | |||
197 | PetscErrorCode ierr; | |||
198 | ||||
199 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 199; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
200 | dm = field->dm; | |||
201 | nc = field->numComponents; | |||
202 | dafield = (DMField_DA *) field->data; | |||
203 | ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),203,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
204 | ierr = VecGetLocalSize(points,&N);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),204,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
205 | if (N % dim) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Point vector size %D not divisible by coordinate dimension %D\n",N,dim)return PetscError(((MPI_Comm)0x44000001),205,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,75,PETSC_ERROR_INITIAL,"Point vector size %D not divisible by coordinate dimension %D\n" ,N,dim); | |||
206 | n = N / dim; | |||
207 | coordRange = &(dafield->coordRange[0]); | |||
208 | ierr = VecGetArrayRead(points,&array);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),208,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
209 | MultilinearEvaluate(dim,coordRange,nc,dafield->cornerCoeffs,dafield->work,n,array,datatype,B,D,H); | |||
210 | ierr = VecRestoreArrayRead(points,&array);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),210,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
211 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
212 | } | |||
213 | ||||
214 | static PetscErrorCode DMFieldEvaluateFE_DA(DMField field, IS cellIS, PetscQuadrature points, PetscDataType datatype, void *B, void *D, void *H) | |||
215 | { | |||
216 | PetscInt c, i, j, k, dim, cellsPer[3] = {0}, first[3] = {0}, whol, half; | |||
217 | PetscReal stepPer[3] = {0.}; | |||
218 | PetscReal cellCoordRange[3][2] = {{0.,1.},{0.,1.},{0.,1.}}; | |||
219 | PetscScalar *cellCoeffs, *work; | |||
220 | DM dm; | |||
221 | DMDALocalInfo info; | |||
222 | PetscInt cStart, cEnd; | |||
223 | PetscInt nq, nc; | |||
224 | const PetscReal *q; | |||
225 | #if defined(PETSC_USE_COMPLEX) | |||
226 | PetscScalar *qs; | |||
227 | #else | |||
228 | const PetscScalar *qs; | |||
229 | #endif | |||
230 | DMField_DA *dafield; | |||
231 | PetscBool isStride; | |||
232 | const PetscInt *cells = NULL((void*)0); | |||
233 | PetscInt sfirst = -1, stride = -1, nCells; | |||
234 | PetscErrorCode ierr; | |||
235 | ||||
236 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 236; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
237 | dafield = (DMField_DA *) field->data; | |||
238 | dm = field->dm; | |||
239 | nc = field->numComponents; | |||
240 | ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),240,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
241 | dim = info.dim; | |||
242 | work = dafield->work; | |||
243 | stepPer[0] = 1./ info.mx; | |||
244 | stepPer[1] = 1./ info.my; | |||
245 | stepPer[2] = 1./ info.mz; | |||
246 | first[0] = info.gxs; | |||
247 | first[1] = info.gys; | |||
248 | first[2] = info.gzs; | |||
249 | cellsPer[0] = info.gxm; | |||
250 | cellsPer[1] = info.gym; | |||
251 | cellsPer[2] = info.gzm; | |||
252 | /* TODO: probably take components into account */ | |||
253 | ierr = PetscQuadratureGetData(points, NULL((void*)0), NULL((void*)0), &nq, &q, NULL((void*)0));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),253,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
254 | #if defined(PETSC_USE_COMPLEX) | |||
255 | ierr = DMGetWorkArray(dm,nq * dim,MPIU_SCALAR((MPI_Datatype)0x4c00080b),&qs);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),255,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
256 | for (i = 0; i < nq * dim; i++) qs[i] = q[i]; | |||
257 | #else | |||
258 | qs = q; | |||
259 | #endif | |||
260 | ierr = DMDAGetHeightStratum(dm,0,&cStart,&cEnd);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),260,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
261 | ierr = DMGetWorkArray(dm,(1 << dim) * nc,MPIU_SCALAR((MPI_Datatype)0x4c00080b),&cellCoeffs);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),261,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
262 | whol = (1 << dim); | |||
263 | half = whol >> 1; | |||
264 | ierr = ISGetLocalSize(cellIS,&nCells);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),264,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
265 | ierr = PetscObjectTypeCompare((PetscObject)cellIS,ISSTRIDE"stride",&isStride);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),265,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
266 | if (isStride) { | |||
267 | ierr = ISStrideGetInfo(cellIS,&sfirst,&stride);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),267,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
268 | } else { | |||
269 | ierr = ISGetIndices(cellIS,&cells);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),269,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
270 | } | |||
271 | for (c = 0; c < nCells; c++) { | |||
272 | PetscInt cell = isStride ? (sfirst + c * stride) : cells[c]; | |||
273 | PetscInt rem = cell; | |||
274 | PetscInt ijk[3] = {0}; | |||
275 | void *cB, *cD, *cH; | |||
276 | ||||
277 | if (datatype == PETSC_SCALARPETSC_DOUBLE) { | |||
278 | cB = B ? &((PetscScalar *)B)[nc * nq * c] : NULL((void*)0); | |||
279 | cD = D ? &((PetscScalar *)D)[nc * nq * dim * c] : NULL((void*)0); | |||
280 | cH = H ? &((PetscScalar *)H)[nc * nq * dim * dim * c] : NULL((void*)0); | |||
281 | } else { | |||
282 | cB = B ? &((PetscReal *)B)[nc * nq * c] : NULL((void*)0); | |||
283 | cD = D ? &((PetscReal *)D)[nc * nq * dim * c] : NULL((void*)0); | |||
284 | cH = H ? &((PetscReal *)H)[nc * nq * dim * dim * c] : NULL((void*)0); | |||
285 | } | |||
286 | if (cell < cStart || cell >= cEnd) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Point %D not a cell [%D,%D), not implemented yet",cell,cStart,cEnd)return PetscError(((MPI_Comm)0x44000001),286,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,63,PETSC_ERROR_INITIAL,"Point %D not a cell [%D,%D), not implemented yet" ,cell,cStart,cEnd); | |||
287 | for (i = 0; i < nc * whol; i++) {work[i] = dafield->cornerCoeffs[i];} | |||
288 | for (j = 0; j < dim; j++) { | |||
289 | PetscReal e, d; | |||
290 | ijk[j] = (rem % cellsPer[j]); | |||
291 | rem /= cellsPer[j]; | |||
292 | ||||
293 | e = 2. * (ijk[j] + first[j] + 0.5) * stepPer[j] - 1.; | |||
294 | d = stepPer[j]; | |||
295 | for (i = 0; i < half; i++) { | |||
296 | for (k = 0; k < nc; k++) { | |||
297 | cellCoeffs[ i * nc + k] = work[ 2 * i * nc + k] * d; | |||
298 | cellCoeffs[(i + half) * nc + k] = work[ 2 * i * nc + k] * e + work[(2 * i + 1) * nc + k]; | |||
299 | } | |||
300 | } | |||
301 | for (i = 0; i < whol * nc; i++) {work[i] = cellCoeffs[i];} | |||
302 | } | |||
303 | MultilinearEvaluate(dim,cellCoordRange,nc,cellCoeffs,dafield->work,nq,qs,datatype,cB,cD,cH); | |||
304 | } | |||
305 | if (!isStride) { | |||
306 | ierr = ISRestoreIndices(cellIS,&cells);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),306,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
307 | } | |||
308 | ierr = DMRestoreWorkArray(dm,(1 << dim) * nc,MPIU_SCALAR((MPI_Datatype)0x4c00080b),&cellCoeffs);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),308,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
309 | #if defined(PETSC_USE_COMPLEX) | |||
310 | ierr = DMRestoreWorkArray(dm,nq * dim,MPIU_SCALAR((MPI_Datatype)0x4c00080b),&qs);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),310,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
311 | #endif | |||
312 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
313 | } | |||
314 | ||||
315 | static PetscErrorCode DMFieldEvaluateFV_DA(DMField field, IS cellIS, PetscDataType datatype, void *B, void *D, void *H) | |||
316 | { | |||
317 | PetscInt c, i, dim, cellsPer[3] = {0}, first[3] = {0}; | |||
| ||||
318 | PetscReal stepPer[3] = {0.}; | |||
319 | DM dm; | |||
320 | DMDALocalInfo info; | |||
321 | PetscInt cStart, cEnd, numCells; | |||
322 | PetscInt nc; | |||
323 | PetscScalar *points; | |||
324 | DMField_DA *dafield; | |||
325 | PetscBool isStride; | |||
326 | const PetscInt *cells = NULL((void*)0); | |||
327 | PetscInt sfirst = -1, stride = -1; | |||
328 | PetscErrorCode ierr; | |||
329 | ||||
330 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 330; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
331 | dafield = (DMField_DA *) field->data; | |||
332 | dm = field->dm; | |||
333 | nc = field->numComponents; | |||
334 | ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),334,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
335 | dim = info.dim; | |||
336 | stepPer[0] = 1./ info.mx; | |||
337 | stepPer[1] = 1./ info.my; | |||
338 | stepPer[2] = 1./ info.mz; | |||
339 | first[0] = info.gxs; | |||
340 | first[1] = info.gys; | |||
341 | first[2] = info.gzs; | |||
342 | cellsPer[0] = info.gxm; | |||
343 | cellsPer[1] = info.gym; | |||
344 | cellsPer[2] = info.gzm; | |||
345 | ierr = DMDAGetHeightStratum(dm,0,&cStart,&cEnd);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),345,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
346 | ierr = ISGetLocalSize(cellIS,&numCells);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),346,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
347 | ierr = DMGetWorkArray(dm,dim * numCells,MPIU_SCALAR((MPI_Datatype)0x4c00080b),&points);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),347,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
348 | ierr = PetscObjectTypeCompare((PetscObject)cellIS,ISSTRIDE"stride",&isStride);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),348,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
349 | if (isStride) { | |||
350 | ierr = ISStrideGetInfo(cellIS,&sfirst,&stride);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),350,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
351 | } else { | |||
352 | ierr = ISGetIndices(cellIS,&cells);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),352,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
353 | } | |||
354 | for (c = 0; c < numCells; c++) { | |||
355 | PetscInt cell = isStride ? (sfirst + c * stride) : cells[c]; | |||
356 | PetscInt rem = cell; | |||
357 | PetscInt ijk[3] = {0}; | |||
358 | ||||
359 | if (cell < cStart || cell >= cEnd) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Point %D not a cell [%D,%D), not implemented yet",cell,cStart,cEnd)return PetscError(((MPI_Comm)0x44000001),359,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,63,PETSC_ERROR_INITIAL,"Point %D not a cell [%D,%D), not implemented yet" ,cell,cStart,cEnd); | |||
360 | for (i = 0; i < dim; i++) { | |||
361 | ijk[i] = (rem % cellsPer[i]); | |||
| ||||
362 | rem /= cellsPer[i]; | |||
363 | points[dim * c + i] = (ijk[i] + first[i] + 0.5) * stepPer[i]; | |||
364 | } | |||
365 | } | |||
366 | if (!isStride) { | |||
367 | ierr = ISRestoreIndices(cellIS,&cells);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),367,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
368 | } | |||
369 | MultilinearEvaluate(dim,dafield->coordRange,nc,dafield->cornerCoeffs,dafield->work,numCells,points,datatype,B,D,H); | |||
370 | ierr = DMRestoreWorkArray(dm,dim * numCells,MPIU_SCALAR((MPI_Datatype)0x4c00080b),&points);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),370,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
371 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
372 | } | |||
373 | ||||
374 | static PetscErrorCode DMFieldGetDegree_DA(DMField field, IS pointIS, PetscInt *minDegree, PetscInt *maxDegree) | |||
375 | { | |||
376 | DM dm; | |||
377 | PetscInt dim, h, imin; | |||
378 | PetscErrorCode ierr; | |||
379 | ||||
380 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 380; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
381 | dm = field->dm; | |||
382 | ierr = ISGetMinMax(pointIS,&imin,NULL((void*)0));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),382,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
383 | ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),383,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
384 | for (h = 0; h <= dim; h++) { | |||
385 | PetscInt hEnd; | |||
386 | ||||
387 | ierr = DMDAGetHeightStratum(dm,h,NULL((void*)0),&hEnd);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),387,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
388 | if (imin < hEnd) break; | |||
389 | } | |||
390 | dim -= h; | |||
391 | if (minDegree) *minDegree = 1; | |||
392 | if (maxDegree) *maxDegree = dim; | |||
393 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
394 | } | |||
395 | ||||
396 | static PetscErrorCode DMFieldCreateDefaultQuadrature_DA(DMField field, IS cellIS, PetscQuadrature *quad) | |||
397 | { | |||
398 | PetscInt h, dim, imax, imin; | |||
399 | DM dm; | |||
400 | PetscErrorCode ierr; | |||
401 | ||||
402 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 402; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
403 | dm = field->dm; | |||
404 | ierr = ISGetMinMax(cellIS,&imin,&imax);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),404,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
405 | ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),405,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
406 | *quad = NULL((void*)0); | |||
407 | for (h = 0; h <= dim; h++) { | |||
408 | PetscInt hStart, hEnd; | |||
409 | ||||
410 | ierr = DMDAGetHeightStratum(dm,h,&hStart,&hEnd);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),410,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
411 | if (imin >= hStart && imax < hEnd) break; | |||
412 | } | |||
413 | dim -= h; | |||
414 | if (dim > 0) { | |||
415 | ierr = PetscDTGaussTensorQuadrature(dim, 1, 1, -1.0, 1.0, quad);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),415,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
416 | } | |||
417 | ||||
418 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
419 | } | |||
420 | ||||
421 | static PetscErrorCode DMFieldInitialize_DA(DMField field) | |||
422 | { | |||
423 | DM dm; | |||
424 | Vec coords = NULL((void*)0); | |||
425 | PetscInt dim, i, j, k; | |||
426 | DMField_DA *dafield = (DMField_DA *) field->data; | |||
427 | PetscErrorCode ierr; | |||
428 | ||||
429 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 429; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
430 | field->ops->destroy = DMFieldDestroy_DA; | |||
431 | field->ops->evaluate = DMFieldEvaluate_DA; | |||
432 | field->ops->evaluateFE = DMFieldEvaluateFE_DA; | |||
433 | field->ops->evaluateFV = DMFieldEvaluateFV_DA; | |||
434 | field->ops->getDegree = DMFieldGetDegree_DA; | |||
435 | field->ops->createDefaultQuadrature = DMFieldCreateDefaultQuadrature_DA; | |||
436 | field->ops->view = DMFieldView_DA; | |||
437 | dm = field->dm; | |||
438 | ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),438,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
439 | if (dm->coordinates) coords = dm->coordinates; | |||
440 | else if (dm->coordinatesLocal) coords = dm->coordinatesLocal; | |||
441 | if (coords) { | |||
442 | PetscInt n; | |||
443 | const PetscScalar *array; | |||
444 | PetscReal mins[3][2] = {{PETSC_MAX_REAL1.7976931348623157e+308,PETSC_MAX_REAL1.7976931348623157e+308},{PETSC_MAX_REAL1.7976931348623157e+308,PETSC_MAX_REAL1.7976931348623157e+308},{PETSC_MAX_REAL1.7976931348623157e+308,PETSC_MAX_REAL1.7976931348623157e+308}}; | |||
445 | ||||
446 | ierr = VecGetLocalSize(coords,&n);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),446,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
447 | n /= dim; | |||
448 | ierr = VecGetArrayRead(coords,&array);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),448,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
449 | for (i = 0, k = 0; i < n; i++) { | |||
450 | for (j = 0; j < dim; j++, k++) { | |||
451 | PetscReal val = PetscRealPart(array[k])(array[k]); | |||
452 | ||||
453 | mins[j][0] = PetscMin(mins[j][0],val)(((mins[j][0])<(val)) ? (mins[j][0]) : (val)); | |||
454 | mins[j][1] = PetscMin(mins[j][1],-val)(((mins[j][1])<(-val)) ? (mins[j][1]) : (-val)); | |||
455 | } | |||
456 | } | |||
457 | ierr = VecRestoreArrayRead(coords,&array);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),457,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
458 | ierr = MPIU_Allreduce((PetscReal *) mins,&(dafield->coordRange[0][0]),2*dim,MPIU_REAL,MPI_MIN,PetscObjectComm((PetscObject)dm))(PetscAllreduceBarrierCheck(PetscObjectComm((PetscObject)dm), 2*dim,458,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ) || ((petsc_allreduce_ct += PetscMPIParallelComm((PetscObjectComm ((PetscObject)dm))),0) || MPI_Allreduce(((PetscReal *) mins), (&(dafield->coordRange[0][0])),(2*dim),(((MPI_Datatype )0x4c00080b)),((MPI_Op)(0x58000002)),(PetscObjectComm((PetscObject )dm)))));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),458,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
459 | for (j = 0; j < dim; j++) { | |||
460 | dafield->coordRange[j][1] = -dafield->coordRange[j][1]; | |||
461 | } | |||
462 | } else { | |||
463 | for (j = 0; j < dim; j++) { | |||
464 | dafield->coordRange[j][0] = 0.; | |||
465 | dafield->coordRange[j][1] = 1.; | |||
466 | } | |||
467 | } | |||
468 | for (j = 0; j < dim; j++) { | |||
469 | PetscReal avg = 0.5 * (dafield->coordRange[j][1] + dafield->coordRange[j][0]); | |||
470 | PetscReal dif = 0.5 * (dafield->coordRange[j][1] - dafield->coordRange[j][0]); | |||
471 | ||||
472 | dafield->coordRange[j][0] = avg; | |||
473 | dafield->coordRange[j][1] = dif; | |||
474 | } | |||
475 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
476 | } | |||
477 | ||||
478 | PETSC_INTERNextern __attribute__((visibility ("hidden"))) PetscErrorCode DMFieldCreate_DA(DMField field) | |||
479 | { | |||
480 | DMField_DA *dafield; | |||
481 | PetscErrorCode ierr; | |||
482 | ||||
483 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 483; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
484 | ierr = PetscNewLog(field,&dafield)(PetscMallocA(1,PETSC_TRUE,484,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,(size_t)(1)*sizeof(**(((&dafield)))),(((&dafield)))) || PetscLogObjectMemory((PetscObject)field,sizeof(**(&dafield ))));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),484,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
485 | field->data = dafield; | |||
486 | ierr = DMFieldInitialize_DA(field);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),486,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
487 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
488 | } | |||
489 | ||||
490 | PetscErrorCode DMFieldCreateDA(DM dm, PetscInt nc, const PetscScalar *cornerValues,DMField *field) | |||
491 | { | |||
492 | DMField b; | |||
493 | DMField_DA *dafield; | |||
494 | PetscInt dim, nv, i, j, k; | |||
495 | PetscInt half; | |||
496 | PetscScalar *cv, *cf, *work; | |||
497 | PetscErrorCode ierr; | |||
498 | ||||
499 | PetscFunctionBegindo { do { ; if (petscstack && (petscstack->currentsize < 64)) { petscstack->function[petscstack->currentsize ] = __func__; petscstack->file[petscstack->currentsize] = "/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ; petscstack->line[petscstack->currentsize] = 499; petscstack ->petscroutine[petscstack->currentsize] = PETSC_TRUE; petscstack ->currentsize++; } if (petscstack) { petscstack->hotdepth += (PETSC_FALSE || petscstack->hotdepth); } ; } while (0) ; ; } while (0); | |||
500 | ierr = DMFieldCreate(dm,nc,DMFIELD_VERTEX,&b);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),500,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
501 | ierr = DMFieldSetType(b,DMFIELDDA"da");CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),501,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
502 | dafield = (DMField_DA *) b->data; | |||
503 | ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),503,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
504 | nv = (1 << dim) * nc; | |||
505 | ierr = PetscMalloc3(nv,&cv,nv,&cf,nv,&work)PetscMallocA(3,PETSC_FALSE,505,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,(size_t)(nv)*sizeof(**(&cv)),(&cv),(size_t)(nv)*sizeof (**(&cf)),(&cf),(size_t)(nv)*sizeof(**(&work)),(& work));CHKERRQ(ierr)do {if (__builtin_expect(!!(ierr),0)) return PetscError(((MPI_Comm )0x44000001),505,__func__,"/sandbox/petsc/petsc.master/src/dm/field/impls/da/dmfieldda.c" ,ierr,PETSC_ERROR_REPEAT," ");} while (0); | |||
506 | for (i = 0; i < nv; i++) cv[i] = cornerValues[i]; | |||
507 | for (i = 0; i < nv; i++) cf[i] = cv[i]; | |||
508 | dafield->cornerVals = cv; | |||
509 | dafield->cornerCoeffs = cf; | |||
510 | dafield->work = work; | |||
511 | half = (1 << (dim - 1)); | |||
512 | for (i = 0; i < dim; i++) { | |||
513 | PetscScalar *w; | |||
514 | ||||
515 | w = work; | |||
516 | for (j = 0; j < half; j++) { | |||
517 | for (k = 0; k < nc; k++) { | |||
518 | w[j * nc + k] = 0.5 * (cf[(2 * j + 1) * nc + k] - cf[(2 * j) * nc + k]); | |||
519 | } | |||
520 | } | |||
521 | w = &work[j * nc]; | |||
522 | for (j = 0; j < half; j++) { | |||
523 | for (k = 0; k < nc; k++) { | |||
524 | w[j * nc + k] = 0.5 * (cf[(2 * j) * nc + k] + cf[(2 * j + 1) * nc + k]); | |||
525 | } | |||
526 | } | |||
527 | for (j = 0; j < nv; j++) {cf[j] = work[j];} | |||
528 | } | |||
529 | *field = b; | |||
530 | PetscFunctionReturn(0)do { do { ; if (petscstack && petscstack->currentsize > 0) { petscstack->currentsize--; petscstack->function [petscstack->currentsize] = 0; petscstack->file[petscstack ->currentsize] = 0; petscstack->line[petscstack->currentsize ] = 0; petscstack->petscroutine[petscstack->currentsize ] = PETSC_FALSE; } if (petscstack) { petscstack->hotdepth = (((petscstack->hotdepth-1)<(0)) ? (0) : (petscstack-> hotdepth-1)); } ; } while (0); return(0);} while (0); | |||
531 | } |