1 /*
2 Include "petsctao.h" so that we can use TAO solvers. Note that this
3 file automatically includes libraries such as:
4 petsc.h - base PETSc routines petscvec.h - vectors
5 petscsys.h - system routines petscmat.h - matrices
6 petscis.h - index sets petscksp.h - Krylov subspace methods
7 petscviewer.h - viewers petscpc.h - preconditioners
8
9 This version tests correlated terms using both vector and listed forms
10 */
11
12 #include <petsctao.h>
13
14 /*
15 Description: These data are the result of a NIST study involving
16 ultrasonic calibration. The response variable is
17 ultrasonic response, and the predictor variable is
18 metal distance.
19
20 Reference: Chwirut, D., NIST (197?).
21 Ultrasonic Reference Block Study.
22 */
23
24 static char help[] = "Finds the nonlinear least-squares solution to the model \n\
25 y = exp[-b1*x]/(b2+b3*x) + e \n";
26
27 #define NOBSERVATIONS 214
28 #define NPARAMETERS 3
29
30 /* User-defined application context */
31 typedef struct {
32 /* Working space */
33 PetscReal t[NOBSERVATIONS]; /* array of independent variables of observation */
34 PetscReal y[NOBSERVATIONS]; /* array of dependent variables */
35 PetscReal j[NOBSERVATIONS][NPARAMETERS]; /* dense jacobian matrix array*/
36 PetscInt idm[NOBSERVATIONS]; /* Matrix indices for jacobian */
37 PetscInt idn[NPARAMETERS];
38 } AppCtx;
39
40 /* User provided Routines */
41 PetscErrorCode InitializeData(AppCtx *user);
42 PetscErrorCode FormStartingPoint(Vec);
43 PetscErrorCode EvaluateFunction(Tao, Vec, Vec, void *);
44 PetscErrorCode EvaluateJacobian(Tao, Vec, Mat, Mat, void *);
45
46 /*--------------------------------------------------------------------*/
main(int argc,char ** argv)47 int main(int argc, char **argv)
48 {
49 PetscInt wtype = 0;
50 Vec x, f; /* solution, function */
51 Vec w; /* weights */
52 Mat J; /* Jacobian matrix */
53 Tao tao; /* Tao solver context */
54 PetscInt i; /* iteration information */
55 PetscReal hist[100], resid[100];
56 PetscInt lits[100];
57 PetscInt w_row[NOBSERVATIONS]; /* explicit weights */
58 PetscInt w_col[NOBSERVATIONS];
59 PetscReal w_vals[NOBSERVATIONS];
60 PetscBool flg;
61 AppCtx user; /* user-defined work context */
62
63 PetscFunctionBeginUser;
64 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
65 PetscCall(PetscOptionsGetInt(NULL, NULL, "-wtype", &wtype, &flg));
66 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "wtype=%" PetscInt_FMT "\n", wtype));
67 /* Allocate vectors */
68 PetscCall(VecCreateSeq(MPI_COMM_SELF, NPARAMETERS, &x));
69 PetscCall(VecCreateSeq(MPI_COMM_SELF, NOBSERVATIONS, &f));
70
71 PetscCall(VecDuplicate(f, &w));
72
73 /* no correlation, but set in different ways */
74 PetscCall(VecSet(w, 1.0));
75 for (i = 0; i < NOBSERVATIONS; i++) {
76 w_row[i] = i;
77 w_col[i] = i;
78 w_vals[i] = 1.0;
79 }
80
81 /* Create the Jacobian matrix. */
82 PetscCall(MatCreateSeqDense(MPI_COMM_SELF, NOBSERVATIONS, NPARAMETERS, NULL, &J));
83
84 for (i = 0; i < NOBSERVATIONS; i++) user.idm[i] = i;
85
86 for (i = 0; i < NPARAMETERS; i++) user.idn[i] = i;
87
88 /* Create TAO solver and set desired solution method */
89 PetscCall(TaoCreate(PETSC_COMM_SELF, &tao));
90 PetscCall(TaoSetType(tao, TAOPOUNDERS));
91
92 /* Set the function and Jacobian routines. */
93 PetscCall(InitializeData(&user));
94 PetscCall(FormStartingPoint(x));
95 PetscCall(TaoSetSolution(tao, x));
96 PetscCall(TaoSetResidualRoutine(tao, f, EvaluateFunction, (void *)&user));
97 if (wtype == 1) {
98 PetscCall(TaoSetResidualWeights(tao, w, 0, NULL, NULL, NULL));
99 } else if (wtype == 2) {
100 PetscCall(TaoSetResidualWeights(tao, NULL, NOBSERVATIONS, w_row, w_col, w_vals));
101 }
102 PetscCall(TaoSetJacobianResidualRoutine(tao, J, J, EvaluateJacobian, (void *)&user));
103 PetscCall(TaoSetTolerances(tao, 1e-5, 0.0, PETSC_CURRENT));
104
105 /* Check for any TAO command line arguments */
106 PetscCall(TaoSetFromOptions(tao));
107
108 PetscCall(TaoSetConvergenceHistory(tao, hist, resid, 0, lits, 100, PETSC_TRUE));
109 /* Perform the Solve */
110 PetscCall(TaoSolve(tao));
111
112 /* Free TAO data structures */
113 PetscCall(TaoDestroy(&tao));
114
115 /* Free PETSc data structures */
116 PetscCall(VecDestroy(&x));
117 PetscCall(VecDestroy(&w));
118 PetscCall(VecDestroy(&f));
119 PetscCall(MatDestroy(&J));
120
121 PetscCall(PetscFinalize());
122 return 0;
123 }
124
125 /*--------------------------------------------------------------------*/
EvaluateFunction(Tao tao,Vec X,Vec F,void * ptr)126 PetscErrorCode EvaluateFunction(Tao tao, Vec X, Vec F, void *ptr)
127 {
128 AppCtx *user = (AppCtx *)ptr;
129 PetscInt i;
130 PetscReal *y = user->y, *f, *t = user->t;
131 const PetscReal *x;
132
133 PetscFunctionBegin;
134 PetscCall(VecGetArrayRead(X, &x));
135 PetscCall(VecGetArray(F, &f));
136
137 for (i = 0; i < NOBSERVATIONS; i++) f[i] = y[i] - PetscExpScalar(-x[0] * t[i]) / (x[1] + x[2] * t[i]);
138 PetscCall(VecRestoreArrayRead(X, &x));
139 PetscCall(VecRestoreArray(F, &f));
140 PetscCall(PetscLogFlops(6 * NOBSERVATIONS));
141 PetscFunctionReturn(PETSC_SUCCESS);
142 }
143
144 /*------------------------------------------------------------*/
145 /* J[i][j] = df[i]/dt[j] */
EvaluateJacobian(Tao tao,Vec X,Mat J,Mat Jpre,void * ptr)146 PetscErrorCode EvaluateJacobian(Tao tao, Vec X, Mat J, Mat Jpre, void *ptr)
147 {
148 AppCtx *user = (AppCtx *)ptr;
149 PetscInt i;
150 PetscReal *t = user->t;
151 const PetscReal *x;
152 PetscReal base;
153
154 PetscFunctionBegin;
155 PetscCall(VecGetArrayRead(X, &x));
156 for (i = 0; i < NOBSERVATIONS; i++) {
157 base = PetscExpScalar(-x[0] * t[i]) / (x[1] + x[2] * t[i]);
158
159 user->j[i][0] = t[i] * base;
160 user->j[i][1] = base / (x[1] + x[2] * t[i]);
161 user->j[i][2] = base * t[i] / (x[1] + x[2] * t[i]);
162 }
163
164 /* Assemble the matrix */
165 PetscCall(MatSetValues(J, NOBSERVATIONS, user->idm, NPARAMETERS, user->idn, (PetscReal *)user->j, INSERT_VALUES));
166 PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
167 PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
168
169 PetscCall(VecRestoreArrayRead(X, &x));
170 PetscCall(PetscLogFlops(NOBSERVATIONS * 13));
171 PetscFunctionReturn(PETSC_SUCCESS);
172 }
173
174 /* ------------------------------------------------------------ */
FormStartingPoint(Vec X)175 PetscErrorCode FormStartingPoint(Vec X)
176 {
177 PetscReal *x;
178
179 PetscFunctionBegin;
180 PetscCall(VecGetArray(X, &x));
181 x[0] = 1.19;
182 x[1] = -1.86;
183 x[2] = 1.08;
184 PetscCall(VecRestoreArray(X, &x));
185 PetscFunctionReturn(PETSC_SUCCESS);
186 }
187
188 /* ---------------------------------------------------------------------- */
InitializeData(AppCtx * user)189 PetscErrorCode InitializeData(AppCtx *user)
190 {
191 PetscReal *t = user->t, *y = user->y;
192 PetscInt i = 0;
193
194 PetscFunctionBegin;
195 y[i] = 92.9000;
196 t[i++] = 0.5000;
197 y[i] = 78.7000;
198 t[i++] = 0.6250;
199 y[i] = 64.2000;
200 t[i++] = 0.7500;
201 y[i] = 64.9000;
202 t[i++] = 0.8750;
203 y[i] = 57.1000;
204 t[i++] = 1.0000;
205 y[i] = 43.3000;
206 t[i++] = 1.2500;
207 y[i] = 31.1000;
208 t[i++] = 1.7500;
209 y[i] = 23.6000;
210 t[i++] = 2.2500;
211 y[i] = 31.0500;
212 t[i++] = 1.7500;
213 y[i] = 23.7750;
214 t[i++] = 2.2500;
215 y[i] = 17.7375;
216 t[i++] = 2.7500;
217 y[i] = 13.8000;
218 t[i++] = 3.2500;
219 y[i] = 11.5875;
220 t[i++] = 3.7500;
221 y[i] = 9.4125;
222 t[i++] = 4.2500;
223 y[i] = 7.7250;
224 t[i++] = 4.7500;
225 y[i] = 7.3500;
226 t[i++] = 5.2500;
227 y[i] = 8.0250;
228 t[i++] = 5.7500;
229 y[i] = 90.6000;
230 t[i++] = 0.5000;
231 y[i] = 76.9000;
232 t[i++] = 0.6250;
233 y[i] = 71.6000;
234 t[i++] = 0.7500;
235 y[i] = 63.6000;
236 t[i++] = 0.8750;
237 y[i] = 54.0000;
238 t[i++] = 1.0000;
239 y[i] = 39.2000;
240 t[i++] = 1.2500;
241 y[i] = 29.3000;
242 t[i++] = 1.7500;
243 y[i] = 21.4000;
244 t[i++] = 2.2500;
245 y[i] = 29.1750;
246 t[i++] = 1.7500;
247 y[i] = 22.1250;
248 t[i++] = 2.2500;
249 y[i] = 17.5125;
250 t[i++] = 2.7500;
251 y[i] = 14.2500;
252 t[i++] = 3.2500;
253 y[i] = 9.4500;
254 t[i++] = 3.7500;
255 y[i] = 9.1500;
256 t[i++] = 4.2500;
257 y[i] = 7.9125;
258 t[i++] = 4.7500;
259 y[i] = 8.4750;
260 t[i++] = 5.2500;
261 y[i] = 6.1125;
262 t[i++] = 5.7500;
263 y[i] = 80.0000;
264 t[i++] = 0.5000;
265 y[i] = 79.0000;
266 t[i++] = 0.6250;
267 y[i] = 63.8000;
268 t[i++] = 0.7500;
269 y[i] = 57.2000;
270 t[i++] = 0.8750;
271 y[i] = 53.2000;
272 t[i++] = 1.0000;
273 y[i] = 42.5000;
274 t[i++] = 1.2500;
275 y[i] = 26.8000;
276 t[i++] = 1.7500;
277 y[i] = 20.4000;
278 t[i++] = 2.2500;
279 y[i] = 26.8500;
280 t[i++] = 1.7500;
281 y[i] = 21.0000;
282 t[i++] = 2.2500;
283 y[i] = 16.4625;
284 t[i++] = 2.7500;
285 y[i] = 12.5250;
286 t[i++] = 3.2500;
287 y[i] = 10.5375;
288 t[i++] = 3.7500;
289 y[i] = 8.5875;
290 t[i++] = 4.2500;
291 y[i] = 7.1250;
292 t[i++] = 4.7500;
293 y[i] = 6.1125;
294 t[i++] = 5.2500;
295 y[i] = 5.9625;
296 t[i++] = 5.7500;
297 y[i] = 74.1000;
298 t[i++] = 0.5000;
299 y[i] = 67.3000;
300 t[i++] = 0.6250;
301 y[i] = 60.8000;
302 t[i++] = 0.7500;
303 y[i] = 55.5000;
304 t[i++] = 0.8750;
305 y[i] = 50.3000;
306 t[i++] = 1.0000;
307 y[i] = 41.0000;
308 t[i++] = 1.2500;
309 y[i] = 29.4000;
310 t[i++] = 1.7500;
311 y[i] = 20.4000;
312 t[i++] = 2.2500;
313 y[i] = 29.3625;
314 t[i++] = 1.7500;
315 y[i] = 21.1500;
316 t[i++] = 2.2500;
317 y[i] = 16.7625;
318 t[i++] = 2.7500;
319 y[i] = 13.2000;
320 t[i++] = 3.2500;
321 y[i] = 10.8750;
322 t[i++] = 3.7500;
323 y[i] = 8.1750;
324 t[i++] = 4.2500;
325 y[i] = 7.3500;
326 t[i++] = 4.7500;
327 y[i] = 5.9625;
328 t[i++] = 5.2500;
329 y[i] = 5.6250;
330 t[i++] = 5.7500;
331 y[i] = 81.5000;
332 t[i++] = .5000;
333 y[i] = 62.4000;
334 t[i++] = .7500;
335 y[i] = 32.5000;
336 t[i++] = 1.5000;
337 y[i] = 12.4100;
338 t[i++] = 3.0000;
339 y[i] = 13.1200;
340 t[i++] = 3.0000;
341 y[i] = 15.5600;
342 t[i++] = 3.0000;
343 y[i] = 5.6300;
344 t[i++] = 6.0000;
345 y[i] = 78.0000;
346 t[i++] = .5000;
347 y[i] = 59.9000;
348 t[i++] = .7500;
349 y[i] = 33.2000;
350 t[i++] = 1.5000;
351 y[i] = 13.8400;
352 t[i++] = 3.0000;
353 y[i] = 12.7500;
354 t[i++] = 3.0000;
355 y[i] = 14.6200;
356 t[i++] = 3.0000;
357 y[i] = 3.9400;
358 t[i++] = 6.0000;
359 y[i] = 76.8000;
360 t[i++] = .5000;
361 y[i] = 61.0000;
362 t[i++] = .7500;
363 y[i] = 32.9000;
364 t[i++] = 1.5000;
365 y[i] = 13.8700;
366 t[i++] = 3.0000;
367 y[i] = 11.8100;
368 t[i++] = 3.0000;
369 y[i] = 13.3100;
370 t[i++] = 3.0000;
371 y[i] = 5.4400;
372 t[i++] = 6.0000;
373 y[i] = 78.0000;
374 t[i++] = .5000;
375 y[i] = 63.5000;
376 t[i++] = .7500;
377 y[i] = 33.8000;
378 t[i++] = 1.5000;
379 y[i] = 12.5600;
380 t[i++] = 3.0000;
381 y[i] = 5.6300;
382 t[i++] = 6.0000;
383 y[i] = 12.7500;
384 t[i++] = 3.0000;
385 y[i] = 13.1200;
386 t[i++] = 3.0000;
387 y[i] = 5.4400;
388 t[i++] = 6.0000;
389 y[i] = 76.8000;
390 t[i++] = .5000;
391 y[i] = 60.0000;
392 t[i++] = .7500;
393 y[i] = 47.8000;
394 t[i++] = 1.0000;
395 y[i] = 32.0000;
396 t[i++] = 1.5000;
397 y[i] = 22.2000;
398 t[i++] = 2.0000;
399 y[i] = 22.5700;
400 t[i++] = 2.0000;
401 y[i] = 18.8200;
402 t[i++] = 2.5000;
403 y[i] = 13.9500;
404 t[i++] = 3.0000;
405 y[i] = 11.2500;
406 t[i++] = 4.0000;
407 y[i] = 9.0000;
408 t[i++] = 5.0000;
409 y[i] = 6.6700;
410 t[i++] = 6.0000;
411 y[i] = 75.8000;
412 t[i++] = .5000;
413 y[i] = 62.0000;
414 t[i++] = .7500;
415 y[i] = 48.8000;
416 t[i++] = 1.0000;
417 y[i] = 35.2000;
418 t[i++] = 1.5000;
419 y[i] = 20.0000;
420 t[i++] = 2.0000;
421 y[i] = 20.3200;
422 t[i++] = 2.0000;
423 y[i] = 19.3100;
424 t[i++] = 2.5000;
425 y[i] = 12.7500;
426 t[i++] = 3.0000;
427 y[i] = 10.4200;
428 t[i++] = 4.0000;
429 y[i] = 7.3100;
430 t[i++] = 5.0000;
431 y[i] = 7.4200;
432 t[i++] = 6.0000;
433 y[i] = 70.5000;
434 t[i++] = .5000;
435 y[i] = 59.5000;
436 t[i++] = .7500;
437 y[i] = 48.5000;
438 t[i++] = 1.0000;
439 y[i] = 35.8000;
440 t[i++] = 1.5000;
441 y[i] = 21.0000;
442 t[i++] = 2.0000;
443 y[i] = 21.6700;
444 t[i++] = 2.0000;
445 y[i] = 21.0000;
446 t[i++] = 2.5000;
447 y[i] = 15.6400;
448 t[i++] = 3.0000;
449 y[i] = 8.1700;
450 t[i++] = 4.0000;
451 y[i] = 8.5500;
452 t[i++] = 5.0000;
453 y[i] = 10.1200;
454 t[i++] = 6.0000;
455 y[i] = 78.0000;
456 t[i++] = .5000;
457 y[i] = 66.0000;
458 t[i++] = .6250;
459 y[i] = 62.0000;
460 t[i++] = .7500;
461 y[i] = 58.0000;
462 t[i++] = .8750;
463 y[i] = 47.7000;
464 t[i++] = 1.0000;
465 y[i] = 37.8000;
466 t[i++] = 1.2500;
467 y[i] = 20.2000;
468 t[i++] = 2.2500;
469 y[i] = 21.0700;
470 t[i++] = 2.2500;
471 y[i] = 13.8700;
472 t[i++] = 2.7500;
473 y[i] = 9.6700;
474 t[i++] = 3.2500;
475 y[i] = 7.7600;
476 t[i++] = 3.7500;
477 y[i] = 5.4400;
478 t[i++] = 4.2500;
479 y[i] = 4.8700;
480 t[i++] = 4.7500;
481 y[i] = 4.0100;
482 t[i++] = 5.2500;
483 y[i] = 3.7500;
484 t[i++] = 5.7500;
485 y[i] = 24.1900;
486 t[i++] = 3.0000;
487 y[i] = 25.7600;
488 t[i++] = 3.0000;
489 y[i] = 18.0700;
490 t[i++] = 3.0000;
491 y[i] = 11.8100;
492 t[i++] = 3.0000;
493 y[i] = 12.0700;
494 t[i++] = 3.0000;
495 y[i] = 16.1200;
496 t[i++] = 3.0000;
497 y[i] = 70.8000;
498 t[i++] = .5000;
499 y[i] = 54.7000;
500 t[i++] = .7500;
501 y[i] = 48.0000;
502 t[i++] = 1.0000;
503 y[i] = 39.8000;
504 t[i++] = 1.5000;
505 y[i] = 29.8000;
506 t[i++] = 2.0000;
507 y[i] = 23.7000;
508 t[i++] = 2.5000;
509 y[i] = 29.6200;
510 t[i++] = 2.0000;
511 y[i] = 23.8100;
512 t[i++] = 2.5000;
513 y[i] = 17.7000;
514 t[i++] = 3.0000;
515 y[i] = 11.5500;
516 t[i++] = 4.0000;
517 y[i] = 12.0700;
518 t[i++] = 5.0000;
519 y[i] = 8.7400;
520 t[i++] = 6.0000;
521 y[i] = 80.7000;
522 t[i++] = .5000;
523 y[i] = 61.3000;
524 t[i++] = .7500;
525 y[i] = 47.5000;
526 t[i++] = 1.0000;
527 y[i] = 29.0000;
528 t[i++] = 1.5000;
529 y[i] = 24.0000;
530 t[i++] = 2.0000;
531 y[i] = 17.7000;
532 t[i++] = 2.5000;
533 y[i] = 24.5600;
534 t[i++] = 2.0000;
535 y[i] = 18.6700;
536 t[i++] = 2.5000;
537 y[i] = 16.2400;
538 t[i++] = 3.0000;
539 y[i] = 8.7400;
540 t[i++] = 4.0000;
541 y[i] = 7.8700;
542 t[i++] = 5.0000;
543 y[i] = 8.5100;
544 t[i++] = 6.0000;
545 y[i] = 66.7000;
546 t[i++] = .5000;
547 y[i] = 59.2000;
548 t[i++] = .7500;
549 y[i] = 40.8000;
550 t[i++] = 1.0000;
551 y[i] = 30.7000;
552 t[i++] = 1.5000;
553 y[i] = 25.7000;
554 t[i++] = 2.0000;
555 y[i] = 16.3000;
556 t[i++] = 2.5000;
557 y[i] = 25.9900;
558 t[i++] = 2.0000;
559 y[i] = 16.9500;
560 t[i++] = 2.5000;
561 y[i] = 13.3500;
562 t[i++] = 3.0000;
563 y[i] = 8.6200;
564 t[i++] = 4.0000;
565 y[i] = 7.2000;
566 t[i++] = 5.0000;
567 y[i] = 6.6400;
568 t[i++] = 6.0000;
569 y[i] = 13.6900;
570 t[i++] = 3.0000;
571 y[i] = 81.0000;
572 t[i++] = .5000;
573 y[i] = 64.5000;
574 t[i++] = .7500;
575 y[i] = 35.5000;
576 t[i++] = 1.5000;
577 y[i] = 13.3100;
578 t[i++] = 3.0000;
579 y[i] = 4.8700;
580 t[i++] = 6.0000;
581 y[i] = 12.9400;
582 t[i++] = 3.0000;
583 y[i] = 5.0600;
584 t[i++] = 6.0000;
585 y[i] = 15.1900;
586 t[i++] = 3.0000;
587 y[i] = 14.6200;
588 t[i++] = 3.0000;
589 y[i] = 15.6400;
590 t[i++] = 3.0000;
591 y[i] = 25.5000;
592 t[i++] = 1.7500;
593 y[i] = 25.9500;
594 t[i++] = 1.7500;
595 y[i] = 81.7000;
596 t[i++] = .5000;
597 y[i] = 61.6000;
598 t[i++] = .7500;
599 y[i] = 29.8000;
600 t[i++] = 1.7500;
601 y[i] = 29.8100;
602 t[i++] = 1.7500;
603 y[i] = 17.1700;
604 t[i++] = 2.7500;
605 y[i] = 10.3900;
606 t[i++] = 3.7500;
607 y[i] = 28.4000;
608 t[i++] = 1.7500;
609 y[i] = 28.6900;
610 t[i++] = 1.7500;
611 y[i] = 81.3000;
612 t[i++] = .5000;
613 y[i] = 60.9000;
614 t[i++] = .7500;
615 y[i] = 16.6500;
616 t[i++] = 2.7500;
617 y[i] = 10.0500;
618 t[i++] = 3.7500;
619 y[i] = 28.9000;
620 t[i++] = 1.7500;
621 y[i] = 28.9500;
622 t[i++] = 1.7500;
623 PetscFunctionReturn(PETSC_SUCCESS);
624 }
625
626 /*TEST
627
628 build:
629 requires: !complex
630
631 test:
632 args: -tao_monitor_short -tao_max_it 100 -tao_type pounders -tao_pounders_delta 0.05 -tao_gatol 1.e-5
633 requires: !single
634 TODO: produces different output for many different systems
635
636 test:
637 suffix: 2
638 args: -tao_monitor_short -tao_max_it 100 -wtype 1 -tao_type pounders -tao_pounders_delta 0.05 -tao_gatol 1.e-5
639 requires: !single
640 TODO: produces different output for many different systems
641
642 test:
643 suffix: 3
644 args: -tao_monitor_short -tao_max_it 100 -wtype 2 -tao_type pounders -tao_pounders_delta 0.05 -tao_gatol 1.e-5
645 requires: !single
646 TODO: produces different output for many different systems
647
648 test:
649 suffix: 4
650 args: -tao_monitor_short -tao_max_it 100 -tao_type pounders -tao_pounders_delta 0.05 -pounders_subsolver_tao_type blmvm -tao_gatol 1.e-5
651 requires: !single
652 TODO: produces different output for many different systems
653
654 TEST*/
655