← Back to team overview

ffc team mailing list archive

Re: [UFL-dev] [HG UFL] Added check for restriction on facet normals and form arguments

 

On Mon, Jun 8, 2009 at 5:45 PM, Anders Logg<logg@xxxxxxxxx> wrote:
> On Mon, Jun 08, 2009 at 05:29:55PM +0200, Martin Sandve Alnæs wrote:
>> On Mon, Jun 8, 2009 at 4:20 PM, Kristian
>> Oelgaard<k.b.oelgaard@xxxxxxxxxx> wrote:
>> > Quoting Harish Narayanan <harish.mlists@xxxxxxxxx>:
>> >
>> >> Anders Logg wrote:
>> >> > On Mon, Jun 08, 2009 at 03:10:58PM +0200, Harish Narayanan wrote:
>> >> >> UFL wrote:
>> >> >>> One or more new changesets pushed to the primary ufl repository.
>> >> >>> A short summary of the last three changesets is included below.
>> >> >>>
>> >> >>> changeset:   897:7c6c19b55ba45ffbec4bfffc68ae0dce9636d20d
>> >> >>> tag:         tip
>> >> >>> user:        "Martin Sandve Alnæs <martinal@xxxxxxxxx>"
>> >> >>> date:        Mon Jun 08 10:29:42 2009 +0200
>> >> >>> files:       sandbox/scratch/restrictiontest.py ufl/algorithms/checks.py
>> >> ufl/algorithms/propagate_restrictions.py
>> >> >>> description:
>> >> >>> Added check for restriction on facet normals and form arguments
>> >> >>> inside the propagate_restrictions algorithm (quick and easy fix).
>> >> >> With this changeset, I now see:
>> >> >>
>> >> >> *** FFC: Form argument must be restricted.
>> >> >> *** FFC: To get more information about this error, rerun FFC with --debug.
>> >> >>
>> >> >>
>> >> >> for any form I try to compile. Must something be added to FFC for it to
>> >> >> work with this UFL change?
>> >> >
>> >> > What do you get with --debug? It should then print out the traceback
>> >> > for the UFL exception.
>> >>
>> >>
>> >> $ ffc --debug -l dolfin Poisson.ufl
>> >> This is FFC, the FEniCS Form Compiler, version 0.6.2.
>> >> For further information, visit http://www.fenics.org/ffc/.
>> >>
>> >> Preprocessing form file: Poisson.ufl --> Poisson.py
>> >>
>> >> Compiler stage 1: Analyzing form
>> >> --------------------------------
>> >> This form is linear in 2 arguments.
>> >>   Name:                               a
>> >>   Rank:                               2
> 5F> >>   Cell:                               <triangle of degree 1>
>> >>   Geometric dimension:                2
>> >>   Topological dimension:              2
>> >>   Number of functions:                0
>> >>   Number of cell integrals:           1
>> >>   Number of exterior facet integrals: 0
>> >>   Number of interior facet integrals: 0
>> >>   Basis functions:                    [v_0, v_1]
>> >>   Functions:                          []
>> >>   Basis function names:               [v0, v1]
>> >>   Function names:                     []
>> >>   Unique elements:                    CG1
>> >>   Unique sub elements:                CG1
>> >>
>> >>   Automatic selection of representation not implemented, defaulting to
>> >> quadrature.
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>   Found element in dof map cache: <CG1 on a <triangle of degree 1>>
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>
>> >> Compiler stage 2: Computing form representation(s)
>> >> --------------------------------------------------
>> >>
>> >>   QR, init, form:
>> >>   { sum_{i_1} (({ A | A_{i_0} = dv_0/dx_i_0 })[i_1]) * (({ A | A_{i_2} =
>> >> dv_1/dx_i_2 })[i_1])  } * dx0
>> >>   Computing quadrature representation
>> >>   Found element in element cache: <CG1 on a <triangle of degree 1>>
>> >>   Derivatives:
>> >> set([SpatialDerivative(BasisFunction(FiniteElement('Lagrange',
>> >> Cell('triangle', 1), 1), 0), MultiIndex((Index(0),), {Index(0): 2})),
>> >> SpatialDerivative(BasisFunction(FiniteElement('Lagrange',
>> >> Cell('triangle', 1), 1), 1), MultiIndex((Index(2),), {Index(2): 2}))])
>> >>   num_derivatives: {FiniteElement('Lagrange', Cell('triangle', 1), 1): 1}
>> >>
>> >>   QR, init, psi_tables:
>> >>   {'cell': {1: {FiniteElement('Lagrange', Cell('triangle', 1), 1):
>> >> {None: [{(0, 0): array([[ 0.33333333],
>> >>          [ 0.33333333],
>> >>          [ 0.33333333]])}, {(0, 1): array([[ -1.00000000e+00],
>> >>          [  7.85046229e-17],
>> >>          [  1.00000000e+00]]), (1, 0): array([[-1.],
>> >>          [ 1.],
>> >>          [ 0.]])}]}}}, 'exterior_facet': {}, 'interior_facet': {}}
>> >>
>> >>   QR, init, quadrature_weights:
>> >>   {'cell': {1: array([ 0.5])}, 'exterior_facet': {}, 'interior_facet': {}}
>> >>
>> >> Compiler stage 3: Optimizing form representation
>> >> ------------------------------------------------
>> >>   Optimization of tensor contraction representation currently broken (to
>> >> be fixed).
>> >>
>> >> Compiler stage 4: Generating code
>> >> ---------------------------------
>> >>   Generating code for finite elements...
>> >>   Removing unused variable: Jinv_11
>> >>   Removing unused variable: Jinv_10
>> >>   Removing unused variable: Jinv_01
>> >>   Removing unused variable: Jinv_00
>> >>   Removing unused variable: Jinv_11
>> >>   Removing unused variable: Jinv_10
>> >>   Removing unused variable: Jinv_01
>> >>   Removing unused variable: Jinv_00
>> >>   Removing unused variable: Jinv_11
>> >>   Removing unused variable: Jinv_10
>> >>   Removing unused variable: Jinv_01
>> >>   Removing unused variable: Jinv_00
>> >>   Removing unused variable: Jinv_11
>> >>   Removing unused variable: Jinv_10
>> >>   Removing unused variable: Jinv_01
>> >>   Removing unused variable: Jinv_00
>> >>   done
>> >>   Generating code for finite dof maps...
>> >>   done
>> >>   Generating code for form...
>> >>   done
>> >>   Generating code using quadrature representation
>> >>
>> >>   QG-utils, psi_tables:
>> >>   {1: {FiniteElement('Lagrange', Cell('triangle', 1), 1): {None: [{(0,
>> >> 0): array([[ 0.33333333],
>> >>          [ 0.33333333],
>> >>          [ 0.33333333]])}, {(0, 1): array([[ -1.00000000e+00],
>> >>          [  7.85046229e-17],
>> >>          [  1.00000000e+00]]), (1, 0): array([[-1.],
>> >>          [ 1.],
>> >>          [ 0.]])}]}}}
>> >>
>> >>   QG-utils, flatten_tables, points:
>> >>   1
>> >>
>> >>   QG-utils, flatten_tables, elem_dict:
>> >>   {FiniteElement('Lagrange', Cell('triangle', 1), 1): {None: [{(0, 0):
>> >> array([[ 0.33333333],
>> >>          [ 0.33333333],
>> >>          [ 0.33333333]])}, {(0, 1): array([[ -1.00000000e+00],
>> >>          [  7.85046229e-17],
>> >>          [  1.00000000e+00]]), (1, 0): array([[-1.],
>> >>          [ 1.],
>> >>          [ 0.]])}]}}
>> >>
>> >>   QG-utils, flatten_tables, elem:
>> >>   <CG1 on a <triangle of degree 1>>
>> >>
>> >>   QG-utils, flatten_tables, facet_tables:
>> >>   {None: [{(0, 0): array([[ 0.33333333],
>> >>          [ 0.33333333],
>> >>          [ 0.33333333]])}, {(0, 1): array([[ -1.00000000e+00],
>> >>          [  7.85046229e-17],
>> >>          [  1.00000000e+00]]), (1, 0): array([[-1.],
>> >>          [ 1.],
>> >>          [ 0.]])}]}
>> >>
>> >>   QG-utils, flatten_tables, derivs:
>> >>   (0, 0)
>> >>
>> >>   QG-utils, flatten_tables, psi_table:
>> >>   [[ 0.33333333]
>> >>    [ 0.33333333]
>> >>    [ 0.33333333]]
>> >>   Table name: FE0
>> >>
>> >>   QG-utils, flatten_tables, derivs:
>> >>   (0, 1)
>> >>
>> >>   QG-utils, flatten_tables, psi_table:
>> >>   [[ -1.00000000e+00]
>> >>    [  7.85046229e-17]
>> >>    [  1.00000000e+00]]
>> >>   Table name: FE0_D01
>> >>
>> >>   QG-utils, flatten_tables, derivs:
>> >>   (1, 0)
>> >>
>> >>   QG-utils, flatten_tables, psi_table:
>> >>   [[-1.]
>> >>    [ 1.]
>> >>    [ 0.]]
>> >>   Table name: FE0_D10
>> >>
>> >>   QG-utils, psi_tables, flat_tables:
>> >>   {'FE0_D10': array([[-1.,  1.,  0.]]), 'FE0_D01': array([[
>> >> -1.00000000e+00,   7.85046229e-17,   1.00000000e+00]]), 'FE0': array([[
>> >> 0.33333333,  0.33333333,  0.33333333]])}
>> >>
>> >>   tables: {'FE0_D10': array([[-1.,  1.,  0.]]), 'FE0_D01': array([[
>> >> -1.00000000e+00,   7.85046229e-17,   1.00000000e+00]]), 'FE0': array([[
>> >> 0.33333333,  0.33333333,  0.33333333]])}
>> >>
>> >>   name_map: {}
>> >>
>> >>   inv_name_map: {'FE0_D10': 'FE0_D10', 'FE0_D01': 'FE0_D01', 'FE0': 'FE0'}
>> >>
>> >>   QG-utils, psi_tables, unique_tables:
>> >>   {'FE0_D10': array([[-1.,  1.,  0.]]), 'FE0_D01': array([[-1.,  0.,
>> >> 1.]]), 'FE0': array([[ 0.33333333,  0.33333333,  0.33333333]])}
>> >>
>> >>   QG-utils, psi_tables, name_map:
>> >>   {'FE0_D10': ('FE0_D10', (), False, False), 'FE0_D01': ('FE0_D01', (),
>> >> False, False), 'FE0': ('FE0', (), False, False)}
>> >>   Generating code for cell integrals using quadrature representation.
>> >>
>> >>   QG, cell_integral, integrals:
>> >>   {1:
>> >>
>> > Integral(IndexSum(Product(Indexed(ComponentTensor(SpatialDerivative(BasisFunction(FiniteElement('Lagrange',
>> >> Cell('triangle', 1), 1), 0), MultiIndex((Index(0),), {Index(0): 2})),
>> >> MultiIndex((Index(0),), {Index(0): 2})), MultiIndex((Index(1),),
>> >> {Index(1): 2})),
>> >>
>> > Indexed(ComponentTensor(SpatialDerivative(BasisFunction(FiniteElement('Lagrange',
>> >> Cell('triangle', 1), 1), 1), MultiIndex((Index(2),), {Index(2): 2})),
>> >> MultiIndex((Index(2),), {Index(2): 2})), MultiIndex((Index(1),),
>> >> {Index(1): 2}))), MultiIndex((Index(1),), {Index(1): 2})),
>> >> Measure('cell', 0, {}))}
>> >>   Looping points: 1
>> >>   integral: { sum_{i_1} (({ A | A_{i_0} = dv_0/dx_i_0 })[i_1]) * (({ A |
>> >> A_{i_2} = dv_1/dx_i_2 })[i_1])  } * dx0<{}>
>> >>
>> >>   Integral tree_format: Integral:
>> >>       domain type: cell
>> >>       domain id: 0
>> >>       integrand:
>> >>           IndexSum
>> >>           (
>> >>               Product
>> >>               (
>> >>                   Indexed
>> >>                   (
>> >>                       ComponentTensor
>> >>                       (
>> >>                           SpatialDerivative
>> >>                           (
>> >>                               BasisFunction(FiniteElement('Lagrange',
>> >> Cell('triangle', 1), 1), 0)
>> >>                               MultiIndex((Index(0),), {Index(0): 2})
>> >>                           )
>> >>                           MultiIndex((Index(0),), {Index(0): 2})
>> >>                       )
>> >>                       MultiIndex((Index(1),), {Index(1): 2})
>> >>                   )
>> >>                   Indexed
>> >>                   (
>> >>                       ComponentTensor
>> >>                       (
>> >>                           SpatialDerivative
>> >>                           (
>> >>                               BasisFunction(FiniteElement('Lagrange',
>> >> Cell('triangle', 1), 1), 1)
>> >>                               MultiIndex((Index(2),), {Index(2): 2})
>> >>                           )
>> >>                           MultiIndex((Index(2),), {Index(2): 2})
>> >>                       )
>> >>                       MultiIndex((Index(1),), {Index(1): 2})
>> >>                   )
>> >>               )
>> >>               MultiIndex((Index(1),), {Index(1): 2})
>> >>           )
>> >>
>> >>   QG, Using Transformer
>> >> Form argument must be restricted.
>> >> Traceback (most recent call last):
>> >>   File "/Users/harish/Work/FEniCS/unstable/build/bin/ffc", line 186, in
>> >> <module>
>> >>     sys.exit(main(sys.argv[1:]))
>> >>   File "/Users/harish/Work/FEniCS/unstable/build/bin/ffc", line 127, in main
>> >>     execfile(script, {})
>> >>   File "Poisson.py", line 34, in <module>
>> >>     compile([a, L, M, element], "Poisson", {'log_level': 10, 'format':
>> >> 'dolfin', 'form_postfix': True, 'quadrature_order': 'auto', 'precision':
>> >> '15', 'cpp optimize': False, 'cache_dir': None, 'split': False,
>> >> 'representation': 'auto', 'optimize': False, 'output_dir': '.'}, globals())
>> >>   File
>> >>
>> > "/Users/harish/Work/FEniCS/unstable/build/lib/python2.5/site-packages/ffc/compiler/compiler.py",
>> >> line 99, in compile
>> >>     form_code = generate_form_code(form_data, representations, prefix,
>> >> format.format, options)
>> >>   File
>> >>
>> > "/Users/harish/Work/FEniCS/unstable/build/lib/python2.5/site-packages/ffc/compiler/compiler.py",
>> >> line 191, in generate_form_code
>> >>     codes.append(code_generator.generate_integrals(representations[i],
>> >> format))
>> >>   File
>> >>
>> > "/unstable/lib/python2.5/site-packages/ffc/compiler/codegeneration/quadrature/quadraturegenerator.py",
>> >> line 74, in generate_integrals
>> >>     code.update(self.generate_cell_integrals(form_representation, format))
>> >>   File
>> >>
>> > "/unstable/lib/python2.5/site-packages/ffc/compiler/codegeneration/quadrature/quadraturegenerator.py",
>> >> line 102, in generate_cell_integrals
>> >>     self.generate_cell_integral(form_representation, transformer,
>> >> integrals, format)
>> >>   File
>> >>
>> > "/unstable/lib/python2.5/site-packages/ffc/compiler/codegeneration/quadrature/quadraturegenerator.py",
>> >> line 165, in generate_cell_integral
>> >>     integrals, Indent, format)
>> >>   File
>> >>
>> > "/unstable/lib/python2.5/site-packages/ffc/compiler/codegeneration/quadrature/quadraturegenerator.py",
>> >> line 352, in __generate_element_tensor
>> >>     generate_code(integral.integrand(), transformer, Indent, format)
>> >>   File
>> >>
>> > "/unstable/lib/python2.5/site-packages/ffc/compiler/codegeneration/quadrature/quadraturetransformer.py",
>> >> line 1132, in generate_code
>> >>     new_integrand = propagate_restrictions(new_integrand)
>> >
>> > I guess the error occurs here. I apply a series of UFL algorithms:
>> >
>> > new_integrand = expand_indices(integrand)
>> > new_integrand = purge_list_tensors(new_integrand)
>> > new_integrand = propagate_restrictions(new_integrand)
>> >
>> > before starting the transformation.
>> >
>> > Should this not be legal/work?
>>
>> I figured propagate_restrictions should only be applied
>> to the integrands of interior facet integrals.
>
> Is propagate_restrictions performing the test?
>
> In that case, it would be natural to have it check that terms on
> interior facets *are* restricted and terms not on interior facet are
> *not* restricted.
>
> Also, would it be a good idea to wrap expand_indices,
> purge_list_tensors and propagate_restrictions into a common function
> to be applied for massaging a form before using it?

No. But propagate_restrictions can probably be merged with
renumber_indices and some other stuff. However, it's best to
wait until things are more streamlined before deciding that.

(purge_list_tensors is currently implemented by calling expand_indices,
which does much more than purge_list_tensors should do)

Martin


References