SyFi  0.3
test_form_argument_code.py
Go to the documentation of this file.
00001 #!/usr/bin/env python
00002 
00003 __authors__ = "Martin Sandve Alnes"
00004 __date__ = "2008-09-04 -- 2012-05-16"
00005 
00006 import unittest
00007 import os, sys, glob, shutil, commands
00008 
00009 import ufl
00010 
00011 from ufl import FiniteElement
00012 from ufl import VectorElement
00013 from ufl import TensorElement
00014 from ufl import MixedElement
00015 
00016 from ufl import Argument
00017 from ufl import TestFunction
00018 from ufl import TrialFunction
00019 
00020 from ufl import Coefficient
00021 from ufl import Constant, VectorConstant, TensorConstant
00022 
00023 from ufl import dx, ds
00024 from ufl.common import product
00025 
00026 import SyFi
00027 import sfc
00028 from sfc.representation.geometryrepresentation import GeometryRepresentation
00029 from sfc.representation.elementrepresentation import ElementRepresentation
00030 import numpy
00031 
00032 from test_tempdir_base import TempDirTestBase
00033 from cell_assembly import assemble_on_cell, cell2volume, num_integrals
00034 
00035 simplices = (ufl.interval, ufl.triangle, ufl.tetrahedron)
00036 quads = (ufl.quadrilateral, ufl.hexahedron)
00037 
00038 #from sfc.representation import FormArgumentRepresentation
00039 from sfc.codegeneration.formargumentcg import FormArgumentCG
00040 
00041 try:
00042     from uflacs.codeutils.format_code_structure import format_code_structure
00043 except:
00044     def format_code_structure(x): return str(x)
00045 
00046 class ElementRepTest(unittest.TestCase):
00047     def test_foo(self):
00048         SyFi.initSyFi(2)
00049 
00050         family = "CG"
00051         ufl_cell = ufl.triangle
00052         degree = 1
00053         ufl_element = FiniteElement(family, ufl_cell, degree)
00054         geomrep = GeometryRepresentation(ufl_cell.domain())
00055         er = ElementRepresentation(ufl_element, geomrep, quad_rule=None, options=None, cache=None)
00056         print
00057         print
00058         print dir(er)
00059         print
00060 
00061         ngrad = 0
00062         code = []
00063         for dof in range(er.local_dimension):
00064             #sub_element_index = er.dof_to_sub_element_index(dof)
00065 
00066             components = [()] # FIXME: permutations of value shape
00067 
00068             if ngrad:
00069                 assert ngrad == 1
00070                 all_derivatives = [(d,) for d in range(nsd)] # FIXME: permutations of ngrad
00071 
00072                 for component in components:
00073                     for directions in all_derivatives:
00074                         expr = er.basis_function_derivative(dof, component, directions)
00075                         name = "FIXME"
00076                         code.append((name, expr.printc()))
00077             else:
00078                 for component in components:
00079                     expr = er.basis_function(dof, component)
00080                     name = "FIXME"
00081                     code.append((name, expr.printc()))
00082 
00083         print "Element code:"
00084         for n,e in code:
00085             print "%s = %s;" % (n, e)
00086 
00087 class FormArgumentCodeGenTest(object):#unittest.TestCase):
00088     def xtest_form_argument_representation_interval(self):
00089         ufl_cell = ufl.interval
00090 
00091         SyFi.initSyFi(ufl_cell.d)
00092 
00093         # FIXME: Must FormArgumentCG really depend on all this stuff?
00094         #formdata, element_reps, geomrep, options = ...
00095         formrep = None #FormRepresentation(formdata, element_reps, geomrep, options)
00096         farep = FormArgumentRepresentation(formrep)
00097         facg = FormArgumentCG(farep, False)
00098         decl = facg.generate_declarations()
00099 
00100         self.assertTrue(decl)
00101         if 1:
00102             print
00103             print ufl_cell
00104             print format_code_structure(decl)
00105             print
00106 
00107     def xtest_form_argument_representation_triangle(self):
00108         ufl_cell = ufl.triangle
00109 
00110         SyFi.initSyFi(ufl_cell.d)
00111 
00112         farep = FormArgumentRepresentation()
00113         facg = FormArgumentCG(farep) # FIXME
00114         decl = facg.generate_declarations()
00115 
00116         self.assertTrue(decl)
00117         if 1:
00118             print
00119             print ufl_cell
00120             print format_code_structure(decl)
00121             print
00122 
00123     def xtest_form_argument_representation_tetrahedron(self):
00124         ufl_cell = ufl.tetrahedron
00125 
00126         SyFi.initSyFi(ufl_cell.d)
00127 
00128         farep = FormArgumentRepresentation()
00129         facg = FormArgumentCG(farep) # FIXME
00130         decl = facg.generate_declarations()
00131 
00132         self.assertTrue(decl)
00133         if 1:
00134             print
00135             print ufl_cell
00136             print format_code_structure(decl)
00137             print
00138 
00139     def xtest_form_argument_representation_quadrilateral(self):
00140         ufl_cell = ufl.quadrilateral
00141     def xtest_form_argument_representation_hexahedron(self):
00142         ufl_cell = ufl.hexahedron
00143 
00144 
00145 class FormArgumentTestMixin(object):
00146 
00147     def assertAlmostEqualNumpy(self, actual, expected, eps=1e-8):
00148         self.assertEqual(actual.shape, expected.shape)
00149         mdiff = expected - actual
00150         vdiff = mdiff.reshape((product(mdiff.shape),))
00151         sdiff = numpy.dot(vdiff,vdiff)
00152         self.assertAlmostEqual(sdiff, 0.0, eps)
00153 
00154     def assertAssembledAlmostEqual(self, ufl_form, ufl_cell, coeffs, expected):
00155         """Test that the integral of 1.0 over a unit cell
00156         equals the length/area/volume of the unit cell."""
00157         ufc_form, module, formdata, prefix = sfc.jit(ufl_form, parameters=self.options)
00158         assembled = assemble_on_cell(ufc_form, ufl_cell, coeffs=coeffs)
00159         if isinstance(assembled, float):
00160             self.assertAlmostEqual(assembled, expected)
00161         else:
00162             self.assertAlmostEqualNumpy(assembled, expected)
00163 
00164 
00165     def _test_constant(self, cell):
00166         "Test that a single constant is handled correctly."
00167         D = cell2volume[cell]
00168         c = Constant(cell)
00169         a = c*dx
00170         coeffs = [1.2]
00171         expected = 1.2 * D
00172         self.assertAssembledAlmostEqual(a, cell, coeffs, expected)
00173 
00174     def xtest_constant_interval(self):
00175         self._test_constant(ufl.interval)
00176     def xtest_constant_triangle(self):
00177         self._test_constant(ufl.triangle)
00178     def xtest_constant_tetrahedron(self):
00179         self._test_constant(ufl.tetrahedron)
00180 
00181     def _test_constant_quadrilateral(self): # TODO: Not supported by dolfin yet
00182         self._test_constant(ufl.quadrilateral)
00183     def _test_constant_hexahedron(self): # TODO: Not supported by dolfin yet
00184         self._test_constant(ufl.hexahedron)
00185 
00186 
00187     def _test_constants(self, cell):
00188         "Test that multiple constants are handled correctly."
00189         D = cell2volume[cell]
00190         c1 = Constant(cell)
00191         c2 = Constant(cell)
00192         c3 = Constant(cell)
00193         a = c1*c2*c3*dx
00194         coeffs = [2,3,5]
00195         expected = coeffs[0]*coeffs[1]*coeffs[2] * D
00196         self.assertAssembledAlmostEqual(a, cell, coeffs, expected)
00197 
00198     def xtest_constants_interval(self):
00199         self._test_constants(ufl.interval)
00200     def xtest_constants_triangle(self):
00201         self._test_constants(ufl.triangle)
00202     def xtest_constants_tetrahedron(self):
00203         self._test_constants(ufl.tetrahedron)
00204 
00205     def _test_constants_quadrilateral(self): # TODO: Not supported by dolfin yet
00206         self._test_constants(ufl.quadrilateral)
00207     def _test_constants_hexahedron(self): # TODO: Not supported by dolfin yet
00208         self._test_constants(ufl.hexahedron)
00209 
00210 
00211     def _test_vector_constant(self, cell):
00212         "Test that a vector valued constant is handled correctly."
00213         D = cell2volume[cell]
00214         c = VectorConstant(cell)
00215         a = sum(c[i] for i in range(cell.d))*dx
00216         coeffs = [(1.2,3.4,5.6)[:cell.d]]
00217         expected = sum(coeffs[0]) * D
00218         self.assertAssembledAlmostEqual(a, cell, coeffs, expected)
00219 
00220     def xtest_vector_constant_interval(self):
00221         self._test_vector_constant(ufl.interval)
00222     def xtest_vector_constant_triangle(self):
00223         self._test_vector_constant(ufl.triangle)
00224     def xtest_vector_constant_tetrahedron(self):
00225         self._test_vector_constant(ufl.tetrahedron)
00226 
00227     def _test_vector_constant_quadrilateral(self): # TODO: Not supported by dolfin yet
00228         self._test_vector_constant(ufl.quadrilateral)
00229     def _test_vector_constant_hexahedron(self): # TODO: Not supported by dolfin yet
00230         self._test_vector_constant(ufl.hexahedron)
00231 
00232 
00233     def _test_tensor_constant(self, cell):
00234         "Test that a tensor valued constant is handled correctly."
00235         D = cell2volume[cell]
00236         c = TensorConstant(cell)
00237         sh = c.shape()
00238         a = sum(c[i,j] for i in range(sh[0]) for j in range(sh[1]))*dx
00239         rows = ( (1.1,1.2,1.3)[:sh[0]],
00240                  (2.1,2.2,2.3)[:sh[0]],
00241                  (3.1,3.2,3.3)[:sh[0]] )[:sh[1]]
00242         coeffs = [rows]
00243         expected = D * sum(rows[i][j] for i in range(sh[0]) for j in range(sh[1]))
00244         self.assertAssembledAlmostEqual(a, cell, coeffs, expected)
00245 
00246     def xtest_tensor_constant_interval(self):
00247         self._test_tensor_constant(ufl.interval)
00248     def xtest_tensor_constant_triangle(self):
00249         self._test_tensor_constant(ufl.triangle)
00250     def xtest_tensor_constant_tetrahedron(self):
00251         self._test_tensor_constant(ufl.tetrahedron)
00252 
00253     def _test_tensor_constant_quadrilateral(self): # TODO: Not supported by dolfin yet
00254         self._test_tensor_constant(ufl.quadrilateral)
00255     def _test_tensor_constant_hexahedron(self): # TODO: Not supported by dolfin yet
00256         self._test_tensor_constant(ufl.hexahedron)
00257 
00258 
00259     def _test_coefficients(self, cell):
00260         self._test_coefficient(cell, "DG", 0)
00261         self._test_coefficient(cell, "DG", 1)
00262         #self._test_coefficient(cell, "DG", 2)
00263         #self._test_coefficient(cell, "DG", 3)
00264 
00265         #self._test_coefficient(cell, "CG", 1)
00266         #self._test_coefficient(cell, "CG", 2)
00267         #self._test_coefficient(cell, "CG", 3)
00268 
00269     def _test_coefficient(self, cell, family="CG", degree=1):
00270         "Test that a single coefficient is handled correctly."
00271         D = cell2volume[cell]
00272         V = FiniteElement(family, cell, degree)
00273         c = Coefficient(V)
00274         a = c*dx
00275         coeffs = [1.2]
00276         expected = 1.2 * D
00277         self.assertAssembledAlmostEqual(a, cell, coeffs, expected)
00278 
00279     # TODO: Make these work with uflacs:
00280     def xtest_coefficients_interval(self):
00281         self._test_coefficients(ufl.interval)
00282     def xtest_coefficients_triangle(self):
00283         self._test_coefficients(ufl.triangle)
00284     def xtest_coefficients_tetrahedron(self):
00285         self._test_coefficients(ufl.tetrahedron)
00286 
00287     def xtest_coefficients_quadrilateral(self): # TODO: Not supported by dolfin yet
00288         self._test_coefficients(ufl.quadrilateral)
00289     def xtest_coefficients_hexahedron(self): # TODO: Not supported by dolfin yet
00290         self._test_coefficients(ufl.hexahedron)
00291 
00292 
00293     def test_code_snippets(self):
00294         print
00295         cell = ufl.triangle
00296 
00297         # Mock up a formrep object
00298         class Mock(object): pass
00299         formrep = Mock()
00300         formrep.geomrep = Mock()
00301         formrep.geomrep.sfc_cell = Mock()
00302         formrep.geomrep.sfc_cell.nsd = cell.d
00303 
00304         facg = FormArgumentCG(formrep, on_facet=False)
00305 
00306         S = FiniteElement("CG", cell, 1)
00307         u = Argument(S, count=1)
00308         f = Coefficient(S, count=2)
00309         print
00310         print facg.gen_argument_name(u, (), (), None)
00311         print facg.gen_argument_name(u, (0,), (), None)
00312         print facg.gen_coefficient_name(f, (), (), None)
00313         print facg.gen_coefficient_name(f, (0,), (), None)
00314         print
00315 
00316         V = VectorElement("CG", cell, 1)
00317         v = Argument(V, count=1)
00318         g = Coefficient(V, count=2)
00319         print
00320         print facg.gen_argument_name(v, (), (0,), None)
00321         print facg.gen_argument_name(v, (0,), (1,), None)
00322         print facg.gen_coefficient_name(g, (), (0,), None)
00323         print facg.gen_coefficient_name(g, (0,), (1,), None)
00324         print
00325 
00326         T = TensorElement("CG", cell, 1)
00327         w = Argument(T, count=1)
00328         h = Coefficient(T, count=2)
00329         print
00330         print facg.gen_argument_name(w, (), (0,1), None)
00331         print facg.gen_argument_name(w, (0,), (1,0), None)
00332         print facg.gen_coefficient_name(h, (), (0,1), None)
00333         print facg.gen_coefficient_name(h, (0,), (1,0), None)
00334         print
00335 
00336         print
00337 
00338 if 1:
00339     class UflacsFormArgumentTest(TempDirTestBase, FormArgumentTestMixin):
00340         def __init__(self, *args, **kwargs):
00341             TempDirTestBase.__init__(self, *args, **kwargs)
00342 
00343         def setUp(self):
00344             TempDirTestBase.setUp(self)
00345             self.options.code.integral.integration_method = "quadrature"
00346             self.set_uflacs_options()
00347 
00348         def set_uflacs_options(self):
00349             try:
00350                 import uflacs
00351                 self.options.code.integral.use_uflacs = True
00352             except:
00353                 print "UFLACS NOT AVAILABLE"
00354 
00355 if 0:
00356     class SfcQuadratureFormArgumentTest(TempDirTestBase, FormArgumentTestMixin):
00357         def __init__(self, *args, **kwargs):
00358             TempDirTestBase.__init__(self, *args, **kwargs)
00359 
00360         def setUp(self):
00361             TempDirTestBase.setUp(self)
00362             self.options.code.integral.integration_method = "quadrature"
00363 
00364 if 0:
00365     class SfcSymbolicFormArgumentTest(TempDirTestBase, FormArgumentTestMixin):
00366         def __init__(self, *args, **kwargs):
00367             TempDirTestBase.__init__(self, *args, **kwargs)
00368 
00369         def setUp(self):
00370             TempDirTestBase.setUp(self)
00371             self.options.code.integral.integration_method = "symbolic"
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Defines