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