1! Program usage: mpiexec -n 1 chwirut1f [-help] [all TAO options] 2! 3! Description: This example demonstrates use of the TAO package to solve a 4! nonlinear least-squares problem on a single processor. We minimize the 5! Chwirut function: 6! sum_{i=0}^{n/2-1} ( alpha*(x_{2i+1}-x_{2i}^2)^2 + (1-x_{2i})^2) 7! 8! The C version of this code is chwirut1.c 9! 10 module chwirut2fmodule 11 use petscmpi ! or mpi or mpi_f08 12 use petsctao 13#include <petsc/finclude/petsctao.h> 14 PetscReal t(0:213) 15 PetscReal y(0:213) 16 PetscInt m,n 17 PetscMPIInt nn 18 PetscMPIInt rank 19 PetscMPIInt size 20 PetscMPIInt idle_tag, die_tag 21 PetscMPIInt zero,one 22 parameter (m=214) 23 parameter (n=3) 24 parameter (nn=n) 25 parameter (idle_tag=2000) 26 parameter (die_tag=3000) 27 parameter (zero=0,one=1) 28 end module chwirut2fmodule 29 30 program main 31 use chwirut2fmodule 32 33! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 34! Variable declarations 35! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 36! 37! See additional variable declarations in the file chwirut2f.h 38 39 PetscErrorCode ierr ! used to check for functions returning nonzeros 40 Vec x ! solution vector 41 Vec f ! vector of functions 42 Tao ta ! Tao context 43 44! Note: Any user-defined Fortran routines (such as FormGradient) 45! MUST be declared as external. 46 47 external FormFunction 48 49! Initialize TAO and PETSc 50 PetscCallA(PetscInitialize(ierr)) 51 PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)) 52 PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)) 53 54! Initialize problem parameters 55 call InitializeData() 56 57 if (rank .eq. 0) then 58! Allocate vectors for the solution and gradient 59 PetscCallA(VecCreateSeq(PETSC_COMM_SELF,n,x,ierr)) 60 PetscCallA(VecCreateSeq(PETSC_COMM_SELF,m,f,ierr)) 61 62! The TAO code begins here 63 64! Create TAO solver 65 PetscCallA(TaoCreate(PETSC_COMM_SELF,ta,ierr)) 66 PetscCallA(TaoSetType(ta,TAOPOUNDERS,ierr)) 67 68! Set routines for function, gradient, and hessian evaluation 69 PetscCallA(TaoSetResidualRoutine(ta,f,FormFunction,0,ierr)) 70 71! Optional: Set initial guess 72 call FormStartingPoint(x) 73 PetscCallA(TaoSetSolution(ta, x, ierr)) 74 75! Check for TAO command line options 76 PetscCallA(TaoSetFromOptions(ta,ierr)) 77! SOLVE THE APPLICATION 78 PetscCallA(TaoSolve(ta,ierr)) 79 80! Free TAO data structures 81 PetscCallA(TaoDestroy(ta,ierr)) 82 83! Free PETSc data structures 84 PetscCallA(VecDestroy(x,ierr)) 85 PetscCallA(VecDestroy(f,ierr)) 86 PetscCallA(StopWorkers(ierr)) 87 88 else 89 PetscCallA(TaskWorker(ierr)) 90 endif 91 92 PetscCallA(PetscFinalize(ierr)) 93 end 94 95! -------------------------------------------------------------------- 96! FormFunction - Evaluates the function f(X) and gradient G(X) 97! 98! Input Parameters: 99! tao - the Tao context 100! X - input vector 101! dummy - not used 102! 103! Output Parameters: 104! f - function vector 105 106 subroutine FormFunction(ta, x, f, dummy, ierr) 107 use chwirut2fmodule 108 109 Tao ta 110 Vec x,f 111 PetscErrorCode ierr 112 113 PetscInt i,checkedin 114 PetscInt finished_tasks 115 PetscMPIInt next_task 116 PetscMPIInt status(MPI_STATUS_SIZE),tag,source 117 PetscInt dummy 118 119 PetscReal, pointer :: f_v(:),x_v(:) 120 PetscReal fval(1) 121 122 ierr = 0 123 124! Get pointers to vector data 125 PetscCall(VecGetArrayRead(x,x_v,ierr)) 126 PetscCall(VecGetArray(f,f_v,ierr)) 127 128! Compute F(X) 129 if (size .eq. 1) then 130 ! Single processor 131 do i=1,m 132 PetscCall(RunSimulation(x_v,i,f_v(i),ierr)) 133 enddo 134 else 135 ! Multiprocessor main 136 next_task = zero 137 finished_tasks = 0 138 checkedin = 0 139 140 do while (finished_tasks .lt. m .or. checkedin .lt. size-1) 141 PetscCallMPI(MPI_Recv(fval,one,MPIU_SCALAR,MPI_ANY_SOURCE,MPI_ANY_TAG,PETSC_COMM_WORLD,status,ierr)) 142 tag = status(MPI_TAG) 143 source = status(MPI_SOURCE) 144 if (tag .eq. IDLE_TAG) then 145 checkedin = checkedin + 1 146 else 147 f_v(tag+1) = fval(1) 148 finished_tasks = finished_tasks + 1 149 endif 150 if (next_task .lt. m) then 151 ! Send task to worker 152 PetscCallMPI(MPI_Send(x_v,nn,MPIU_SCALAR,source,next_task,PETSC_COMM_WORLD,ierr)) 153 next_task = next_task + one 154 else 155 ! Send idle message to worker 156 PetscCallMPI(MPI_Send(x_v,nn,MPIU_SCALAR,source,IDLE_TAG,PETSC_COMM_WORLD,ierr)) 157 end if 158 enddo 159 endif 160 161! Restore vectors 162 PetscCall(VecRestoreArrayRead(x,x_v,ierr)) 163 PetscCall(VecRestoreArray(F,f_v,ierr)) 164 end 165 166 subroutine FormStartingPoint(x) 167 use chwirut2fmodule 168 169 Vec x 170 PetscReal, pointer :: x_v(:) 171 PetscErrorCode ierr 172 173 PetscCall(VecGetArray(x,x_v,ierr)) 174 x_v(1) = 0.15 175 x_v(2) = 0.008 176 x_v(3) = 0.01 177 PetscCall(VecRestoreArray(x,x_v,ierr)) 178 end 179 180 subroutine InitializeData() 181 use chwirut2fmodule 182 183 PetscInt i 184 i=0 185 y(i) = 92.9000; t(i) = 0.5000; i=i+1 186 y(i) = 78.7000; t(i) = 0.6250; i=i+1 187 y(i) = 64.2000; t(i) = 0.7500; i=i+1 188 y(i) = 64.9000; t(i) = 0.8750; i=i+1 189 y(i) = 57.1000; t(i) = 1.0000; i=i+1 190 y(i) = 43.3000; t(i) = 1.2500; i=i+1 191 y(i) = 31.1000; t(i) = 1.7500; i=i+1 192 y(i) = 23.6000; t(i) = 2.2500; i=i+1 193 y(i) = 31.0500; t(i) = 1.7500; i=i+1 194 y(i) = 23.7750; t(i) = 2.2500; i=i+1 195 y(i) = 17.7375; t(i) = 2.7500; i=i+1 196 y(i) = 13.8000; t(i) = 3.2500; i=i+1 197 y(i) = 11.5875; t(i) = 3.7500; i=i+1 198 y(i) = 9.4125; t(i) = 4.2500; i=i+1 199 y(i) = 7.7250; t(i) = 4.7500; i=i+1 200 y(i) = 7.3500; t(i) = 5.2500; i=i+1 201 y(i) = 8.0250; t(i) = 5.7500; i=i+1 202 y(i) = 90.6000; t(i) = 0.5000; i=i+1 203 y(i) = 76.9000; t(i) = 0.6250; i=i+1 204 y(i) = 71.6000; t(i) = 0.7500; i=i+1 205 y(i) = 63.6000; t(i) = 0.8750; i=i+1 206 y(i) = 54.0000; t(i) = 1.0000; i=i+1 207 y(i) = 39.2000; t(i) = 1.2500; i=i+1 208 y(i) = 29.3000; t(i) = 1.7500; i=i+1 209 y(i) = 21.4000; t(i) = 2.2500; i=i+1 210 y(i) = 29.1750; t(i) = 1.7500; i=i+1 211 y(i) = 22.1250; t(i) = 2.2500; i=i+1 212 y(i) = 17.5125; t(i) = 2.7500; i=i+1 213 y(i) = 14.2500; t(i) = 3.2500; i=i+1 214 y(i) = 9.4500; t(i) = 3.7500; i=i+1 215 y(i) = 9.1500; t(i) = 4.2500; i=i+1 216 y(i) = 7.9125; t(i) = 4.7500; i=i+1 217 y(i) = 8.4750; t(i) = 5.2500; i=i+1 218 y(i) = 6.1125; t(i) = 5.7500; i=i+1 219 y(i) = 80.0000; t(i) = 0.5000; i=i+1 220 y(i) = 79.0000; t(i) = 0.6250; i=i+1 221 y(i) = 63.8000; t(i) = 0.7500; i=i+1 222 y(i) = 57.2000; t(i) = 0.8750; i=i+1 223 y(i) = 53.2000; t(i) = 1.0000; i=i+1 224 y(i) = 42.5000; t(i) = 1.2500; i=i+1 225 y(i) = 26.8000; t(i) = 1.7500; i=i+1 226 y(i) = 20.4000; t(i) = 2.2500; i=i+1 227 y(i) = 26.8500; t(i) = 1.7500; i=i+1 228 y(i) = 21.0000; t(i) = 2.2500; i=i+1 229 y(i) = 16.4625; t(i) = 2.7500; i=i+1 230 y(i) = 12.5250; t(i) = 3.2500; i=i+1 231 y(i) = 10.5375; t(i) = 3.7500; i=i+1 232 y(i) = 8.5875; t(i) = 4.2500; i=i+1 233 y(i) = 7.1250; t(i) = 4.7500; i=i+1 234 y(i) = 6.1125; t(i) = 5.2500; i=i+1 235 y(i) = 5.9625; t(i) = 5.7500; i=i+1 236 y(i) = 74.1000; t(i) = 0.5000; i=i+1 237 y(i) = 67.3000; t(i) = 0.6250; i=i+1 238 y(i) = 60.8000; t(i) = 0.7500; i=i+1 239 y(i) = 55.5000; t(i) = 0.8750; i=i+1 240 y(i) = 50.3000; t(i) = 1.0000; i=i+1 241 y(i) = 41.0000; t(i) = 1.2500; i=i+1 242 y(i) = 29.4000; t(i) = 1.7500; i=i+1 243 y(i) = 20.4000; t(i) = 2.2500; i=i+1 244 y(i) = 29.3625; t(i) = 1.7500; i=i+1 245 y(i) = 21.1500; t(i) = 2.2500; i=i+1 246 y(i) = 16.7625; t(i) = 2.7500; i=i+1 247 y(i) = 13.2000; t(i) = 3.2500; i=i+1 248 y(i) = 10.8750; t(i) = 3.7500; i=i+1 249 y(i) = 8.1750; t(i) = 4.2500; i=i+1 250 y(i) = 7.3500; t(i) = 4.7500; i=i+1 251 y(i) = 5.9625; t(i) = 5.2500; i=i+1 252 y(i) = 5.6250; t(i) = 5.7500; i=i+1 253 y(i) = 81.5000; t(i) = .5000; i=i+1 254 y(i) = 62.4000; t(i) = .7500; i=i+1 255 y(i) = 32.5000; t(i) = 1.5000; i=i+1 256 y(i) = 12.4100; t(i) = 3.0000; i=i+1 257 y(i) = 13.1200; t(i) = 3.0000; i=i+1 258 y(i) = 15.5600; t(i) = 3.0000; i=i+1 259 y(i) = 5.6300; t(i) = 6.0000; i=i+1 260 y(i) = 78.0000; t(i) = .5000; i=i+1 261 y(i) = 59.9000; t(i) = .7500; i=i+1 262 y(i) = 33.2000; t(i) = 1.5000; i=i+1 263 y(i) = 13.8400; t(i) = 3.0000; i=i+1 264 y(i) = 12.7500; t(i) = 3.0000; i=i+1 265 y(i) = 14.6200; t(i) = 3.0000; i=i+1 266 y(i) = 3.9400; t(i) = 6.0000; i=i+1 267 y(i) = 76.8000; t(i) = .5000; i=i+1 268 y(i) = 61.0000; t(i) = .7500; i=i+1 269 y(i) = 32.9000; t(i) = 1.5000; i=i+1 270 y(i) = 13.8700; t(i) = 3.0000; i=i+1 271 y(i) = 11.8100; t(i) = 3.0000; i=i+1 272 y(i) = 13.3100; t(i) = 3.0000; i=i+1 273 y(i) = 5.4400; t(i) = 6.0000; i=i+1 274 y(i) = 78.0000; t(i) = .5000; i=i+1 275 y(i) = 63.5000; t(i) = .7500; i=i+1 276 y(i) = 33.8000; t(i) = 1.5000; i=i+1 277 y(i) = 12.5600; t(i) = 3.0000; i=i+1 278 y(i) = 5.6300; t(i) = 6.0000; i=i+1 279 y(i) = 12.7500; t(i) = 3.0000; i=i+1 280 y(i) = 13.1200; t(i) = 3.0000; i=i+1 281 y(i) = 5.4400; t(i) = 6.0000; i=i+1 282 y(i) = 76.8000; t(i) = .5000; i=i+1 283 y(i) = 60.0000; t(i) = .7500; i=i+1 284 y(i) = 47.8000; t(i) = 1.0000; i=i+1 285 y(i) = 32.0000; t(i) = 1.5000; i=i+1 286 y(i) = 22.2000; t(i) = 2.0000; i=i+1 287 y(i) = 22.5700; t(i) = 2.0000; i=i+1 288 y(i) = 18.8200; t(i) = 2.5000; i=i+1 289 y(i) = 13.9500; t(i) = 3.0000; i=i+1 290 y(i) = 11.2500; t(i) = 4.0000; i=i+1 291 y(i) = 9.0000; t(i) = 5.0000; i=i+1 292 y(i) = 6.6700; t(i) = 6.0000; i=i+1 293 y(i) = 75.8000; t(i) = .5000; i=i+1 294 y(i) = 62.0000; t(i) = .7500; i=i+1 295 y(i) = 48.8000; t(i) = 1.0000; i=i+1 296 y(i) = 35.2000; t(i) = 1.5000; i=i+1 297 y(i) = 20.0000; t(i) = 2.0000; i=i+1 298 y(i) = 20.3200; t(i) = 2.0000; i=i+1 299 y(i) = 19.3100; t(i) = 2.5000; i=i+1 300 y(i) = 12.7500; t(i) = 3.0000; i=i+1 301 y(i) = 10.4200; t(i) = 4.0000; i=i+1 302 y(i) = 7.3100; t(i) = 5.0000; i=i+1 303 y(i) = 7.4200; t(i) = 6.0000; i=i+1 304 y(i) = 70.5000; t(i) = .5000; i=i+1 305 y(i) = 59.5000; t(i) = .7500; i=i+1 306 y(i) = 48.5000; t(i) = 1.0000; i=i+1 307 y(i) = 35.8000; t(i) = 1.5000; i=i+1 308 y(i) = 21.0000; t(i) = 2.0000; i=i+1 309 y(i) = 21.6700; t(i) = 2.0000; i=i+1 310 y(i) = 21.0000; t(i) = 2.5000; i=i+1 311 y(i) = 15.6400; t(i) = 3.0000; i=i+1 312 y(i) = 8.1700; t(i) = 4.0000; i=i+1 313 y(i) = 8.5500; t(i) = 5.0000; i=i+1 314 y(i) = 10.1200; t(i) = 6.0000; i=i+1 315 y(i) = 78.0000; t(i) = .5000; i=i+1 316 y(i) = 66.0000; t(i) = .6250; i=i+1 317 y(i) = 62.0000; t(i) = .7500; i=i+1 318 y(i) = 58.0000; t(i) = .8750; i=i+1 319 y(i) = 47.7000; t(i) = 1.0000; i=i+1 320 y(i) = 37.8000; t(i) = 1.2500; i=i+1 321 y(i) = 20.2000; t(i) = 2.2500; i=i+1 322 y(i) = 21.0700; t(i) = 2.2500; i=i+1 323 y(i) = 13.8700; t(i) = 2.7500; i=i+1 324 y(i) = 9.6700; t(i) = 3.2500; i=i+1 325 y(i) = 7.7600; t(i) = 3.7500; i=i+1 326 y(i) = 5.4400; t(i) = 4.2500; i=i+1 327 y(i) = 4.8700; t(i) = 4.7500; i=i+1 328 y(i) = 4.0100; t(i) = 5.2500; i=i+1 329 y(i) = 3.7500; t(i) = 5.7500; i=i+1 330 y(i) = 24.1900; t(i) = 3.0000; i=i+1 331 y(i) = 25.7600; t(i) = 3.0000; i=i+1 332 y(i) = 18.0700; t(i) = 3.0000; i=i+1 333 y(i) = 11.8100; t(i) = 3.0000; i=i+1 334 y(i) = 12.0700; t(i) = 3.0000; i=i+1 335 y(i) = 16.1200; t(i) = 3.0000; i=i+1 336 y(i) = 70.8000; t(i) = .5000; i=i+1 337 y(i) = 54.7000; t(i) = .7500; i=i+1 338 y(i) = 48.0000; t(i) = 1.0000; i=i+1 339 y(i) = 39.8000; t(i) = 1.5000; i=i+1 340 y(i) = 29.8000; t(i) = 2.0000; i=i+1 341 y(i) = 23.7000; t(i) = 2.5000; i=i+1 342 y(i) = 29.6200; t(i) = 2.0000; i=i+1 343 y(i) = 23.8100; t(i) = 2.5000; i=i+1 344 y(i) = 17.7000; t(i) = 3.0000; i=i+1 345 y(i) = 11.5500; t(i) = 4.0000; i=i+1 346 y(i) = 12.0700; t(i) = 5.0000; i=i+1 347 y(i) = 8.7400; t(i) = 6.0000; i=i+1 348 y(i) = 80.7000; t(i) = .5000; i=i+1 349 y(i) = 61.3000; t(i) = .7500; i=i+1 350 y(i) = 47.5000; t(i) = 1.0000; i=i+1 351 y(i) = 29.0000; t(i) = 1.5000; i=i+1 352 y(i) = 24.0000; t(i) = 2.0000; i=i+1 353 y(i) = 17.7000; t(i) = 2.5000; i=i+1 354 y(i) = 24.5600; t(i) = 2.0000; i=i+1 355 y(i) = 18.6700; t(i) = 2.5000; i=i+1 356 y(i) = 16.2400; t(i) = 3.0000; i=i+1 357 y(i) = 8.7400; t(i) = 4.0000; i=i+1 358 y(i) = 7.8700; t(i) = 5.0000; i=i+1 359 y(i) = 8.5100; t(i) = 6.0000; i=i+1 360 y(i) = 66.7000; t(i) = .5000; i=i+1 361 y(i) = 59.2000; t(i) = .7500; i=i+1 362 y(i) = 40.8000; t(i) = 1.0000; i=i+1 363 y(i) = 30.7000; t(i) = 1.5000; i=i+1 364 y(i) = 25.7000; t(i) = 2.0000; i=i+1 365 y(i) = 16.3000; t(i) = 2.5000; i=i+1 366 y(i) = 25.9900; t(i) = 2.0000; i=i+1 367 y(i) = 16.9500; t(i) = 2.5000; i=i+1 368 y(i) = 13.3500; t(i) = 3.0000; i=i+1 369 y(i) = 8.6200; t(i) = 4.0000; i=i+1 370 y(i) = 7.2000; t(i) = 5.0000; i=i+1 371 y(i) = 6.6400; t(i) = 6.0000; i=i+1 372 y(i) = 13.6900; t(i) = 3.0000; i=i+1 373 y(i) = 81.0000; t(i) = .5000; i=i+1 374 y(i) = 64.5000; t(i) = .7500; i=i+1 375 y(i) = 35.5000; t(i) = 1.5000; i=i+1 376 y(i) = 13.3100; t(i) = 3.0000; i=i+1 377 y(i) = 4.8700; t(i) = 6.0000; i=i+1 378 y(i) = 12.9400; t(i) = 3.0000; i=i+1 379 y(i) = 5.0600; t(i) = 6.0000; i=i+1 380 y(i) = 15.1900; t(i) = 3.0000; i=i+1 381 y(i) = 14.6200; t(i) = 3.0000; i=i+1 382 y(i) = 15.6400; t(i) = 3.0000; i=i+1 383 y(i) = 25.5000; t(i) = 1.7500; i=i+1 384 y(i) = 25.9500; t(i) = 1.7500; i=i+1 385 y(i) = 81.7000; t(i) = .5000; i=i+1 386 y(i) = 61.6000; t(i) = .7500; i=i+1 387 y(i) = 29.8000; t(i) = 1.7500; i=i+1 388 y(i) = 29.8100; t(i) = 1.7500; i=i+1 389 y(i) = 17.1700; t(i) = 2.7500; i=i+1 390 y(i) = 10.3900; t(i) = 3.7500; i=i+1 391 y(i) = 28.4000; t(i) = 1.7500; i=i+1 392 y(i) = 28.6900; t(i) = 1.7500; i=i+1 393 y(i) = 81.3000; t(i) = .5000; i=i+1 394 y(i) = 60.9000; t(i) = .7500; i=i+1 395 y(i) = 16.6500; t(i) = 2.7500; i=i+1 396 y(i) = 10.0500; t(i) = 3.7500; i=i+1 397 y(i) = 28.9000; t(i) = 1.7500; i=i+1 398 y(i) = 28.9500; t(i) = 1.7500; i=i+1 399 400 end 401 402 subroutine TaskWorker(ierr) 403 use chwirut2fmodule 404 405 PetscErrorCode ierr 406 PetscReal x(n),f(1) 407 PetscMPIInt tag 408 PetscInt index 409 PetscMPIInt status(MPI_STATUS_SIZE) 410 411 tag = IDLE_TAG 412 f = 0.0 413 ! Send check-in message to rank-0 414 PetscCallMPI(MPI_Send(f,one,MPIU_SCALAR,zero,IDLE_TAG,PETSC_COMM_WORLD,ierr)) 415 do while (tag .ne. DIE_TAG) 416 PetscCallMPI(MPI_Recv(x,nn,MPIU_SCALAR,zero,MPI_ANY_TAG,PETSC_COMM_WORLD,status,ierr)) 417 tag = status(MPI_TAG) 418 if (tag .eq. IDLE_TAG) then 419 PetscCallMPI(MPI_Send(f,one,MPIU_SCALAR,zero,IDLE_TAG,PETSC_COMM_WORLD,ierr)) 420 else if (tag .ne. DIE_TAG) then 421 index = tag 422 ! Compute local part of residual 423 PetscCall(RunSimulation(x,index,f(1),ierr)) 424 425 ! Return residual to rank-0 426 PetscCallMPI(MPI_Send(f,one,MPIU_SCALAR,zero,tag,PETSC_COMM_WORLD,ierr)) 427 end if 428 enddo 429 ierr = 0 430 end 431 432 subroutine RunSimulation(x,i,f,ierr) 433 use chwirut2fmodule 434 435 PetscReal x(n),f 436 PetscInt i 437 PetscErrorCode ierr 438 f = y(i) - exp(-x(1)*t(i))/(x(2)+x(3)*t(i)) 439 ierr = 0 440 end 441 442 subroutine StopWorkers(ierr) 443 use chwirut2fmodule 444 445 integer checkedin 446 PetscMPIInt status(MPI_STATUS_SIZE) 447 PetscMPIInt source 448 PetscReal f(1),x(n) 449 PetscErrorCode ierr 450 PetscInt i 451 452 checkedin=0 453 do while (checkedin .lt. size-1) 454 PetscCallMPI(MPI_Recv(f,one,MPIU_SCALAR,MPI_ANY_SOURCE,MPI_ANY_TAG,PETSC_COMM_WORLD,status,ierr)) 455 checkedin=checkedin+1 456 source = status(MPI_SOURCE) 457 do i=1,n 458 x(i) = 0.0 459 enddo 460 PetscCallMPI(MPI_Send(x,nn,MPIU_SCALAR,source,DIE_TAG,PETSC_COMM_WORLD,ierr)) 461 enddo 462 ierr = 0 463 end 464 465!/*TEST 466! 467! build: 468! requires: !complex 469! 470! test: 471! nsize: 3 472! args: -tao_monitor_short -tao_max_it 100 -tao_type pounders -tao_gatol 1.e-5 473! requires: !single 474! 475! 476!TEST*/ 477