| Thread Previous • Date Previous • Date Next • Thread Next |
On 13/04/10 21:00, Anders Logg wrote:
On Tue, Apr 13, 2010 at 08:19:33PM +0800, Garth N. Wells wrote:On 13/04/10 17:59, Anders Logg wrote:On Tue, Apr 13, 2010 at 11:31:38AM +0200, Mehdi Nikbakht wrote:On Tue, 2010-04-13 at 10:43 +0200, Anders Logg wrote:On Tue, Apr 13, 2010 at 10:31:51AM +0200, Mehdi Nikbakht wrote:On Tue, 2010-04-13 at 09:45 +0200, Anders Logg wrote:On Tue, Apr 13, 2010 at 07:17:28AM +0800, Garth N. Wells wrote:On 12/04/10 23:35, Anders Logg wrote:On Mon, Apr 12, 2010 at 10:20:13PM +0800, Garth N. Wells wrote:On 12/04/10 21:49, Anders Logg wrote:On Mon, Apr 12, 2010 at 09:34:38PM +0800, Garth N. Wells wrote:On 12/04/10 21:29, Anders Logg wrote:On Mon, Apr 12, 2010 at 09:21:32PM +0800, Garth N. Wells wrote:On 12/04/10 21:19, Garth N. Wells wrote:On 12/04/10 20:47, Anders Logg wrote:We are doing some work where we need to do run-time quadrature over arbitrary polyhedra.We (Mehdi and I) do this already (using UFC), so I don't see why a new function is required. Can you explain why evaluate_tensor is not enough?I meant 'tabulate_tensor'.Which function do you call for evaluating the integrand?We evaluate it inside ufc::tabulate_tensor. We construct our forms with an extra argument, say an object "CutCellIntegrator", which can provide quadrature schemes which depend on the considered cell.That would require a special purpose code generator.What's wrong with that? FFC won't (and shouldn't) be able to do everything. Just adding a function to UFC won't make FFC do what we do now. We reuse FFC (import modules) and add special purpose extensions.Exactly, it won't make FFC do what we need, but we could *use* FFC to do what we need (without adding a special-purpose code generator).Having evaluate_integrand would allow more flexibility for users to implement their own special quadrature scheme.We make "CutCellIntegrator" an abstract base class, so the user has *complete* freedom to define the quadrature scheme and the generated code does not depend on the scheme, since the scheme may depend on things like how the cell 'cut' is represented.Then it sounds to me like that generated code is not at all special, but instead general purpose and should be added to UFC/FFC. And the most general interface would (I think) be an interface for evaluating the integrand at a given point. We already have the same for basis functions (evaluate_basis_function) so it is a natural extension.I still don't see the need for 'evaluate_integrand' unless you plan to call it directly from the assembler side (i.e. DOLFIN). Is that the case? Perhaps you can give me a concrete example of how you plan to use it.Yes, that's the plan. In pseudo-code, this is what we want to do: for polyhedron in intersection.cut_cells: quadrature_rule = QuadratureRule(polyhedron) AK = 0 for (x, w) in quadrature_rule: AK += w * evaluate_integrand(x) A += AK All data representing the geometry, the polyhedra, mappings from those polyhedra to the original cells etc is in the Intersection class (in the sandbox): intersection = Intersection(mesh0, mesh1) Eventually, we might want to move part of the functionality into either DOLFIN or FFC, but having access to and evaluate_integrand function makes it possible to experiment (from the C++ side) without the need for building complex abstractions at this point.As far as I understood you want to compute the integration points for polyhedrons inside Dolfin and evaluate_integrands will just compute that integrand in that specific integration points. If it is the case how would you determine what is the order of quadrature rule that you want to use?The quadrature rule would be a simple option for the user to set. Currently we only have one general rule implemented which is barycentric quadrature.Then if you use higher order elements, you need to update your quadrature rule manually. I think it would be nice to compute your own quadrature rule inside tabualte_tensor by using the standard quadrature rule.Yes, but the problem is that the computation of the quadrature rule is nontrivial and it might be better to do it from C++ than as part of the generated code, especially when the code depends on external libraries like CGAL. See BarycenterQuadrature.cpp.What we plan is that the generated code would provide a 1D scheme (or possible just a polynomial order) and the C++ quadrature object could apply it on a sub-triangulations.I don't understand. How would the 1D scheme connect to a triangulation?
Construct a rule for triangles from the 1D scheme (like FIAT does) and apply it on each sub-triangle.
Btw, we have considered subtriangulating the polyhedra to get a quadrature rule, but decided to just to barycentric quadrature for now to save time, both implementation and run time...
We sub-triangulate. I'm hoping that some of the work can be moved to CGAL.
Since to evaluate integrands you need to tabulate basis functions and their derivatives on arbitrary points. How do you want to tabulate basis functions and their derivatives inside evaluate_integrands?That's a good point. We would need to evaluate the basis functions and their derivatives at a given arbitrary point which is not known at compile-time.Yes, you need to use the tabulate_basis* functions implemented by Kristian. Then I am not the only one who is using them. Good for Kristian. ;)ok, good. Then there is no principal problem of generating code for evaluate_integrand (if it is allowed to call tabaulate_basis). I don't see what the problem is of adding evaluate_integrand. It is a natural extension (we have evaluate_basis already), it would be "simple" to implement (Kristian can correct me if I'm wrong) and it makes generated code useful for assembly over cut meshes (without the need for writing a special-purpose code generator).It's not a problem - just good to scrutinise heavily additions toAgree.UFC to avoid bloat. I'm happy for it to be added. We should make clear that it's for testing/special cases since the performance for computing element tensors will be poor compared to tabulate_tensor.ok. So we have the following proposed additions to UFC: evaluate_integrand tabulate_tensor(x, w) evaluate_basis (two versions) + quite a few other suggestions as blueprints https://blueprints.launchpad.net/ufc Perhaps it's a good time to think about UFC 2.0 and take care of all the blueprints.
I think that we're settled on evaluate_basis, so this can go into UFC-dev. The other needs more discussion (and more work), and I won't have time for a while.
Garth Or perhaps we should wait until after DOLFIN 1.0?
I have limited resources at the moment so we might want to postpone this. -- Anders
| Thread Previous • Date Previous • Date Next • Thread Next |
