syfi team mailing list archive

[Kent-Andre Mardal <kent-and@xxxxxxxxx>]

fre, 02.05.2008 kl. 14.25 +0200, skrev Martin Sandve Alnæs:
> 2008/5/2 Jed Brown <jed@xxxxxxxx>:
> > On Fri 2008-05-02 12:34, Martin Sandve Alnæs wrote:
> >> Let me try to summarize what's already there.
> >>
> >> In FFC, I think the action of a form can be accomplished by simply changing
> >>   u = TrialFunction(element)
> >> to
> >>   u = Function(element)
> >> such that u becomes a coefficient, the vector to multiply with.
> >
> > Oh, right.  We can use a linear form to compute the action of the bilinear form
> > without getting matrix entries.  So this part is already implemented and well
> > tested.  Great.
> >
> >> In SFC, you can write
> >>   b = Action(a)
> >> if a is a Form. This will simply do the element matrix-vector
> >> product symbolically prior to compilation.
> >> (This hasn't been used or tested very well, but it is short and simple)
> >
> > How do you think this would compare to using a linear form?  For high order
> > elements, the symbolic expression could be so huge the compiler has trouble
> > doing the necessary CSE.
> 
> It shouldn't be very different, but I think those problems will be the
> same any way.
> 
> Frankly, I don't think SFC will cope well with p=10 anyway.
> I plan to fix some bottlenecks to make it more robust.
> However, there are some things I need UFL to be able to do,
> so it may take a little time.
> 

I have just printed out the basis functions on a unit-cube p=10. 
It looks like: 

fe.N(0)=(1.0/9.3329542656000000e+10)*
( 602595.0*(y*(y*y*y)*(y*y*y))*pow( y-1.0,3.0)+87498.0*y*pow( y-1.0,9.0)+4536.0*pow( y-1.0,10.0)+87498.0*(((y*y)*(y*y))*((y*y)*(y*y)))*pow( y-1.0,2.0)+1996867.0*(y*y*y)*pow( y-1.0,7.0)+602595.0*(y*y)*pow( y-1.0,8.0)+3558504.0*(y*(y*y)*(y*y))*pow( y-1.0,5.0)+4536.0*(y*((y*y)*(y*y))*((y*y)*(y*y)))*( y-1.0)+3558504.0*((y*y)*(y*y))*pow( y-1.0,6.0)+1996867.0*((y*y*y)*(y*y*y))*pow( y-1.0,4.0))*( 602595.0*pow( x-1.0,3.0)*(x*(x*x*x)*(x*x*x))+4536.0*pow( x-1.0,10.0)+602595.0*pow( x-1.0,8.0)*(x*x)+1996867.0*pow( x-1.0,7.0)*(x*x*x)+4536.0*( x-1.0)*(x*((x*x)*(x*x))*((x*x)*(x*x)))+3558504.0*pow( x-1.0,6.0)*((x*x)*(x*x))+1996867.0*pow( x-1.0,4.0)*((x*x*x)*(x*x*x))+87498.0*pow( x-1.0,9.0)*x+3558504.0*pow( x-1.0,5.0)*(x*(x*x)*(x*x))+87498.0*pow( x-1.0,2.0)*(((x*x)*(x*x))*((x*x)*(x*x))))*( 602595.0*(z*(z*z*z)*(z*z*z))*pow( z-1.0,3.0)+87498.0*z*pow( z-1.0,9.0)+87498.0*(((z*z)*(z*z))*((z*z)*(z*z)))*pow( z-1.0,2.0)+1996867.0*(z*z*z)*pow( z-1.0,7.0)+602595.0*(z*z)*pow( z-1.0,8.0)+3558504.0*(z*(z*z)*(z*z))*pow( z-1.0,5.0)+4536.0*(z*((z*z)*(z*z))*((z*z)*(z*z)))*( z-1.0)+3558504.0*((z*z)*(z*z))*pow( z-1.0,6.0)+1996867.0*((z*z*z)*(z*z*z))*pow( z-1.0,4.0)+4536.0*pow( z-1.0,10.0))

The expression is a bit large and its more than a thousand of them so
CSE probably is not a good idea. However, one can run CSE on a group of
functions to start with and then go trough all in a systematic fashion.
SFC will not do this automatically, but it should still be feasible.
Doing this in the "proper way" could also merit a paper. It should
also be possible to pre-compute the quantities in quadrature points
etc.  

SyFi is sufficiently fast that it is possible to experiment with
different algorithms here. Computing the element functions, computing
one entry in the element matrix etc typically takes less than 1s, even
with this analytical approach. 

> 
> >> SFC supports hexahedral elements, FFC does not.
> >
> > It sounds like I can get what I want if I write a form file, specify high order
> > hex elements, just have a linear form, and specify quadrature.  Can I expect
> > that applying this form will be O(p^4)?
> 
> Note that SFC doesn't yet have form files in the same sense as FFC does.
> I'm working on the unified form language (UFL), but it's not finished.
> There are some SFC demos in the newest syfi under examples/solvers/.
> 
> 
> >> Adding mappings to the form compilers shouldn't be too hard with quadrature.
> >>
> >> UFC doesn't cover high order geometries, varying cell types in a mesh
> >> or varying element orders. This is probably the largest obstacle.
> >
> > Indeed.  We can do some proof of concept work without these things, but they are
> > clearly desirable in the long term.  It would greatly complicate the spec,
> > though.  Would it be better to create a new HpUFC with these features?  Of
> > course, if hp-adaptivity is made easy, everyone will be doing it :-)
> 
> That is a good idea, and it can be merged back in UFC for UFC 2.0.
> 
> --
> Martin

We would _really_ like to have hp in UFC, but it might be that the best
way to start is a HpUFC. This would alow experimentation. I think this
could be a very good project. 


Kent