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 test_chwirut1.c 9! 10 11! 12! ---------------------------------------------------------------------- 13! 14module chwirut1fmodule 15#include <petsc/finclude/petsctao.h> 16 use petsctao 17 implicit none 18 19 PetscReal t(0:213) 20 PetscReal y(0:213) 21 PetscInt m, n 22end module chwirut1fmodule 23 24program main 25 use chwirut1fmodule 26! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 27! Variable declarations 28! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 29! 30! See additional variable declarations in the file chwirut1f.h 31 32 PetscErrorCode ierr ! used to check for functions returning nonzeros 33 Vec x ! solution vector 34 Vec f ! vector of functions 35 Tao ta ! Tao context 36 PetscInt nhist 37 PetscMPIInt size, rank ! number of processes running 38 PetscReal hist(100) ! objective value history 39 PetscReal resid(100)! residual history 40 PetscReal cnorm(100)! cnorm history 41 PetscInt lits(100) ! #ksp history 42 PetscInt oh 43 TaoConvergedReason reason 44 45! Note: Any user-defined Fortran routines (such as FormGradient) 46! MUST be declared as external. 47 48 external FormFunction 49 50! Initialize TAO and PETSc 51 PetscCallA(PetscInitialize(ierr)) 52 53 PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, size, ierr)) 54 PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr)) 55 PetscCheckA(size == 1, PETSC_COMM_SELF, PETSC_ERR_WRONG_MPI_SIZE, 'This is a uniprocessor example only') 56 57! Initialize problem parameters 58 m = 214 59 n = 3 60 61! Allocate vectors for the solution and gradient 62 PetscCallA(VecCreateSeq(PETSC_COMM_SELF, n, x, ierr)) 63 PetscCallA(VecCreateSeq(PETSC_COMM_SELF, m, f, ierr)) 64 65! The TAO code begins here 66 67! Create TAO solver 68 PetscCallA(TaoCreate(PETSC_COMM_SELF, ta, ierr)) 69 PetscCallA(TaoSetType(ta, TAOPOUNDERS, ierr)) 70! Set routines for function, gradient, and hessian evaluation 71 72 PetscCallA(TaoSetResidualRoutine(ta, f, FormFunction, 0, ierr)) 73 74! Optional: Set initial guess 75 call InitializeData() 76 call FormStartingPoint(x) 77 PetscCallA(TaoSetSolution(ta, x, ierr)) 78 79! Check for TAO command line options 80 PetscCallA(TaoSetFromOptions(ta, ierr)) 81 oh = 100 82 PetscCallA(TaoSetConvergenceHistory(ta, hist, resid, cnorm, lits, oh, PETSC_TRUE, ierr)) 83! SOLVE THE APPLICATION 84 PetscCallA(TaoSolve(ta, ierr)) 85 PetscCallA(TaoGetConvergenceHistory(ta, nhist, ierr)) 86 PetscCallA(TaoGetConvergedReason(ta, reason, ierr)) 87 if (reason%v <= 0) then 88 print *, 'Tao failed.' 89 print *, 'Try a different TAO method, adjust some parameters,' 90 print *, 'or check the function evaluation routines.' 91 end if 92 93! Free TAO data structures 94 PetscCallA(TaoDestroy(ta, ierr)) 95 96! Free PETSc data structures 97 PetscCallA(VecDestroy(x, ierr)) 98 PetscCallA(VecDestroy(f, ierr)) 99 100 PetscCallA(PetscFinalize(ierr)) 101 102end 103 104! -------------------------------------------------------------------- 105! FormFunction - Evaluates the function f(X) and gradient G(X) 106! 107! Input Parameters: 108! tao - the Tao context 109! X - input vector 110! dummy - not used 111! 112! Output Parameters: 113! f - function vector 114 115subroutine FormFunction(ta, x, f, dummy, ierr) 116 use chwirut1fmodule 117 118 Tao ta 119 Vec x, f 120 PetscErrorCode ierr 121 PetscInt dummy 122 123 PetscInt i 124 PetscScalar, pointer, dimension(:) :: x_v, f_v 125 126 ierr = 0 127 128! Get pointers to vector data 129 PetscCall(VecGetArray(x, x_v, ierr)) 130 PetscCall(VecGetArray(f, f_v, ierr)) 131 132! Compute F(X) 133 do i = 0, m - 1 134 f_v(i + 1) = y(i) - exp(-x_v(1)*t(i))/(x_v(2) + x_v(3)*t(i)) 135 end do 136 137! Restore vectors 138 PetscCall(VecRestoreArray(X, x_v, ierr)) 139 PetscCall(VecRestoreArray(F, f_v, ierr)) 140 141end 142 143subroutine FormStartingPoint(x) 144 use chwirut1fmodule 145 146 Vec x 147 PetscScalar, pointer, dimension(:) :: x_v 148 PetscErrorCode ierr 149 150 PetscCall(VecGetArray(x, x_v, ierr)) 151 x_v(1) = 0.15 152 x_v(2) = 0.008 153 x_v(3) = 0.01 154 PetscCall(VecRestoreArray(x, x_v, ierr)) 155end 156 157subroutine InitializeData() 158 use chwirut1fmodule 159 160 integer i 161 i = 0 162 y(i) = 92.9000; t(i) = 0.5000; i = i + 1 163 y(i) = 78.7000; t(i) = 0.6250; i = i + 1 164 y(i) = 64.2000; t(i) = 0.7500; i = i + 1 165 y(i) = 64.9000; t(i) = 0.8750; i = i + 1 166 y(i) = 57.1000; t(i) = 1.0000; i = i + 1 167 y(i) = 43.3000; t(i) = 1.2500; i = i + 1 168 y(i) = 31.1000; t(i) = 1.7500; i = i + 1 169 y(i) = 23.6000; t(i) = 2.2500; i = i + 1 170 y(i) = 31.0500; t(i) = 1.7500; i = i + 1 171 y(i) = 23.7750; t(i) = 2.2500; i = i + 1 172 y(i) = 17.7375; t(i) = 2.7500; i = i + 1 173 y(i) = 13.8000; t(i) = 3.2500; i = i + 1 174 y(i) = 11.5875; t(i) = 3.7500; i = i + 1 175 y(i) = 9.4125; t(i) = 4.2500; i = i + 1 176 y(i) = 7.7250; t(i) = 4.7500; i = i + 1 177 y(i) = 7.3500; t(i) = 5.2500; i = i + 1 178 y(i) = 8.0250; t(i) = 5.7500; i = i + 1 179 y(i) = 90.6000; t(i) = 0.5000; i = i + 1 180 y(i) = 76.9000; t(i) = 0.6250; i = i + 1 181 y(i) = 71.6000; t(i) = 0.7500; i = i + 1 182 y(i) = 63.6000; t(i) = 0.8750; i = i + 1 183 y(i) = 54.0000; t(i) = 1.0000; i = i + 1 184 y(i) = 39.2000; t(i) = 1.2500; i = i + 1 185 y(i) = 29.3000; t(i) = 1.7500; i = i + 1 186 y(i) = 21.4000; t(i) = 2.2500; i = i + 1 187 y(i) = 29.1750; t(i) = 1.7500; i = i + 1 188 y(i) = 22.1250; t(i) = 2.2500; i = i + 1 189 y(i) = 17.5125; t(i) = 2.7500; i = i + 1 190 y(i) = 14.2500; t(i) = 3.2500; i = i + 1 191 y(i) = 9.4500; t(i) = 3.7500; i = i + 1 192 y(i) = 9.1500; t(i) = 4.2500; i = i + 1 193 y(i) = 7.9125; t(i) = 4.7500; i = i + 1 194 y(i) = 8.4750; t(i) = 5.2500; i = i + 1 195 y(i) = 6.1125; t(i) = 5.7500; i = i + 1 196 y(i) = 80.0000; t(i) = 0.5000; i = i + 1 197 y(i) = 79.0000; t(i) = 0.6250; i = i + 1 198 y(i) = 63.8000; t(i) = 0.7500; i = i + 1 199 y(i) = 57.2000; t(i) = 0.8750; i = i + 1 200 y(i) = 53.2000; t(i) = 1.0000; i = i + 1 201 y(i) = 42.5000; t(i) = 1.2500; i = i + 1 202 y(i) = 26.8000; t(i) = 1.7500; i = i + 1 203 y(i) = 20.4000; t(i) = 2.2500; i = i + 1 204 y(i) = 26.8500; t(i) = 1.7500; i = i + 1 205 y(i) = 21.0000; t(i) = 2.2500; i = i + 1 206 y(i) = 16.4625; t(i) = 2.7500; i = i + 1 207 y(i) = 12.5250; t(i) = 3.2500; i = i + 1 208 y(i) = 10.5375; t(i) = 3.7500; i = i + 1 209 y(i) = 8.5875; t(i) = 4.2500; i = i + 1 210 y(i) = 7.1250; t(i) = 4.7500; i = i + 1 211 y(i) = 6.1125; t(i) = 5.2500; i = i + 1 212 y(i) = 5.9625; t(i) = 5.7500; i = i + 1 213 y(i) = 74.1000; t(i) = 0.5000; i = i + 1 214 y(i) = 67.3000; t(i) = 0.6250; i = i + 1 215 y(i) = 60.8000; t(i) = 0.7500; i = i + 1 216 y(i) = 55.5000; t(i) = 0.8750; i = i + 1 217 y(i) = 50.3000; t(i) = 1.0000; i = i + 1 218 y(i) = 41.0000; t(i) = 1.2500; i = i + 1 219 y(i) = 29.4000; t(i) = 1.7500; i = i + 1 220 y(i) = 20.4000; t(i) = 2.2500; i = i + 1 221 y(i) = 29.3625; t(i) = 1.7500; i = i + 1 222 y(i) = 21.1500; t(i) = 2.2500; i = i + 1 223 y(i) = 16.7625; t(i) = 2.7500; i = i + 1 224 y(i) = 13.2000; t(i) = 3.2500; i = i + 1 225 y(i) = 10.8750; t(i) = 3.7500; i = i + 1 226 y(i) = 8.1750; t(i) = 4.2500; i = i + 1 227 y(i) = 7.3500; t(i) = 4.7500; i = i + 1 228 y(i) = 5.9625; t(i) = 5.2500; i = i + 1 229 y(i) = 5.6250; t(i) = 5.7500; i = i + 1 230 y(i) = 81.5000; t(i) = .5000; i = i + 1 231 y(i) = 62.4000; t(i) = .7500; i = i + 1 232 y(i) = 32.5000; t(i) = 1.5000; i = i + 1 233 y(i) = 12.4100; t(i) = 3.0000; i = i + 1 234 y(i) = 13.1200; t(i) = 3.0000; i = i + 1 235 y(i) = 15.5600; t(i) = 3.0000; i = i + 1 236 y(i) = 5.6300; t(i) = 6.0000; i = i + 1 237 y(i) = 78.0000; t(i) = .5000; i = i + 1 238 y(i) = 59.9000; t(i) = .7500; i = i + 1 239 y(i) = 33.2000; t(i) = 1.5000; i = i + 1 240 y(i) = 13.8400; t(i) = 3.0000; i = i + 1 241 y(i) = 12.7500; t(i) = 3.0000; i = i + 1 242 y(i) = 14.6200; t(i) = 3.0000; i = i + 1 243 y(i) = 3.9400; t(i) = 6.0000; i = i + 1 244 y(i) = 76.8000; t(i) = .5000; i = i + 1 245 y(i) = 61.0000; t(i) = .7500; i = i + 1 246 y(i) = 32.9000; t(i) = 1.5000; i = i + 1 247 y(i) = 13.8700; t(i) = 3.0000; i = i + 1 248 y(i) = 11.8100; t(i) = 3.0000; i = i + 1 249 y(i) = 13.3100; t(i) = 3.0000; i = i + 1 250 y(i) = 5.4400; t(i) = 6.0000; i = i + 1 251 y(i) = 78.0000; t(i) = .5000; i = i + 1 252 y(i) = 63.5000; t(i) = .7500; i = i + 1 253 y(i) = 33.8000; t(i) = 1.5000; i = i + 1 254 y(i) = 12.5600; t(i) = 3.0000; i = i + 1 255 y(i) = 5.6300; t(i) = 6.0000; i = i + 1 256 y(i) = 12.7500; t(i) = 3.0000; i = i + 1 257 y(i) = 13.1200; t(i) = 3.0000; i = i + 1 258 y(i) = 5.4400; t(i) = 6.0000; i = i + 1 259 y(i) = 76.8000; t(i) = .5000; i = i + 1 260 y(i) = 60.0000; t(i) = .7500; i = i + 1 261 y(i) = 47.8000; t(i) = 1.0000; i = i + 1 262 y(i) = 32.0000; t(i) = 1.5000; i = i + 1 263 y(i) = 22.2000; t(i) = 2.0000; i = i + 1 264 y(i) = 22.5700; t(i) = 2.0000; i = i + 1 265 y(i) = 18.8200; t(i) = 2.5000; i = i + 1 266 y(i) = 13.9500; t(i) = 3.0000; i = i + 1 267 y(i) = 11.2500; t(i) = 4.0000; i = i + 1 268 y(i) = 9.0000; t(i) = 5.0000; i = i + 1 269 y(i) = 6.6700; t(i) = 6.0000; i = i + 1 270 y(i) = 75.8000; t(i) = .5000; i = i + 1 271 y(i) = 62.0000; t(i) = .7500; i = i + 1 272 y(i) = 48.8000; t(i) = 1.0000; i = i + 1 273 y(i) = 35.2000; t(i) = 1.5000; i = i + 1 274 y(i) = 20.0000; t(i) = 2.0000; i = i + 1 275 y(i) = 20.3200; t(i) = 2.0000; i = i + 1 276 y(i) = 19.3100; t(i) = 2.5000; i = i + 1 277 y(i) = 12.7500; t(i) = 3.0000; i = i + 1 278 y(i) = 10.4200; t(i) = 4.0000; i = i + 1 279 y(i) = 7.3100; t(i) = 5.0000; i = i + 1 280 y(i) = 7.4200; t(i) = 6.0000; i = i + 1 281 y(i) = 70.5000; t(i) = .5000; i = i + 1 282 y(i) = 59.5000; t(i) = .7500; i = i + 1 283 y(i) = 48.5000; t(i) = 1.0000; i = i + 1 284 y(i) = 35.8000; t(i) = 1.5000; i = i + 1 285 y(i) = 21.0000; t(i) = 2.0000; i = i + 1 286 y(i) = 21.6700; t(i) = 2.0000; i = i + 1 287 y(i) = 21.0000; t(i) = 2.5000; i = i + 1 288 y(i) = 15.6400; t(i) = 3.0000; i = i + 1 289 y(i) = 8.1700; t(i) = 4.0000; i = i + 1 290 y(i) = 8.5500; t(i) = 5.0000; i = i + 1 291 y(i) = 10.1200; t(i) = 6.0000; i = i + 1 292 y(i) = 78.0000; t(i) = .5000; i = i + 1 293 y(i) = 66.0000; t(i) = .6250; i = i + 1 294 y(i) = 62.0000; t(i) = .7500; i = i + 1 295 y(i) = 58.0000; t(i) = .8750; i = i + 1 296 y(i) = 47.7000; t(i) = 1.0000; i = i + 1 297 y(i) = 37.8000; t(i) = 1.2500; i = i + 1 298 y(i) = 20.2000; t(i) = 2.2500; i = i + 1 299 y(i) = 21.0700; t(i) = 2.2500; i = i + 1 300 y(i) = 13.8700; t(i) = 2.7500; i = i + 1 301 y(i) = 9.6700; t(i) = 3.2500; i = i + 1 302 y(i) = 7.7600; t(i) = 3.7500; i = i + 1 303 y(i) = 5.4400; t(i) = 4.2500; i = i + 1 304 y(i) = 4.8700; t(i) = 4.7500; i = i + 1 305 y(i) = 4.0100; t(i) = 5.2500; i = i + 1 306 y(i) = 3.7500; t(i) = 5.7500; i = i + 1 307 y(i) = 24.1900; t(i) = 3.0000; i = i + 1 308 y(i) = 25.7600; t(i) = 3.0000; i = i + 1 309 y(i) = 18.0700; t(i) = 3.0000; i = i + 1 310 y(i) = 11.8100; t(i) = 3.0000; i = i + 1 311 y(i) = 12.0700; t(i) = 3.0000; i = i + 1 312 y(i) = 16.1200; t(i) = 3.0000; i = i + 1 313 y(i) = 70.8000; t(i) = .5000; i = i + 1 314 y(i) = 54.7000; t(i) = .7500; i = i + 1 315 y(i) = 48.0000; t(i) = 1.0000; i = i + 1 316 y(i) = 39.8000; t(i) = 1.5000; i = i + 1 317 y(i) = 29.8000; t(i) = 2.0000; i = i + 1 318 y(i) = 23.7000; t(i) = 2.5000; i = i + 1 319 y(i) = 29.6200; t(i) = 2.0000; i = i + 1 320 y(i) = 23.8100; t(i) = 2.5000; i = i + 1 321 y(i) = 17.7000; t(i) = 3.0000; i = i + 1 322 y(i) = 11.5500; t(i) = 4.0000; i = i + 1 323 y(i) = 12.0700; t(i) = 5.0000; i = i + 1 324 y(i) = 8.7400; t(i) = 6.0000; i = i + 1 325 y(i) = 80.7000; t(i) = .5000; i = i + 1 326 y(i) = 61.3000; t(i) = .7500; i = i + 1 327 y(i) = 47.5000; t(i) = 1.0000; i = i + 1 328 y(i) = 29.0000; t(i) = 1.5000; i = i + 1 329 y(i) = 24.0000; t(i) = 2.0000; i = i + 1 330 y(i) = 17.7000; t(i) = 2.5000; i = i + 1 331 y(i) = 24.5600; t(i) = 2.0000; i = i + 1 332 y(i) = 18.6700; t(i) = 2.5000; i = i + 1 333 y(i) = 16.2400; t(i) = 3.0000; i = i + 1 334 y(i) = 8.7400; t(i) = 4.0000; i = i + 1 335 y(i) = 7.8700; t(i) = 5.0000; i = i + 1 336 y(i) = 8.5100; t(i) = 6.0000; i = i + 1 337 y(i) = 66.7000; t(i) = .5000; i = i + 1 338 y(i) = 59.2000; t(i) = .7500; i = i + 1 339 y(i) = 40.8000; t(i) = 1.0000; i = i + 1 340 y(i) = 30.7000; t(i) = 1.5000; i = i + 1 341 y(i) = 25.7000; t(i) = 2.0000; i = i + 1 342 y(i) = 16.3000; t(i) = 2.5000; i = i + 1 343 y(i) = 25.9900; t(i) = 2.0000; i = i + 1 344 y(i) = 16.9500; t(i) = 2.5000; i = i + 1 345 y(i) = 13.3500; t(i) = 3.0000; i = i + 1 346 y(i) = 8.6200; t(i) = 4.0000; i = i + 1 347 y(i) = 7.2000; t(i) = 5.0000; i = i + 1 348 y(i) = 6.6400; t(i) = 6.0000; i = i + 1 349 y(i) = 13.6900; t(i) = 3.0000; i = i + 1 350 y(i) = 81.0000; t(i) = .5000; i = i + 1 351 y(i) = 64.5000; t(i) = .7500; i = i + 1 352 y(i) = 35.5000; t(i) = 1.5000; i = i + 1 353 y(i) = 13.3100; t(i) = 3.0000; i = i + 1 354 y(i) = 4.8700; t(i) = 6.0000; i = i + 1 355 y(i) = 12.9400; t(i) = 3.0000; i = i + 1 356 y(i) = 5.0600; t(i) = 6.0000; i = i + 1 357 y(i) = 15.1900; t(i) = 3.0000; i = i + 1 358 y(i) = 14.6200; t(i) = 3.0000; i = i + 1 359 y(i) = 15.6400; t(i) = 3.0000; i = i + 1 360 y(i) = 25.5000; t(i) = 1.7500; i = i + 1 361 y(i) = 25.9500; t(i) = 1.7500; i = i + 1 362 y(i) = 81.7000; t(i) = .5000; i = i + 1 363 y(i) = 61.6000; t(i) = .7500; i = i + 1 364 y(i) = 29.8000; t(i) = 1.7500; i = i + 1 365 y(i) = 29.8100; t(i) = 1.7500; i = i + 1 366 y(i) = 17.1700; t(i) = 2.7500; i = i + 1 367 y(i) = 10.3900; t(i) = 3.7500; i = i + 1 368 y(i) = 28.4000; t(i) = 1.7500; i = i + 1 369 y(i) = 28.6900; t(i) = 1.7500; i = i + 1 370 y(i) = 81.3000; t(i) = .5000; i = i + 1 371 y(i) = 60.9000; t(i) = .7500; i = i + 1 372 y(i) = 16.6500; t(i) = 2.7500; i = i + 1 373 y(i) = 10.0500; t(i) = 3.7500; i = i + 1 374 y(i) = 28.9000; t(i) = 1.7500; i = i + 1 375 y(i) = 28.9500; t(i) = 1.7500; i = i + 1 376 377end 378 379!/*TEST 380! 381! build: 382! requires: !complex 383! 384! test: 385! args: -tao_monitor_short -tao_max_it 100 -tao_type pounders -tao_gatol 1.e-5 386! requires: !single 387! 388!TEST*/ 389