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