On Wed, 27 Aug 2014 12:13:11 +0200 Kristian Ølgaard <[email protected]> wrote:
> On 27 August 2014 11:36, Garth N. Wells <[email protected]> wrote: > > > > > > > On Wed, 27 Aug, 2014 at 10:21 AM, Kristian Ølgaard > > <[email protected]> wrote: > > > >> On 27 August 2014 10:44, Garth N. Wells <[email protected]> wrote: > >> > >>> > >>> > >>> On Wed, 27 Aug, 2014 at 9:15 AM, Martin Sandve Alnæs > >>> <[email protected]> wrote: > >>> > >>>> On 27 August 2014 08:25, Kristian Ølgaard > >>>> <[email protected]> wrote: > >>>> > >>>>> On 26 August 2014 15:59, Martin Sandve Alnæs > >>>>> <[email protected]> wrote: > >>>>> > >>>>>> On 26 August 2014 15:21, Kristian Ølgaard > >>>>>> <[email protected]> wrote: > >>>>>> > >>>>>>> > >>>>>>> > >>>>>>> ---------- Forwarded message ---------- > >>>>>>> From: Kristian Ølgaard <[email protected]> > >>>>>>> Date: 26 August 2014 15:20 > >>>>>>> Subject: Re: [FEniCS] `Expression`s and their silent > >>>>>>> interpolation To: Jan Blechta <[email protected]> > >>>>>>> > >>>>>>> > >>>>>>> On 26 August 2014 14:18, Jan Blechta > >>>>>>> <[email protected]> wrote: > >>>>>>> > >>>>>>>> On Tue, 26 Aug 2014 09:50:23 +0100 > >>>>>>>> "Garth N. Wells" <[email protected]> wrote: > >>>>>>>> > >>>>>>>> > To summarise this thread, it seems we need to introduce the > >>>>>>>> concept > >>>>>>>> > of an 'Expression' that can be evaluated at arbitrary > >>>>>>>> > points. It should not be a Quadrature{Element/Function} > >>>>>>>> > because the proposed object could be used in different > >>>>>>>> > forms with different evaluation > >>>>>>>> > >>>>>>> > >>>>>> Agree. Also it should not have any notion of degree. Pointwise > >>>>>> is pointwise. > >>>>>> > >>>>>> > >>>>>> > points. The follow-on on issue is then how a 'point-wise' > >>>>>>>> expression > >>>>>>>> > should be treated in forms. We could estimate the > >>>>>>>> > quadrature > >>>>>>>> scheme > >>>>>>>> > when test/trial functions are present, and in the case of > >>>>>>>> functionals > >>>>>>>> > throw an error if the user doesn't supply the quadrature > >>>>>>>> > degree. > >>>>>>>> > >>>>>>>> There's no principal difference regarding rank of the form. > >>>>>>>> Consider > >>>>>>>> > >>>>>>>> f = PointwiseExpression(eval_formula) > >>>>>>>> u, v = TrialFunction(V), TestFunction(V) > >>>>>>>> a = f*u*v*dx > >>>>>>>> L = f*v*dx > >>>>>>>> F = f*dx > >>>>>>>> > >>>>>>>> Still, you need to know what is the polynomial degree of f > >>>>>>>> to have exact quadrature of any of these forms. Ignoring > >>>>>>>> non-zero degree of f > >>>>>>>> (which seems to me you do suggest for a and L) means that > >>>>>>>> you're underintegrating any of those three forms. This is > >>>>>>>> analogical to integrating F with scheme of order zero. I > >>>>>>>> don't see any good reason why having distinct behaviour > >>>>>>>> based on rank of the respective form. > >>>>>>>> > >>>>>>> > >>>>>> Agree. > >>>>>>> For PointwiseExpression, one should define EITHER the > >>>>>>> polynomial degree that the user would like the use for the > >>>>>>> approximation (of e.g., 'sin(x)') OR the (degree of) > >>>>>>> quadrature rule for the measure. The latter should take > >>>>>>> precedence if both are defined, just as it does currently. > >>>>>>> > >>>>>> > >>>>>> Please, no. Isn't that basically the situation we're trying to > >>>>>> get away from? A pointwise expression doesn't have a degree > >>>>>> and it's not a good abstraction to assign one to it. The rules > >>>>>> become complex which makes the source code hard to follow, the > >>>>>> documentation poor, and confuses the users and developers > >>>>>> alike. > >>>>>> > >>>>> > >>>>> Perhaps I got a little confused. Is PointwiseExpression > >>>>> supposed to replace Expression or will they co-exist? > >>>>> > >>>>> If PointwiseExpression will replace Expression, my suggestion to > >>>>> forcing the user to supply the degree (or element) was intended > >>>>> to solve Nico's original problem. > >>>>> > >>>> > >>>> Good question. There are two sides to that: > >>>> > >>>> From the UFL/FFC/UFC/Assembler side there needs to be two > >>>> distinct concepts: If a function is pointwise evaluated FFC must > >>>> generate function calls instead of using basis tables and dofs. > >>>> Whether this is implemented by adding a PointwiseExpression or a > >>>> "Pointwise" family to UFL is an implementation detail that > >>>> affects other work in progress so lets not go there in this > >>>> thread at least. This is basically the new functionality we're > >>>> discussing here, and I don't think anyone disagrees with this, > >>>> apart from annotating the PointwiseExpression with a "virtual > >>>> degree" for integration degree purposes. > >>>> > >>>> > >> Agree, let's discuss this later. > >> > >> From the DOLFIN user interface side, there are several ways to > >> go, and > >>>> this is where most opinions in this thread differ. > >>>> > >>>> The main interface point is whether to introduce new notation > >>>> PointwiseExpression("x[0]") and deprecate Expression("x[0]"), or > >>>> to change Expression("x[0]") to mean pointwise and remove the > >>>> implicit interpolation. > >>>> > >>>> Deprecating Expression("x[0]") breaks all demos and lots of user > >>>> programs. > >>>> > >>>> Changing the behaviour of Expression("x[0]") changes the > >>>> numerical results of lots of programs. > >>>> > >>> > >>> I'd be reluctant to introduce yet another object in the form of > >>> PointwiseExpression. How about passing a string or element > >>> argument to Expression? For now, the default could be the present > >>> behaviour (interpolation) with a warning that the default will > >>> change in release 1.foo to pointwise? > >>> > >> > >> Do we need two arguments to do that? > >> > >> Expression('x[0]*x[0]', pointwise=False, element=None) > >> > > > > We need to allow passing of different types, but it could be just > > one argument and the Expression constructor can do whatever it > > needs to do depending on the type. > > > OK. > > > > > > > > This should then default to the current behaviour and interpolate > > the > >> expression using a linear element. > >> > > > > > > For now, yes. The default isn't a linear element. I think there is > > some magic to pick a suitable order. > > > Could be, I don't recall it exactly. I don't think there is any parser of C++/Python code guessing the degree. It is not anywhere in DOLFIN and UFL does not see Expression::eval code. Jan > > > > > > > > One can then supply the element to get the 'exact' integration > >> > >> Expression('x[0]*x[0]', pointwise=False, element=quadratic_element) > >> > > > > No need to pass pointwise. Just an element would do. > > > Sure, that was just to show the complete argument list, but passing > one argument should be enough. > > > > > > > > and test the implementation of the pointwise feature by > >> > >> Expression('x[0]*x[0]', pointwise=True) > >> > >> which, depending on implementation details, will have a virtual > >> degree equal to zero. > >> > > > > If one did > > > > f = Expression('x[0]*x[0]', pointwise=True) > > assemeble(f*dx) > > > > it should throw an error. Integrating f would require something like > > > > assemeble(f*dx(degree=2)) > > > Yes. > > Kristian > > > > > > > > > > In the future we can either set pointwise=True by default, or > > remove this > >> argument such that failing to provide the element implies pointwise > >> evaluation. > >> > > > > Yes. > > > > Garth > > > > > > Kristian > >> > >> > >> > >> > >> > >>> Garth > >>> > >>> > >>> > >>> Martin > >>>> > >>>> > >>>> These are two distinct issues: > >>>>>> 1) We need a "PointwiseExpression" with no degree and no hidden > >>>>>> interpolation under the hood. This expression is evaluated in > >>>>>> quadrature points - this is a clean concept and easy to > >>>>>> understand. > >>>>>> > >>>>> > >>>>> If PointwiseExpression and Expression will co-exist, I agree to > >>>>> this definition. If it makes implementation cleaner for > >>>>> quadrature degree estimation, we can set the degree equal to > >>>>> zero, but hidden from users. > >>>>> > >>>>> > >>>>> 2) Degree estimation is not exact and some people are confused > >>>>> by > >>>>>> that. But it is not exact today, never was claimed to be, and > >>>>>> never will be. If that's not acceptable, we can just as well > >>>>>> disable it completely. Disabling it where it isn't exact will > >>>>>> break a _lot_ of programs. What we _can_ do without breaking > >>>>>> programs or making the interface more cumbersome than today, > >>>>>> is to make it more obvious how to control the integration > >>>>>> degree, and to document it better. > >>>>>> > >>>>> > >>>>> I don't think anyone can disagree with this. > >>>>> > >>>>> Kristian > >>>>> > >>>>> > >>>>>> Martin > >>>>>> > >>>>> > >>>>> > >>>> > >>> > >> > > _______________________________________________ fenics mailing list [email protected] http://fenicsproject.org/mailman/listinfo/fenics
