dolfin team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2008/8/21  <kent-and@xxxxxxxxx>:
>> 2008/8/20 Anders Logg <logg@xxxxxxxxx>:
>>> On Wed, Aug 20, 2008 at 08:43:08AM -0600, Ostien, Jakob T wrote:
>>>> Hi,
>>>>
>>>> I need to be able to assemble an active set of integration points.
>>>> Essentially, to determine the set I loop over cells and then loop again
>>>> over the integration points in the cell and determine if that
>>>> integration point is active with some criteria.  Then I'd like to be
>>>> able to assemble that set.
>>>>
>>>> This is a problem because currently the element_tensor does not break
>>>> down into integration points, and I need to take derivatives, so the
>>>> QuadratureElement in FFC is also ruled out.  I suppose I could
>>>> calculate derivatives and then project on the QuadratureElement, but
>>>> that seems sort of unclean.
>>>>
>>>> It seems like the recent discussion about integrating at a point (on
>>>> the UFC list) might help me out here.
>>>>
>>>> Any other thoughts on how I might go about this?
>>>>
>>>> Jake
>>>
>>> I don't know yet. I'm thinking about how much of this should go into
>>> the compiler and how much should be left to the user.
>>>
>>> We're working on similar things (and I know Garth is too). I think a
>>> common denominator would at least be to be able to evaluate a form at
>>> a given point.
>>>
>>> I like Kent's suggestion earlier about adding a new UFC class to go
>>> along with the other three integral classes, something like
>>>
>>>  point_integral
>>>
>>> with tabulate_tensor taking the coordinates for a point as additional
>>> input.
>>>
>>> --
>>> Anders
>>
>> I liked your version with adding a new function to each *_integral class
>> better.
>> The reason is that each *_integral object will correspond to one
>> foo*dx(i) or bar*ds(j)
>> in the form definition, and each *_integral object will be able to
>> compute both an
>> integral over some cell and the corresponding integrand.
>>
>> Although a separate point_integral might make sense for something I
>> don't know about,
>> matching them to corresponding cell_integral, exterior_facet_integral
>> and interior_facet_integral
>> objects will become messy.
>>
>> --
>> Martin
>
> One may consider point evaluation also as a form, that was my point.
>
> There is an advantage with adding a new function to the classes in the sense
> that the operators and functions make sense in certain places.
> Eg. jump can be evaluated on points in the interior of a facet whereas it
> should not be used on the interior of a cell (note that the jump function
> does not
> have the same meaning in eg. a vertex). This coincide (of course) with the
> the standard usage of dx, ds, and dS. Although with dx, ds, dS one knows that
> eg jump is not evaluated at vertices because of the choosen quadrature
> points.

I see.

> An advantage of having a separate form object and a separate dx-like
> symbol in
> the form language, say dp or delta, is that one may then define forms that
> are not
> based on the dx, ds, and dS forms. Put in another way, one can eg. define
> functions
> in this way;
>
> f = u*sin(x)*dp
>
> (although dp does not look nice and I don't have any good suggestion).
>
> Kent

Maybe we should add both ways to UFC?
Seems to me they cover two different concepts,
evaluating the integrand of an existing form at a point,
and evaluating a form that is defined only at a point.

// Corresponding to "f = u*sin(x)*dp":
class point_integral
{
    virtual void tabulate_tensor(...) const = 0;
};

// Corresponding to "a = u*v*dx":
class cell_integral
{
    virtual void tabulate_tensor(...) const = 0;
    virtual void tabulate_integrand(...) const = 0;
};
+ same for other integrals.

-- 
Martin