Actual source code: ex56.c
petsc-3.6.4 2016-04-12
1: static char help[] = "3D, tri-linear quadrilateral (Q1), displacement finite element formulation\n\
2: of linear elasticity. E=1.0, nu=0.25.\n\
3: Unit square domain with Dirichelet boundary condition on the y=0 side only.\n\
4: Load of 1.0 in x + 2y direction on all nodes (not a true uniform load).\n\
5: -ne <size> : number of (square) quadrilateral elements in each dimension\n\
6: -alpha <v> : scaling of material coeficient in embedded circle\n\n";
8: #include <petscksp.h>
10: static PetscBool log_stages = PETSC_TRUE;
11: static PetscErrorCode MaybeLogStagePush(PetscLogStage stage) { return log_stages ? PetscLogStagePush(stage) : 0; }
12: static PetscErrorCode MaybeLogStagePop() { return log_stages ? PetscLogStagePop() : 0; }
13: PetscErrorCode elem_3d_elast_v_25(PetscScalar *);
17: int main(int argc,char **args)
18: {
19: Mat Amat;
21: PetscInt m,nn,M,Istart,Iend,i,j,k,ii,jj,kk,ic,ne=4,id;
22: PetscReal x,y,z,h,*coords,soft_alpha=1.e-3;
23: PetscBool two_solves=PETSC_FALSE,test_nonzero_cols=PETSC_FALSE,use_nearnullspace=PETSC_FALSE;
24: Vec xx,bb;
25: KSP ksp;
26: MPI_Comm comm;
27: PetscMPIInt npe,mype;
28: PC pc;
29: PetscScalar DD[24][24],DD2[24][24];
30: PetscLogStage stage[6];
31: PetscScalar DD1[24][24];
32: PCType type;
34: PetscInitialize(&argc,&args,(char*)0,help);
35: comm = PETSC_COMM_WORLD;
36: MPI_Comm_rank(comm, &mype);
37: MPI_Comm_size(comm, &npe);
39: PetscOptionsBegin(comm,NULL,"3D bilinear Q1 elasticity options","");
40: {
41: char nestring[256];
42: PetscSNPrintf(nestring,sizeof nestring,"number of elements in each direction, ne+1 must be a multiple of %D (sizes^{1/3})",(PetscInt)(PetscPowReal((PetscReal)npe,1./3.) + .5));
43: PetscOptionsInt("-ne",nestring,"",ne,&ne,NULL);
44: PetscOptionsBool("-log_stages","Log stages of solve separately","",log_stages,&log_stages,NULL);
45: PetscOptionsReal("-alpha","material coefficient inside circle","",soft_alpha,&soft_alpha,NULL);
46: PetscOptionsBool("-two_solves","solve additional variant of the problem","",two_solves,&two_solves,NULL);
47: PetscOptionsBool("-test_nonzero_cols","nonzero test","",test_nonzero_cols,&test_nonzero_cols,NULL);
48: PetscOptionsBool("-use_mat_nearnullspace","MatNearNullSpace API test","",use_nearnullspace,&use_nearnullspace,NULL);
49: }
50: PetscOptionsEnd();
52: if (log_stages) {
53: PetscLogStageRegister("Setup", &stage[0]);
54: PetscLogStageRegister("Solve", &stage[1]);
55: PetscLogStageRegister("2nd Setup", &stage[2]);
56: PetscLogStageRegister("2nd Solve", &stage[3]);
57: PetscLogStageRegister("3rd Setup", &stage[4]);
58: PetscLogStageRegister("3rd Solve", &stage[5]);
59: } else {
60: for (i=0; i<(PetscInt)(sizeof(stage)/sizeof(stage[0])); i++) stage[i] = -1;
61: }
63: h = 1./ne; nn = ne+1;
64: /* ne*ne; number of global elements */
65: M = 3*nn*nn*nn; /* global number of equations */
66: if (npe==2) {
67: if (mype==1) m=0;
68: else m = nn*nn*nn;
69: npe = 1;
70: } else {
71: m = nn*nn*nn/npe;
72: if (mype==npe-1) m = nn*nn*nn - (npe-1)*m;
73: }
74: m *= 3; /* number of equations local*/
75: /* Setup solver, get PC type and pc */
76: KSPCreate(PETSC_COMM_WORLD,&ksp);
77: KSPSetType(ksp, KSPCG);
78: KSPSetComputeSingularValues(ksp, PETSC_TRUE);
79: KSPGetPC(ksp, &pc);
80: PCSetType(pc, PCGAMG); /* default */
81: KSPSetFromOptions(ksp);
82: PCGetType(pc, &type);
84: {
85: /* configureation */
86: const PetscInt NP = (PetscInt)(PetscPowReal((PetscReal)npe,1./3.) + .5);
87: const PetscInt ipx = mype%NP, ipy = (mype%(NP*NP))/NP, ipz = mype/(NP*NP);
88: const PetscInt Ni0 = ipx*(nn/NP), Nj0 = ipy*(nn/NP), Nk0 = ipz*(nn/NP);
89: const PetscInt Ni1 = Ni0 + (m>0 ? (nn/NP) : 0), Nj1 = Nj0 + (nn/NP), Nk1 = Nk0 + (nn/NP);
90: const PetscInt NN = nn/NP, id0 = ipz*nn*nn*NN + ipy*nn*NN*NN + ipx*NN*NN*NN;
91: PetscInt *d_nnz, *o_nnz,osz[4]={0,9,15,19},nbc;
92: PetscScalar vv[24], v2[24];
93: if (npe!=NP*NP*NP) SETERRQ1(comm,PETSC_ERR_ARG_WRONG, "npe=%d: npe^{1/3} must be integer",npe);
94: if (nn!=NP*(nn/NP)) SETERRQ1(comm,PETSC_ERR_ARG_WRONG, "-ne %d: (ne+1)%(npe^{1/3}) must equal zero",ne);
96: /* count nnz */
97: PetscMalloc1(m+1, &d_nnz);
98: PetscMalloc1(m+1, &o_nnz);
99: for (i=Ni0,ic=0; i<Ni1; i++) {
100: for (j=Nj0; j<Nj1; j++) {
101: for (k=Nk0; k<Nk1; k++) {
102: nbc = 0;
103: if (i==Ni0 || i==Ni1-1) nbc++;
104: if (j==Nj0 || j==Nj1-1) nbc++;
105: if (k==Nk0 || k==Nk1-1) nbc++;
106: for (jj=0; jj<3; jj++,ic++) {
107: d_nnz[ic] = 3*(27-osz[nbc]);
108: o_nnz[ic] = 3*osz[nbc];
109: }
110: }
111: }
112: }
113: if (ic != m) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"ic %D does not equal m %D",ic,m);
115: /* create stiffness matrix */
116: MatCreate(comm,&Amat);
117: MatSetSizes(Amat,m,m,M,M);
118: MatSetBlockSize(Amat,3);
119: MatSetType(Amat,MATAIJ);
120: MatSeqAIJSetPreallocation(Amat,0,d_nnz);
121: MatMPIAIJSetPreallocation(Amat,0,d_nnz,0,o_nnz);
123: PetscFree(d_nnz);
124: PetscFree(o_nnz);
126: MatGetOwnershipRange(Amat,&Istart,&Iend);
128: if (m != Iend - Istart) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"m %D does not equal Iend %D - Istart %D",m,Iend,Istart);
129: /* Generate vectors */
130: VecCreate(comm,&xx);
131: VecSetSizes(xx,m,M);
132: VecSetBlockSize(xx,3);
133: VecSetFromOptions(xx);
134: VecDuplicate(xx,&bb);
135: VecSet(bb,.0);
136: /* generate element matrices */
137: {
138: PetscBool hasData = PETSC_TRUE;
139: if (!hasData) {
140: PetscPrintf(PETSC_COMM_WORLD,"\t No data is provided\n");
141: for (i=0; i<24; i++) {
142: for (j=0; j<24; j++) {
143: if (i==j) DD1[i][j] = 1.0;
144: else DD1[i][j] = -.25;
145: }
146: }
147: } else {
148: /* Get array data */
149: elem_3d_elast_v_25((PetscScalar*)DD1);
150: }
151:
152: /* BC version of element */
153: for (i=0; i<24; i++) {
154: for (j=0; j<24; j++) {
155: if (i<12 || (j < 12 && !test_nonzero_cols)) {
156: if (i==j) DD2[i][j] = 0.1*DD1[i][j];
157: else DD2[i][j] = 0.0;
158: } else DD2[i][j] = DD1[i][j];
159: }
160: }
161: /* element residual/load vector */
162: for (i=0; i<24; i++) {
163: if (i%3==0) vv[i] = h*h;
164: else if (i%3==1) vv[i] = 2.0*h*h;
165: else vv[i] = .0;
166: }
167: for (i=0; i<24; i++) {
168: if (i%3==0 && i>=12) v2[i] = h*h;
169: else if (i%3==1 && i>=12) v2[i] = 2.0*h*h;
170: else v2[i] = .0;
171: }
172: }
174: PetscMalloc1(m+1, &coords);
175: coords[m] = -99.0;
177: /* forms the element stiffness and coordinates */
178: for (i=Ni0,ic=0,ii=0; i<Ni1; i++,ii++) {
179: for (j=Nj0,jj=0; j<Nj1; j++,jj++) {
180: for (k=Nk0,kk=0; k<Nk1; k++,kk++,ic++) {
182: /* coords */
183: x = coords[3*ic] = h*(PetscReal)i;
184: y = coords[3*ic+1] = h*(PetscReal)j;
185: z = coords[3*ic+2] = h*(PetscReal)k;
186: /* matrix */
187: id = id0 + ii + NN*jj + NN*NN*kk;
189: if (i<ne && j<ne && k<ne) {
190: /* radius */
191: PetscReal radius = PetscSqrtReal((x-.5+h/2)*(x-.5+h/2)+(y-.5+h/2)*(y-.5+h/2)+(z-.5+h/2)*(z-.5+h/2));
192: PetscReal alpha = 1.0;
193: PetscInt jx,ix,idx[8] = { id, id+1, id+NN+1, id+NN,
194: id + NN*NN, id+1 + NN*NN,
195: id+NN+1 + NN*NN, id+NN + NN*NN };
197: /* correct indices */
198: if (i==Ni1-1 && Ni1!=nn) {
199: idx[1] += NN*(NN*NN-1);
200: idx[2] += NN*(NN*NN-1);
201: idx[5] += NN*(NN*NN-1);
202: idx[6] += NN*(NN*NN-1);
203: }
204: if (j==Nj1-1 && Nj1!=nn) {
205: idx[2] += NN*NN*(nn-1);
206: idx[3] += NN*NN*(nn-1);
207: idx[6] += NN*NN*(nn-1);
208: idx[7] += NN*NN*(nn-1);
209: }
210: if (k==Nk1-1 && Nk1!=nn) {
211: idx[4] += NN*(nn*nn-NN*NN);
212: idx[5] += NN*(nn*nn-NN*NN);
213: idx[6] += NN*(nn*nn-NN*NN);
214: idx[7] += NN*(nn*nn-NN*NN);
215: }
217: if (radius < 0.25) alpha = soft_alpha;
219: for (ix=0; ix<24; ix++) {
220: for (jx=0;jx<24;jx++) DD[ix][jx] = alpha*DD1[ix][jx];
221: }
222: if (k>0) {
223: MatSetValuesBlocked(Amat,8,idx,8,idx,(const PetscScalar*)DD,ADD_VALUES);
224: VecSetValuesBlocked(bb,8,idx,(const PetscScalar*)vv,ADD_VALUES);
225: } else {
226: /* a BC */
227: for (ix=0;ix<24;ix++) {
228: for (jx=0;jx<24;jx++) DD[ix][jx] = alpha*DD2[ix][jx];
229: }
230: MatSetValuesBlocked(Amat,8,idx,8,idx,(const PetscScalar*)DD,ADD_VALUES);
231: VecSetValuesBlocked(bb,8,idx,(const PetscScalar*)v2,ADD_VALUES);
232: }
233: }
234: }
235: }
237: }
238: MatAssemblyBegin(Amat,MAT_FINAL_ASSEMBLY);
239: MatAssemblyEnd(Amat,MAT_FINAL_ASSEMBLY);
240: VecAssemblyBegin(bb);
241: VecAssemblyEnd(bb);
242: }
244: if (!PETSC_TRUE) {
245: PetscViewer viewer;
246: PetscViewerASCIIOpen(comm, "Amat.m", &viewer);
247: PetscViewerSetFormat(viewer, PETSC_VIEWER_ASCII_MATLAB);
248: MatView(Amat,viewer);
249: PetscViewerDestroy(&viewer);
250: }
252: /* finish KSP/PC setup */
253: KSPSetOperators(ksp, Amat, Amat);
254: if (use_nearnullspace) {
255: MatNullSpace matnull;
256: Vec vec_coords;
257: PetscScalar *c;
259: VecCreate(MPI_COMM_WORLD,&vec_coords);
260: VecSetBlockSize(vec_coords,3);
261: VecSetSizes(vec_coords,m,PETSC_DECIDE);
262: VecSetUp(vec_coords);
263: VecGetArray(vec_coords,&c);
264: for (i=0; i<m; i++) c[i] = coords[i]; /* Copy since Scalar type might be Complex */
265: VecRestoreArray(vec_coords,&c);
266: MatNullSpaceCreateRigidBody(vec_coords,&matnull);
267: MatSetNearNullSpace(Amat,matnull);
268: MatNullSpaceDestroy(&matnull);
269: VecDestroy(&vec_coords);
270: } else {
271: PCSetCoordinates(pc, 3, m/3, coords);
272: }
274: MaybeLogStagePush(stage[0]);
276: /* PC setup basically */
277: KSPSetUp(ksp);
279: MaybeLogStagePop();
280: MaybeLogStagePush(stage[1]);
282: /* test BCs */
283: if (test_nonzero_cols) {
284: VecZeroEntries(xx);
285: if (mype==0) VecSetValue(xx,0,1.0,INSERT_VALUES);
286: VecAssemblyBegin(xx);
287: VecAssemblyEnd(xx);
288: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
289: }
291: /* 1st solve */
292: KSPSolve(ksp, bb, xx);
294: MaybeLogStagePop();
296: /* 2nd solve */
297: if (two_solves) {
298: PetscReal emax, emin;
299: PetscReal norm,norm2;
300: Vec res;
302: MaybeLogStagePush(stage[2]);
303: /* PC setup basically */
304: MatScale(Amat, 100000.0);
305: KSPSetOperators(ksp, Amat, Amat);
306: KSPSetUp(ksp);
308: MaybeLogStagePop();
309: MaybeLogStagePush(stage[3]);
310: KSPSolve(ksp, bb, xx);
311: KSPComputeExtremeSingularValues(ksp, &emax, &emin);
313: MaybeLogStagePop();
314: MaybeLogStagePush(stage[4]);
316: /* 3rd solve */
317: MatScale(Amat, 100000.0);
318: KSPSetOperators(ksp, Amat, Amat);
319: KSPSetUp(ksp);
321: MaybeLogStagePop();
322: MaybeLogStagePush(stage[5]);
324: KSPSolve(ksp, bb, xx);
326: MaybeLogStagePop();
329: VecNorm(bb, NORM_2, &norm2);
331: VecDuplicate(xx, &res);
332: MatMult(Amat, xx, res);
333: VecAXPY(bb, -1.0, res);
334: VecDestroy(&res);
335: VecNorm(bb, NORM_2, &norm);
336: PetscPrintf(PETSC_COMM_WORLD,"[%d]%s |b-Ax|/|b|=%e, |b|=%e, emax=%e\n",0,__FUNCT__,(double)(norm/norm2),(double)norm2,(double)emax);
337: /*PetscViewerASCIIOpen(comm, "residual.m", &viewer);
338: PetscViewerSetFormat(viewer, PETSC_VIEWER_ASCII_MATLAB);
339: VecView(bb,viewer);
340: PetscViewerDestroy(&viewer);*/
343: /* PetscViewerASCIIOpen(comm, "rhs.m", &viewer); */
344: /* PetscViewerSetFormat(viewer, PETSC_VIEWER_ASCII_MATLAB); */
345: /* */
346: /* VecView(bb,viewer); */
347: /* PetscViewerDestroy(&viewer); */
349: /* PetscViewerASCIIOpen(comm, "solution.m", &viewer); */
350: /* PetscViewerSetFormat(viewer, PETSC_VIEWER_ASCII_MATLAB); */
351: /* */
352: /* VecView(xx, viewer); */
353: /* PetscViewerDestroy(&viewer); */
354: }
356: /* Free work space */
357: KSPDestroy(&ksp);
358: VecDestroy(&xx);
359: VecDestroy(&bb);
360: MatDestroy(&Amat);
361: PetscFree(coords);
363: PetscFinalize();
364: return 0;
365: }
367: /* Data was previously provided in the file data/elem_3d_elast_v_25.tx */
370: PetscErrorCode elem_3d_elast_v_25(PetscScalar *dd)
371: {
373: PetscScalar DD[] = {
374: 0.18981481481481474 ,
375: 5.27777777777777568E-002,
376: 5.27777777777777568E-002,
377: -5.64814814814814659E-002,
378: -1.38888888888889072E-002,
379: -1.38888888888889089E-002,
380: -8.24074074074073876E-002,
381: -5.27777777777777429E-002,
382: 1.38888888888888725E-002,
383: 4.90740740740740339E-002,
384: 1.38888888888889124E-002,
385: 4.72222222222222071E-002,
386: 4.90740740740740339E-002,
387: 4.72222222222221932E-002,
388: 1.38888888888888968E-002,
389: -8.24074074074073876E-002,
390: 1.38888888888888673E-002,
391: -5.27777777777777429E-002,
392: -7.87037037037036785E-002,
393: -4.72222222222221932E-002,
394: -4.72222222222222071E-002,
395: 1.20370370370370180E-002,
396: -1.38888888888888742E-002,
397: -1.38888888888888829E-002,
398: 5.27777777777777568E-002,
399: 0.18981481481481474 ,
400: 5.27777777777777568E-002,
401: 1.38888888888889124E-002,
402: 4.90740740740740269E-002,
403: 4.72222222222221932E-002,
404: -5.27777777777777637E-002,
405: -8.24074074074073876E-002,
406: 1.38888888888888725E-002,
407: -1.38888888888889037E-002,
408: -5.64814814814814728E-002,
409: -1.38888888888888985E-002,
410: 4.72222222222221932E-002,
411: 4.90740740740740478E-002,
412: 1.38888888888888968E-002,
413: -1.38888888888888673E-002,
414: 1.20370370370370058E-002,
415: -1.38888888888888742E-002,
416: -4.72222222222221932E-002,
417: -7.87037037037036785E-002,
418: -4.72222222222222002E-002,
419: 1.38888888888888742E-002,
420: -8.24074074074073598E-002,
421: -5.27777777777777568E-002,
422: 5.27777777777777568E-002,
423: 5.27777777777777568E-002,
424: 0.18981481481481474 ,
425: 1.38888888888889055E-002,
426: 4.72222222222222002E-002,
427: 4.90740740740740269E-002,
428: -1.38888888888888829E-002,
429: -1.38888888888888829E-002,
430: 1.20370370370370180E-002,
431: 4.72222222222222002E-002,
432: 1.38888888888888985E-002,
433: 4.90740740740740339E-002,
434: -1.38888888888888985E-002,
435: -1.38888888888888968E-002,
436: -5.64814814814814520E-002,
437: -5.27777777777777568E-002,
438: 1.38888888888888777E-002,
439: -8.24074074074073876E-002,
440: -4.72222222222222002E-002,
441: -4.72222222222221932E-002,
442: -7.87037037037036646E-002,
443: 1.38888888888888794E-002,
444: -5.27777777777777568E-002,
445: -8.24074074074073598E-002,
446: -5.64814814814814659E-002,
447: 1.38888888888889124E-002,
448: 1.38888888888889055E-002,
449: 0.18981481481481474 ,
450: -5.27777777777777568E-002,
451: -5.27777777777777499E-002,
452: 4.90740740740740269E-002,
453: -1.38888888888889072E-002,
454: -4.72222222222221932E-002,
455: -8.24074074074073876E-002,
456: 5.27777777777777568E-002,
457: -1.38888888888888812E-002,
458: -8.24074074074073876E-002,
459: -1.38888888888888742E-002,
460: 5.27777777777777499E-002,
461: 4.90740740740740269E-002,
462: -4.72222222222221863E-002,
463: -1.38888888888889089E-002,
464: 1.20370370370370162E-002,
465: 1.38888888888888673E-002,
466: 1.38888888888888742E-002,
467: -7.87037037037036785E-002,
468: 4.72222222222222002E-002,
469: 4.72222222222222071E-002,
470: -1.38888888888889072E-002,
471: 4.90740740740740269E-002,
472: 4.72222222222222002E-002,
473: -5.27777777777777568E-002,
474: 0.18981481481481480 ,
475: 5.27777777777777568E-002,
476: 1.38888888888889020E-002,
477: -5.64814814814814728E-002,
478: -1.38888888888888951E-002,
479: 5.27777777777777637E-002,
480: -8.24074074074073876E-002,
481: 1.38888888888888881E-002,
482: 1.38888888888888742E-002,
483: 1.20370370370370232E-002,
484: -1.38888888888888812E-002,
485: -4.72222222222221863E-002,
486: 4.90740740740740339E-002,
487: 1.38888888888888933E-002,
488: -1.38888888888888812E-002,
489: -8.24074074074073876E-002,
490: -5.27777777777777568E-002,
491: 4.72222222222222071E-002,
492: -7.87037037037036924E-002,
493: -4.72222222222222140E-002,
494: -1.38888888888889089E-002,
495: 4.72222222222221932E-002,
496: 4.90740740740740269E-002,
497: -5.27777777777777499E-002,
498: 5.27777777777777568E-002,
499: 0.18981481481481477 ,
500: -4.72222222222222071E-002,
501: 1.38888888888888968E-002,
502: 4.90740740740740131E-002,
503: 1.38888888888888812E-002,
504: -1.38888888888888708E-002,
505: 1.20370370370370267E-002,
506: 5.27777777777777568E-002,
507: 1.38888888888888812E-002,
508: -8.24074074074073876E-002,
509: 1.38888888888889124E-002,
510: -1.38888888888889055E-002,
511: -5.64814814814814589E-002,
512: -1.38888888888888812E-002,
513: -5.27777777777777568E-002,
514: -8.24074074074073737E-002,
515: 4.72222222222222002E-002,
516: -4.72222222222222002E-002,
517: -7.87037037037036924E-002,
518: -8.24074074074073876E-002,
519: -5.27777777777777637E-002,
520: -1.38888888888888829E-002,
521: 4.90740740740740269E-002,
522: 1.38888888888889020E-002,
523: -4.72222222222222071E-002,
524: 0.18981481481481480 ,
525: 5.27777777777777637E-002,
526: -5.27777777777777637E-002,
527: -5.64814814814814728E-002,
528: -1.38888888888889037E-002,
529: 1.38888888888888951E-002,
530: -7.87037037037036785E-002,
531: -4.72222222222222002E-002,
532: 4.72222222222221932E-002,
533: 1.20370370370370128E-002,
534: -1.38888888888888725E-002,
535: 1.38888888888888812E-002,
536: 4.90740740740740408E-002,
537: 4.72222222222222002E-002,
538: -1.38888888888888951E-002,
539: -8.24074074074073876E-002,
540: 1.38888888888888812E-002,
541: 5.27777777777777637E-002,
542: -5.27777777777777429E-002,
543: -8.24074074074073876E-002,
544: -1.38888888888888829E-002,
545: -1.38888888888889072E-002,
546: -5.64814814814814728E-002,
547: 1.38888888888888968E-002,
548: 5.27777777777777637E-002,
549: 0.18981481481481480 ,
550: -5.27777777777777568E-002,
551: 1.38888888888888916E-002,
552: 4.90740740740740339E-002,
553: -4.72222222222222210E-002,
554: -4.72222222222221932E-002,
555: -7.87037037037036924E-002,
556: 4.72222222222222002E-002,
557: 1.38888888888888742E-002,
558: -8.24074074074073876E-002,
559: 5.27777777777777429E-002,
560: 4.72222222222222002E-002,
561: 4.90740740740740269E-002,
562: -1.38888888888888951E-002,
563: -1.38888888888888846E-002,
564: 1.20370370370370267E-002,
565: 1.38888888888888916E-002,
566: 1.38888888888888725E-002,
567: 1.38888888888888725E-002,
568: 1.20370370370370180E-002,
569: -4.72222222222221932E-002,
570: -1.38888888888888951E-002,
571: 4.90740740740740131E-002,
572: -5.27777777777777637E-002,
573: -5.27777777777777568E-002,
574: 0.18981481481481480 ,
575: -1.38888888888888968E-002,
576: -4.72222222222221932E-002,
577: 4.90740740740740339E-002,
578: 4.72222222222221932E-002,
579: 4.72222222222222071E-002,
580: -7.87037037037036646E-002,
581: -1.38888888888888742E-002,
582: 5.27777777777777499E-002,
583: -8.24074074074073737E-002,
584: 1.38888888888888933E-002,
585: 1.38888888888889020E-002,
586: -5.64814814814814589E-002,
587: 5.27777777777777568E-002,
588: -1.38888888888888794E-002,
589: -8.24074074074073876E-002,
590: 4.90740740740740339E-002,
591: -1.38888888888889037E-002,
592: 4.72222222222222002E-002,
593: -8.24074074074073876E-002,
594: 5.27777777777777637E-002,
595: 1.38888888888888812E-002,
596: -5.64814814814814728E-002,
597: 1.38888888888888916E-002,
598: -1.38888888888888968E-002,
599: 0.18981481481481480 ,
600: -5.27777777777777499E-002,
601: 5.27777777777777707E-002,
602: 1.20370370370370180E-002,
603: 1.38888888888888812E-002,
604: -1.38888888888888812E-002,
605: -7.87037037037036785E-002,
606: 4.72222222222222002E-002,
607: -4.72222222222222071E-002,
608: -8.24074074074073876E-002,
609: -1.38888888888888742E-002,
610: -5.27777777777777568E-002,
611: 4.90740740740740616E-002,
612: -4.72222222222222002E-002,
613: 1.38888888888888846E-002,
614: 1.38888888888889124E-002,
615: -5.64814814814814728E-002,
616: 1.38888888888888985E-002,
617: 5.27777777777777568E-002,
618: -8.24074074074073876E-002,
619: -1.38888888888888708E-002,
620: -1.38888888888889037E-002,
621: 4.90740740740740339E-002,
622: -4.72222222222221932E-002,
623: -5.27777777777777499E-002,
624: 0.18981481481481480 ,
625: -5.27777777777777568E-002,
626: -1.38888888888888673E-002,
627: -8.24074074074073598E-002,
628: 5.27777777777777429E-002,
629: 4.72222222222222002E-002,
630: -7.87037037037036785E-002,
631: 4.72222222222222002E-002,
632: 1.38888888888888708E-002,
633: 1.20370370370370128E-002,
634: 1.38888888888888760E-002,
635: -4.72222222222222002E-002,
636: 4.90740740740740478E-002,
637: -1.38888888888888951E-002,
638: 4.72222222222222071E-002,
639: -1.38888888888888985E-002,
640: 4.90740740740740339E-002,
641: -1.38888888888888812E-002,
642: 1.38888888888888881E-002,
643: 1.20370370370370267E-002,
644: 1.38888888888888951E-002,
645: -4.72222222222222210E-002,
646: 4.90740740740740339E-002,
647: 5.27777777777777707E-002,
648: -5.27777777777777568E-002,
649: 0.18981481481481477 ,
650: 1.38888888888888829E-002,
651: 5.27777777777777707E-002,
652: -8.24074074074073598E-002,
653: -4.72222222222222140E-002,
654: 4.72222222222222140E-002,
655: -7.87037037037036646E-002,
656: -5.27777777777777707E-002,
657: -1.38888888888888829E-002,
658: -8.24074074074073876E-002,
659: -1.38888888888888881E-002,
660: 1.38888888888888881E-002,
661: -5.64814814814814589E-002,
662: 4.90740740740740339E-002,
663: 4.72222222222221932E-002,
664: -1.38888888888888985E-002,
665: -8.24074074074073876E-002,
666: 1.38888888888888742E-002,
667: 5.27777777777777568E-002,
668: -7.87037037037036785E-002,
669: -4.72222222222221932E-002,
670: 4.72222222222221932E-002,
671: 1.20370370370370180E-002,
672: -1.38888888888888673E-002,
673: 1.38888888888888829E-002,
674: 0.18981481481481469 ,
675: 5.27777777777777429E-002,
676: -5.27777777777777429E-002,
677: -5.64814814814814659E-002,
678: -1.38888888888889055E-002,
679: 1.38888888888889055E-002,
680: -8.24074074074074153E-002,
681: -5.27777777777777429E-002,
682: -1.38888888888888760E-002,
683: 4.90740740740740408E-002,
684: 1.38888888888888968E-002,
685: -4.72222222222222071E-002,
686: 4.72222222222221932E-002,
687: 4.90740740740740478E-002,
688: -1.38888888888888968E-002,
689: -1.38888888888888742E-002,
690: 1.20370370370370232E-002,
691: 1.38888888888888812E-002,
692: -4.72222222222222002E-002,
693: -7.87037037037036924E-002,
694: 4.72222222222222071E-002,
695: 1.38888888888888812E-002,
696: -8.24074074074073598E-002,
697: 5.27777777777777707E-002,
698: 5.27777777777777429E-002,
699: 0.18981481481481477 ,
700: -5.27777777777777499E-002,
701: 1.38888888888889107E-002,
702: 4.90740740740740478E-002,
703: -4.72222222222221932E-002,
704: -5.27777777777777568E-002,
705: -8.24074074074074153E-002,
706: -1.38888888888888812E-002,
707: -1.38888888888888846E-002,
708: -5.64814814814814659E-002,
709: 1.38888888888888812E-002,
710: 1.38888888888888968E-002,
711: 1.38888888888888968E-002,
712: -5.64814814814814520E-002,
713: 5.27777777777777499E-002,
714: -1.38888888888888812E-002,
715: -8.24074074074073876E-002,
716: 4.72222222222221932E-002,
717: 4.72222222222222002E-002,
718: -7.87037037037036646E-002,
719: -1.38888888888888812E-002,
720: 5.27777777777777429E-002,
721: -8.24074074074073598E-002,
722: -5.27777777777777429E-002,
723: -5.27777777777777499E-002,
724: 0.18981481481481474 ,
725: -1.38888888888888985E-002,
726: -4.72222222222221863E-002,
727: 4.90740740740740339E-002,
728: 1.38888888888888829E-002,
729: 1.38888888888888777E-002,
730: 1.20370370370370249E-002,
731: -4.72222222222222002E-002,
732: -1.38888888888888933E-002,
733: 4.90740740740740339E-002,
734: -8.24074074074073876E-002,
735: -1.38888888888888673E-002,
736: -5.27777777777777568E-002,
737: 4.90740740740740269E-002,
738: -4.72222222222221863E-002,
739: 1.38888888888889124E-002,
740: 1.20370370370370128E-002,
741: 1.38888888888888742E-002,
742: -1.38888888888888742E-002,
743: -7.87037037037036785E-002,
744: 4.72222222222222002E-002,
745: -4.72222222222222140E-002,
746: -5.64814814814814659E-002,
747: 1.38888888888889107E-002,
748: -1.38888888888888985E-002,
749: 0.18981481481481474 ,
750: -5.27777777777777499E-002,
751: 5.27777777777777499E-002,
752: 4.90740740740740339E-002,
753: -1.38888888888889055E-002,
754: 4.72222222222221932E-002,
755: -8.24074074074074153E-002,
756: 5.27777777777777499E-002,
757: 1.38888888888888829E-002,
758: 1.38888888888888673E-002,
759: 1.20370370370370058E-002,
760: 1.38888888888888777E-002,
761: -4.72222222222221863E-002,
762: 4.90740740740740339E-002,
763: -1.38888888888889055E-002,
764: -1.38888888888888725E-002,
765: -8.24074074074073876E-002,
766: 5.27777777777777499E-002,
767: 4.72222222222222002E-002,
768: -7.87037037037036785E-002,
769: 4.72222222222222140E-002,
770: -1.38888888888889055E-002,
771: 4.90740740740740478E-002,
772: -4.72222222222221863E-002,
773: -5.27777777777777499E-002,
774: 0.18981481481481469 ,
775: -5.27777777777777499E-002,
776: 1.38888888888889072E-002,
777: -5.64814814814814659E-002,
778: 1.38888888888889003E-002,
779: 5.27777777777777429E-002,
780: -8.24074074074074153E-002,
781: -1.38888888888888812E-002,
782: -5.27777777777777429E-002,
783: -1.38888888888888742E-002,
784: -8.24074074074073876E-002,
785: -1.38888888888889089E-002,
786: 1.38888888888888933E-002,
787: -5.64814814814814589E-002,
788: 1.38888888888888812E-002,
789: 5.27777777777777429E-002,
790: -8.24074074074073737E-002,
791: -4.72222222222222071E-002,
792: 4.72222222222222002E-002,
793: -7.87037037037036646E-002,
794: 1.38888888888889055E-002,
795: -4.72222222222221932E-002,
796: 4.90740740740740339E-002,
797: 5.27777777777777499E-002,
798: -5.27777777777777499E-002,
799: 0.18981481481481474 ,
800: 4.72222222222222002E-002,
801: -1.38888888888888985E-002,
802: 4.90740740740740339E-002,
803: -1.38888888888888846E-002,
804: 1.38888888888888812E-002,
805: 1.20370370370370284E-002,
806: -7.87037037037036785E-002,
807: -4.72222222222221932E-002,
808: -4.72222222222222002E-002,
809: 1.20370370370370162E-002,
810: -1.38888888888888812E-002,
811: -1.38888888888888812E-002,
812: 4.90740740740740408E-002,
813: 4.72222222222222002E-002,
814: 1.38888888888888933E-002,
815: -8.24074074074073876E-002,
816: 1.38888888888888708E-002,
817: -5.27777777777777707E-002,
818: -8.24074074074074153E-002,
819: -5.27777777777777568E-002,
820: 1.38888888888888829E-002,
821: 4.90740740740740339E-002,
822: 1.38888888888889072E-002,
823: 4.72222222222222002E-002,
824: 0.18981481481481477 ,
825: 5.27777777777777429E-002,
826: 5.27777777777777568E-002,
827: -5.64814814814814659E-002,
828: -1.38888888888888846E-002,
829: -1.38888888888888881E-002,
830: -4.72222222222221932E-002,
831: -7.87037037037036785E-002,
832: -4.72222222222221932E-002,
833: 1.38888888888888673E-002,
834: -8.24074074074073876E-002,
835: -5.27777777777777568E-002,
836: 4.72222222222222002E-002,
837: 4.90740740740740269E-002,
838: 1.38888888888889020E-002,
839: -1.38888888888888742E-002,
840: 1.20370370370370128E-002,
841: -1.38888888888888829E-002,
842: -5.27777777777777429E-002,
843: -8.24074074074074153E-002,
844: 1.38888888888888777E-002,
845: -1.38888888888889055E-002,
846: -5.64814814814814659E-002,
847: -1.38888888888888985E-002,
848: 5.27777777777777429E-002,
849: 0.18981481481481469 ,
850: 5.27777777777777429E-002,
851: 1.38888888888888933E-002,
852: 4.90740740740740339E-002,
853: 4.72222222222222071E-002,
854: -4.72222222222222071E-002,
855: -4.72222222222222002E-002,
856: -7.87037037037036646E-002,
857: 1.38888888888888742E-002,
858: -5.27777777777777568E-002,
859: -8.24074074074073737E-002,
860: -1.38888888888888951E-002,
861: -1.38888888888888951E-002,
862: -5.64814814814814589E-002,
863: -5.27777777777777568E-002,
864: 1.38888888888888760E-002,
865: -8.24074074074073876E-002,
866: -1.38888888888888760E-002,
867: -1.38888888888888812E-002,
868: 1.20370370370370249E-002,
869: 4.72222222222221932E-002,
870: 1.38888888888889003E-002,
871: 4.90740740740740339E-002,
872: 5.27777777777777568E-002,
873: 5.27777777777777429E-002,
874: 0.18981481481481474 ,
875: 1.38888888888888933E-002,
876: 4.72222222222222071E-002,
877: 4.90740740740740339E-002,
878: 1.20370370370370180E-002,
879: 1.38888888888888742E-002,
880: 1.38888888888888794E-002,
881: -7.87037037037036785E-002,
882: 4.72222222222222071E-002,
883: 4.72222222222222002E-002,
884: -8.24074074074073876E-002,
885: -1.38888888888888846E-002,
886: 5.27777777777777568E-002,
887: 4.90740740740740616E-002,
888: -4.72222222222222002E-002,
889: -1.38888888888888881E-002,
890: 4.90740740740740408E-002,
891: -1.38888888888888846E-002,
892: -4.72222222222222002E-002,
893: -8.24074074074074153E-002,
894: 5.27777777777777429E-002,
895: -1.38888888888888846E-002,
896: -5.64814814814814659E-002,
897: 1.38888888888888933E-002,
898: 1.38888888888888933E-002,
899: 0.18981481481481477 ,
900: -5.27777777777777568E-002,
901: -5.27777777777777637E-002,
902: -1.38888888888888742E-002,
903: -8.24074074074073598E-002,
904: -5.27777777777777568E-002,
905: 4.72222222222222002E-002,
906: -7.87037037037036924E-002,
907: -4.72222222222222002E-002,
908: 1.38888888888888812E-002,
909: 1.20370370370370267E-002,
910: -1.38888888888888794E-002,
911: -4.72222222222222002E-002,
912: 4.90740740740740478E-002,
913: 1.38888888888888881E-002,
914: 1.38888888888888968E-002,
915: -5.64814814814814659E-002,
916: -1.38888888888888933E-002,
917: 5.27777777777777499E-002,
918: -8.24074074074074153E-002,
919: 1.38888888888888812E-002,
920: -1.38888888888888846E-002,
921: 4.90740740740740339E-002,
922: 4.72222222222222071E-002,
923: -5.27777777777777568E-002,
924: 0.18981481481481477 ,
925: 5.27777777777777637E-002,
926: -1.38888888888888829E-002,
927: -5.27777777777777568E-002,
928: -8.24074074074073598E-002,
929: 4.72222222222222071E-002,
930: -4.72222222222222140E-002,
931: -7.87037037037036924E-002,
932: 5.27777777777777637E-002,
933: 1.38888888888888916E-002,
934: -8.24074074074073876E-002,
935: 1.38888888888888846E-002,
936: -1.38888888888888951E-002,
937: -5.64814814814814589E-002,
938: -4.72222222222222071E-002,
939: 1.38888888888888812E-002,
940: 4.90740740740740339E-002,
941: 1.38888888888888829E-002,
942: -1.38888888888888812E-002,
943: 1.20370370370370284E-002,
944: -1.38888888888888881E-002,
945: 4.72222222222222071E-002,
946: 4.90740740740740339E-002,
947: -5.27777777777777637E-002,
948: 5.27777777777777637E-002,
949: 0.18981481481481477 ,
950: };
953: PetscMemcpy(dd,DD,sizeof(PetscScalar)*576);
954: return(0);
955: }