Actual source code: ex6.c
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: x_t = (y - ws) y_t = (ws/2H)*(Pm - Pmax*sin(x))
7: Boundary conditions: -bc_type 0 => Zero dirichlet boundary
8: -bc_type 1 => Steady state boundary condition
9: Steady state boundary condition found by setting p_t = 0
10: */
12: #include <petscdm.h>
13: #include <petscdmda.h>
14: #include <petscts.h>
16: /*
17: User-defined data structures and routines
18: */
19: typedef struct {
20: PetscScalar ws; /* Synchronous speed */
21: PetscScalar H; /* Inertia constant */
22: PetscScalar D; /* Damping constant */
23: PetscScalar Pmax; /* Maximum power output of generator */
24: PetscScalar PM_min; /* Mean mechanical power input */
25: PetscScalar lambda; /* correlation time */
26: PetscScalar q; /* noise strength */
27: PetscScalar mux; /* Initial average angle */
28: PetscScalar sigmax; /* Standard deviation of initial angle */
29: PetscScalar muy; /* Average speed */
30: PetscScalar sigmay; /* standard deviation of initial speed */
31: PetscScalar rho; /* Cross-correlation coefficient */
32: PetscScalar t0; /* Initial time */
33: PetscScalar tmax; /* Final time */
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: PetscInt bc; /* Boundary conditions */
41: PetscScalar disper_coe; /* Dispersion coefficient */
42: DM da;
43: } AppCtx;
45: PetscErrorCode Parameter_settings(AppCtx*);
46: PetscErrorCode ini_bou(Vec,AppCtx*);
47: PetscErrorCode IFunction(TS,PetscReal,Vec,Vec,Vec,void*);
48: PetscErrorCode IJacobian(TS,PetscReal,Vec,Vec,PetscReal,Mat,Mat,void*);
49: PetscErrorCode PostStep(TS);
51: int main(int argc, char **argv)
52: {
53: Vec x; /* Solution vector */
54: TS ts; /* Time-stepping context */
55: AppCtx user; /* Application context */
56: Mat J;
57: PetscViewer viewer;
59: PetscInitialize(&argc,&argv,"petscopt_ex6", help);
60: /* Get physics and time parameters */
61: Parameter_settings(&user);
62: /* Create a 2D DA with dof = 1 */
63: DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,4,4,PETSC_DECIDE,PETSC_DECIDE,1,1,NULL,NULL,&user.da);
64: DMSetFromOptions(user.da);
65: DMSetUp(user.da);
66: /* Set x and y coordinates */
67: DMDASetUniformCoordinates(user.da,user.xmin,user.xmax,user.ymin,user.ymax,0.0,1.0);
69: /* Get global vector x from DM */
70: DMCreateGlobalVector(user.da,&x);
72: ini_bou(x,&user);
73: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"ini_x",FILE_MODE_WRITE,&viewer);
74: VecView(x,viewer);
75: PetscViewerDestroy(&viewer);
77: /* Get Jacobian matrix structure from the da */
78: DMSetMatType(user.da,MATAIJ);
79: DMCreateMatrix(user.da,&J);
81: TSCreate(PETSC_COMM_WORLD,&ts);
82: TSSetProblemType(ts,TS_NONLINEAR);
83: TSSetIFunction(ts,NULL,IFunction,&user);
84: TSSetIJacobian(ts,J,J,IJacobian,&user);
85: TSSetApplicationContext(ts,&user);
86: TSSetMaxTime(ts,user.tmax);
87: TSSetExactFinalTime(ts,TS_EXACTFINALTIME_MATCHSTEP);
88: TSSetTime(ts,user.t0);
89: TSSetTimeStep(ts,.005);
90: TSSetFromOptions(ts);
91: TSSetPostStep(ts,PostStep);
92: TSSolve(ts,x);
94: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"fin_x",FILE_MODE_WRITE,&viewer);
95: VecView(x,viewer);
96: PetscViewerDestroy(&viewer);
98: VecDestroy(&x);
99: MatDestroy(&J);
100: DMDestroy(&user.da);
101: TSDestroy(&ts);
102: PetscFinalize();
103: return 0;
104: }
106: PetscErrorCode PostStep(TS ts)
107: {
108: Vec X;
109: AppCtx *user;
110: PetscScalar sum;
111: PetscReal t;
113: TSGetApplicationContext(ts,&user);
114: TSGetTime(ts,&t);
115: TSGetSolution(ts,&X);
116: VecSum(X,&sum);
117: PetscPrintf(PETSC_COMM_WORLD,"sum(p)*dw*dtheta at t = %3.2f = %3.6f\n",(double)t,(double)sum*user->dx*user->dy);
118: return 0;
119: }
121: PetscErrorCode ini_bou(Vec X,AppCtx* user)
122: {
123: DM cda;
124: DMDACoor2d **coors;
125: PetscScalar **p;
126: Vec gc;
127: PetscInt i,j;
128: PetscInt xs,ys,xm,ym,M,N;
129: PetscScalar xi,yi;
130: PetscScalar sigmax=user->sigmax,sigmay=user->sigmay;
131: PetscScalar rho =user->rho;
132: PetscScalar mux =user->mux,muy=user->muy;
133: PetscMPIInt rank;
136: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
137: DMDAGetInfo(user->da,NULL,&M,&N,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
138: user->dx = (user->xmax - user->xmin)/(M-1); user->dy = (user->ymax - user->ymin)/(N-1);
139: DMGetCoordinateDM(user->da,&cda);
140: DMGetCoordinates(user->da,&gc);
141: DMDAVecGetArray(cda,gc,&coors);
142: DMDAVecGetArray(user->da,X,&p);
143: DMDAGetCorners(cda,&xs,&ys,0,&xm,&ym,0);
144: for (i=xs; i < xs+xm; i++) {
145: for (j=ys; j < ys+ym; j++) {
146: xi = coors[j][i].x; yi = coors[j][i].y;
147: if (i == 0 || j == 0 || i == M-1 || j == N-1) p[j][i] = 0.0;
148: else 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)));
149: }
150: }
151: /* p[N/2+N%2][M/2+M%2] = 1/(user->dx*user->dy); */
153: DMDAVecRestoreArray(cda,gc,&coors);
154: DMDAVecRestoreArray(user->da,X,&p);
155: return 0;
156: }
158: /* First advection term */
159: PetscErrorCode adv1(PetscScalar **p,PetscScalar y,PetscInt i,PetscInt j,PetscInt M,PetscScalar *p1,AppCtx *user)
160: {
161: PetscScalar f;
162: /* PetscScalar v1,v2,v3,v4,v5,s1,s2,s3; */
163: /* if (i > 2 && i < M-2) {
164: v1 = (y-user->ws)*(p[j][i-2] - p[j][i-3])/user->dx;
165: v2 = (y-user->ws)*(p[j][i-1] - p[j][i-2])/user->dx;
166: v3 = (y-user->ws)*(p[j][i] - p[j][i-1])/user->dx;
167: v4 = (y-user->ws)*(p[j][i+1] - p[j][i])/user->dx;
168: v5 = (y-user->ws)*(p[j][i+1] - p[j][i+2])/user->dx;
170: s1 = v1/3.0 - (7.0/6.0)*v2 + (11.0/6.0)*v3;
171: s2 =-v2/6.0 + (5.0/6.0)*v3 + (1.0/3.0)*v4;
172: s3 = v3/3.0 + (5.0/6.0)*v4 - (1.0/6.0)*v5;
174: *p1 = 0.1*s1 + 0.6*s2 + 0.3*s3;
175: } else *p1 = 0.0; */
176: f = (y - user->ws);
177: *p1 = f*(p[j][i+1] - p[j][i-1])/(2*user->dx);
178: return 0;
179: }
181: /* Second advection term */
182: PetscErrorCode adv2(PetscScalar **p,PetscScalar x,PetscInt i,PetscInt j,PetscInt N,PetscScalar *p2,AppCtx *user)
183: {
184: PetscScalar f;
185: /* PetscScalar v1,v2,v3,v4,v5,s1,s2,s3; */
186: /* if (j > 2 && j < N-2) {
187: v1 = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*sin(x))*(p[j-2][i] - p[j-3][i])/user->dy;
188: v2 = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*sin(x))*(p[j-1][i] - p[j-2][i])/user->dy;
189: v3 = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*sin(x))*(p[j][i] - p[j-1][i])/user->dy;
190: v4 = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*sin(x))*(p[j+1][i] - p[j][i])/user->dy;
191: v5 = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*sin(x))*(p[j+2][i] - p[j+1][i])/user->dy;
193: s1 = v1/3.0 - (7.0/6.0)*v2 + (11.0/6.0)*v3;
194: s2 =-v2/6.0 + (5.0/6.0)*v3 + (1.0/3.0)*v4;
195: s3 = v3/3.0 + (5.0/6.0)*v4 - (1.0/6.0)*v5;
197: *p2 = 0.1*s1 + 0.6*s2 + 0.3*s3;
198: } else *p2 = 0.0; */
199: f = (user->ws/(2*user->H))*(user->PM_min - user->Pmax*PetscSinScalar(x));
200: *p2 = f*(p[j+1][i] - p[j-1][i])/(2*user->dy);
201: return 0;
202: }
204: /* Diffusion term */
205: PetscErrorCode diffuse(PetscScalar **p,PetscInt i,PetscInt j,PetscReal t,PetscScalar *p_diff,AppCtx * user)
206: {
209: *p_diff = user->disper_coe*((p[j-1][i] - 2*p[j][i] + p[j+1][i])/(user->dy*user->dy));
210: return 0;
211: }
213: PetscErrorCode BoundaryConditions(PetscScalar **p,DMDACoor2d **coors,PetscInt i,PetscInt j,PetscInt M, PetscInt N,PetscScalar **f,AppCtx *user)
214: {
215: PetscScalar fwc,fthetac;
216: PetscScalar w=coors[j][i].y,theta=coors[j][i].x;
219: if (user->bc == 0) { /* Natural boundary condition */
220: f[j][i] = p[j][i];
221: } else { /* Steady state boundary condition */
222: fthetac = user->ws/(2*user->H)*(user->PM_min - user->Pmax*PetscSinScalar(theta));
223: fwc = (w*w/2.0 - user->ws*w);
224: if (i == 0 && j == 0) { /* left bottom corner */
225: f[j][i] = fwc*(p[j][i+1] - p[j][i])/user->dx + fthetac*p[j][i] - user->disper_coe*(p[j+1][i] - p[j][i])/user->dy;
226: } else if (i == 0 && j == N-1) { /* right bottom corner */
227: f[j][i] = fwc*(p[j][i+1] - p[j][i])/user->dx + fthetac*p[j][i] - user->disper_coe*(p[j][i] - p[j-1][i])/user->dy;
228: } else if (i == M-1 && j == 0) { /* left top corner */
229: f[j][i] = fwc*(p[j][i] - p[j][i-1])/user->dx + fthetac*p[j][i] - user->disper_coe*(p[j+1][i] - p[j][i])/user->dy;
230: } else if (i == M-1 && j == N-1) { /* right top corner */
231: f[j][i] = fwc*(p[j][i] - p[j][i-1])/user->dx + fthetac*p[j][i] - user->disper_coe*(p[j][i] - p[j-1][i])/user->dy;
232: } else if (i == 0) { /* Bottom edge */
233: f[j][i] = fwc*(p[j][i+1] - p[j][i])/(user->dx) + fthetac*p[j][i] - user->disper_coe*(p[j+1][i] - p[j-1][i])/(2*user->dy);
234: } else if (i == M-1) { /* Top edge */
235: f[j][i] = fwc*(p[j][i] - p[j][i-1])/(user->dx) + fthetac*p[j][i] - user->disper_coe*(p[j+1][i] - p[j-1][i])/(2*user->dy);
236: } else if (j == 0) { /* Left edge */
237: f[j][i] = fwc*(p[j][i+1] - p[j][i-1])/(2*user->dx) + fthetac*p[j][i] - user->disper_coe*(p[j+1][i] - p[j][i])/(user->dy);
238: } else if (j == N-1) { /* Right edge */
239: f[j][i] = fwc*(p[j][i+1] - p[j][i-1])/(2*user->dx) + fthetac*p[j][i] - user->disper_coe*(p[j][i] - p[j-1][i])/(user->dy);
240: }
241: }
242: return 0;
243: }
245: PetscErrorCode IFunction(TS ts,PetscReal t,Vec X,Vec Xdot,Vec F,void *ctx)
246: {
247: AppCtx *user=(AppCtx*)ctx;
248: DM cda;
249: DMDACoor2d **coors;
250: PetscScalar **p,**f,**pdot;
251: PetscInt i,j;
252: PetscInt xs,ys,xm,ym,M,N;
253: Vec localX,gc,localXdot;
254: PetscScalar p_adv1,p_adv2,p_diff;
257: DMDAGetInfo(user->da,NULL,&M,&N,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
258: DMGetCoordinateDM(user->da,&cda);
259: DMDAGetCorners(cda,&xs,&ys,0,&xm,&ym,0);
261: DMGetLocalVector(user->da,&localX);
262: DMGetLocalVector(user->da,&localXdot);
264: DMGlobalToLocalBegin(user->da,X,INSERT_VALUES,localX);
265: DMGlobalToLocalEnd(user->da,X,INSERT_VALUES,localX);
266: DMGlobalToLocalBegin(user->da,Xdot,INSERT_VALUES,localXdot);
267: DMGlobalToLocalEnd(user->da,Xdot,INSERT_VALUES,localXdot);
269: DMGetCoordinatesLocal(user->da,&gc);
271: DMDAVecGetArrayRead(cda,gc,&coors);
272: DMDAVecGetArrayRead(user->da,localX,&p);
273: DMDAVecGetArrayRead(user->da,localXdot,&pdot);
274: DMDAVecGetArray(user->da,F,&f);
276: user->disper_coe = PetscPowScalar((user->lambda*user->ws)/(2*user->H),2)*user->q*(1.0-PetscExpScalar(-t/user->lambda));
277: for (i=xs; i < xs+xm; i++) {
278: for (j=ys; j < ys+ym; j++) {
279: if (i == 0 || j == 0 || i == M-1 || j == N-1) {
280: BoundaryConditions(p,coors,i,j,M,N,f,user);
281: } else {
282: adv1(p,coors[j][i].y,i,j,M,&p_adv1,user);
283: adv2(p,coors[j][i].x,i,j,N,&p_adv2,user);
284: diffuse(p,i,j,t,&p_diff,user);
285: f[j][i] = -p_adv1 - p_adv2 + p_diff - pdot[j][i];
286: }
287: }
288: }
289: DMDAVecRestoreArrayRead(user->da,localX,&p);
290: DMDAVecRestoreArrayRead(user->da,localX,&pdot);
291: DMRestoreLocalVector(user->da,&localX);
292: DMRestoreLocalVector(user->da,&localXdot);
293: DMDAVecRestoreArray(user->da,F,&f);
294: DMDAVecRestoreArrayRead(cda,gc,&coors);
296: return 0;
297: }
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;
312: DMDAGetInfo(user->da,NULL,&M,&N,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
313: DMGetCoordinateDM(user->da,&cda);
314: DMDAGetCorners(cda,&xs,&ys,0,&xm,&ym,0);
316: DMGetCoordinatesLocal(user->da,&gc);
317: 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; yi = coors[j][i].y;
322: row.i = i; row.j = j;
323: if (i == 0 || j == 0 || i == M-1 || j == N-1) {
324: if (user->bc == 0) {
325: col[nc].i = i; col[nc].j = j; val[nc++] = 1.0;
326: } else {
327: PetscScalar fthetac,fwc;
328: fthetac = user->ws/(2*user->H)*(user->PM_min - user->Pmax*PetscSinScalar(xi));
329: fwc = (yi*yi/2.0 - user->ws*yi);
330: if (i==0 && j==0) {
331: col[nc].i = i+1; col[nc].j = j; val[nc++] = fwc/user->dx;
332: col[nc].i = i; col[nc].j = j+1; val[nc++] = -user->disper_coe/user->dy;
333: col[nc].i = i; col[nc].j = j; val[nc++] = -fwc/user->dx + fthetac + user->disper_coe/user->dy;
334: } else if (i==0 && j == N-1) {
335: col[nc].i = i+1; col[nc].j = j; val[nc++] = fwc/user->dx;
336: col[nc].i = i; col[nc].j = j-1; val[nc++] = user->disper_coe/user->dy;
337: col[nc].i = i; col[nc].j = j; val[nc++] = -fwc/user->dx + fthetac - user->disper_coe/user->dy;
338: } else if (i== M-1 && j == 0) {
339: col[nc].i = i-1; col[nc].j = j; val[nc++] = -fwc/user->dx;
340: col[nc].i = i; col[nc].j = j+1; val[nc++] = -user->disper_coe/user->dy;
341: col[nc].i = i; col[nc].j = j; val[nc++] = fwc/user->dx + fthetac + user->disper_coe/user->dy;
342: } else if (i == M-1 && j == N-1) {
343: col[nc].i = i-1; col[nc].j = j; val[nc++] = -fwc/user->dx;
344: col[nc].i = i; col[nc].j = j-1; val[nc++] = user->disper_coe/user->dy;
345: col[nc].i = i; col[nc].j = j; val[nc++] = fwc/user->dx + fthetac - user->disper_coe/user->dy;
346: } else if (i==0) {
347: col[nc].i = i+1; col[nc].j = j; val[nc++] = fwc/user->dx;
348: col[nc].i = i; col[nc].j = j+1; val[nc++] = -user->disper_coe/(2*user->dy);
349: col[nc].i = i; col[nc].j = j-1; val[nc++] = user->disper_coe/(2*user->dy);
350: col[nc].i = i; col[nc].j = j; val[nc++] = -fwc/user->dx + fthetac;
351: } else if (i == M-1) {
352: col[nc].i = i-1; col[nc].j = j; val[nc++] = -fwc/user->dx;
353: col[nc].i = i; col[nc].j = j+1; val[nc++] = -user->disper_coe/(2*user->dy);
354: col[nc].i = i; col[nc].j = j-1; val[nc++] = user->disper_coe/(2*user->dy);
355: col[nc].i = i; col[nc].j = j; val[nc++] = fwc/user->dx + fthetac;
356: } else if (j==0) {
357: col[nc].i = i+1; col[nc].j = j; val[nc++] = fwc/(2*user->dx);
358: col[nc].i = i-1; col[nc].j = j; val[nc++] = -fwc/(2*user->dx);
359: col[nc].i = i; col[nc].j = j+1; val[nc++] = -user->disper_coe/user->dy;
360: col[nc].i = i; col[nc].j = j; val[nc++] = user->disper_coe/user->dy + fthetac;
361: } else if (j == N-1) {
362: col[nc].i = i+1; col[nc].j = j; val[nc++] = fwc/(2*user->dx);
363: col[nc].i = i-1; col[nc].j = j; val[nc++] = -fwc/(2*user->dx);
364: col[nc].i = i; col[nc].j = j-1; val[nc++] = user->disper_coe/user->dy;
365: col[nc].i = i; col[nc].j = j; val[nc++] = -user->disper_coe/user->dy + fthetac;
366: }
367: }
368: } else {
369: c1 = (yi-user->ws)/(2*user->dx);
370: c3 = (user->ws/(2.0*user->H))*(user->PM_min - user->Pmax*PetscSinScalar(xi))/(2*user->dy);
371: c5 = (PetscPowScalar((user->lambda*user->ws)/(2*user->H),2)*user->q*(1.0-PetscExpScalar(-t/user->lambda)))/(user->dy*user->dy);
372: col[nc].i = i-1; col[nc].j = j; val[nc++] = c1;
373: col[nc].i = i+1; col[nc].j = j; val[nc++] = -c1;
374: col[nc].i = i; col[nc].j = j-1; val[nc++] = c3 + c5;
375: col[nc].i = i; col[nc].j = j+1; val[nc++] = -c3 + c5;
376: col[nc].i = i; col[nc].j = j; val[nc++] = -2*c5 -a;
377: }
378: MatSetValuesStencil(Jpre,1,&row,nc,col,val,INSERT_VALUES);
379: }
380: }
381: DMDAVecRestoreArrayRead(cda,gc,&coors);
383: MatAssemblyBegin(Jpre,MAT_FINAL_ASSEMBLY);
384: MatAssemblyEnd(Jpre,MAT_FINAL_ASSEMBLY);
385: if (J != Jpre) {
386: MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);
387: MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);
388: }
389: return 0;
390: }
392: PetscErrorCode Parameter_settings(AppCtx *user)
393: {
394: PetscBool flg;
398: /* Set default parameters */
399: user->ws = 1.0;
400: user->H = 5.0; user->Pmax = 2.1;
401: user->PM_min = 1.0; user->lambda = 0.1;
402: user->q = 1.0; user->mux = PetscAsinScalar(user->PM_min/user->Pmax);
403: user->sigmax = 0.1;
404: user->sigmay = 0.1; user->rho = 0.0;
405: user->t0 = 0.0; user->tmax = 2.0;
406: user->xmin = -1.0; user->xmax = 10.0;
407: user->ymin = -1.0; user->ymax = 10.0;
408: user->bc = 0;
410: PetscOptionsGetScalar(NULL,NULL,"-ws",&user->ws,&flg);
411: PetscOptionsGetScalar(NULL,NULL,"-Inertia",&user->H,&flg);
412: PetscOptionsGetScalar(NULL,NULL,"-Pmax",&user->Pmax,&flg);
413: PetscOptionsGetScalar(NULL,NULL,"-PM_min",&user->PM_min,&flg);
414: PetscOptionsGetScalar(NULL,NULL,"-lambda",&user->lambda,&flg);
415: PetscOptionsGetScalar(NULL,NULL,"-q",&user->q,&flg);
416: PetscOptionsGetScalar(NULL,NULL,"-mux",&user->mux,&flg);
417: PetscOptionsGetScalar(NULL,NULL,"-sigmax",&user->sigmax,&flg);
418: PetscOptionsGetScalar(NULL,NULL,"-muy",&user->muy,&flg);
419: PetscOptionsGetScalar(NULL,NULL,"-sigmay",&user->sigmay,&flg);
420: PetscOptionsGetScalar(NULL,NULL,"-rho",&user->rho,&flg);
421: PetscOptionsGetScalar(NULL,NULL,"-t0",&user->t0,&flg);
422: PetscOptionsGetScalar(NULL,NULL,"-tmax",&user->tmax,&flg);
423: PetscOptionsGetScalar(NULL,NULL,"-xmin",&user->xmin,&flg);
424: PetscOptionsGetScalar(NULL,NULL,"-xmax",&user->xmax,&flg);
425: PetscOptionsGetScalar(NULL,NULL,"-ymin",&user->ymin,&flg);
426: PetscOptionsGetScalar(NULL,NULL,"-ymax",&user->ymax,&flg);
427: PetscOptionsGetInt(NULL,NULL,"-bc_type",&user->bc,&flg);
428: user->muy = user->ws;
429: return 0;
430: }
432: /*TEST
434: build:
435: requires: !complex
437: test:
438: args: -nox -ts_max_steps 2
439: localrunfiles: petscopt_ex6
440: timeoutfactor: 4
442: TEST*/