Actual source code: ex56.c

petsc-3.8.4 2018-03-24
Report Typos and Errors
  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 *);

 15: int main(int argc,char **args)
 16: {
 17:   Mat            Amat;
 19:   PetscInt       m,nn,M,Istart,Iend,i,j,k,ii,jj,kk,ic,ne=4,id;
 20:   PetscReal      x,y,z,h,*coords,soft_alpha=1.e-3;
 21:   PetscBool      two_solves=PETSC_FALSE,test_nonzero_cols=PETSC_FALSE,use_nearnullspace=PETSC_FALSE,test_late_bs=PETSC_FALSE;
 22:   Vec            xx,bb;
 23:   KSP            ksp;
 24:   MPI_Comm       comm;
 25:   PetscMPIInt    npe,mype;
 26:   PetscScalar    DD[24][24],DD2[24][24];
 27:   PetscLogStage  stage[6];
 28:   PetscScalar    DD1[24][24];

 30:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
 31:   comm = PETSC_COMM_WORLD;
 32:   MPI_Comm_rank(comm, &mype);
 33:   MPI_Comm_size(comm, &npe);

 35:   PetscOptionsBegin(comm,NULL,"3D bilinear Q1 elasticity options","");
 36:   {
 37:     char nestring[256];
 38:     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));
 39:     PetscOptionsInt("-ne",nestring,"",ne,&ne,NULL);
 40:     PetscOptionsBool("-log_stages","Log stages of solve separately","",log_stages,&log_stages,NULL);
 41:     PetscOptionsReal("-alpha","material coefficient inside circle","",soft_alpha,&soft_alpha,NULL);
 42:     PetscOptionsBool("-two_solves","solve additional variant of the problem","",two_solves,&two_solves,NULL);
 43:     PetscOptionsBool("-test_nonzero_cols","nonzero test","",test_nonzero_cols,&test_nonzero_cols,NULL);
 44:     PetscOptionsBool("-use_mat_nearnullspace","MatNearNullSpace API test","",use_nearnullspace,&use_nearnullspace,NULL);
 45:     PetscOptionsBool("-test_late_bs","","",test_late_bs,&test_late_bs,NULL);
 46:   }
 47:   PetscOptionsEnd();

 49:   if (log_stages) {
 50:     PetscLogStageRegister("Setup", &stage[0]);
 51:     PetscLogStageRegister("Solve", &stage[1]);
 52:     PetscLogStageRegister("2nd Setup", &stage[2]);
 53:     PetscLogStageRegister("2nd Solve", &stage[3]);
 54:     PetscLogStageRegister("3rd Setup", &stage[4]);
 55:     PetscLogStageRegister("3rd Solve", &stage[5]);
 56:   } else {
 57:     for (i=0; i<(PetscInt)(sizeof(stage)/sizeof(stage[0])); i++) stage[i] = -1;
 58:   }

 60:   h = 1./ne; nn = ne+1;
 61:   /* ne*ne; number of global elements */
 62:   M = 3*nn*nn*nn; /* global number of equations */
 63:   if (npe==2) {
 64:     if (mype==1) m=0;
 65:     else m = nn*nn*nn;
 66:     npe = 1;
 67:   } else {
 68:     m = nn*nn*nn/npe;
 69:     if (mype==npe-1) m = nn*nn*nn - (npe-1)*m;
 70:   }
 71:   m *= 3; /* number of equations local*/
 72:   /* Setup solver */
 73:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 74:   KSPSetComputeSingularValues(ksp, PETSC_TRUE);
 75:   KSPSetFromOptions(ksp);

 77:   {
 78:     /* configuration */
 79:     const PetscInt NP = (PetscInt)(PetscPowReal((PetscReal)npe,1./3.) + .5);
 80:     const PetscInt ipx = mype%NP, ipy = (mype%(NP*NP))/NP, ipz = mype/(NP*NP);
 81:     const PetscInt Ni0 = ipx*(nn/NP), Nj0 = ipy*(nn/NP), Nk0 = ipz*(nn/NP);
 82:     const PetscInt Ni1 = Ni0 + (m>0 ? (nn/NP) : 0), Nj1 = Nj0 + (nn/NP), Nk1 = Nk0 + (nn/NP);
 83:     const PetscInt NN  = nn/NP, id0 = ipz*nn*nn*NN + ipy*nn*NN*NN + ipx*NN*NN*NN;
 84:     PetscInt       *d_nnz, *o_nnz,osz[4]={0,9,15,19},nbc;
 85:     PetscScalar    vv[24], v2[24];
 86:     if (npe!=NP*NP*NP) SETERRQ1(comm,PETSC_ERR_ARG_WRONG, "npe=%d: npe^{1/3} must be integer",npe);
 87:     if (nn!=NP*(nn/NP)) SETERRQ1(comm,PETSC_ERR_ARG_WRONG, "-ne %d: (ne+1)%(npe^{1/3}) must equal zero",ne);

 89:     /* count nnz */
 90:     PetscMalloc1(m+1, &d_nnz);
 91:     PetscMalloc1(m+1, &o_nnz);
 92:     for (i=Ni0,ic=0; i<Ni1; i++) {
 93:       for (j=Nj0; j<Nj1; j++) {
 94:         for (k=Nk0; k<Nk1; k++) {
 95:           nbc = 0;
 96:           if (i==Ni0 || i==Ni1-1) nbc++;
 97:           if (j==Nj0 || j==Nj1-1) nbc++;
 98:           if (k==Nk0 || k==Nk1-1) nbc++;
 99:           for (jj=0; jj<3; jj++,ic++) {
100:             d_nnz[ic] = 3*(27-osz[nbc]);
101:             o_nnz[ic] = 3*osz[nbc];
102:           }
103:         }
104:       }
105:     }
106:     if (ic != m) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"ic %D does not equal m %D",ic,m);

108:     /* create stiffness matrix */
109:     MatCreate(comm,&Amat);
110:     MatSetSizes(Amat,m,m,M,M);
111:     if (!test_late_bs) {
112:       MatSetBlockSize(Amat,3);
113:     }
114:     MatSetType(Amat,MATAIJ);
115:     MatSeqAIJSetPreallocation(Amat,0,d_nnz);
116:     MatMPIAIJSetPreallocation(Amat,0,d_nnz,0,o_nnz);

118:     PetscFree(d_nnz);
119:     PetscFree(o_nnz);

121:     MatGetOwnershipRange(Amat,&Istart,&Iend);

123:     if (m != Iend - Istart) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"m %D does not equal Iend %D - Istart %D",m,Iend,Istart);
124:     /* Generate vectors */
125:     VecCreate(comm,&xx);
126:     VecSetSizes(xx,m,M);
127:     if (!test_late_bs) {
128:       VecSetBlockSize(xx,3);
129:     }
130:     VecSetFromOptions(xx);
131:     VecDuplicate(xx,&bb);
132:     VecSet(bb,.0);
133:     /* generate element matrices */
134:     {
135:       PetscBool hasData = PETSC_TRUE;
136:       if (!hasData) {
137:         PetscPrintf(PETSC_COMM_WORLD,"\t No data is provided\n");
138:         for (i=0; i<24; i++) {
139:           for (j=0; j<24; j++) {
140:             if (i==j) DD1[i][j] = 1.0;
141:             else DD1[i][j] = -.25;
142:           }
143:         }
144:       } else {
145:         /* Get array data */
146:         elem_3d_elast_v_25((PetscScalar*)DD1);
147:       }

149:       /* BC version of element */
150:       for (i=0; i<24; i++) {
151:         for (j=0; j<24; j++) {
152:           if (i<12 || (j < 12 && !test_nonzero_cols)) {
153:             if (i==j) DD2[i][j] = 0.1*DD1[i][j];
154:             else DD2[i][j] = 0.0;
155:           } else DD2[i][j] = DD1[i][j];
156:         }
157:       }
158:       /* element residual/load vector */
159:       for (i=0; i<24; i++) {
160:         if (i%3==0) vv[i] = h*h;
161:         else if (i%3==1) vv[i] = 2.0*h*h;
162:         else vv[i] = .0;
163:       }
164:       for (i=0; i<24; i++) {
165:         if (i%3==0 && i>=12) v2[i] = h*h;
166:         else if (i%3==1 && i>=12) v2[i] = 2.0*h*h;
167:         else v2[i] = .0;
168:       }
169:     }

171:     PetscMalloc1(m+1, &coords);
172:     coords[m] = -99.0;

174:     /* forms the element stiffness and coordinates */
175:     for (i=Ni0,ic=0,ii=0; i<Ni1; i++,ii++) {
176:       for (j=Nj0,jj=0; j<Nj1; j++,jj++) {
177:         for (k=Nk0,kk=0; k<Nk1; k++,kk++,ic++) {

179:           /* coords */
180:           x = coords[3*ic] = h*(PetscReal)i;
181:           y = coords[3*ic+1] = h*(PetscReal)j;
182:           z = coords[3*ic+2] = h*(PetscReal)k;
183:           /* matrix */
184:           id = id0 + ii + NN*jj + NN*NN*kk;

186:           if (i<ne && j<ne && k<ne) {
187:             /* radius */
188:             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));
189:             PetscReal alpha = 1.0;
190:             PetscInt  jx,ix,idx[8],idx3[24];
191:             idx[0] = id;
192:             idx[1] = id+1;
193:             idx[2] = id+NN+1;
194:             idx[3] = id+NN;
195:             idx[4] = id + NN*NN;
196:             idx[5] = id+1 + NN*NN;
197:             idx[6] = id+NN+1 + NN*NN;
198:             idx[7] = id+NN + NN*NN;

200:             /* correct indices */
201:             if (i==Ni1-1 && Ni1!=nn) {
202:               idx[1] += NN*(NN*NN-1);
203:               idx[2] += NN*(NN*NN-1);
204:               idx[5] += NN*(NN*NN-1);
205:               idx[6] += NN*(NN*NN-1);
206:             }
207:             if (j==Nj1-1 && Nj1!=nn) {
208:               idx[2] += NN*NN*(nn-1);
209:               idx[3] += NN*NN*(nn-1);
210:               idx[6] += NN*NN*(nn-1);
211:               idx[7] += NN*NN*(nn-1);
212:             }
213:             if (k==Nk1-1 && Nk1!=nn) {
214:               idx[4] += NN*(nn*nn-NN*NN);
215:               idx[5] += NN*(nn*nn-NN*NN);
216:               idx[6] += NN*(nn*nn-NN*NN);
217:               idx[7] += NN*(nn*nn-NN*NN);
218:             }

220:             if (radius < 0.25) alpha = soft_alpha;

222:             for (ix=0; ix<24; ix++) {
223:               for (jx=0;jx<24;jx++) DD[ix][jx] = alpha*DD1[ix][jx];
224:             }
225:             if (k>0) {
226:               if (!test_late_bs) {
227:                 MatSetValuesBlocked(Amat,8,idx,8,idx,(const PetscScalar*)DD,ADD_VALUES);
228:                 VecSetValuesBlocked(bb,8,idx,(const PetscScalar*)vv,ADD_VALUES);
229:               } else {
230:                 for (ix=0; ix<8; ix++) { idx3[3*ix] = 3*idx[ix]; idx3[3*ix+1] = 3*idx[ix]+1; idx3[3*ix+2] = 3*idx[ix]+2; }
231:                 MatSetValues(Amat,24,idx3,24,idx3,(const PetscScalar*)DD,ADD_VALUES);
232:                 VecSetValues(bb,24,idx3,(const PetscScalar*)vv,ADD_VALUES);
233:               }
234:             } else {
235:               /* a BC */
236:               for (ix=0;ix<24;ix++) {
237:                 for (jx=0;jx<24;jx++) DD[ix][jx] = alpha*DD2[ix][jx];
238:               }
239:               if (!test_late_bs) {
240:                 MatSetValuesBlocked(Amat,8,idx,8,idx,(const PetscScalar*)DD,ADD_VALUES);
241:                 VecSetValuesBlocked(bb,8,idx,(const PetscScalar*)v2,ADD_VALUES);
242:               } else {
243:                 for (ix=0; ix<8; ix++) { idx3[3*ix] = 3*idx[ix]; idx3[3*ix+1] = 3*idx[ix]+1; idx3[3*ix+2] = 3*idx[ix]+2; }
244:                 MatSetValues(Amat,24,idx3,24,idx3,(const PetscScalar*)DD,ADD_VALUES);
245:                 VecSetValues(bb,24,idx3,(const PetscScalar*)v2,ADD_VALUES);
246:               }
247:             }
248:           }
249:         }
250:       }

252:     }
253:     MatAssemblyBegin(Amat,MAT_FINAL_ASSEMBLY);
254:     MatAssemblyEnd(Amat,MAT_FINAL_ASSEMBLY);
255:     VecAssemblyBegin(bb);
256:     VecAssemblyEnd(bb);
257:   }

259:   if (test_late_bs) {
260:     VecSetBlockSize(xx,3);
261:     VecSetBlockSize(bb,3);
262:     MatSetBlockSize(Amat,3);
263:   }

265:   if (!PETSC_TRUE) {
266:     PetscViewer viewer;
267:     PetscViewerASCIIOpen(comm, "Amat.m", &viewer);
268:     PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_MATLAB);
269:     MatView(Amat,viewer);
270:     PetscViewerPopFormat(viewer);
271:     PetscViewerDestroy(&viewer);
272:   }

274:   /* finish KSP/PC setup */
275:   KSPSetOperators(ksp, Amat, Amat);
276:   if (use_nearnullspace) {
277:     MatNullSpace matnull;
278:     Vec          vec_coords;
279:     PetscScalar  *c;

281:     VecCreate(MPI_COMM_WORLD,&vec_coords);
282:     VecSetBlockSize(vec_coords,3);
283:     VecSetSizes(vec_coords,m,PETSC_DECIDE);
284:     VecSetUp(vec_coords);
285:     VecGetArray(vec_coords,&c);
286:     for (i=0; i<m; i++) c[i] = coords[i]; /* Copy since Scalar type might be Complex */
287:     VecRestoreArray(vec_coords,&c);
288:     MatNullSpaceCreateRigidBody(vec_coords,&matnull);
289:     MatSetNearNullSpace(Amat,matnull);
290:     MatNullSpaceDestroy(&matnull);
291:     VecDestroy(&vec_coords);
292:   } else {
293:     PC             pc;
294:     KSPGetPC(ksp,&pc);
295:     PCSetCoordinates(pc, 3, m/3, coords);
296:   }

298:   MaybeLogStagePush(stage[0]);

300:   /* PC setup basically */
301:   KSPSetUp(ksp);

303:   MaybeLogStagePop();
304:   MaybeLogStagePush(stage[1]);

306:   /* test BCs */
307:   if (test_nonzero_cols) {
308:     VecZeroEntries(xx);
309:     if (mype==0) VecSetValue(xx,0,1.0,INSERT_VALUES);
310:     VecAssemblyBegin(xx);
311:     VecAssemblyEnd(xx);
312:     KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
313:   }

315:   /* 1st solve */
316:   KSPSolve(ksp, bb, xx);

318:   MaybeLogStagePop();

320:   /* 2nd solve */
321:   if (two_solves) {
322:     PetscReal emax, emin;
323:     PetscReal norm,norm2;
324:     Vec       res;

326:     MaybeLogStagePush(stage[2]);
327:     /* PC setup basically */
328:     MatScale(Amat, 100000.0);
329:     KSPSetOperators(ksp, Amat, Amat);
330:     KSPSetUp(ksp);

332:     MaybeLogStagePop();
333:     MaybeLogStagePush(stage[3]);
334:     KSPSolve(ksp, bb, xx);
335:     KSPComputeExtremeSingularValues(ksp, &emax, &emin);

337:     MaybeLogStagePop();
338:     MaybeLogStagePush(stage[4]);

340:     MaybeLogStagePop();
341:     MaybeLogStagePush(stage[5]);

343:     /* 3rd solve */
344:     KSPSolve(ksp, bb, xx);

346:     MaybeLogStagePop();


349:     VecNorm(bb, NORM_2, &norm2);

351:     VecDuplicate(xx, &res);
352:     MatMult(Amat, xx, res);
353:     VecAXPY(bb, -1.0, res);
354:     VecDestroy(&res);
355:     VecNorm(bb, NORM_2, &norm);
356:     PetscPrintf(PETSC_COMM_WORLD,"[%d]%s |b-Ax|/|b|=%e, |b|=%e, emax=%e\n",0,PETSC_FUNCTION_NAME,(double)(norm/norm2),(double)norm2,(double)emax);
357:   }

359:   /* Free work space */
360:   KSPDestroy(&ksp);
361:   VecDestroy(&xx);
362:   VecDestroy(&bb);
363:   MatDestroy(&Amat);
364:   PetscFree(coords);

366:   PetscFinalize();
367:   return ierr;
368: }

370: /* Data was previously provided in the file data/elem_3d_elast_v_25.tx */
371: PetscErrorCode elem_3d_elast_v_25(PetscScalar *dd)
372: {
374:   PetscScalar    DD[] = {
375:   0.18981481481481474     ,
376:   5.27777777777777568E-002,
377:   5.27777777777777568E-002,
378:  -5.64814814814814659E-002,
379:  -1.38888888888889072E-002,
380:  -1.38888888888889089E-002,
381:  -8.24074074074073876E-002,
382:  -5.27777777777777429E-002,
383:   1.38888888888888725E-002,
384:   4.90740740740740339E-002,
385:   1.38888888888889124E-002,
386:   4.72222222222222071E-002,
387:   4.90740740740740339E-002,
388:   4.72222222222221932E-002,
389:   1.38888888888888968E-002,
390:  -8.24074074074073876E-002,
391:   1.38888888888888673E-002,
392:  -5.27777777777777429E-002,
393:  -7.87037037037036785E-002,
394:  -4.72222222222221932E-002,
395:  -4.72222222222222071E-002,
396:   1.20370370370370180E-002,
397:  -1.38888888888888742E-002,
398:  -1.38888888888888829E-002,
399:   5.27777777777777568E-002,
400:   0.18981481481481474     ,
401:   5.27777777777777568E-002,
402:   1.38888888888889124E-002,
403:   4.90740740740740269E-002,
404:   4.72222222222221932E-002,
405:  -5.27777777777777637E-002,
406:  -8.24074074074073876E-002,
407:   1.38888888888888725E-002,
408:  -1.38888888888889037E-002,
409:  -5.64814814814814728E-002,
410:  -1.38888888888888985E-002,
411:   4.72222222222221932E-002,
412:   4.90740740740740478E-002,
413:   1.38888888888888968E-002,
414:  -1.38888888888888673E-002,
415:   1.20370370370370058E-002,
416:  -1.38888888888888742E-002,
417:  -4.72222222222221932E-002,
418:  -7.87037037037036785E-002,
419:  -4.72222222222222002E-002,
420:   1.38888888888888742E-002,
421:  -8.24074074074073598E-002,
422:  -5.27777777777777568E-002,
423:   5.27777777777777568E-002,
424:   5.27777777777777568E-002,
425:   0.18981481481481474     ,
426:   1.38888888888889055E-002,
427:   4.72222222222222002E-002,
428:   4.90740740740740269E-002,
429:  -1.38888888888888829E-002,
430:  -1.38888888888888829E-002,
431:   1.20370370370370180E-002,
432:   4.72222222222222002E-002,
433:   1.38888888888888985E-002,
434:   4.90740740740740339E-002,
435:  -1.38888888888888985E-002,
436:  -1.38888888888888968E-002,
437:  -5.64814814814814520E-002,
438:  -5.27777777777777568E-002,
439:   1.38888888888888777E-002,
440:  -8.24074074074073876E-002,
441:  -4.72222222222222002E-002,
442:  -4.72222222222221932E-002,
443:  -7.87037037037036646E-002,
444:   1.38888888888888794E-002,
445:  -5.27777777777777568E-002,
446:  -8.24074074074073598E-002,
447:  -5.64814814814814659E-002,
448:   1.38888888888889124E-002,
449:   1.38888888888889055E-002,
450:   0.18981481481481474     ,
451:  -5.27777777777777568E-002,
452:  -5.27777777777777499E-002,
453:   4.90740740740740269E-002,
454:  -1.38888888888889072E-002,
455:  -4.72222222222221932E-002,
456:  -8.24074074074073876E-002,
457:   5.27777777777777568E-002,
458:  -1.38888888888888812E-002,
459:  -8.24074074074073876E-002,
460:  -1.38888888888888742E-002,
461:   5.27777777777777499E-002,
462:   4.90740740740740269E-002,
463:  -4.72222222222221863E-002,
464:  -1.38888888888889089E-002,
465:   1.20370370370370162E-002,
466:   1.38888888888888673E-002,
467:   1.38888888888888742E-002,
468:  -7.87037037037036785E-002,
469:   4.72222222222222002E-002,
470:   4.72222222222222071E-002,
471:  -1.38888888888889072E-002,
472:   4.90740740740740269E-002,
473:   4.72222222222222002E-002,
474:  -5.27777777777777568E-002,
475:   0.18981481481481480     ,
476:   5.27777777777777568E-002,
477:   1.38888888888889020E-002,
478:  -5.64814814814814728E-002,
479:  -1.38888888888888951E-002,
480:   5.27777777777777637E-002,
481:  -8.24074074074073876E-002,
482:   1.38888888888888881E-002,
483:   1.38888888888888742E-002,
484:   1.20370370370370232E-002,
485:  -1.38888888888888812E-002,
486:  -4.72222222222221863E-002,
487:   4.90740740740740339E-002,
488:   1.38888888888888933E-002,
489:  -1.38888888888888812E-002,
490:  -8.24074074074073876E-002,
491:  -5.27777777777777568E-002,
492:   4.72222222222222071E-002,
493:  -7.87037037037036924E-002,
494:  -4.72222222222222140E-002,
495:  -1.38888888888889089E-002,
496:   4.72222222222221932E-002,
497:   4.90740740740740269E-002,
498:  -5.27777777777777499E-002,
499:   5.27777777777777568E-002,
500:   0.18981481481481477     ,
501:  -4.72222222222222071E-002,
502:   1.38888888888888968E-002,
503:   4.90740740740740131E-002,
504:   1.38888888888888812E-002,
505:  -1.38888888888888708E-002,
506:   1.20370370370370267E-002,
507:   5.27777777777777568E-002,
508:   1.38888888888888812E-002,
509:  -8.24074074074073876E-002,
510:   1.38888888888889124E-002,
511:  -1.38888888888889055E-002,
512:  -5.64814814814814589E-002,
513:  -1.38888888888888812E-002,
514:  -5.27777777777777568E-002,
515:  -8.24074074074073737E-002,
516:   4.72222222222222002E-002,
517:  -4.72222222222222002E-002,
518:  -7.87037037037036924E-002,
519:  -8.24074074074073876E-002,
520:  -5.27777777777777637E-002,
521:  -1.38888888888888829E-002,
522:   4.90740740740740269E-002,
523:   1.38888888888889020E-002,
524:  -4.72222222222222071E-002,
525:   0.18981481481481480     ,
526:   5.27777777777777637E-002,
527:  -5.27777777777777637E-002,
528:  -5.64814814814814728E-002,
529:  -1.38888888888889037E-002,
530:   1.38888888888888951E-002,
531:  -7.87037037037036785E-002,
532:  -4.72222222222222002E-002,
533:   4.72222222222221932E-002,
534:   1.20370370370370128E-002,
535:  -1.38888888888888725E-002,
536:   1.38888888888888812E-002,
537:   4.90740740740740408E-002,
538:   4.72222222222222002E-002,
539:  -1.38888888888888951E-002,
540:  -8.24074074074073876E-002,
541:   1.38888888888888812E-002,
542:   5.27777777777777637E-002,
543:  -5.27777777777777429E-002,
544:  -8.24074074074073876E-002,
545:  -1.38888888888888829E-002,
546:  -1.38888888888889072E-002,
547:  -5.64814814814814728E-002,
548:   1.38888888888888968E-002,
549:   5.27777777777777637E-002,
550:   0.18981481481481480     ,
551:  -5.27777777777777568E-002,
552:   1.38888888888888916E-002,
553:   4.90740740740740339E-002,
554:  -4.72222222222222210E-002,
555:  -4.72222222222221932E-002,
556:  -7.87037037037036924E-002,
557:   4.72222222222222002E-002,
558:   1.38888888888888742E-002,
559:  -8.24074074074073876E-002,
560:   5.27777777777777429E-002,
561:   4.72222222222222002E-002,
562:   4.90740740740740269E-002,
563:  -1.38888888888888951E-002,
564:  -1.38888888888888846E-002,
565:   1.20370370370370267E-002,
566:   1.38888888888888916E-002,
567:   1.38888888888888725E-002,
568:   1.38888888888888725E-002,
569:   1.20370370370370180E-002,
570:  -4.72222222222221932E-002,
571:  -1.38888888888888951E-002,
572:   4.90740740740740131E-002,
573:  -5.27777777777777637E-002,
574:  -5.27777777777777568E-002,
575:   0.18981481481481480     ,
576:  -1.38888888888888968E-002,
577:  -4.72222222222221932E-002,
578:   4.90740740740740339E-002,
579:   4.72222222222221932E-002,
580:   4.72222222222222071E-002,
581:  -7.87037037037036646E-002,
582:  -1.38888888888888742E-002,
583:   5.27777777777777499E-002,
584:  -8.24074074074073737E-002,
585:   1.38888888888888933E-002,
586:   1.38888888888889020E-002,
587:  -5.64814814814814589E-002,
588:   5.27777777777777568E-002,
589:  -1.38888888888888794E-002,
590:  -8.24074074074073876E-002,
591:   4.90740740740740339E-002,
592:  -1.38888888888889037E-002,
593:   4.72222222222222002E-002,
594:  -8.24074074074073876E-002,
595:   5.27777777777777637E-002,
596:   1.38888888888888812E-002,
597:  -5.64814814814814728E-002,
598:   1.38888888888888916E-002,
599:  -1.38888888888888968E-002,
600:   0.18981481481481480     ,
601:  -5.27777777777777499E-002,
602:   5.27777777777777707E-002,
603:   1.20370370370370180E-002,
604:   1.38888888888888812E-002,
605:  -1.38888888888888812E-002,
606:  -7.87037037037036785E-002,
607:   4.72222222222222002E-002,
608:  -4.72222222222222071E-002,
609:  -8.24074074074073876E-002,
610:  -1.38888888888888742E-002,
611:  -5.27777777777777568E-002,
612:   4.90740740740740616E-002,
613:  -4.72222222222222002E-002,
614:   1.38888888888888846E-002,
615:   1.38888888888889124E-002,
616:  -5.64814814814814728E-002,
617:   1.38888888888888985E-002,
618:   5.27777777777777568E-002,
619:  -8.24074074074073876E-002,
620:  -1.38888888888888708E-002,
621:  -1.38888888888889037E-002,
622:   4.90740740740740339E-002,
623:  -4.72222222222221932E-002,
624:  -5.27777777777777499E-002,
625:   0.18981481481481480     ,
626:  -5.27777777777777568E-002,
627:  -1.38888888888888673E-002,
628:  -8.24074074074073598E-002,
629:   5.27777777777777429E-002,
630:   4.72222222222222002E-002,
631:  -7.87037037037036785E-002,
632:   4.72222222222222002E-002,
633:   1.38888888888888708E-002,
634:   1.20370370370370128E-002,
635:   1.38888888888888760E-002,
636:  -4.72222222222222002E-002,
637:   4.90740740740740478E-002,
638:  -1.38888888888888951E-002,
639:   4.72222222222222071E-002,
640:  -1.38888888888888985E-002,
641:   4.90740740740740339E-002,
642:  -1.38888888888888812E-002,
643:   1.38888888888888881E-002,
644:   1.20370370370370267E-002,
645:   1.38888888888888951E-002,
646:  -4.72222222222222210E-002,
647:   4.90740740740740339E-002,
648:   5.27777777777777707E-002,
649:  -5.27777777777777568E-002,
650:   0.18981481481481477     ,
651:   1.38888888888888829E-002,
652:   5.27777777777777707E-002,
653:  -8.24074074074073598E-002,
654:  -4.72222222222222140E-002,
655:   4.72222222222222140E-002,
656:  -7.87037037037036646E-002,
657:  -5.27777777777777707E-002,
658:  -1.38888888888888829E-002,
659:  -8.24074074074073876E-002,
660:  -1.38888888888888881E-002,
661:   1.38888888888888881E-002,
662:  -5.64814814814814589E-002,
663:   4.90740740740740339E-002,
664:   4.72222222222221932E-002,
665:  -1.38888888888888985E-002,
666:  -8.24074074074073876E-002,
667:   1.38888888888888742E-002,
668:   5.27777777777777568E-002,
669:  -7.87037037037036785E-002,
670:  -4.72222222222221932E-002,
671:   4.72222222222221932E-002,
672:   1.20370370370370180E-002,
673:  -1.38888888888888673E-002,
674:   1.38888888888888829E-002,
675:   0.18981481481481469     ,
676:   5.27777777777777429E-002,
677:  -5.27777777777777429E-002,
678:  -5.64814814814814659E-002,
679:  -1.38888888888889055E-002,
680:   1.38888888888889055E-002,
681:  -8.24074074074074153E-002,
682:  -5.27777777777777429E-002,
683:  -1.38888888888888760E-002,
684:   4.90740740740740408E-002,
685:   1.38888888888888968E-002,
686:  -4.72222222222222071E-002,
687:   4.72222222222221932E-002,
688:   4.90740740740740478E-002,
689:  -1.38888888888888968E-002,
690:  -1.38888888888888742E-002,
691:   1.20370370370370232E-002,
692:   1.38888888888888812E-002,
693:  -4.72222222222222002E-002,
694:  -7.87037037037036924E-002,
695:   4.72222222222222071E-002,
696:   1.38888888888888812E-002,
697:  -8.24074074074073598E-002,
698:   5.27777777777777707E-002,
699:   5.27777777777777429E-002,
700:   0.18981481481481477     ,
701:  -5.27777777777777499E-002,
702:   1.38888888888889107E-002,
703:   4.90740740740740478E-002,
704:  -4.72222222222221932E-002,
705:  -5.27777777777777568E-002,
706:  -8.24074074074074153E-002,
707:  -1.38888888888888812E-002,
708:  -1.38888888888888846E-002,
709:  -5.64814814814814659E-002,
710:   1.38888888888888812E-002,
711:   1.38888888888888968E-002,
712:   1.38888888888888968E-002,
713:  -5.64814814814814520E-002,
714:   5.27777777777777499E-002,
715:  -1.38888888888888812E-002,
716:  -8.24074074074073876E-002,
717:   4.72222222222221932E-002,
718:   4.72222222222222002E-002,
719:  -7.87037037037036646E-002,
720:  -1.38888888888888812E-002,
721:   5.27777777777777429E-002,
722:  -8.24074074074073598E-002,
723:  -5.27777777777777429E-002,
724:  -5.27777777777777499E-002,
725:   0.18981481481481474     ,
726:  -1.38888888888888985E-002,
727:  -4.72222222222221863E-002,
728:   4.90740740740740339E-002,
729:   1.38888888888888829E-002,
730:   1.38888888888888777E-002,
731:   1.20370370370370249E-002,
732:  -4.72222222222222002E-002,
733:  -1.38888888888888933E-002,
734:   4.90740740740740339E-002,
735:  -8.24074074074073876E-002,
736:  -1.38888888888888673E-002,
737:  -5.27777777777777568E-002,
738:   4.90740740740740269E-002,
739:  -4.72222222222221863E-002,
740:   1.38888888888889124E-002,
741:   1.20370370370370128E-002,
742:   1.38888888888888742E-002,
743:  -1.38888888888888742E-002,
744:  -7.87037037037036785E-002,
745:   4.72222222222222002E-002,
746:  -4.72222222222222140E-002,
747:  -5.64814814814814659E-002,
748:   1.38888888888889107E-002,
749:  -1.38888888888888985E-002,
750:   0.18981481481481474     ,
751:  -5.27777777777777499E-002,
752:   5.27777777777777499E-002,
753:   4.90740740740740339E-002,
754:  -1.38888888888889055E-002,
755:   4.72222222222221932E-002,
756:  -8.24074074074074153E-002,
757:   5.27777777777777499E-002,
758:   1.38888888888888829E-002,
759:   1.38888888888888673E-002,
760:   1.20370370370370058E-002,
761:   1.38888888888888777E-002,
762:  -4.72222222222221863E-002,
763:   4.90740740740740339E-002,
764:  -1.38888888888889055E-002,
765:  -1.38888888888888725E-002,
766:  -8.24074074074073876E-002,
767:   5.27777777777777499E-002,
768:   4.72222222222222002E-002,
769:  -7.87037037037036785E-002,
770:   4.72222222222222140E-002,
771:  -1.38888888888889055E-002,
772:   4.90740740740740478E-002,
773:  -4.72222222222221863E-002,
774:  -5.27777777777777499E-002,
775:   0.18981481481481469     ,
776:  -5.27777777777777499E-002,
777:   1.38888888888889072E-002,
778:  -5.64814814814814659E-002,
779:   1.38888888888889003E-002,
780:   5.27777777777777429E-002,
781:  -8.24074074074074153E-002,
782:  -1.38888888888888812E-002,
783:  -5.27777777777777429E-002,
784:  -1.38888888888888742E-002,
785:  -8.24074074074073876E-002,
786:  -1.38888888888889089E-002,
787:   1.38888888888888933E-002,
788:  -5.64814814814814589E-002,
789:   1.38888888888888812E-002,
790:   5.27777777777777429E-002,
791:  -8.24074074074073737E-002,
792:  -4.72222222222222071E-002,
793:   4.72222222222222002E-002,
794:  -7.87037037037036646E-002,
795:   1.38888888888889055E-002,
796:  -4.72222222222221932E-002,
797:   4.90740740740740339E-002,
798:   5.27777777777777499E-002,
799:  -5.27777777777777499E-002,
800:   0.18981481481481474     ,
801:   4.72222222222222002E-002,
802:  -1.38888888888888985E-002,
803:   4.90740740740740339E-002,
804:  -1.38888888888888846E-002,
805:   1.38888888888888812E-002,
806:   1.20370370370370284E-002,
807:  -7.87037037037036785E-002,
808:  -4.72222222222221932E-002,
809:  -4.72222222222222002E-002,
810:   1.20370370370370162E-002,
811:  -1.38888888888888812E-002,
812:  -1.38888888888888812E-002,
813:   4.90740740740740408E-002,
814:   4.72222222222222002E-002,
815:   1.38888888888888933E-002,
816:  -8.24074074074073876E-002,
817:   1.38888888888888708E-002,
818:  -5.27777777777777707E-002,
819:  -8.24074074074074153E-002,
820:  -5.27777777777777568E-002,
821:   1.38888888888888829E-002,
822:   4.90740740740740339E-002,
823:   1.38888888888889072E-002,
824:   4.72222222222222002E-002,
825:   0.18981481481481477     ,
826:   5.27777777777777429E-002,
827:   5.27777777777777568E-002,
828:  -5.64814814814814659E-002,
829:  -1.38888888888888846E-002,
830:  -1.38888888888888881E-002,
831:  -4.72222222222221932E-002,
832:  -7.87037037037036785E-002,
833:  -4.72222222222221932E-002,
834:   1.38888888888888673E-002,
835:  -8.24074074074073876E-002,
836:  -5.27777777777777568E-002,
837:   4.72222222222222002E-002,
838:   4.90740740740740269E-002,
839:   1.38888888888889020E-002,
840:  -1.38888888888888742E-002,
841:   1.20370370370370128E-002,
842:  -1.38888888888888829E-002,
843:  -5.27777777777777429E-002,
844:  -8.24074074074074153E-002,
845:   1.38888888888888777E-002,
846:  -1.38888888888889055E-002,
847:  -5.64814814814814659E-002,
848:  -1.38888888888888985E-002,
849:   5.27777777777777429E-002,
850:   0.18981481481481469     ,
851:   5.27777777777777429E-002,
852:   1.38888888888888933E-002,
853:   4.90740740740740339E-002,
854:   4.72222222222222071E-002,
855:  -4.72222222222222071E-002,
856:  -4.72222222222222002E-002,
857:  -7.87037037037036646E-002,
858:   1.38888888888888742E-002,
859:  -5.27777777777777568E-002,
860:  -8.24074074074073737E-002,
861:  -1.38888888888888951E-002,
862:  -1.38888888888888951E-002,
863:  -5.64814814814814589E-002,
864:  -5.27777777777777568E-002,
865:   1.38888888888888760E-002,
866:  -8.24074074074073876E-002,
867:  -1.38888888888888760E-002,
868:  -1.38888888888888812E-002,
869:   1.20370370370370249E-002,
870:   4.72222222222221932E-002,
871:   1.38888888888889003E-002,
872:   4.90740740740740339E-002,
873:   5.27777777777777568E-002,
874:   5.27777777777777429E-002,
875:   0.18981481481481474     ,
876:   1.38888888888888933E-002,
877:   4.72222222222222071E-002,
878:   4.90740740740740339E-002,
879:   1.20370370370370180E-002,
880:   1.38888888888888742E-002,
881:   1.38888888888888794E-002,
882:  -7.87037037037036785E-002,
883:   4.72222222222222071E-002,
884:   4.72222222222222002E-002,
885:  -8.24074074074073876E-002,
886:  -1.38888888888888846E-002,
887:   5.27777777777777568E-002,
888:   4.90740740740740616E-002,
889:  -4.72222222222222002E-002,
890:  -1.38888888888888881E-002,
891:   4.90740740740740408E-002,
892:  -1.38888888888888846E-002,
893:  -4.72222222222222002E-002,
894:  -8.24074074074074153E-002,
895:   5.27777777777777429E-002,
896:  -1.38888888888888846E-002,
897:  -5.64814814814814659E-002,
898:   1.38888888888888933E-002,
899:   1.38888888888888933E-002,
900:   0.18981481481481477     ,
901:  -5.27777777777777568E-002,
902:  -5.27777777777777637E-002,
903:  -1.38888888888888742E-002,
904:  -8.24074074074073598E-002,
905:  -5.27777777777777568E-002,
906:   4.72222222222222002E-002,
907:  -7.87037037037036924E-002,
908:  -4.72222222222222002E-002,
909:   1.38888888888888812E-002,
910:   1.20370370370370267E-002,
911:  -1.38888888888888794E-002,
912:  -4.72222222222222002E-002,
913:   4.90740740740740478E-002,
914:   1.38888888888888881E-002,
915:   1.38888888888888968E-002,
916:  -5.64814814814814659E-002,
917:  -1.38888888888888933E-002,
918:   5.27777777777777499E-002,
919:  -8.24074074074074153E-002,
920:   1.38888888888888812E-002,
921:  -1.38888888888888846E-002,
922:   4.90740740740740339E-002,
923:   4.72222222222222071E-002,
924:  -5.27777777777777568E-002,
925:   0.18981481481481477     ,
926:   5.27777777777777637E-002,
927:  -1.38888888888888829E-002,
928:  -5.27777777777777568E-002,
929:  -8.24074074074073598E-002,
930:   4.72222222222222071E-002,
931:  -4.72222222222222140E-002,
932:  -7.87037037037036924E-002,
933:   5.27777777777777637E-002,
934:   1.38888888888888916E-002,
935:  -8.24074074074073876E-002,
936:   1.38888888888888846E-002,
937:  -1.38888888888888951E-002,
938:  -5.64814814814814589E-002,
939:  -4.72222222222222071E-002,
940:   1.38888888888888812E-002,
941:   4.90740740740740339E-002,
942:   1.38888888888888829E-002,
943:  -1.38888888888888812E-002,
944:   1.20370370370370284E-002,
945:  -1.38888888888888881E-002,
946:   4.72222222222222071E-002,
947:   4.90740740740740339E-002,
948:  -5.27777777777777637E-002,
949:   5.27777777777777637E-002,
950:   0.18981481481481477     ,
951:   };

954:   PetscMemcpy(dd,DD,sizeof(PetscScalar)*576);
955:   return(0);
956: }