<html><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; "><div><br></div> Any ideas?<div><br><div><br><div>Begin forwarded message:</div><br class="Apple-interchange-newline"><blockquote type="cite"><div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px; "><font face="Helvetica" size="3" color="#000000" style="font: 12.0px Helvetica; color: #000000"><b>From: </b></font><font face="Helvetica" size="3" style="font: 12.0px Helvetica">Dennis Zeki Kaya <<a href="mailto:Dennis.Kaya.3319@student.uu.se">Dennis.Kaya.3319@student.uu.se</a>></font></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px; "><font face="Helvetica" size="3" color="#000000" style="font: 12.0px Helvetica; color: #000000"><b>Date: </b></font><font face="Helvetica" size="3" style="font: 12.0px Helvetica">December 4, 2009 5:06:35 AM CST</font></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px; "><font face="Helvetica" size="3" color="#000000" style="font: 12.0px Helvetica; color: #000000"><b>To: </b></font><font face="Helvetica" size="3" style="font: 12.0px Helvetica">petsc-maint Maintenance <<a href="mailto:petsc-maint@mcs.anl.gov">petsc-maint@mcs.anl.gov</a>></font></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px; "><font face="Helvetica" size="3" color="#000000" style="font: 12.0px Helvetica; color: #000000"><b>Subject: </b></font><font face="Helvetica" size="3" style="font: 12.0px Helvetica"><b>Re: [petsc-maint #38493] TSSetIJacobian()</b></font></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px; "><font face="Helvetica" size="3" color="#000000" style="font: 12.0px Helvetica; color: #000000"><b>Reply-To: </b></font><font face="Helvetica" size="3" style="font: 12.0px Helvetica"><a href="mailto:petsc-maint@mcs.anl.gov">petsc-maint@mcs.anl.gov</a>, Dennis Zeki Kaya <<a href="mailto:Dennis.Kaya.3319@student.uu.se">Dennis.Kaya.3319@student.uu.se</a>></font></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px; min-height: 14px; "><br></div> </div>Hello,<br><br>I'm using the TSTHETA solver to solve my DAE problem since I have a <br>nonsingular Jacobian. I'm following the procedure in ex8.c in <br>petsc-dev/src/ts/examples/tutorials as close as I can. I have attached <br>a outline of my Petsc_Solve() routine to this email. The only big <br>difference in my implementation to that used in ex8.c is that I create <br>a global Jacobian matrix with my routine Petsc_JacCreate().<br><br><br><br>Quoting Barry Smith <<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>>:<br><br><blockquote type="cite"><br></blockquote><blockquote type="cite"> Did you set the type to Backward Euler? The default Euler doesn't<br></blockquote><blockquote type="cite">use Jacobians.<br></blockquote><blockquote type="cite"><br></blockquote><blockquote type="cite"> Barry<br></blockquote><blockquote type="cite"><br></blockquote><blockquote type="cite">On Dec 3, 2009, at 6:28 PM, Dennis Zeki Kaya wrote:<br></blockquote><blockquote type="cite"><br></blockquote><blockquote type="cite"><blockquote type="cite">Hello,<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">I have a problem with a routine that evaluates a Jacobian matrix. The<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">problem is that PETSc never enters the RHSJacobian() routine. I know<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">that PETSc enters the RHSFunction() rotine threw a print statement.<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">The output is something like:<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">At t[0] = 0.00e+00 u= 0.00e+00 at the center<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">...<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">At t[10] = 1.00e-02 u= 0.00e+00 at the center<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">RHSFunction!<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">...<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">I have taken bits and pieces for monitoring the solution procedure<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">from the examples provided. Since the Jacobian is never computed the<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">solution stays zero at each time-step, indifferent of parameters like<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">the size of time-step etc.<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">This is the calls to the setting of the functions:<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">/* Set user-provided RHSFunction and RHSJacobian */<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">ierr = TSSetIFunction(ts, RHSFunction, ctx);CHKERRQ(ierr);<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">ierr = TSSetIJacobian(ts, J, J, RHSJacobian, ctx);CHKERRQ(ierr);<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">I can attach the whole Petsc_Solve() routine with the RHSFunction()<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">and the RHSJacobian() routines if needed. I just thought that there<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">maybe is an easy answer to my question regarding the time-stepper not<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">entering the RHSJacobian() routine.<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">Thank you in advance,<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite">Dennis Kaya<br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><blockquote type="cite"><blockquote type="cite"><br></blockquote></blockquote><br><br><br>PetscErrorCode Petsc_Solve(/* Parameters for setting the DAE userdefined context */) {<br><br> PetscErrorCode ierr;<br> TS ts; /* nonlinear solver */<br> PetscInt time_steps=100, steps, iout, NOUT=1;<br> PetscReal ftime;<br> Mat J; /* Jacobian matrix */<br> Vec u; /* solution */<br> PC pc; /* preconditioner */<br> DataCtx *ctx; /* user-defined DAE-data context */<br> PetscReal dt;<br> PetscTruth flg;<br> const TSType tstype;<br> PetscViewer viewer;<br> PetscViewer<span class="Apple-tab-span" style="white-space:pre"> </span> viewfile;<br> char pcinfo[120], tsinfo[120];<br><br> SNES ts_snes;<br> KSP ksp;<br><br> PetscFunctionBegin;<br><br> /* Create Jacobian matrix */<br> ierr = Petsc_JacCreate(&J, ctx);CHKERRQ(ierr);<br><br> /* Set the solution vector to have the same parallel layout as the Jacobian matrix. */<br> ierr = MatGetVecs(J, &u, PETSC_NULL);CHKERRQ(ierr);<br><br> /* Create time-step context */<br> ierr = TSCreate(PETSC_COMM_WORLD, &ts);CHKERRQ(ierr);<br> ierr = TSSetProblemType(ts, TS_NONLINEAR);CHKERRQ(ierr);<br> ierr = TSSetType(ts, TSTHETA);CHKERRQ(ierr);<br><br> /* Set user-provided RHSFunction and RHSJacobian */<br> ierr = TSSetIFunction(ts, RHSFunction, ctx);CHKERRQ(ierr);<br> ierr = TSSetIJacobian(ts, J, J, RHSJacobian, ctx);CHKERRQ(ierr);<br> ierr = TSSetDuration(ts, time_steps, 0.3);CHKERRQ(ierr);<br><br> ierr = TSMonitorSet(ts, Monitor, ctx, PETSC_NULL);CHKERRQ(ierr); <br><br> /* Set up the preconditioner type to be used */<br> ierr = TSGetSNES(ts, &ts_snes);CHKERRQ(ierr);<br> ierr = SNESGetKSP(ts_snes, &ksp);CHKERRQ(ierr);<br> ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);<br> /* Set preconditioner to Jacobi, i.e diagonal scaling preconditioning. */<br> ierr = PCSetType(pc, PCJACOBI);CHKERRQ(ierr);<br><br> dt = 0.001;<br> /* Set the initial solution vector to zero. */<br> ierr = Initial(0.0, u, ctx);CHKERRQ(ierr);<br> ierr = TSSetInitialTimeStep(ts, 0.0, dt);CHKERRQ(ierr); <br> ierr = TSSetSolution(ts, u);CHKERRQ(ierr);<br><br> /* Test TSSetPostStep() */<br> ierr = PetscOptionsHasName(PETSC_NULL,"-test_PostStep",&flg);CHKERRQ(ierr);<br> if (flg) {<br> ierr = TSSetPostStep(ts,PostStep);CHKERRQ(ierr);<br> }<br><br> /* Set up the time-stepping context from the the options given */<br> ierr = TSSetFromOptions(ts);CHKERRQ(ierr);<br><br> ierr = PetscOptionsGetInt(PETSC_NULL,"-NOUT",&NOUT,PETSC_NULL);CHKERRQ(ierr);<br> for (iout=1; iout<=NOUT; iout++){<br> ierr = TSSetDuration(ts,time_steps,iout*1.0/NOUT);CHKERRQ(ierr);<br> ierr = TSStep(ts,&steps,&ftime);CHKERRQ(ierr);<br> ierr = TSSetInitialTimeStep(ts,ftime,dt);CHKERRQ(ierr);<br> }<br><br> ierr = PetscPrintf(PETSC_COMM_WORLD,"steps %D, ftime %G\n",steps,ftime);CHKERRQ(ierr);<br><br> ierr = TSGetType(ts, &tstype);CHKERRQ(ierr);<br> ierr = PetscViewerStringOpen(PETSC_COMM_WORLD, tsinfo, 120, &viewer);CHKERRQ(ierr);<br> ierr = TSView(ts, viewer);CHKERRQ(ierr);<br> ierr = PetscViewerDestroy(viewer);CHKERRQ(ierr);<br> ierr = PetscViewerStringOpen(PETSC_COMM_WORLD, pcinfo, 120, &viewer);CHKERRQ(ierr);<br> ierr = PCView(pc, viewer);CHKERRQ(ierr);<br> ierr = PetscPrintf(PETSC_COMM_WORLD, "%d Procs,%s TSType, %s Preconditioner\n",<br><span class="Apple-tab-span" style="white-space:pre"> </span><span class="Apple-tab-span" style="white-space:pre"> </span> size, tsinfo, pcinfo);CHKERRQ(ierr);<br> ierr = PetscViewerDestroy(viewer);CHKERRQ(ierr);<br><br> ierr = PetscOptionsHasName(PETSC_NULL,"-matlab_view",&flg);CHKERRQ(ierr);<br> if (flg) { /* print solution into a MATLAB file */<br> ierr = TSGetSolution(ts, &u);CHKERRQ(ierr);<br> ierr = PetscViewerASCIIOpen(PETSC_COMM_WORLD,"out.m", &viewfile);CHKERRQ(ierr);<br> ierr = PetscViewerSetFormat(viewfile, PETSC_VIEWER_ASCII_MATLAB);CHKERRQ(ierr);<br> ierr = VecView(u, viewfile);CHKERRQ(ierr);<br> ierr = PetscViewerDestroy(viewfile);CHKERRQ(ierr);<br> }<br><br> /* Free memory used by this routine. */<br> ierr = TSDestroy(ts);CHKERRQ(ierr);<br> ierr = VecDestroy(u);CHKERRQ(ierr);<br> ierr = PetscFree(ctx);CHKERRQ(ierr);<br><br> PetscFunctionReturn(0);<br>}/*.Petsc_Solve.*/<br></blockquote></div><br></div></body></html>