1 static char help[] = "Time-dependent PDE in 2d for calculating joint PDF. \n"; 2 /* 3 p_t = -x_t*p_x -y_t*p_y + f(t)*p_yy 4 xmin < x < xmax, ymin < y < ymax; 5 6 Boundary conditions: 7 Zero dirichlet in y using ghosted values 8 Periodic in x 9 10 Note that x_t and y_t in the above are given functions of x and y; they are not derivatives of x and y. 11 x_t = (y - ws) 12 y_t = (ws/2H)*(Pm - Pmax*sin(x) - D*(w - ws)) 13 14 In this example, we can see the effect of a fault, that zeroes the electrical power output 15 Pmax*sin(x), on the PDF. The fault on/off times can be controlled by options -tf and -tcl respectively. 16 17 */ 18 19 #include <petscdm.h> 20 #include <petscdmda.h> 21 #include <petscts.h> 22 23 /* 24 User-defined data structures and routines 25 */ 26 typedef struct { 27 PetscScalar ws; /* Synchronous speed */ 28 PetscScalar H; /* Inertia constant */ 29 PetscScalar D; /* Damping constant */ 30 PetscScalar Pmax,Pmax_s; /* Maximum power output of generator */ 31 PetscScalar PM_min; /* Mean mechanical power input */ 32 PetscScalar lambda; /* correlation time */ 33 PetscScalar q; /* noise strength */ 34 PetscScalar mux; /* Initial average angle */ 35 PetscScalar sigmax; /* Standard deviation of initial angle */ 36 PetscScalar muy; /* Average speed */ 37 PetscScalar sigmay; /* standard deviation of initial speed */ 38 PetscScalar rho; /* Cross-correlation coefficient */ 39 PetscScalar xmin; /* left boundary of angle */ 40 PetscScalar xmax; /* right boundary of angle */ 41 PetscScalar ymin; /* bottom boundary of speed */ 42 PetscScalar ymax; /* top boundary of speed */ 43 PetscScalar dx; /* x step size */ 44 PetscScalar dy; /* y step size */ 45 PetscScalar disper_coe; /* Dispersion coefficient */ 46 DM da; 47 PetscInt st_width; /* Stencil width */ 48 DMBoundaryType bx; /* x boundary type */ 49 DMBoundaryType by; /* y boundary type */ 50 PetscReal tf,tcl; /* Fault incidence and clearing times */ 51 } AppCtx; 52 53 PetscErrorCode Parameter_settings(AppCtx*); 54 PetscErrorCode ini_bou(Vec,AppCtx*); 55 PetscErrorCode IFunction(TS,PetscReal,Vec,Vec,Vec,void*); 56 PetscErrorCode IJacobian(TS,PetscReal,Vec,Vec,PetscReal,Mat,Mat,void*); 57 PetscErrorCode PostStep(TS); 58 59 int main(int argc, char **argv) 60 { 61 PetscErrorCode ierr; 62 Vec x; /* Solution vector */ 63 TS ts; /* Time-stepping context */ 64 AppCtx user; /* Application context */ 65 PetscViewer viewer; 66 67 ierr = PetscInitialize(&argc,&argv,"petscopt_ex8", help);if (ierr) return ierr; 68 69 /* Get physics and time parameters */ 70 ierr = Parameter_settings(&user);CHKERRQ(ierr); 71 /* Create a 2D DA with dof = 1 */ 72 ierr = DMDACreate2d(PETSC_COMM_WORLD,user.bx,user.by,DMDA_STENCIL_STAR,4,4,PETSC_DECIDE,PETSC_DECIDE,1,user.st_width,NULL,NULL,&user.da);CHKERRQ(ierr); 73 ierr = DMSetFromOptions(user.da);CHKERRQ(ierr); 74 ierr = DMSetUp(user.da);CHKERRQ(ierr); 75 /* Set x and y coordinates */ 76 ierr = DMDASetUniformCoordinates(user.da,user.xmin,user.xmax,user.ymin,user.ymax,0,0);CHKERRQ(ierr); 77 ierr = DMDASetCoordinateName(user.da,0,"X - the angle");CHKERRQ(ierr); 78 ierr = DMDASetCoordinateName(user.da,1,"Y - the speed");CHKERRQ(ierr); 79 80 /* Get global vector x from DM */ 81 ierr = DMCreateGlobalVector(user.da,&x);CHKERRQ(ierr); 82 83 ierr = ini_bou(x,&user);CHKERRQ(ierr); 84 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"ini_x",FILE_MODE_WRITE,&viewer);CHKERRQ(ierr); 85 ierr = VecView(x,viewer);CHKERRQ(ierr); 86 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 87 88 ierr = TSCreate(PETSC_COMM_WORLD,&ts);CHKERRQ(ierr); 89 ierr = TSSetDM(ts,user.da);CHKERRQ(ierr); 90 ierr = TSSetProblemType(ts,TS_NONLINEAR);CHKERRQ(ierr); 91 ierr = TSSetType(ts,TSARKIMEX);CHKERRQ(ierr); 92 ierr = TSSetIFunction(ts,NULL,IFunction,&user);CHKERRQ(ierr); 93 /* ierr = TSSetIJacobian(ts,NULL,NULL,IJacobian,&user);CHKERRQ(ierr); */ 94 ierr = TSSetApplicationContext(ts,&user);CHKERRQ(ierr); 95 ierr = TSSetTimeStep(ts,.005);CHKERRQ(ierr); 96 ierr = TSSetFromOptions(ts);CHKERRQ(ierr); 97 ierr = TSSetPostStep(ts,PostStep);CHKERRQ(ierr); 98 ierr = TSSolve(ts,x);CHKERRQ(ierr); 99 100 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"fin_x",FILE_MODE_WRITE,&viewer);CHKERRQ(ierr); 101 ierr = VecView(x,viewer);CHKERRQ(ierr); 102 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 103 104 ierr = VecDestroy(&x);CHKERRQ(ierr); 105 ierr = DMDestroy(&user.da);CHKERRQ(ierr); 106 ierr = TSDestroy(&ts);CHKERRQ(ierr); 107 ierr = PetscFinalize(); 108 return ierr; 109 } 110 111 PetscErrorCode PostStep(TS ts) 112 { 113 PetscErrorCode ierr; 114 Vec X; 115 AppCtx *user; 116 PetscReal t; 117 PetscScalar asum; 118 119 PetscFunctionBegin; 120 ierr = TSGetApplicationContext(ts,&user);CHKERRQ(ierr); 121 ierr = TSGetTime(ts,&t);CHKERRQ(ierr); 122 ierr = TSGetSolution(ts,&X);CHKERRQ(ierr); 123 /* 124 if (t >= .2) { 125 ierr = TSGetSolution(ts,&X);CHKERRQ(ierr); 126 ierr = VecView(X,PETSC_VIEWER_BINARY_WORLD);CHKERRQ(ierr); 127 exit(0); 128 results in initial conditions after fault in binaryoutput 129 }*/ 130 131 if ((t > user->tf) && (t < user->tcl)) user->Pmax = 0.0; /* A short-circuit that drives the electrical power output (Pmax*sin(delta)) to zero */ 132 else user->Pmax = user->Pmax_s; 133 134 ierr = VecSum(X,&asum);CHKERRQ(ierr); 135 ierr = PetscPrintf(PETSC_COMM_WORLD,"sum(p) at t = %f = %f\n",(double)t,(double)(asum));CHKERRQ(ierr); 136 PetscFunctionReturn(0); 137 } 138 139 PetscErrorCode ini_bou(Vec X,AppCtx* user) 140 { 141 PetscErrorCode ierr; 142 DM cda; 143 DMDACoor2d **coors; 144 PetscScalar **p; 145 Vec gc; 146 PetscInt M,N,Ir,J; 147 PetscMPIInt rank; 148 149 PetscFunctionBeginUser; 150 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr); 151 ierr = DMDAGetInfo(user->da,NULL,&M,&N,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);CHKERRQ(ierr); 152 user->dx = (user->xmax - user->xmin)/(M-1); user->dy = (user->ymax - user->ymin)/(N-1); 153 ierr = DMGetCoordinateDM(user->da,&cda);CHKERRQ(ierr); 154 ierr = DMGetCoordinates(user->da,&gc);CHKERRQ(ierr); 155 ierr = DMDAVecGetArrayRead(cda,gc,&coors);CHKERRQ(ierr); 156 ierr = DMDAVecGetArray(user->da,X,&p);CHKERRQ(ierr); 157 158 /* Point mass at (mux,muy) */ 159 ierr = PetscPrintf(PETSC_COMM_WORLD,"Original user->mux = %f, user->muy = %f\n",user->mux,user->muy);CHKERRQ(ierr); 160 ierr = DMDAGetLogicalCoordinate(user->da,user->mux,user->muy,0.0,&Ir,&J,NULL,&user->mux,&user->muy,NULL);CHKERRQ(ierr); 161 user->PM_min = user->Pmax*PetscSinScalar(user->mux); 162 ierr = PetscPrintf(PETSC_COMM_WORLD,"Corrected user->mux = %f, user->muy = %f user->PM_min = %f,user->dx = %f\n",user->mux,user->muy,user->PM_min,user->dx);CHKERRQ(ierr); 163 if (Ir > -1 && J > -1) { 164 p[J][Ir] = 1.0; 165 } 166 167 ierr = DMDAVecRestoreArrayRead(cda,gc,&coors);CHKERRQ(ierr); 168 ierr = DMDAVecRestoreArray(user->da,X,&p);CHKERRQ(ierr); 169 PetscFunctionReturn(0); 170 } 171 172 /* First advection term */ 173 PetscErrorCode adv1(PetscScalar **p,PetscScalar y,PetscInt i,PetscInt j,PetscInt M,PetscScalar *p1,AppCtx *user) 174 { 175 PetscScalar f,fpos,fneg; 176 PetscFunctionBegin; 177 f = (y - user->ws); 178 fpos = PetscMax(f,0); 179 fneg = PetscMin(f,0); 180 if (user->st_width == 1) { 181 *p1 = fpos*(p[j][i] - p[j][i-1])/user->dx + fneg*(p[j][i+1] - p[j][i])/user->dx; 182 } else if (user->st_width == 2) { 183 *p1 = fpos*(3*p[j][i] - 4*p[j][i-1] + p[j][i-2])/(2*user->dx) + fneg*(-p[j][i+2] + 4*p[j][i+1] - 3*p[j][i])/(2*user->dx); 184 } else if (user->st_width == 3) { 185 *p1 = fpos*(2*p[j][i+1] + 3*p[j][i] - 6*p[j][i-1] + p[j][i-2])/(6*user->dx) + fneg*(-p[j][i+2] + 6*p[j][i+1] - 3*p[j][i] - 2*p[j][i-1])/(6*user->dx); 186 } 187 /* *p1 = f*(p[j][i+1] - p[j][i-1])/user->dx;*/ 188 PetscFunctionReturn(0); 189 } 190 191 /* Second advection term */ 192 PetscErrorCode adv2(PetscScalar **p,PetscScalar x,PetscScalar y,PetscInt i,PetscInt j,PetscInt N,PetscScalar *p2,AppCtx *user) 193 { 194 PetscScalar f,fpos,fneg; 195 PetscFunctionBegin; 196 f = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*PetscSinScalar(x) - user->D*(y - user->ws)); 197 fpos = PetscMax(f,0); 198 fneg = PetscMin(f,0); 199 if (user->st_width == 1) { 200 *p2 = fpos*(p[j][i] - p[j-1][i])/user->dy + fneg*(p[j+1][i] - p[j][i])/user->dy; 201 } else if (user->st_width ==2) { 202 *p2 = fpos*(3*p[j][i] - 4*p[j-1][i] + p[j-2][i])/(2*user->dy) + fneg*(-p[j+2][i] + 4*p[j+1][i] - 3*p[j][i])/(2*user->dy); 203 } else if (user->st_width == 3) { 204 *p2 = fpos*(2*p[j+1][i] + 3*p[j][i] - 6*p[j-1][i] + p[j-2][i])/(6*user->dy) + fneg*(-p[j+2][i] + 6*p[j+1][i] - 3*p[j][i] - 2*p[j-1][i])/(6*user->dy); 205 } 206 207 /* *p2 = f*(p[j+1][i] - p[j-1][i])/user->dy;*/ 208 PetscFunctionReturn(0); 209 } 210 211 /* Diffusion term */ 212 PetscErrorCode diffuse(PetscScalar **p,PetscInt i,PetscInt j,PetscReal t,PetscScalar *p_diff,AppCtx * user) 213 { 214 PetscFunctionBeginUser; 215 if (user->st_width == 1) { 216 *p_diff = user->disper_coe*((p[j-1][i] - 2*p[j][i] + p[j+1][i])/(user->dy*user->dy)); 217 } else if (user->st_width == 2) { 218 *p_diff = user->disper_coe*((-p[j-2][i] + 16*p[j-1][i] - 30*p[j][i] + 16*p[j+1][i] - p[j+2][i])/(12.0*user->dy*user->dy)); 219 } else if (user->st_width == 3) { 220 *p_diff = user->disper_coe*((2*p[j-3][i] - 27*p[j-2][i] + 270*p[j-1][i] - 490*p[j][i] + 270*p[j+1][i] - 27*p[j+2][i] + 2*p[j+3][i])/(180.0*user->dy*user->dy)); 221 } 222 PetscFunctionReturn(0); 223 } 224 225 PetscErrorCode IFunction(TS ts,PetscReal t,Vec X,Vec Xdot,Vec F,void *ctx) 226 { 227 PetscErrorCode ierr; 228 AppCtx *user=(AppCtx*)ctx; 229 DM cda; 230 DMDACoor2d **coors; 231 PetscScalar **p,**f,**pdot; 232 PetscInt i,j; 233 PetscInt xs,ys,xm,ym,M,N; 234 Vec localX,gc,localXdot; 235 PetscScalar p_adv1 = 0.0,p_adv2 = 0.0,p_diff = 0; 236 PetscScalar diffuse1,gamma; 237 238 PetscFunctionBeginUser; 239 ierr = DMDAGetInfo(user->da,NULL,&M,&N,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);CHKERRQ(ierr); 240 ierr = DMGetCoordinateDM(user->da,&cda);CHKERRQ(ierr); 241 ierr = DMDAGetCorners(cda,&xs,&ys,0,&xm,&ym,0);CHKERRQ(ierr); 242 243 ierr = DMGetLocalVector(user->da,&localX);CHKERRQ(ierr); 244 ierr = DMGetLocalVector(user->da,&localXdot);CHKERRQ(ierr); 245 246 ierr = DMGlobalToLocalBegin(user->da,X,INSERT_VALUES,localX);CHKERRQ(ierr); 247 ierr = DMGlobalToLocalEnd(user->da,X,INSERT_VALUES,localX);CHKERRQ(ierr); 248 ierr = DMGlobalToLocalBegin(user->da,Xdot,INSERT_VALUES,localXdot);CHKERRQ(ierr); 249 ierr = DMGlobalToLocalEnd(user->da,Xdot,INSERT_VALUES,localXdot);CHKERRQ(ierr); 250 251 ierr = DMGetCoordinatesLocal(user->da,&gc);CHKERRQ(ierr); 252 253 ierr = DMDAVecGetArrayRead(cda,gc,&coors);CHKERRQ(ierr); 254 ierr = DMDAVecGetArrayRead(user->da,localX,&p);CHKERRQ(ierr); 255 ierr = DMDAVecGetArrayRead(user->da,localXdot,&pdot);CHKERRQ(ierr); 256 ierr = DMDAVecGetArray(user->da,F,&f);CHKERRQ(ierr); 257 258 gamma = user->D*user->ws/(2*user->H); 259 diffuse1 = user->lambda*user->lambda*user->q/(user->lambda*gamma+1)*(1.0 - PetscExpScalar(-t*(gamma+1.0)/user->lambda)); 260 user->disper_coe = user->ws*user->ws/(4*user->H*user->H)*diffuse1; 261 262 for (i=xs; i < xs+xm; i++) { 263 for (j=ys; j < ys+ym; j++) { 264 ierr = adv1(p,coors[j][i].y,i,j,M,&p_adv1,user);CHKERRQ(ierr); 265 ierr = adv2(p,coors[j][i].x,coors[j][i].y,i,j,N,&p_adv2,user);CHKERRQ(ierr); 266 ierr = diffuse(p,i,j,t,&p_diff,user);CHKERRQ(ierr); 267 f[j][i] = -p_adv1 - p_adv2 + p_diff - pdot[j][i]; 268 } 269 } 270 ierr = DMDAVecRestoreArrayRead(user->da,localX,&p);CHKERRQ(ierr); 271 ierr = DMDAVecRestoreArrayRead(user->da,localX,&pdot);CHKERRQ(ierr); 272 ierr = DMRestoreLocalVector(user->da,&localX);CHKERRQ(ierr); 273 ierr = DMRestoreLocalVector(user->da,&localXdot);CHKERRQ(ierr); 274 ierr = DMDAVecRestoreArray(user->da,F,&f);CHKERRQ(ierr); 275 ierr = DMDAVecRestoreArrayRead(cda,gc,&coors);CHKERRQ(ierr); 276 277 PetscFunctionReturn(0); 278 } 279 280 PetscErrorCode IJacobian(TS ts,PetscReal t,Vec X,Vec Xdot,PetscReal a,Mat J,Mat Jpre,void *ctx) 281 { 282 PetscErrorCode ierr; 283 AppCtx *user=(AppCtx*)ctx; 284 DM cda; 285 DMDACoor2d **coors; 286 PetscInt i,j; 287 PetscInt xs,ys,xm,ym,M,N; 288 Vec gc; 289 PetscScalar val[5],xi,yi; 290 MatStencil row,col[5]; 291 PetscScalar c1,c3,c5,c1pos,c1neg,c3pos,c3neg; 292 293 PetscFunctionBeginUser; 294 ierr = DMDAGetInfo(user->da,NULL,&M,&N,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);CHKERRQ(ierr); 295 ierr = DMGetCoordinateDM(user->da,&cda);CHKERRQ(ierr); 296 ierr = DMDAGetCorners(cda,&xs,&ys,0,&xm,&ym,0);CHKERRQ(ierr); 297 298 ierr = DMGetCoordinatesLocal(user->da,&gc);CHKERRQ(ierr); 299 ierr = DMDAVecGetArrayRead(cda,gc,&coors);CHKERRQ(ierr); 300 for (i=xs; i < xs+xm; i++) { 301 for (j=ys; j < ys+ym; j++) { 302 PetscInt nc = 0; 303 xi = coors[j][i].x; yi = coors[j][i].y; 304 row.i = i; row.j = j; 305 c1 = (yi-user->ws)/user->dx; 306 c1pos = PetscMax(c1,0); 307 c1neg = PetscMin(c1,0); 308 c3 = (user->ws/(2.0*user->H))*(user->PM_min - user->Pmax*PetscSinScalar(xi) - user->D*(yi - user->ws))/user->dy; 309 c3pos = PetscMax(c3,0); 310 c3neg = PetscMin(c3,0); 311 c5 = (PetscPowScalar((user->lambda*user->ws)/(2*user->H),2)*user->q*(1.0-PetscExpScalar(-t/user->lambda)))/(user->dy*user->dy); 312 col[nc].i = i-1; col[nc].j = j; val[nc++] = c1pos; 313 col[nc].i = i+1; col[nc].j = j; val[nc++] = -c1neg; 314 col[nc].i = i; col[nc].j = j-1; val[nc++] = c3pos + c5; 315 col[nc].i = i; col[nc].j = j+1; val[nc++] = -c3neg + c5; 316 col[nc].i = i; col[nc].j = j; val[nc++] = -c1pos + c1neg -c3pos + c3neg -2*c5 -a; 317 ierr = MatSetValuesStencil(Jpre,1,&row,nc,col,val,INSERT_VALUES);CHKERRQ(ierr); 318 } 319 } 320 ierr = DMDAVecRestoreArrayRead(cda,gc,&coors);CHKERRQ(ierr); 321 322 ierr = MatAssemblyBegin(Jpre,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 323 ierr = MatAssemblyEnd(Jpre,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 324 if (J != Jpre) { 325 ierr = MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 326 ierr = MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 327 } 328 PetscFunctionReturn(0); 329 } 330 331 PetscErrorCode Parameter_settings(AppCtx *user) 332 { 333 PetscErrorCode ierr; 334 PetscBool flg; 335 336 PetscFunctionBeginUser; 337 338 /* Set default parameters */ 339 user->ws = 1.0; 340 user->H = 5.0; 341 user->D = 0.0; 342 user->Pmax = user->Pmax_s = 2.1; 343 user->PM_min = 1.0; 344 user->lambda = 0.1; 345 user->q = 1.0; 346 user->mux = PetscAsinScalar(user->PM_min/user->Pmax); 347 user->sigmax = 0.1; 348 user->sigmay = 0.1; 349 user->rho = 0.0; 350 user->xmin = -PETSC_PI; 351 user->xmax = PETSC_PI; 352 user->bx = DM_BOUNDARY_PERIODIC; 353 user->by = DM_BOUNDARY_GHOSTED; 354 user->tf = user->tcl = -1; 355 user->ymin = -2.0; 356 user->ymax = 2.0; 357 user->st_width = 1; 358 359 ierr = PetscOptionsGetScalar(NULL,NULL,"-ws",&user->ws,&flg);CHKERRQ(ierr); 360 ierr = PetscOptionsGetScalar(NULL,NULL,"-Inertia",&user->H,&flg);CHKERRQ(ierr); 361 ierr = PetscOptionsGetScalar(NULL,NULL,"-D",&user->D,&flg);CHKERRQ(ierr); 362 ierr = PetscOptionsGetScalar(NULL,NULL,"-Pmax",&user->Pmax,&flg);CHKERRQ(ierr); 363 ierr = PetscOptionsGetScalar(NULL,NULL,"-PM_min",&user->PM_min,&flg);CHKERRQ(ierr); 364 ierr = PetscOptionsGetScalar(NULL,NULL,"-lambda",&user->lambda,&flg);CHKERRQ(ierr); 365 ierr = PetscOptionsGetScalar(NULL,NULL,"-q",&user->q,&flg);CHKERRQ(ierr); 366 ierr = PetscOptionsGetScalar(NULL,NULL,"-mux",&user->mux,&flg);CHKERRQ(ierr); 367 ierr = PetscOptionsGetScalar(NULL,NULL,"-muy",&user->muy,&flg);CHKERRQ(ierr); 368 if (flg == 0) { 369 user->muy = user->ws; 370 } 371 ierr = PetscOptionsGetScalar(NULL,NULL,"-xmin",&user->xmin,&flg);CHKERRQ(ierr); 372 ierr = PetscOptionsGetScalar(NULL,NULL,"-xmax",&user->xmax,&flg);CHKERRQ(ierr); 373 ierr = PetscOptionsGetScalar(NULL,NULL,"-ymin",&user->ymin,&flg);CHKERRQ(ierr); 374 ierr = PetscOptionsGetScalar(NULL,NULL,"-ymax",&user->ymax,&flg);CHKERRQ(ierr); 375 ierr = PetscOptionsGetInt(NULL,NULL,"-stencil_width",&user->st_width,&flg);CHKERRQ(ierr); 376 ierr = PetscOptionsGetEnum(NULL,NULL,"-bx",DMBoundaryTypes,(PetscEnum*)&user->bx,&flg);CHKERRQ(ierr); 377 ierr = PetscOptionsGetEnum(NULL,NULL,"-by",DMBoundaryTypes,(PetscEnum*)&user->by,&flg);CHKERRQ(ierr); 378 ierr = PetscOptionsGetReal(NULL,NULL,"-tf",&user->tf,&flg);CHKERRQ(ierr); 379 ierr = PetscOptionsGetReal(NULL,NULL,"-tcl",&user->tcl,&flg);CHKERRQ(ierr); 380 PetscFunctionReturn(0); 381 } 382 383 /*TEST 384 385 build: 386 requires: !complex x 387 388 test: 389 args: -ts_max_steps 1 390 localrunfiles: petscopt_ex8 391 timeoutfactor: 3 392 393 TEST*/ 394