<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Nov 9, 2016 at 12:32 AM, Julian Andrej <span dir="ltr"><<a href="mailto:juan@tf.uni-kiel.de" target="_blank">juan@tf.uni-kiel.de</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">I tested a few other things and it turned out that the function added<br>
as a dirichlet boundary condition via PetscDSAddBoundary doesn't<br>
receive the time values correctly (PetscReal time is always 0.0).<br></blockquote><div><br></div><div>I will check out your example and fix this. However, I wonder if there is a bug in your implementation:</div><div><br></div><div><div>void f0_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,</div><div> const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],</div><div> const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],</div><div> PetscReal t, const PetscReal x[], PetscScalar f0[])</div><div>{</div><div> f0[0] = u_t[0] - 2.0;</div><div>}</div></div><div><br></div><div>I get</div><div><br></div><div> du/dt - \Delta u = -2</div><div><br></div><div>so that if u = x^2 + y^2 + 2t, then</div><div><br></div><div> 2 - (2 + 2) = -2</div><div><br></div><div>so you would want u_t + 2.</div><div><br></div><div>If you give the Laplacian the opposite sign, this is an unstable formulation.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
I also see no effect adding DMTSSetIJacobianLocal vs. not setting the<br>
jacobian function. The analytic solution i used in my attached example<br>
is the same as in the src/ts/examples/tutorials/<wbr>ex32.c.<br>
<br>
Am i doing anything obviously wrong here? My next step would be to try<br>
if assembling the matrices myself and writing a custom<br>
IFunction/IJacobian which assembles the different parts of the matrix<br>
like M and J for the stiff ODE with nontrivial mass matrix (see manual<br>
section of TS -> F = M u' - f) but i think this should be obsolete<br>
right?<br>
<br>
Appreciate your help.<br>
<br>
Regards<br>
<span class="gmail-HOEnZb"><font color="#888888">Julian<br>
</font></span><div class="gmail-HOEnZb"><div class="gmail-h5"><br>
On Mon, Nov 7, 2016 at 10:49 AM, Julian Andrej <<a href="mailto:juan@tf.uni-kiel.de">juan@tf.uni-kiel.de</a>> wrote:<br>
> Hello,<br>
><br>
> i'm using PetscFE in combination with SNES and a hand written backward<br>
> euler solver right now and like to hand that over to TS. I have been<br>
> fiddling quite a while getting TS to work with PetscFE and have<br>
> encountered a few unclarities.<br>
><br>
> I have looked at examples which are outdated i guess<br>
> (src/ts/examples/tutorials/<wbr>ex32.c) and i am confused on how i have to<br>
> formulate the discretization of the jacobian and the residual. It<br>
> seems obvious to use the DMTSSetIFunctionLocal and<br>
> DMTSSetIJacobianLocal functions because there are prepared wrappers<br>
> from DMPlexTSComputeIFunctionFEM and DMPlexTSComputeIJacobianFEM.<br>
><br>
> I need a Mass matrix for my problem, is that formed from the function<br>
> space information or do i have to form it myself? Is there any working<br>
> example which uses PetscFE and TS to form the DMTSSetIFunctionLocal<br>
> and DMTSSetIJacobianLocal? (grepping the tree tells me "no"<br>
> unfortunately).<br>
><br>
> Regards<br>
> Julian<br>
</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>