1 2 static char help[] = "Time-dependent PDE in 2d for calculating joint PDF. \n"; 3 /* 4 p_t = -x_t*p_x -y_t*p_y + f(t)*p_yy 5 xmin < x < xmax, ymin < y < ymax; 6 7 Boundary conditions Neumman using mirror values 8 9 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. 10 x_t = (y - ws) y_t = (ws/2H)*(Pm - Pmax*sin(x)) 11 12 */ 13 14 #include <petscdm.h> 15 #include <petscdmda.h> 16 #include <petscts.h> 17 18 /* 19 User-defined data structures and routines 20 */ 21 typedef struct { 22 PetscScalar ws; /* Synchronous speed */ 23 PetscScalar H; /* Inertia constant */ 24 PetscScalar D; /* Damping constant */ 25 PetscScalar Pmax; /* Maximum power output of generator */ 26 PetscScalar PM_min; /* Mean mechanical power input */ 27 PetscScalar lambda; /* correlation time */ 28 PetscScalar q; /* noise strength */ 29 PetscScalar mux; /* Initial average angle */ 30 PetscScalar sigmax; /* Standard deviation of initial angle */ 31 PetscScalar muy; /* Average speed */ 32 PetscScalar sigmay; /* standard deviation of initial speed */ 33 PetscScalar rho; /* Cross-correlation coefficient */ 34 PetscScalar xmin; /* left boundary of angle */ 35 PetscScalar xmax; /* right boundary of angle */ 36 PetscScalar ymin; /* bottom boundary of speed */ 37 PetscScalar ymax; /* top boundary of speed */ 38 PetscScalar dx; /* x step size */ 39 PetscScalar dy; /* y step size */ 40 PetscScalar disper_coe; /* Dispersion coefficient */ 41 DM da; 42 PetscInt st_width; /* Stencil width */ 43 DMBoundaryType bx; /* x boundary type */ 44 DMBoundaryType by; /* y boundary type */ 45 PetscBool nonoiseinitial; 46 } AppCtx; 47 48 PetscErrorCode Parameter_settings(AppCtx *); 49 PetscErrorCode ini_bou(Vec, AppCtx *); 50 PetscErrorCode IFunction(TS, PetscReal, Vec, Vec, Vec, void *); 51 PetscErrorCode IJacobian(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *); 52 PetscErrorCode PostStep(TS); 53 54 int main(int argc, char **argv) 55 { 56 Vec x; /* Solution vector */ 57 TS ts; /* Time-stepping context */ 58 AppCtx user; /* Application context */ 59 PetscViewer viewer; 60 61 PetscFunctionBeginUser; 62 PetscCall(PetscInitialize(&argc, &argv, "petscopt_ex7", help)); 63 64 /* Get physics and time parameters */ 65 PetscCall(Parameter_settings(&user)); 66 /* Create a 2D DA with dof = 1 */ 67 PetscCall(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)); 68 PetscCall(DMSetFromOptions(user.da)); 69 PetscCall(DMSetUp(user.da)); 70 /* Set x and y coordinates */ 71 PetscCall(DMDASetUniformCoordinates(user.da, user.xmin, user.xmax, user.ymin, user.ymax, 0, 0)); 72 PetscCall(DMDASetCoordinateName(user.da, 0, "X - the angle")); 73 PetscCall(DMDASetCoordinateName(user.da, 1, "Y - the speed")); 74 75 /* Get global vector x from DM */ 76 PetscCall(DMCreateGlobalVector(user.da, &x)); 77 78 PetscCall(ini_bou(x, &user)); 79 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "ini_x", FILE_MODE_WRITE, &viewer)); 80 PetscCall(VecView(x, viewer)); 81 PetscCall(PetscViewerDestroy(&viewer)); 82 83 PetscCall(TSCreate(PETSC_COMM_WORLD, &ts)); 84 PetscCall(TSSetDM(ts, user.da)); 85 PetscCall(TSSetProblemType(ts, TS_NONLINEAR)); 86 PetscCall(TSSetType(ts, TSARKIMEX)); 87 PetscCall(TSSetIFunction(ts, NULL, IFunction, &user)); 88 /* PetscCall(TSSetIJacobian(ts,NULL,NULL,IJacobian,&user)); */ 89 PetscCall(TSSetApplicationContext(ts, &user)); 90 PetscCall(TSSetTimeStep(ts, .005)); 91 PetscCall(TSSetFromOptions(ts)); 92 PetscCall(TSSetPostStep(ts, PostStep)); 93 PetscCall(TSSolve(ts, x)); 94 95 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "fin_x", FILE_MODE_WRITE, &viewer)); 96 PetscCall(VecView(x, viewer)); 97 PetscCall(PetscViewerDestroy(&viewer)); 98 99 PetscCall(VecDestroy(&x)); 100 PetscCall(DMDestroy(&user.da)); 101 PetscCall(TSDestroy(&ts)); 102 PetscCall(PetscFinalize()); 103 return 0; 104 } 105 106 PetscErrorCode PostStep(TS ts) 107 { 108 Vec X, gc; 109 AppCtx *user; 110 PetscScalar sum = 0, asum; 111 PetscReal t, **p; 112 DMDACoor2d **coors; 113 DM cda; 114 PetscInt i, j, xs, ys, xm, ym; 115 116 PetscFunctionBegin; 117 PetscCall(TSGetApplicationContext(ts, &user)); 118 PetscCall(TSGetTime(ts, &t)); 119 PetscCall(TSGetSolution(ts, &X)); 120 121 PetscCall(DMGetCoordinateDM(user->da, &cda)); 122 PetscCall(DMDAGetCorners(cda, &xs, &ys, 0, &xm, &ym, 0)); 123 PetscCall(DMGetCoordinates(user->da, &gc)); 124 PetscCall(DMDAVecGetArrayRead(cda, gc, &coors)); 125 PetscCall(DMDAVecGetArrayRead(user->da, X, &p)); 126 for (i = xs; i < xs + xm; i++) { 127 for (j = ys; j < ys + ym; j++) { 128 if (coors[j][i].y < 5) sum += p[j][i]; 129 } 130 } 131 PetscCall(DMDAVecRestoreArrayRead(cda, gc, &coors)); 132 PetscCall(DMDAVecRestoreArrayRead(user->da, X, &p)); 133 PetscCallMPI(MPI_Allreduce(&sum, &asum, 1, MPIU_SCALAR, MPIU_SUM, PetscObjectComm((PetscObject)ts))); 134 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "sum(p)*dw*dtheta at t = %f = %f\n", (double)t, (double)(asum))); 135 if (sum < 1.0e-2) { 136 PetscCall(TSSetConvergedReason(ts, TS_CONVERGED_USER)); 137 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Exiting TS as the integral of PDF is almost zero\n")); 138 } 139 PetscFunctionReturn(PETSC_SUCCESS); 140 } 141 142 PetscErrorCode ini_bou(Vec X, AppCtx *user) 143 { 144 DM cda; 145 DMDACoor2d **coors; 146 PetscScalar **p; 147 Vec gc; 148 PetscInt i, j; 149 PetscInt xs, ys, xm, ym, M, N; 150 PetscScalar xi, yi; 151 PetscScalar sigmax = user->sigmax, sigmay = user->sigmay; 152 PetscScalar rho = user->rho; 153 PetscScalar muy = user->muy, mux; 154 PetscMPIInt rank; 155 PetscScalar sum; 156 157 PetscFunctionBeginUser; 158 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 159 PetscCall(DMDAGetInfo(user->da, NULL, &M, &N, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL)); 160 user->dx = (user->xmax - user->xmin) / (M - 1); 161 user->dy = (user->ymax - user->ymin) / (N - 1); 162 PetscCall(DMGetCoordinateDM(user->da, &cda)); 163 PetscCall(DMGetCoordinates(user->da, &gc)); 164 PetscCall(DMDAVecGetArray(cda, gc, &coors)); 165 PetscCall(DMDAVecGetArray(user->da, X, &p)); 166 PetscCall(DMDAGetCorners(cda, &xs, &ys, 0, &xm, &ym, 0)); 167 168 /* mux and muy need to be grid points in the x and y-direction otherwise the solution goes unstable 169 muy is set by choosing the y domain, no. of grid points along y-direction so that muy is a grid point 170 in the y-direction. We only modify mux here 171 */ 172 mux = user->mux = coors[0][M / 2 + 10].x; /* For -pi < x < pi, this should be some angle between 0 and pi/2 */ 173 if (user->nonoiseinitial) { 174 for (i = xs; i < xs + xm; i++) { 175 for (j = ys; j < ys + ym; j++) { 176 xi = coors[j][i].x; 177 yi = coors[j][i].y; 178 if ((xi == mux) && (yi == muy)) p[j][i] = 1.0; 179 } 180 } 181 } else { 182 /* Change PM_min accordingly */ 183 user->PM_min = user->Pmax * PetscSinScalar(mux); 184 for (i = xs; i < xs + xm; i++) { 185 for (j = ys; j < ys + ym; j++) { 186 xi = coors[j][i].x; 187 yi = coors[j][i].y; 188 p[j][i] = (0.5 / (PETSC_PI * sigmax * sigmay * PetscSqrtScalar(1.0 - rho * rho))) * PetscExpScalar(-0.5 / (1 - rho * rho) * (PetscPowScalar((xi - mux) / sigmax, 2) + PetscPowScalar((yi - muy) / sigmay, 2) - 2 * rho * (xi - mux) * (yi - muy) / (sigmax * sigmay))); 189 } 190 } 191 } 192 PetscCall(DMDAVecRestoreArray(cda, gc, &coors)); 193 PetscCall(DMDAVecRestoreArray(user->da, X, &p)); 194 PetscCall(VecSum(X, &sum)); 195 PetscCall(VecScale(X, 1.0 / sum)); 196 PetscFunctionReturn(PETSC_SUCCESS); 197 } 198 199 /* First advection term */ 200 PetscErrorCode adv1(PetscScalar **p, PetscScalar y, PetscInt i, PetscInt j, PetscInt M, PetscScalar *p1, AppCtx *user) 201 { 202 PetscScalar f, fpos, fneg; 203 PetscFunctionBegin; 204 f = (y - user->ws); 205 fpos = PetscMax(f, 0); 206 fneg = PetscMin(f, 0); 207 if (user->st_width == 1) { 208 *p1 = fpos * (p[j][i] - p[j][i - 1]) / user->dx + fneg * (p[j][i + 1] - p[j][i]) / user->dx; 209 } else if (user->st_width == 2) { 210 *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); 211 } else if (user->st_width == 3) { 212 *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); 213 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "No support for wider stencils"); 214 PetscFunctionReturn(PETSC_SUCCESS); 215 } 216 217 /* Second advection term */ 218 PetscErrorCode adv2(PetscScalar **p, PetscScalar x, PetscInt i, PetscInt j, PetscInt N, PetscScalar *p2, AppCtx *user) 219 { 220 PetscScalar f, fpos, fneg; 221 PetscFunctionBegin; 222 f = (user->ws / (2 * user->H)) * (user->PM_min - user->Pmax * PetscSinScalar(x)); 223 fpos = PetscMax(f, 0); 224 fneg = PetscMin(f, 0); 225 if (user->st_width == 1) { 226 *p2 = fpos * (p[j][i] - p[j - 1][i]) / user->dy + fneg * (p[j + 1][i] - p[j][i]) / user->dy; 227 } else if (user->st_width == 2) { 228 *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); 229 } else if (user->st_width == 3) { 230 *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); 231 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "No support for wider stencils"); 232 PetscFunctionReturn(PETSC_SUCCESS); 233 } 234 235 /* Diffusion term */ 236 PetscErrorCode diffuse(PetscScalar **p, PetscInt i, PetscInt j, PetscReal t, PetscScalar *p_diff, AppCtx *user) 237 { 238 PetscFunctionBeginUser; 239 if (user->st_width == 1) { 240 *p_diff = user->disper_coe * ((p[j - 1][i] - 2 * p[j][i] + p[j + 1][i]) / (user->dy * user->dy)); 241 } else if (user->st_width == 2) { 242 *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)); 243 } else if (user->st_width == 3) { 244 *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)); 245 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "No support for wider stencils"); 246 PetscFunctionReturn(PETSC_SUCCESS); 247 } 248 249 PetscErrorCode IFunction(TS ts, PetscReal t, Vec X, Vec Xdot, Vec F, void *ctx) 250 { 251 AppCtx *user = (AppCtx *)ctx; 252 DM cda; 253 DMDACoor2d **coors; 254 PetscScalar **p, **f, **pdot; 255 PetscInt i, j; 256 PetscInt xs, ys, xm, ym, M, N; 257 Vec localX, gc, localXdot; 258 PetscScalar p_adv1 = 0.0, p_adv2 = 0.0, p_diff; /* initialize to prevent incorrect compiler warnings */ 259 260 PetscFunctionBeginUser; 261 PetscCall(DMDAGetInfo(user->da, NULL, &M, &N, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL)); 262 PetscCall(DMGetCoordinateDM(user->da, &cda)); 263 PetscCall(DMDAGetCorners(cda, &xs, &ys, 0, &xm, &ym, 0)); 264 265 PetscCall(DMGetLocalVector(user->da, &localX)); 266 PetscCall(DMGetLocalVector(user->da, &localXdot)); 267 268 PetscCall(DMGlobalToLocalBegin(user->da, X, INSERT_VALUES, localX)); 269 PetscCall(DMGlobalToLocalEnd(user->da, X, INSERT_VALUES, localX)); 270 PetscCall(DMGlobalToLocalBegin(user->da, Xdot, INSERT_VALUES, localXdot)); 271 PetscCall(DMGlobalToLocalEnd(user->da, Xdot, INSERT_VALUES, localXdot)); 272 273 PetscCall(DMGetCoordinatesLocal(user->da, &gc)); 274 275 PetscCall(DMDAVecGetArrayRead(cda, gc, &coors)); 276 PetscCall(DMDAVecGetArrayRead(user->da, localX, &p)); 277 PetscCall(DMDAVecGetArrayRead(user->da, localXdot, &pdot)); 278 PetscCall(DMDAVecGetArray(user->da, F, &f)); 279 280 user->disper_coe = PetscPowScalar((user->lambda * user->ws) / (2 * user->H), 2) * user->q * (1.0 - PetscExpScalar(-t / user->lambda)); 281 for (i = xs; i < xs + xm; i++) { 282 for (j = ys; j < ys + ym; j++) { 283 PetscCall(adv1(p, coors[j][i].y, i, j, M, &p_adv1, user)); 284 PetscCall(adv2(p, coors[j][i].x, i, j, N, &p_adv2, user)); 285 PetscCall(diffuse(p, i, j, t, &p_diff, user)); 286 f[j][i] = -p_adv1 - p_adv2 + p_diff - pdot[j][i]; 287 } 288 } 289 PetscCall(DMDAVecRestoreArrayRead(user->da, localX, &p)); 290 PetscCall(DMDAVecRestoreArrayRead(user->da, localX, &pdot)); 291 PetscCall(DMRestoreLocalVector(user->da, &localX)); 292 PetscCall(DMRestoreLocalVector(user->da, &localXdot)); 293 PetscCall(DMDAVecRestoreArray(user->da, F, &f)); 294 PetscCall(DMDAVecRestoreArrayRead(cda, gc, &coors)); 295 296 PetscFunctionReturn(PETSC_SUCCESS); 297 } 298 299 PetscErrorCode IJacobian(TS ts, PetscReal t, Vec X, Vec Xdot, PetscReal a, Mat J, Mat Jpre, void *ctx) 300 { 301 AppCtx *user = (AppCtx *)ctx; 302 DM cda; 303 DMDACoor2d **coors; 304 PetscInt i, j; 305 PetscInt xs, ys, xm, ym, M, N; 306 Vec gc; 307 PetscScalar val[5], xi, yi; 308 MatStencil row, col[5]; 309 PetscScalar c1, c3, c5, c1pos, c1neg, c3pos, c3neg; 310 311 PetscFunctionBeginUser; 312 PetscCall(DMDAGetInfo(user->da, NULL, &M, &N, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL)); 313 PetscCall(DMGetCoordinateDM(user->da, &cda)); 314 PetscCall(DMDAGetCorners(cda, &xs, &ys, 0, &xm, &ym, 0)); 315 316 PetscCall(DMGetCoordinatesLocal(user->da, &gc)); 317 PetscCall(DMDAVecGetArrayRead(cda, gc, &coors)); 318 for (i = xs; i < xs + xm; i++) { 319 for (j = ys; j < ys + ym; j++) { 320 PetscInt nc = 0; 321 xi = coors[j][i].x; 322 yi = coors[j][i].y; 323 row.i = i; 324 row.j = j; 325 c1 = (yi - user->ws) / user->dx; 326 c1pos = PetscMax(c1, 0); 327 c1neg = PetscMin(c1, 0); 328 c3 = (user->ws / (2.0 * user->H)) * (user->PM_min - user->Pmax * PetscSinScalar(xi)) / user->dy; 329 c3pos = PetscMax(c3, 0); 330 c3neg = PetscMin(c3, 0); 331 c5 = (PetscPowScalar((user->lambda * user->ws) / (2 * user->H), 2) * user->q * (1.0 - PetscExpScalar(-t / user->lambda))) / (user->dy * user->dy); 332 col[nc].i = i - 1; 333 col[nc].j = j; 334 val[nc++] = c1pos; 335 col[nc].i = i + 1; 336 col[nc].j = j; 337 val[nc++] = -c1neg; 338 col[nc].i = i; 339 col[nc].j = j - 1; 340 val[nc++] = c3pos + c5; 341 col[nc].i = i; 342 col[nc].j = j + 1; 343 val[nc++] = -c3neg + c5; 344 col[nc].i = i; 345 col[nc].j = j; 346 val[nc++] = -c1pos + c1neg - c3pos + c3neg - 2 * c5 - a; 347 PetscCall(MatSetValuesStencil(Jpre, 1, &row, nc, col, val, INSERT_VALUES)); 348 } 349 } 350 PetscCall(DMDAVecRestoreArrayRead(cda, gc, &coors)); 351 352 PetscCall(MatAssemblyBegin(Jpre, MAT_FINAL_ASSEMBLY)); 353 PetscCall(MatAssemblyEnd(Jpre, MAT_FINAL_ASSEMBLY)); 354 if (J != Jpre) { 355 PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY)); 356 PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY)); 357 } 358 PetscFunctionReturn(PETSC_SUCCESS); 359 } 360 361 PetscErrorCode Parameter_settings(AppCtx *user) 362 { 363 PetscBool flg; 364 365 PetscFunctionBeginUser; 366 367 /* Set default parameters */ 368 user->ws = 1.0; 369 user->H = 5.0; 370 user->Pmax = 2.1; 371 user->PM_min = 1.0; 372 user->lambda = 0.1; 373 user->q = 1.0; 374 user->mux = PetscAsinScalar(user->PM_min / user->Pmax); 375 user->sigmax = 0.1; 376 user->sigmay = 0.1; 377 user->rho = 0.0; 378 user->xmin = -PETSC_PI; 379 user->xmax = PETSC_PI; 380 user->bx = DM_BOUNDARY_PERIODIC; 381 user->by = DM_BOUNDARY_MIRROR; 382 user->nonoiseinitial = PETSC_FALSE; 383 384 /* 385 ymin of -3 seems to let the unstable solution move up and leave a zero in its wake 386 with an ymin of -1 the wake is never exactly zero 387 */ 388 user->ymin = -3.0; 389 user->ymax = 10.0; 390 user->st_width = 1; 391 392 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-ws", &user->ws, &flg)); 393 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-Inertia", &user->H, &flg)); 394 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-Pmax", &user->Pmax, &flg)); 395 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-PM_min", &user->PM_min, &flg)); 396 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-lambda", &user->lambda, &flg)); 397 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-q", &user->q, &flg)); 398 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-mux", &user->mux, &flg)); 399 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-sigmax", &user->sigmax, &flg)); 400 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-muy", &user->muy, &flg)); 401 if (flg == 0) user->muy = user->ws; 402 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-sigmay", &user->sigmay, &flg)); 403 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-rho", &user->rho, &flg)); 404 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-xmin", &user->xmin, &flg)); 405 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-xmax", &user->xmax, &flg)); 406 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-ymin", &user->ymin, &flg)); 407 PetscCall(PetscOptionsGetScalar(NULL, NULL, "-ymax", &user->ymax, &flg)); 408 PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_width", &user->st_width, &flg)); 409 PetscCall(PetscOptionsGetEnum(NULL, NULL, "-bx", DMBoundaryTypes, (PetscEnum *)&user->bx, &flg)); 410 PetscCall(PetscOptionsGetEnum(NULL, NULL, "-by", DMBoundaryTypes, (PetscEnum *)&user->by, &flg)); 411 PetscCall(PetscOptionsGetBool(NULL, NULL, "-nonoiseinitial", &user->nonoiseinitial, &flg)); 412 413 PetscFunctionReturn(PETSC_SUCCESS); 414 } 415 416 /*TEST 417 418 build: 419 requires: !complex !single 420 421 test: 422 args: -ts_max_steps 2 423 localrunfiles: petscopt_ex7 424 425 test: 426 suffix: 2 427 args: -ts_max_steps 2 -snes_mf_operator 428 output_file: output/ex7_1.out 429 localrunfiles: petscopt_ex7 430 timeoutfactor: 2 431 432 test: 433 suffix: 3 434 args: -ts_max_steps 2 -snes_mf -pc_type none 435 output_file: output/ex7_1.out 436 localrunfiles: petscopt_ex7 437 timeoutfactor: 2 438 439 TEST*/ 440