1! Program usage: mpiexec -n 1 rosenbrock1f [-help] [all TAO options] 2! 3! Description: This example demonstrates use of the TAO package to solve an 4! unconstrained minimization problem on a single processor. We minimize the 5! extended Rosenbrock 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 rosenbrock1.c 9! 10 11! 12 13! ---------------------------------------------------------------------- 14! 15#include "petsc/finclude/petsctao.h" 16module rosenbrock1fmodule 17 use petsctao 18 implicit none 19contains 20! -------------------------------------------------------------------- 21! FormFunctionGradient - Evaluates the function f(X) and gradient G(X) 22! 23! Input Parameters: 24! tao - the Tao context 25! X - input vector 26! dummy - not used 27! 28! Output Parameters: 29! G - vector containing the newly evaluated gradient 30! f - function value 31 32 subroutine FormFunctionGradient(ta, X, f, G, dummy, ierr) 33 type(tTao) ta 34 type(tVec) X, G 35 PetscReal f 36 PetscErrorCode ierr 37 PetscInt dummy 38 39 PetscReal ff, t1, t2 40 PetscInt i, nn 41 PetscReal, pointer :: g_v(:), x_v(:) 42 PetscReal alpha 43 PetscInt n 44 common/params/alpha, n 45 46 ierr = 0 47 nn = n/2 48 ff = 0 49 50! Get pointers to vector data 51 PetscCall(VecGetArrayRead(X, x_v, ierr)) 52 PetscCall(VecGetArray(G, g_v, ierr)) 53 54! Compute G(X) 55 do i = 0, nn - 1 56 t1 = x_v(1 + 2*i + 1) - x_v(1 + 2*i)*x_v(1 + 2*i) 57 t2 = 1.0 - x_v(1 + 2*i) 58 ff = ff + alpha*t1*t1 + t2*t2 59 g_v(1 + 2*i) = -4*alpha*t1*x_v(1 + 2*i) - 2.0*t2 60 g_v(1 + 2*i + 1) = 2.0*alpha*t1 61 end do 62 63! Restore vectors 64 PetscCall(VecRestoreArrayRead(X, x_v, ierr)) 65 PetscCall(VecRestoreArray(G, g_v, ierr)) 66 67 f = ff 68 PetscCall(PetscLogFlops(15.0d0*nn, ierr)) 69 70 end 71 72! 73! --------------------------------------------------------------------- 74! 75! FormHessian - Evaluates Hessian matrix. 76! 77! Input Parameters: 78! tao - the Tao context 79! X - input vector 80! dummy - optional user-defined context, as set by SNESSetHessian() 81! (not used here) 82! 83! Output Parameters: 84! H - Hessian matrix 85! PrecH - optionally different matrix used to compute the preconditioner (not used here) 86! ierr - error code 87! 88! Note: Providing the Hessian may not be necessary. Only some solvers 89! require this matrix. 90 91 subroutine FormHessian(ta, X, H, PrecH, dummy, ierr) 92 93! Input/output variables: 94 type(tTao) ta 95 type(tVec) X 96 type(tMat) H, PrecH 97 PetscErrorCode ierr 98 PetscInt dummy 99 100 PetscReal v(0:1, 0:1) 101 PetscBool assembled 102 103! PETSc's VecGetArray acts differently in Fortran than it does in C. 104! Calling VecGetArray((Vec) X, (PetscReal) x_array(0:1), (PetscOffset) x_index, ierr) 105! will return an array of doubles referenced by x_array offset by x_index. 106! i.e., to reference the kth element of X, use x_array(k + x_index). 107! Notice that by declaring the arrays with range (0:1), we are using the C 0-indexing practice. 108 PetscReal, pointer :: x_v(:) 109 PetscInt i, nn, ind(0:1), i2 110 PetscReal alpha 111 PetscInt n 112 common/params/alpha, n 113 114 ierr = 0 115 nn = n/2 116 i2 = 2 117 118! Zero existing matrix entries 119 PetscCall(MatAssembled(H, assembled, ierr)) 120 if (assembled .eqv. PETSC_TRUE) PetscCall(MatZeroEntries(H, ierr)) 121 122! Get a pointer to vector data 123 124 PetscCall(VecGetArrayRead(X, x_v, ierr)) 125 126! Compute Hessian entries 127 128 do i = 0, nn - 1 129 v(1, 1) = 2.0*alpha 130 v(0, 0) = -4.0*alpha*(x_v(1 + 2*i + 1) - 3*x_v(1 + 2*i)*x_v(1 + 2*i)) + 2 131 v(1, 0) = -4.0*alpha*x_v(1 + 2*i) 132 v(0, 1) = v(1, 0) 133 ind(0) = 2*i 134 ind(1) = 2*i + 1 135 PetscCall(MatSetValues(H, i2, ind, i2, ind, reshape(v, [i2*i2]), INSERT_VALUES, ierr)) 136 end do 137 138! Restore vector 139 140 PetscCall(VecRestoreArrayRead(X, x_v, ierr)) 141 142! Assemble matrix 143 144 PetscCall(MatAssemblyBegin(H, MAT_FINAL_ASSEMBLY, ierr)) 145 PetscCall(MatAssemblyEnd(H, MAT_FINAL_ASSEMBLY, ierr)) 146 147 PetscCall(PetscLogFlops(9.0d0*nn, ierr)) 148 149 end 150end module rosenbrock1fmodule 151 152program rosenbrock1f 153 use petsctao 154 use rosenbrock1fmodule 155 implicit none 156 157! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 158! Variable declarations 159! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 160! 161! See additional variable declarations in the file rosenbrock1f.h 162 163 PetscErrorCode ierr ! used to check for functions returning nonzeros 164 type(tVec) x ! solution vector 165 type(tMat) H ! hessian matrix 166 type(tTao) ta ! TAO_SOVER context 167 PetscBool flg 168 PetscInt i2, i1 169 PetscMPIInt size 170 PetscReal zero 171 PetscReal alpha 172 PetscInt n 173 common/params/alpha, n 174 175 zero = 0.0d0 176 i2 = 2 177 i1 = 1 178 179! Initialize TAO and PETSc 180 PetscCallA(PetscInitialize(ierr)) 181 182 PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, size, ierr)) 183 PetscCheckA(size == 1, PETSC_COMM_SELF, PETSC_ERR_WRONG_MPI_SIZE, 'This is a uniprocessor example only') 184 185! Initialize problem parameters 186 n = 2 187 alpha = 99.0d0 188 189! Check for command line arguments to override defaults 190 PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr)) 191 PetscCallA(PetscOptionsGetReal(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-alpha', alpha, flg, ierr)) 192 193! Allocate vectors for the solution and gradient 194 PetscCallA(VecCreateSeq(PETSC_COMM_SELF, n, x, ierr)) 195 196! Allocate storage space for Hessian 197 PetscCallA(MatCreateSeqBAIJ(PETSC_COMM_SELF, i2, n, n, i1, PETSC_NULL_INTEGER_ARRAY, H, ierr)) 198 199 PetscCallA(MatSetOption(H, MAT_SYMMETRIC, PETSC_TRUE, ierr)) 200 201! The TAO code begins here 202 203! Create TAO solver 204 PetscCallA(TaoCreate(PETSC_COMM_SELF, ta, ierr)) 205 PetscCallA(TaoSetType(ta, TAOLMVM, ierr)) 206 207! Set routines for function, gradient, and hessian evaluation 208 PetscCallA(TaoSetObjectiveAndGradient(ta, PETSC_NULL_VEC, FormFunctionGradient, 0, ierr)) 209 PetscCallA(TaoSetHessian(ta, H, H, FormHessian, 0, ierr)) 210 211! Optional: Set initial guess 212 PetscCallA(VecSet(x, zero, ierr)) 213 PetscCallA(TaoSetSolution(ta, x, ierr)) 214 215! Check for TAO command line options 216 PetscCallA(TaoSetFromOptions(ta, ierr)) 217 218! SOLVE THE APPLICATION 219 PetscCallA(TaoSolve(ta, ierr)) 220 221! TaoView() prints ierr about the TAO solver; the option 222! -tao_view 223! can alternatively be used to activate this at runtime. 224! PetscCallA(TaoView(ta,PETSC_VIEWER_STDOUT_SELF,ierr)) 225 226! Free TAO data structures 227 PetscCallA(TaoDestroy(ta, ierr)) 228 229! Free PETSc data structures 230 PetscCallA(VecDestroy(x, ierr)) 231 PetscCallA(MatDestroy(H, ierr)) 232 233 PetscCallA(PetscFinalize(ierr)) 234end program rosenbrock1f 235! 236!/*TEST 237! 238! build: 239! requires: !complex 240! 241! test: 242! args: -tao_monitor_short -tao_type ntr -tao_gatol 1.e-5 243! requires: !single 244! 245!TEST*/ 246