Mani, I apologize for my slow reply.

Mani Chandra <[email protected]> writes:

> Hi,
>
> Suppose the following is the residual function that TS needs:
>
> void residualFunction(TS ts, PetscScalar t, Vec X, Vec dX_dt, Vec F, void
> *ptr)
>
> and this returns the following finite volume residual at each grid point
>
> dU_dt + gradF + sourceTerms = 0
>
> where dU_dt are the time derivatives of the conserved variables.
>
> The value of dU_dt needs to be computed from the values of the time
> derivatives of the primitive variables given in dX_dt, i.e. dU_dt is an
> analytic function of dX_dt. 

I recommend using this analytic relationship, i.e.,

  dU_dt = dU/dX dX_dt

This can be done by algorithmic differentiation or differencing, if you
prefer.

> But is it possible to compute dU_dt numerically using the vector X and
> it's value at a previous time step like the following?
>
> dU_dt = (U(X) - U(X_old) )/dt -- (1)
>
> (instead of analytically computing dU_dt which is a function of the vectors
> X and dX_dt which is hard to do if the function is very complicated)
>
> So, is computing dU_dt using (1) permissible using TS? The jacobian will be
> correctly assembled for implicit methods?

No, your equation is a very inaccurate time derivative and will spoil
the accuracy and stability properties of the high-order time
integrators.  You could add U as an algebraic constraint:

  U - eq_state(X) = 0

if you want.  This equation is trivial to solve, but it requires that
you use methods capable of solving DAEs.

Attachment: pgp3YYkds2we4.pgp
Description: PGP signature

Reply via email to