[petsc-users] advection-diffusion pointwise function for FEM

[Justin Chang]

Never mind I think I figured it out. Will let you know if I think something
is wrong. Thanks anyway

On Wed, May 6, 2015 at 8:06 AM, Matthew Knepley <knepley at gmail.com> wrote:

> On Wed, May 6, 2015 at 2:07 AM, Justin Chang <jychang48 at gmail.com> wrote:
>
>> Matt,
>>
>> Yes I would actually use a[] in this case. Back then, I think I was
>> trying to take the derivative of the auxiliary which i had declared to be
>> piece-wise constant.
>>
>> The auxiliary a[] is the gradient of pressure (times permeability over
>> viscosity) obtained from solving the Darcy equation as a Laplacian. This
>> part is solved separately, projected as cell-wise velocity, and inputted as
>> an auxiliary for the advection-diffusion problem.
>>
>> Speaking of which, have you had a chance to look into
>> DMPlexProjectFieldLocal(...) not projecting dirichlet BC constraints? I
>> realize that the "small example" I gave you wasn't small at all, so I can
>> send you a smaller and simpler one if needed.
>>
>
> For some reason, I thought that was taken care of. Please resend. Sorry
> about that.
>
>   Thanks,
>
>      Matt
>
>
>> Thanks,
>> Justin
>>
>> On Tue, May 5, 2015 at 3:34 PM, Matthew Knepley <knepley at gmail.com>
>> wrote:
>>
>>> On Sat, Apr 4, 2015 at 11:54 AM, Justin Chang <jchang27 at uh.edu> wrote:
>>>
>>>> Hi everyone,
>>>>
>>>> I want to include advection into my diffusion FEM code (I am assuming a
>>>> small Peclet number so stability isn't an issue ATM). That is I want to
>>>> incorporate the second term as a pointwise function:
>>>>
>>>> du/dt + v * grad[c] - div[grad[c]] = f
>>>>
>>>> Where v is the velocity (obtained from the auxiliary term a_x[]). For
>>>> the residual, would it be of the following form:
>>>>
>>>
>>> 1) I would think you would use a[] instead. What is your velocity the
>>> gradient of?
>>>
>>>
>>>> f0_u(const PetscScalar u[], const PetscScalar u_t[], const PetscScalar
>>>> u_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar
>>>> a_x[], const PetscReal x[], PetscScalar f0[]) {
>>>>
>>>>   PetscInt d;
>>>>   f[0] = u_t[0];
>>>>   for (d = 0; d < spatialDim; ++d) f0[0] += a_x[d]*u_x[d];
>>>>
>>>> }
>>>>
>>>> What about the jacobian? My guess would be to use g1(...) but what
>>>> would the inside of this function be?
>>>>
>>>
>>> Yes it would be g1. The indices for the output are f,g,dg. I am guessing
>>> that c is a scalar, so f = {0}, g = {0}, dg = {0, 1}
>>> for 2D, so g1 would have two terms,
>>>
>>>   g1[0] = a[0];
>>>   g1[1] = a[1];
>>>
>>>   Thanks,
>>>
>>>     Matt
>>>
>>>
>>>> Thanks,
>>>>
>>>>
>>>> --
>>>> Justin Chang
>>>> PhD Candidate, Civil Engineering - Computational Sciences
>>>> University of Houston, Department of Civil and Environmental Engineering
>>>> Houston, TX 77004
>>>> (512) 963-3262
>>>>
>>>
>>>
>>>
>>> --
>>> What most experimenters take for granted before they begin their
>>> experiments is infinitely more interesting than any results to which their
>>> experiments lead.
>>> -- Norbert Wiener
>>>
>>
>>
>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which their
> experiments lead.
> -- Norbert Wiener
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150506/2d75a652/attachment.html>