<span id="mailbox-conversation"><div>I’m sorry, I’m not trying to be difficult, but I’m not following.</div>
<div><br></div>
<div>The manual states (for my special case):</div>
<div><div>
<div>
<div>
<ul><li>
<div id="mb-reply">u ̇ = A(t)u. Use
</div>
<div id="mb-reply"><br></div>
<div id="mb-reply">TSSetProblemType(ts,TS LINEAR);
TSSetRHSFunction(ts,NULL,TSComputeRHSFunctionLinear,NULL);
TSSetRHSJacobian(ts,A,A,YourComputeRHSJacobian,&appctx);
</div>
<div id="mb-reply"><br></div>
<div>where YourComputeRHSJacobian() is a function you provide that computes A as a function of
time. Or use ...</div>
</li>
</ul></div>
</div>
</div></div>
<div>My `func` does this. It is 7 lines:</div>
<div><br></div>
<div>
<div> context* c = static_cast<context*>( G_u );</div>
<div> PetscScalar e = c->E( t_ );</div>
<div> MatCopy( c->D, A, SAME_NONZERO_PATTERN );</div>
<div> MatShift( A, e );</div>
<div> MatDiagonalSet( A, c->H0, INSERT_VALUES);</div>
<div> MatShift( A, std::complex<double>( 0, -1 ) );</div>
<div> return 0;</div>
</div>
<div><br></div>
<div>SHOULD `func` touch U? If so, what should `func` do to U? I thought that the RHSJacobian function was only meant to create A, since dG/du = A(t) (for this special case).</div>
<div><br></div>
<div>-Andrew</div></span><div class="mailbox_signature"><br></div>
<br><br><div class="gmail_quote"><p>On Fri, Mar 20, 2015 at 3:26 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br></p><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;"><div><div dir="ltr"><div class="gmail_extra">
<div class="gmail_quote">On Fri, Mar 20, 2015 at 3:09 PM, Andrew Spott <span dir="ltr"><<a href="mailto:ansp6066@colorado.edu">ansp6066@colorado.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<u></u>
<div>
<span><div>So, it doesn’t seem that zeroing the given vector in the function passed to TSSetRHSJacobian is the problem. When I do that, it just zeros out the solution.</div></span>
</div>
</blockquote>
<div><br></div>
<div>I would think you would zero the residual vector (if you add to it to construct the residual, as in FEM methods), not the solution.</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><span>
<div>The function that is passed to TSSetRHSJacobian has only one responsibility — to create the jacobian — correct? In my case this is A(t). The solution vector is given for when you are solving nonlinear problems (A(t) also depends on U(t)). In my case, I don’t even look at the solution vector (because my A(t) doesn’t depend on it).</div></span></div></blockquote>
<div><br></div>
<div>Are you initializing the Jacobian to 0 first?</div>
<div><br></div>
<div> Thanks,</div>
<div><br></div>
<div> Matt</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div>
<span>
<div>Is this the case? or is there some other responsibility of said function?</div>
<span class="HOEnZb"><font color="#888888">
<div><br></div>
<div>-Andrew</div></font></span><span class="">
<div><br></div>
<div>>Ah ha!</div>
<div>></div>
<div>>The function passed to TSSetRHSJacobian needs to zero the solution vector?</div>
<div>></div>
<div>>As a point, this isn’t mentioned in any documentation that I can find.</div>
<div>></div>
<div>>-Andrew</div>
<br></span><div><div class="h5"><span style="display:inline">On Friday, Mar 20, 2015 at 2:17 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>>, wrote:<br></span></div></div></span><div><div class="h5">
<span><blockquote class="gmail_quote"><div>
<p dir="ltr">This sounds like a problem in your calculation function where a Vec or Mat does not get reset to 0, but it does in your by hand code.</p>
<p dir="ltr"> Matt</p>
<div class="gmail_quote">On Mar 20, 2015 2:52 PM, "Andrew Spott" <<a href="mailto:ansp6066@colorado.edu">ansp6066@colorado.edu</a>> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<u></u>
<div>
<span><div>I have a fairly simple problem that I’m trying to timestep:</div>
<div><br></div>
<div>u’ = A(t) u</div>
<div><br></div>
<div>I’m using the crank-nicholson method, which I understand (for this problem) to be:</div>
<div><br></div>
<div>u(t + h) = u(t) + h/2[A(t+h)*u(t+h) + A(t)*u(t)]</div>
<div>or</div>
<div>[1 - h/2 * A(t+1)] u(t+1) = [1 + h/2 * A(t)] u(t)</div>
<div><br></div>
<div>When I attempt to timestep using PETSc, the norm of `u` blows up. When I do it directly (using the above), the norm of `u` doesn’t blow up.</div>
<div><br></div>
<div>It is important to note that the solution generated after the first step is identical for both, but the second step for Petsc has a norm of ~2, while for the directly calculated version it is ~1. The third step for petsc has a norm of ~4, while the directly calculated version it is still ~1.</div>
<div><br></div>
<div>I’m not sure what I’m doing wrong.</div>
<div><br></div>
<div>PETSc code is taken out of the manual and is pretty simple:</div>
<div><br></div>
<div>
<div> TSCreate( comm, &ts );</div>
<div> TSSetProblemType( ts, TS_LINEAR);</div>
<div> TSSetType( ts, TSCN );</div>
<div> TSSetInitialTimeStep( ts, 0, 0.01 );</div>
<div> TSSetDuration( ts, 5, 0.03 );</div>
<div> TSSetFromOptions( ts );</div>
<div> TSSetRHSFunction( ts, NULL, TSComputeRHSFunctionLinear, NULL );</div>
<div> TSSetRHSJacobian( ts, A, A, func, &cntx );</div>
<div> TSSolve( ts, psi0 );</div>
<div><br></div>
<div>`func` just constructs A(t) at the time given. The same code for calculating A(t) is used in both calculations, along with the same initial vector psi0, and the same time steps.</div>
<div><br></div>
<div>Let me know what other information is needed. I’m not sure what could be the problem. `func` doesn’t touch U at all (should it?).</div>
<div><br></div>
<div>-Andrew</div>
</div></span><div><br></div>
</div>
</blockquote>
</div>
</div></blockquote></span>
</div></div>
</div>
</blockquote>
</div>
<br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div>
</div></div></div></blockquote></div><br>