1c4762a1bSJed Brown /* 2c4762a1bSJed Brown Include "petsctao.h" so that we can use TAO solvers. Note that this 3c4762a1bSJed Brown file automatically includes libraries such as: 4c4762a1bSJed Brown petsc.h - base PETSc routines petscvec.h - vectors 5a5b23f4aSJose E. Roman petscsys.h - system routines petscmat.h - matrices 6c4762a1bSJed Brown petscis.h - index sets petscksp.h - Krylov subspace methods 7c4762a1bSJed Brown petscviewer.h - viewers petscpc.h - preconditioners 8c4762a1bSJed Brown 9c4762a1bSJed Brown */ 10c4762a1bSJed Brown 11c4762a1bSJed Brown #include <petsctao.h> 12c4762a1bSJed Brown 13c4762a1bSJed Brown /* 14c4762a1bSJed Brown Description: These data are the result of a NIST study involving 15c4762a1bSJed Brown ultrasonic calibration. The response variable is 16c4762a1bSJed Brown ultrasonic response, and the predictor variable is 17c4762a1bSJed Brown metal distance. 18c4762a1bSJed Brown 19c4762a1bSJed Brown Reference: Chwirut, D., NIST (197?). 20c4762a1bSJed Brown Ultrasonic Reference Block Study. 21c4762a1bSJed Brown */ 22c4762a1bSJed Brown 23c4762a1bSJed Brown static char help[]="Finds the nonlinear least-squares solution to the model \n\ 24c4762a1bSJed Brown y = exp[-b1*x]/(b2+b3*x) + e \n"; 25c4762a1bSJed Brown 26c4762a1bSJed Brown /* T 27c4762a1bSJed Brown Concepts: TAO^Solving a system of nonlinear equations, nonlinear least squares 28c4762a1bSJed Brown Routines: TaoCreate(); 29c4762a1bSJed Brown Routines: TaoSetType(); 30c4762a1bSJed Brown Routines: TaoSetResidualRoutine(); 31c4762a1bSJed Brown Routines: TaoSetMonitor(); 32a82e8c82SStefano Zampini Routines: TaoSetSolution(); 33c4762a1bSJed Brown Routines: TaoSetFromOptions(); 34c4762a1bSJed Brown Routines: TaoSolve(); 35c4762a1bSJed Brown Routines: TaoDestroy(); 36c4762a1bSJed Brown Processors: n 37c4762a1bSJed Brown T*/ 38c4762a1bSJed Brown 39c4762a1bSJed Brown #define NOBSERVATIONS 214 40c4762a1bSJed Brown #define NPARAMETERS 3 41c4762a1bSJed Brown 42c4762a1bSJed Brown #define DIE_TAG 2000 43c4762a1bSJed Brown #define IDLE_TAG 1000 44c4762a1bSJed Brown 45c4762a1bSJed Brown /* User-defined application context */ 46c4762a1bSJed Brown typedef struct { 47c4762a1bSJed Brown /* Working space */ 48c4762a1bSJed Brown PetscReal t[NOBSERVATIONS]; /* array of independent variables of observation */ 49c4762a1bSJed Brown PetscReal y[NOBSERVATIONS]; /* array of dependent variables */ 50c4762a1bSJed Brown PetscMPIInt size,rank; 51c4762a1bSJed Brown } AppCtx; 52c4762a1bSJed Brown 53c4762a1bSJed Brown /* User provided Routines */ 54c4762a1bSJed Brown PetscErrorCode InitializeData(AppCtx *user); 55c4762a1bSJed Brown PetscErrorCode FormStartingPoint(Vec); 56c4762a1bSJed Brown PetscErrorCode EvaluateFunction(Tao, Vec, Vec, void *); 57c4762a1bSJed Brown PetscErrorCode TaskWorker(AppCtx *user); 58c4762a1bSJed Brown PetscErrorCode StopWorkers(AppCtx *user); 59c4762a1bSJed Brown PetscErrorCode RunSimulation(PetscReal *x, PetscInt i, PetscReal*f, AppCtx *user); 60c4762a1bSJed Brown 61c4762a1bSJed Brown /*--------------------------------------------------------------------*/ 62c4762a1bSJed Brown int main(int argc,char **argv) 63c4762a1bSJed Brown { 64c4762a1bSJed Brown Vec x, f; /* solution, function */ 65c4762a1bSJed Brown Tao tao; /* Tao solver context */ 66c4762a1bSJed Brown AppCtx user; /* user-defined work context */ 67c4762a1bSJed Brown 68c4762a1bSJed Brown /* Initialize TAO and PETSc */ 69*b122ec5aSJacob Faibussowitsch CHKERRQ(PetscInitialize(&argc,&argv,(char *)0,help)); 70c4762a1bSJed Brown MPI_Comm_size(MPI_COMM_WORLD,&user.size); 71c4762a1bSJed Brown MPI_Comm_rank(MPI_COMM_WORLD,&user.rank); 725f80ce2aSJacob Faibussowitsch CHKERRQ(InitializeData(&user)); 73c4762a1bSJed Brown 74c4762a1bSJed Brown /* Run optimization on rank 0 */ 75c4762a1bSJed Brown if (user.rank == 0) { 76c4762a1bSJed Brown /* Allocate vectors */ 775f80ce2aSJacob Faibussowitsch CHKERRQ(VecCreateSeq(PETSC_COMM_SELF,NPARAMETERS,&x)); 785f80ce2aSJacob Faibussowitsch CHKERRQ(VecCreateSeq(PETSC_COMM_SELF,NOBSERVATIONS,&f)); 79c4762a1bSJed Brown 80c4762a1bSJed Brown /* TAO code begins here */ 81c4762a1bSJed Brown 82c4762a1bSJed Brown /* Create TAO solver and set desired solution method */ 835f80ce2aSJacob Faibussowitsch CHKERRQ(TaoCreate(PETSC_COMM_SELF,&tao)); 845f80ce2aSJacob Faibussowitsch CHKERRQ(TaoSetType(tao,TAOPOUNDERS)); 85c4762a1bSJed Brown 86c4762a1bSJed Brown /* Set the function and Jacobian routines. */ 875f80ce2aSJacob Faibussowitsch CHKERRQ(FormStartingPoint(x)); 885f80ce2aSJacob Faibussowitsch CHKERRQ(TaoSetSolution(tao,x)); 895f80ce2aSJacob Faibussowitsch CHKERRQ(TaoSetResidualRoutine(tao,f,EvaluateFunction,(void*)&user)); 90c4762a1bSJed Brown 91c4762a1bSJed Brown /* Check for any TAO command line arguments */ 925f80ce2aSJacob Faibussowitsch CHKERRQ(TaoSetFromOptions(tao)); 93c4762a1bSJed Brown 94c4762a1bSJed Brown /* Perform the Solve */ 955f80ce2aSJacob Faibussowitsch CHKERRQ(TaoSolve(tao)); 96c4762a1bSJed Brown 97c4762a1bSJed Brown /* Free TAO data structures */ 985f80ce2aSJacob Faibussowitsch CHKERRQ(TaoDestroy(&tao)); 99c4762a1bSJed Brown 100c4762a1bSJed Brown /* Free PETSc data structures */ 1015f80ce2aSJacob Faibussowitsch CHKERRQ(VecDestroy(&x)); 1025f80ce2aSJacob Faibussowitsch CHKERRQ(VecDestroy(&f)); 103c4762a1bSJed Brown StopWorkers(&user); 104c4762a1bSJed Brown } else { 105c4762a1bSJed Brown TaskWorker(&user); 106c4762a1bSJed Brown } 107*b122ec5aSJacob Faibussowitsch CHKERRQ(PetscFinalize()); 108*b122ec5aSJacob Faibussowitsch return 0; 109c4762a1bSJed Brown } 110c4762a1bSJed Brown 111c4762a1bSJed Brown /*--------------------------------------------------------------------*/ 112c4762a1bSJed Brown PetscErrorCode EvaluateFunction(Tao tao, Vec X, Vec F, void *ptr) 113c4762a1bSJed Brown { 114c4762a1bSJed Brown AppCtx *user = (AppCtx *)ptr; 115c4762a1bSJed Brown PetscInt i; 116c4762a1bSJed Brown PetscReal *x,*f; 117c4762a1bSJed Brown 118c4762a1bSJed Brown PetscFunctionBegin; 1195f80ce2aSJacob Faibussowitsch CHKERRQ(VecGetArray(X,&x)); 1205f80ce2aSJacob Faibussowitsch CHKERRQ(VecGetArray(F,&f)); 121c4762a1bSJed Brown if (user->size == 1) { 122c4762a1bSJed Brown /* Single processor */ 123c4762a1bSJed Brown for (i=0;i<NOBSERVATIONS;i++) { 1245f80ce2aSJacob Faibussowitsch CHKERRQ(RunSimulation(x,i,&f[i],user)); 125c4762a1bSJed Brown } 126c4762a1bSJed Brown } else { 1279dddd249SSatish Balay /* Multiprocessor main */ 128c4762a1bSJed Brown PetscMPIInt tag; 129c4762a1bSJed Brown PetscInt finishedtasks,next_task,checkedin; 130c4762a1bSJed Brown PetscReal f_i=0.0; 131c4762a1bSJed Brown MPI_Status status; 132c4762a1bSJed Brown 133c4762a1bSJed Brown next_task=0; 134c4762a1bSJed Brown finishedtasks=0; 135c4762a1bSJed Brown checkedin=0; 136c4762a1bSJed Brown 137c4762a1bSJed Brown while (finishedtasks < NOBSERVATIONS || checkedin < user->size-1) { 1385f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Recv(&f_i,1,MPIU_REAL,MPI_ANY_SOURCE,MPI_ANY_TAG,PETSC_COMM_WORLD,&status)); 139c4762a1bSJed Brown if (status.MPI_TAG == IDLE_TAG) { 140c4762a1bSJed Brown checkedin++; 141c4762a1bSJed Brown } else { 142c4762a1bSJed Brown 143c4762a1bSJed Brown tag = status.MPI_TAG; 144c4762a1bSJed Brown f[tag] = (PetscReal)f_i; 145c4762a1bSJed Brown finishedtasks++; 146c4762a1bSJed Brown } 147c4762a1bSJed Brown 148c4762a1bSJed Brown if (next_task<NOBSERVATIONS) { 1495f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Send(x,NPARAMETERS,MPIU_REAL,status.MPI_SOURCE,next_task,PETSC_COMM_WORLD)); 150c4762a1bSJed Brown next_task++; 151c4762a1bSJed Brown 152c4762a1bSJed Brown } else { 153c4762a1bSJed Brown /* Send idle message */ 1545f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Send(x,NPARAMETERS,MPIU_REAL,status.MPI_SOURCE,IDLE_TAG,PETSC_COMM_WORLD)); 155c4762a1bSJed Brown } 156c4762a1bSJed Brown } 157c4762a1bSJed Brown } 1585f80ce2aSJacob Faibussowitsch CHKERRQ(VecRestoreArray(X,&x)); 1595f80ce2aSJacob Faibussowitsch CHKERRQ(VecRestoreArray(F,&f)); 160c4762a1bSJed Brown PetscLogFlops(6*NOBSERVATIONS); 161c4762a1bSJed Brown PetscFunctionReturn(0); 162c4762a1bSJed Brown } 163c4762a1bSJed Brown 164c4762a1bSJed Brown /* ------------------------------------------------------------ */ 165c4762a1bSJed Brown PetscErrorCode FormStartingPoint(Vec X) 166c4762a1bSJed Brown { 167c4762a1bSJed Brown PetscReal *x; 168c4762a1bSJed Brown 169c4762a1bSJed Brown PetscFunctionBegin; 1705f80ce2aSJacob Faibussowitsch CHKERRQ(VecGetArray(X,&x)); 171c4762a1bSJed Brown x[0] = 0.15; 172c4762a1bSJed Brown x[1] = 0.008; 173c4762a1bSJed Brown x[2] = 0.010; 1745f80ce2aSJacob Faibussowitsch CHKERRQ(VecRestoreArray(X,&x)); 175c4762a1bSJed Brown PetscFunctionReturn(0); 176c4762a1bSJed Brown } 177c4762a1bSJed Brown 178c4762a1bSJed Brown /* ---------------------------------------------------------------------- */ 179c4762a1bSJed Brown PetscErrorCode InitializeData(AppCtx *user) 180c4762a1bSJed Brown { 181c4762a1bSJed Brown PetscReal *t=user->t,*y=user->y; 182c4762a1bSJed Brown PetscInt i=0; 183c4762a1bSJed Brown 184c4762a1bSJed Brown PetscFunctionBegin; 185c4762a1bSJed Brown y[i] = 92.9000; t[i++] = 0.5000; 186c4762a1bSJed Brown y[i] = 78.7000; t[i++] = 0.6250; 187c4762a1bSJed Brown y[i] = 64.2000; t[i++] = 0.7500; 188c4762a1bSJed Brown y[i] = 64.9000; t[i++] = 0.8750; 189c4762a1bSJed Brown y[i] = 57.1000; t[i++] = 1.0000; 190c4762a1bSJed Brown y[i] = 43.3000; t[i++] = 1.2500; 191c4762a1bSJed Brown y[i] = 31.1000; t[i++] = 1.7500; 192c4762a1bSJed Brown y[i] = 23.6000; t[i++] = 2.2500; 193c4762a1bSJed Brown y[i] = 31.0500; t[i++] = 1.7500; 194c4762a1bSJed Brown y[i] = 23.7750; t[i++] = 2.2500; 195c4762a1bSJed Brown y[i] = 17.7375; t[i++] = 2.7500; 196c4762a1bSJed Brown y[i] = 13.8000; t[i++] = 3.2500; 197c4762a1bSJed Brown y[i] = 11.5875; t[i++] = 3.7500; 198c4762a1bSJed Brown y[i] = 9.4125; t[i++] = 4.2500; 199c4762a1bSJed Brown y[i] = 7.7250; t[i++] = 4.7500; 200c4762a1bSJed Brown y[i] = 7.3500; t[i++] = 5.2500; 201c4762a1bSJed Brown y[i] = 8.0250; t[i++] = 5.7500; 202c4762a1bSJed Brown y[i] = 90.6000; t[i++] = 0.5000; 203c4762a1bSJed Brown y[i] = 76.9000; t[i++] = 0.6250; 204c4762a1bSJed Brown y[i] = 71.6000; t[i++] = 0.7500; 205c4762a1bSJed Brown y[i] = 63.6000; t[i++] = 0.8750; 206c4762a1bSJed Brown y[i] = 54.0000; t[i++] = 1.0000; 207c4762a1bSJed Brown y[i] = 39.2000; t[i++] = 1.2500; 208c4762a1bSJed Brown y[i] = 29.3000; t[i++] = 1.7500; 209c4762a1bSJed Brown y[i] = 21.4000; t[i++] = 2.2500; 210c4762a1bSJed Brown y[i] = 29.1750; t[i++] = 1.7500; 211c4762a1bSJed Brown y[i] = 22.1250; t[i++] = 2.2500; 212c4762a1bSJed Brown y[i] = 17.5125; t[i++] = 2.7500; 213c4762a1bSJed Brown y[i] = 14.2500; t[i++] = 3.2500; 214c4762a1bSJed Brown y[i] = 9.4500; t[i++] = 3.7500; 215c4762a1bSJed Brown y[i] = 9.1500; t[i++] = 4.2500; 216c4762a1bSJed Brown y[i] = 7.9125; t[i++] = 4.7500; 217c4762a1bSJed Brown y[i] = 8.4750; t[i++] = 5.2500; 218c4762a1bSJed Brown y[i] = 6.1125; t[i++] = 5.7500; 219c4762a1bSJed Brown y[i] = 80.0000; t[i++] = 0.5000; 220c4762a1bSJed Brown y[i] = 79.0000; t[i++] = 0.6250; 221c4762a1bSJed Brown y[i] = 63.8000; t[i++] = 0.7500; 222c4762a1bSJed Brown y[i] = 57.2000; t[i++] = 0.8750; 223c4762a1bSJed Brown y[i] = 53.2000; t[i++] = 1.0000; 224c4762a1bSJed Brown y[i] = 42.5000; t[i++] = 1.2500; 225c4762a1bSJed Brown y[i] = 26.8000; t[i++] = 1.7500; 226c4762a1bSJed Brown y[i] = 20.4000; t[i++] = 2.2500; 227c4762a1bSJed Brown y[i] = 26.8500; t[i++] = 1.7500; 228c4762a1bSJed Brown y[i] = 21.0000; t[i++] = 2.2500; 229c4762a1bSJed Brown y[i] = 16.4625; t[i++] = 2.7500; 230c4762a1bSJed Brown y[i] = 12.5250; t[i++] = 3.2500; 231c4762a1bSJed Brown y[i] = 10.5375; t[i++] = 3.7500; 232c4762a1bSJed Brown y[i] = 8.5875; t[i++] = 4.2500; 233c4762a1bSJed Brown y[i] = 7.1250; t[i++] = 4.7500; 234c4762a1bSJed Brown y[i] = 6.1125; t[i++] = 5.2500; 235c4762a1bSJed Brown y[i] = 5.9625; t[i++] = 5.7500; 236c4762a1bSJed Brown y[i] = 74.1000; t[i++] = 0.5000; 237c4762a1bSJed Brown y[i] = 67.3000; t[i++] = 0.6250; 238c4762a1bSJed Brown y[i] = 60.8000; t[i++] = 0.7500; 239c4762a1bSJed Brown y[i] = 55.5000; t[i++] = 0.8750; 240c4762a1bSJed Brown y[i] = 50.3000; t[i++] = 1.0000; 241c4762a1bSJed Brown y[i] = 41.0000; t[i++] = 1.2500; 242c4762a1bSJed Brown y[i] = 29.4000; t[i++] = 1.7500; 243c4762a1bSJed Brown y[i] = 20.4000; t[i++] = 2.2500; 244c4762a1bSJed Brown y[i] = 29.3625; t[i++] = 1.7500; 245c4762a1bSJed Brown y[i] = 21.1500; t[i++] = 2.2500; 246c4762a1bSJed Brown y[i] = 16.7625; t[i++] = 2.7500; 247c4762a1bSJed Brown y[i] = 13.2000; t[i++] = 3.2500; 248c4762a1bSJed Brown y[i] = 10.8750; t[i++] = 3.7500; 249c4762a1bSJed Brown y[i] = 8.1750; t[i++] = 4.2500; 250c4762a1bSJed Brown y[i] = 7.3500; t[i++] = 4.7500; 251c4762a1bSJed Brown y[i] = 5.9625; t[i++] = 5.2500; 252c4762a1bSJed Brown y[i] = 5.6250; t[i++] = 5.7500; 253c4762a1bSJed Brown y[i] = 81.5000; t[i++] = .5000; 254c4762a1bSJed Brown y[i] = 62.4000; t[i++] = .7500; 255c4762a1bSJed Brown y[i] = 32.5000; t[i++] = 1.5000; 256c4762a1bSJed Brown y[i] = 12.4100; t[i++] = 3.0000; 257c4762a1bSJed Brown y[i] = 13.1200; t[i++] = 3.0000; 258c4762a1bSJed Brown y[i] = 15.5600; t[i++] = 3.0000; 259c4762a1bSJed Brown y[i] = 5.6300; t[i++] = 6.0000; 260c4762a1bSJed Brown y[i] = 78.0000; t[i++] = .5000; 261c4762a1bSJed Brown y[i] = 59.9000; t[i++] = .7500; 262c4762a1bSJed Brown y[i] = 33.2000; t[i++] = 1.5000; 263c4762a1bSJed Brown y[i] = 13.8400; t[i++] = 3.0000; 264c4762a1bSJed Brown y[i] = 12.7500; t[i++] = 3.0000; 265c4762a1bSJed Brown y[i] = 14.6200; t[i++] = 3.0000; 266c4762a1bSJed Brown y[i] = 3.9400; t[i++] = 6.0000; 267c4762a1bSJed Brown y[i] = 76.8000; t[i++] = .5000; 268c4762a1bSJed Brown y[i] = 61.0000; t[i++] = .7500; 269c4762a1bSJed Brown y[i] = 32.9000; t[i++] = 1.5000; 270c4762a1bSJed Brown y[i] = 13.8700; t[i++] = 3.0000; 271c4762a1bSJed Brown y[i] = 11.8100; t[i++] = 3.0000; 272c4762a1bSJed Brown y[i] = 13.3100; t[i++] = 3.0000; 273c4762a1bSJed Brown y[i] = 5.4400; t[i++] = 6.0000; 274c4762a1bSJed Brown y[i] = 78.0000; t[i++] = .5000; 275c4762a1bSJed Brown y[i] = 63.5000; t[i++] = .7500; 276c4762a1bSJed Brown y[i] = 33.8000; t[i++] = 1.5000; 277c4762a1bSJed Brown y[i] = 12.5600; t[i++] = 3.0000; 278c4762a1bSJed Brown y[i] = 5.6300; t[i++] = 6.0000; 279c4762a1bSJed Brown y[i] = 12.7500; t[i++] = 3.0000; 280c4762a1bSJed Brown y[i] = 13.1200; t[i++] = 3.0000; 281c4762a1bSJed Brown y[i] = 5.4400; t[i++] = 6.0000; 282c4762a1bSJed Brown y[i] = 76.8000; t[i++] = .5000; 283c4762a1bSJed Brown y[i] = 60.0000; t[i++] = .7500; 284c4762a1bSJed Brown y[i] = 47.8000; t[i++] = 1.0000; 285c4762a1bSJed Brown y[i] = 32.0000; t[i++] = 1.5000; 286c4762a1bSJed Brown y[i] = 22.2000; t[i++] = 2.0000; 287c4762a1bSJed Brown y[i] = 22.5700; t[i++] = 2.0000; 288c4762a1bSJed Brown y[i] = 18.8200; t[i++] = 2.5000; 289c4762a1bSJed Brown y[i] = 13.9500; t[i++] = 3.0000; 290c4762a1bSJed Brown y[i] = 11.2500; t[i++] = 4.0000; 291c4762a1bSJed Brown y[i] = 9.0000; t[i++] = 5.0000; 292c4762a1bSJed Brown y[i] = 6.6700; t[i++] = 6.0000; 293c4762a1bSJed Brown y[i] = 75.8000; t[i++] = .5000; 294c4762a1bSJed Brown y[i] = 62.0000; t[i++] = .7500; 295c4762a1bSJed Brown y[i] = 48.8000; t[i++] = 1.0000; 296c4762a1bSJed Brown y[i] = 35.2000; t[i++] = 1.5000; 297c4762a1bSJed Brown y[i] = 20.0000; t[i++] = 2.0000; 298c4762a1bSJed Brown y[i] = 20.3200; t[i++] = 2.0000; 299c4762a1bSJed Brown y[i] = 19.3100; t[i++] = 2.5000; 300c4762a1bSJed Brown y[i] = 12.7500; t[i++] = 3.0000; 301c4762a1bSJed Brown y[i] = 10.4200; t[i++] = 4.0000; 302c4762a1bSJed Brown y[i] = 7.3100; t[i++] = 5.0000; 303c4762a1bSJed Brown y[i] = 7.4200; t[i++] = 6.0000; 304c4762a1bSJed Brown y[i] = 70.5000; t[i++] = .5000; 305c4762a1bSJed Brown y[i] = 59.5000; t[i++] = .7500; 306c4762a1bSJed Brown y[i] = 48.5000; t[i++] = 1.0000; 307c4762a1bSJed Brown y[i] = 35.8000; t[i++] = 1.5000; 308c4762a1bSJed Brown y[i] = 21.0000; t[i++] = 2.0000; 309c4762a1bSJed Brown y[i] = 21.6700; t[i++] = 2.0000; 310c4762a1bSJed Brown y[i] = 21.0000; t[i++] = 2.5000; 311c4762a1bSJed Brown y[i] = 15.6400; t[i++] = 3.0000; 312c4762a1bSJed Brown y[i] = 8.1700; t[i++] = 4.0000; 313c4762a1bSJed Brown y[i] = 8.5500; t[i++] = 5.0000; 314c4762a1bSJed Brown y[i] = 10.1200; t[i++] = 6.0000; 315c4762a1bSJed Brown y[i] = 78.0000; t[i++] = .5000; 316c4762a1bSJed Brown y[i] = 66.0000; t[i++] = .6250; 317c4762a1bSJed Brown y[i] = 62.0000; t[i++] = .7500; 318c4762a1bSJed Brown y[i] = 58.0000; t[i++] = .8750; 319c4762a1bSJed Brown y[i] = 47.7000; t[i++] = 1.0000; 320c4762a1bSJed Brown y[i] = 37.8000; t[i++] = 1.2500; 321c4762a1bSJed Brown y[i] = 20.2000; t[i++] = 2.2500; 322c4762a1bSJed Brown y[i] = 21.0700; t[i++] = 2.2500; 323c4762a1bSJed Brown y[i] = 13.8700; t[i++] = 2.7500; 324c4762a1bSJed Brown y[i] = 9.6700; t[i++] = 3.2500; 325c4762a1bSJed Brown y[i] = 7.7600; t[i++] = 3.7500; 326c4762a1bSJed Brown y[i] = 5.4400; t[i++] = 4.2500; 327c4762a1bSJed Brown y[i] = 4.8700; t[i++] = 4.7500; 328c4762a1bSJed Brown y[i] = 4.0100; t[i++] = 5.2500; 329c4762a1bSJed Brown y[i] = 3.7500; t[i++] = 5.7500; 330c4762a1bSJed Brown y[i] = 24.1900; t[i++] = 3.0000; 331c4762a1bSJed Brown y[i] = 25.7600; t[i++] = 3.0000; 332c4762a1bSJed Brown y[i] = 18.0700; t[i++] = 3.0000; 333c4762a1bSJed Brown y[i] = 11.8100; t[i++] = 3.0000; 334c4762a1bSJed Brown y[i] = 12.0700; t[i++] = 3.0000; 335c4762a1bSJed Brown y[i] = 16.1200; t[i++] = 3.0000; 336c4762a1bSJed Brown y[i] = 70.8000; t[i++] = .5000; 337c4762a1bSJed Brown y[i] = 54.7000; t[i++] = .7500; 338c4762a1bSJed Brown y[i] = 48.0000; t[i++] = 1.0000; 339c4762a1bSJed Brown y[i] = 39.8000; t[i++] = 1.5000; 340c4762a1bSJed Brown y[i] = 29.8000; t[i++] = 2.0000; 341c4762a1bSJed Brown y[i] = 23.7000; t[i++] = 2.5000; 342c4762a1bSJed Brown y[i] = 29.6200; t[i++] = 2.0000; 343c4762a1bSJed Brown y[i] = 23.8100; t[i++] = 2.5000; 344c4762a1bSJed Brown y[i] = 17.7000; t[i++] = 3.0000; 345c4762a1bSJed Brown y[i] = 11.5500; t[i++] = 4.0000; 346c4762a1bSJed Brown y[i] = 12.0700; t[i++] = 5.0000; 347c4762a1bSJed Brown y[i] = 8.7400; t[i++] = 6.0000; 348c4762a1bSJed Brown y[i] = 80.7000; t[i++] = .5000; 349c4762a1bSJed Brown y[i] = 61.3000; t[i++] = .7500; 350c4762a1bSJed Brown y[i] = 47.5000; t[i++] = 1.0000; 351c4762a1bSJed Brown y[i] = 29.0000; t[i++] = 1.5000; 352c4762a1bSJed Brown y[i] = 24.0000; t[i++] = 2.0000; 353c4762a1bSJed Brown y[i] = 17.7000; t[i++] = 2.5000; 354c4762a1bSJed Brown y[i] = 24.5600; t[i++] = 2.0000; 355c4762a1bSJed Brown y[i] = 18.6700; t[i++] = 2.5000; 356c4762a1bSJed Brown y[i] = 16.2400; t[i++] = 3.0000; 357c4762a1bSJed Brown y[i] = 8.7400; t[i++] = 4.0000; 358c4762a1bSJed Brown y[i] = 7.8700; t[i++] = 5.0000; 359c4762a1bSJed Brown y[i] = 8.5100; t[i++] = 6.0000; 360c4762a1bSJed Brown y[i] = 66.7000; t[i++] = .5000; 361c4762a1bSJed Brown y[i] = 59.2000; t[i++] = .7500; 362c4762a1bSJed Brown y[i] = 40.8000; t[i++] = 1.0000; 363c4762a1bSJed Brown y[i] = 30.7000; t[i++] = 1.5000; 364c4762a1bSJed Brown y[i] = 25.7000; t[i++] = 2.0000; 365c4762a1bSJed Brown y[i] = 16.3000; t[i++] = 2.5000; 366c4762a1bSJed Brown y[i] = 25.9900; t[i++] = 2.0000; 367c4762a1bSJed Brown y[i] = 16.9500; t[i++] = 2.5000; 368c4762a1bSJed Brown y[i] = 13.3500; t[i++] = 3.0000; 369c4762a1bSJed Brown y[i] = 8.6200; t[i++] = 4.0000; 370c4762a1bSJed Brown y[i] = 7.2000; t[i++] = 5.0000; 371c4762a1bSJed Brown y[i] = 6.6400; t[i++] = 6.0000; 372c4762a1bSJed Brown y[i] = 13.6900; t[i++] = 3.0000; 373c4762a1bSJed Brown y[i] = 81.0000; t[i++] = .5000; 374c4762a1bSJed Brown y[i] = 64.5000; t[i++] = .7500; 375c4762a1bSJed Brown y[i] = 35.5000; t[i++] = 1.5000; 376c4762a1bSJed Brown y[i] = 13.3100; t[i++] = 3.0000; 377c4762a1bSJed Brown y[i] = 4.8700; t[i++] = 6.0000; 378c4762a1bSJed Brown y[i] = 12.9400; t[i++] = 3.0000; 379c4762a1bSJed Brown y[i] = 5.0600; t[i++] = 6.0000; 380c4762a1bSJed Brown y[i] = 15.1900; t[i++] = 3.0000; 381c4762a1bSJed Brown y[i] = 14.6200; t[i++] = 3.0000; 382c4762a1bSJed Brown y[i] = 15.6400; t[i++] = 3.0000; 383c4762a1bSJed Brown y[i] = 25.5000; t[i++] = 1.7500; 384c4762a1bSJed Brown y[i] = 25.9500; t[i++] = 1.7500; 385c4762a1bSJed Brown y[i] = 81.7000; t[i++] = .5000; 386c4762a1bSJed Brown y[i] = 61.6000; t[i++] = .7500; 387c4762a1bSJed Brown y[i] = 29.8000; t[i++] = 1.7500; 388c4762a1bSJed Brown y[i] = 29.8100; t[i++] = 1.7500; 389c4762a1bSJed Brown y[i] = 17.1700; t[i++] = 2.7500; 390c4762a1bSJed Brown y[i] = 10.3900; t[i++] = 3.7500; 391c4762a1bSJed Brown y[i] = 28.4000; t[i++] = 1.7500; 392c4762a1bSJed Brown y[i] = 28.6900; t[i++] = 1.7500; 393c4762a1bSJed Brown y[i] = 81.3000; t[i++] = .5000; 394c4762a1bSJed Brown y[i] = 60.9000; t[i++] = .7500; 395c4762a1bSJed Brown y[i] = 16.6500; t[i++] = 2.7500; 396c4762a1bSJed Brown y[i] = 10.0500; t[i++] = 3.7500; 397c4762a1bSJed Brown y[i] = 28.9000; t[i++] = 1.7500; 398c4762a1bSJed Brown y[i] = 28.9500; t[i++] = 1.7500; 399c4762a1bSJed Brown PetscFunctionReturn(0); 400c4762a1bSJed Brown } 401c4762a1bSJed Brown 402c4762a1bSJed Brown PetscErrorCode TaskWorker(AppCtx *user) 403c4762a1bSJed Brown { 404c4762a1bSJed Brown PetscReal x[NPARAMETERS],f = 0.0; 405c4762a1bSJed Brown PetscMPIInt tag=IDLE_TAG; 406c4762a1bSJed Brown PetscInt index; 407c4762a1bSJed Brown MPI_Status status; 408c4762a1bSJed Brown 409c4762a1bSJed Brown PetscFunctionBegin; 4109dddd249SSatish Balay /* Send check-in message to rank-0 */ 411c4762a1bSJed Brown 4125f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Send(&f,1,MPIU_REAL,0,IDLE_TAG,PETSC_COMM_WORLD)); 413c4762a1bSJed Brown while (tag != DIE_TAG) { 4145f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Recv(x,NPARAMETERS,MPIU_REAL,0,MPI_ANY_TAG,PETSC_COMM_WORLD,&status)); 415c4762a1bSJed Brown tag = status.MPI_TAG; 416c4762a1bSJed Brown if (tag == IDLE_TAG) { 4175f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Send(&f,1,MPIU_REAL,0,IDLE_TAG,PETSC_COMM_WORLD)); 418c4762a1bSJed Brown } else if (tag != DIE_TAG) { 419c4762a1bSJed Brown index = (PetscInt)tag; 4205f80ce2aSJacob Faibussowitsch CHKERRQ(RunSimulation(x,index,&f,user)); 4215f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Send(&f,1,MPIU_REAL,0,tag,PETSC_COMM_WORLD)); 422c4762a1bSJed Brown } 423c4762a1bSJed Brown } 424c4762a1bSJed Brown PetscFunctionReturn(0); 425c4762a1bSJed Brown } 426c4762a1bSJed Brown 427c4762a1bSJed Brown PetscErrorCode RunSimulation(PetscReal *x, PetscInt i, PetscReal*f, AppCtx *user) 428c4762a1bSJed Brown { 429c4762a1bSJed Brown PetscReal *t = user->t; 430c4762a1bSJed Brown PetscReal *y = user->y; 431c4762a1bSJed Brown #if defined(PETSC_USE_REAL_SINGLE) 432e1dfdf8eSBarry Smith *f = y[i] - exp(-x[0]*t[i])/(x[1] + x[2]*t[i]); /* expf() for single-precision breaks this example on Freebsd, Valgrind errors on Linux */ 433c4762a1bSJed Brown #else 434c4762a1bSJed Brown *f = y[i] - PetscExpScalar(-x[0]*t[i])/(x[1] + x[2]*t[i]); 435c4762a1bSJed Brown #endif 436c4762a1bSJed Brown return(0); 437c4762a1bSJed Brown } 438c4762a1bSJed Brown 439c4762a1bSJed Brown PetscErrorCode StopWorkers(AppCtx *user) 440c4762a1bSJed Brown { 441c4762a1bSJed Brown PetscInt checkedin; 442c4762a1bSJed Brown MPI_Status status; 443c4762a1bSJed Brown PetscReal f,x[NPARAMETERS]; 444c4762a1bSJed Brown 445c4762a1bSJed Brown PetscFunctionBegin; 446c4762a1bSJed Brown checkedin=0; 447c4762a1bSJed Brown while (checkedin < user->size-1) { 4485f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Recv(&f,1,MPIU_REAL,MPI_ANY_SOURCE,MPI_ANY_TAG,PETSC_COMM_WORLD,&status)); 449c4762a1bSJed Brown checkedin++; 4505f80ce2aSJacob Faibussowitsch CHKERRQ(PetscArrayzero(x,NPARAMETERS)); 4515f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Send(x,NPARAMETERS,MPIU_REAL,status.MPI_SOURCE,DIE_TAG,PETSC_COMM_WORLD)); 452c4762a1bSJed Brown } 453c4762a1bSJed Brown PetscFunctionReturn(0); 454c4762a1bSJed Brown } 455c4762a1bSJed Brown 456c4762a1bSJed Brown /*TEST 457c4762a1bSJed Brown 458c4762a1bSJed Brown build: 459c4762a1bSJed Brown requires: !complex 460c4762a1bSJed Brown 461c4762a1bSJed Brown test: 462c4762a1bSJed Brown nsize: 3 463c4762a1bSJed Brown requires: !single 464c4762a1bSJed Brown args: -tao_smonitor -tao_max_it 100 -tao_type pounders -tao_gatol 1.e-5 465c4762a1bSJed Brown 466c4762a1bSJed Brown TEST*/ 467