Actual source code: rosenbrock1f.F90
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: !
11: !
13: ! ----------------------------------------------------------------------
14: !
15: #include "petsc/finclude/petsctao.h"
16: module rosenbrock1fmodule
17: use petsctao
18: implicit none
19: contains
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
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
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
46: ierr = 0
47: nn = n/2
48: ff = 0
50: ! Get pointers to vector data
51: PetscCall(VecGetArrayRead(X, x_v, ierr))
52: PetscCall(VecGetArray(G, g_v, ierr))
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
63: ! Restore vectors
64: PetscCall(VecRestoreArrayRead(X, x_v, ierr))
65: PetscCall(VecRestoreArray(G, g_v, ierr))
67: f = ff
68: PetscCall(PetscLogFlops(15.0d0*nn, ierr))
70: end
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.
91: subroutine FormHessian(ta, X, H, PrecH, dummy, ierr)
93: ! Input/output variables:
94: type(tTao) ta
95: type(tVec) X
96: type(tMat) H, PrecH
97: PetscErrorCode ierr
98: PetscInt dummy
100: PetscReal v(0:1, 0:1)
101: PetscBool assembled
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
114: ierr = 0
115: nn = n/2
116: i2 = 2
118: ! Zero existing matrix entries
119: PetscCall(MatAssembled(H, assembled, ierr))
120: if (assembled .eqv. PETSC_TRUE) PetscCall(MatZeroEntries(H, ierr))
122: ! Get a pointer to vector data
124: PetscCall(VecGetArrayRead(X, x_v, ierr))
126: ! Compute Hessian entries
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
138: ! Restore vector
140: PetscCall(VecRestoreArrayRead(X, x_v, ierr))
142: ! Assemble matrix
144: PetscCall(MatAssemblyBegin(H, MAT_FINAL_ASSEMBLY, ierr))
145: PetscCall(MatAssemblyEnd(H, MAT_FINAL_ASSEMBLY, ierr))
147: PetscCall(PetscLogFlops(9.0d0*nn, ierr))
149: end
150: end module rosenbrock1fmodule
152: program rosenbrock1f
153: use petsctao
154: use rosenbrock1fmodule
155: implicit none
157: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
158: ! Variable declarations
159: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
160: !
161: ! See additional variable declarations in the file rosenbrock1f.h
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
175: zero = 0.0d0
176: i2 = 2
177: i1 = 1
179: ! Initialize TAO and PETSc
180: PetscCallA(PetscInitialize(ierr))
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')
185: ! Initialize problem parameters
186: n = 2
187: alpha = 99.0d0
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))
193: ! Allocate vectors for the solution and gradient
194: PetscCallA(VecCreateSeq(PETSC_COMM_SELF, n, x, ierr))
196: ! Allocate storage space for Hessian
197: PetscCallA(MatCreateSeqBAIJ(PETSC_COMM_SELF, i2, n, n, i1, PETSC_NULL_INTEGER_ARRAY, H, ierr))
199: PetscCallA(MatSetOption(H, MAT_SYMMETRIC, PETSC_TRUE, ierr))
201: ! The TAO code begins here
203: ! Create TAO solver
204: PetscCallA(TaoCreate(PETSC_COMM_SELF, ta, ierr))
205: PetscCallA(TaoSetType(ta, TAOLMVM, ierr))
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))
211: ! Optional: Set initial guess
212: PetscCallA(VecSet(x, zero, ierr))
213: PetscCallA(TaoSetSolution(ta, x, ierr))
215: ! Check for TAO command line options
216: PetscCallA(TaoSetFromOptions(ta, ierr))
218: ! SOLVE THE APPLICATION
219: PetscCallA(TaoSolve(ta, ierr))
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))
226: ! Free TAO data structures
227: PetscCallA(TaoDestroy(ta, ierr))
229: ! Free PETSc data structures
230: PetscCallA(VecDestroy(x, ierr))
231: PetscCallA(MatDestroy(H, ierr))
233: PetscCallA(PetscFinalize(ierr))
234: end 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*/