Package sfc :: Package codegeneration :: Module exteriorfacetintegralcg
[hide private]
[frames] | no frames]

Source Code for Module sfc.codegeneration.exteriorfacetintegralcg

  1  # -*- coding: utf-8 -*- 
  2  """ 
  3  This module contains code generation tools for the ufc::exterior_facet_integral class. 
  4  """ 
  5   
  6  # Copyright (C) 2008-2009 Martin Sandve Alnes and Simula Resarch Laboratory 
  7  # 
  8  # This file is part of SyFi. 
  9  # 
 10  # SyFi is free software: you can redistribute it and/or modify 
 11  # it under the terms of the GNU General Public License as published by 
 12  # the Free Software Foundation, either version 2 of the License, or 
 13  # (at your option) any later version. 
 14  # 
 15  # SyFi is distributed in the hope that it will be useful, 
 16  # but WITHOUT ANY WARRANTY; without even the implied warranty of 
 17  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 
 18  # GNU General Public License for more details. 
 19  # 
 20  # You should have received a copy of the GNU General Public License 
 21  # along with SyFi. If not, see <http://www.gnu.org/licenses/>. 
 22  # 
 23  # First added:  2008-08-13 
 24  # Last changed: 2009-03-09 
 25   
 26  from itertools import tee 
 27  import swiginac 
 28  from sfc.common.output import sfc_debug, sfc_error, sfc_assert, sfc_warning 
 29  from sfc.quadrature import gen_quadrature_rule_definition, gen_quadrature_rule_definitions 
 30  from sfc.codegeneration.codeformatting import indent, CodeFormatter, gen_const_token_definitions, gen_token_prints 
 31  from sfc.codegeneration.codeformatting import gen_symbol_declarations, gen_token_declarations, gen_token_assignments, gen_token_additions 
 32  from sfc.codegeneration.integralcg import IntegralCG 
 33   
34 -class ExteriorFacetIntegralCG(IntegralCG):
35 - def __init__(self, itgrep):
36 IntegralCG.__init__(self, itgrep) 37 # Generator trick: 38 # Iterators over members are used in gen_members and 39 # gen_constructor, avoid double iteration with tee: 40 41 # FIXME: Update this for multiple integral data, and for precomputation without constructor code 42 for integral in self.itgrep.quadrature_integrals: 43 sfc_assert(len(self.itgrep.quadrature_integrals) <= 1, "Not implemented!") 44 data = self.itgrep.integral_data[integral.measure()] 45 46 a, b = tee(self.itgrep.iter_member_quad_tokens(data)) 47 self._iter_member_quad_tokens1 = a 48 self._iter_member_quad_tokens2 = b
49
50 - def gen_members(self):
51 sfc_debug("entering ExteriorFacetIntegral.gen_members") 52 code = "" 53 54 # FIXME: Update this for multiple integral data, and for precomputation without constructor code 55 for integral in self.itgrep.quadrature_integrals: 56 sfc_assert(len(self.itgrep.quadrature_integrals) <= 1, "Not implemented!") 57 data = self.itgrep.integral_data[integral.measure()] 58 59 member_quad_tokens = self._iter_member_quad_tokens1 #self.itgrep.iter_member_quad_tokens() 60 member_quad_token_names = (str(s[0]).split(".")[-1] for s in member_quad_tokens) 61 decl = gen_symbol_declarations(member_quad_token_names) 62 if decl: 63 code += "\n\n" 64 code += "struct member_quad_tokens\n{\n" 65 code += indent(decl) 66 code += "\n};\n" 67 code += "member_quad_tokens qt[%d][%d];\n" % (self.itgrep.formrep.cell.num_facets, data.facet_quad_rules[0].num_points) 68 69 sfc_debug("leaving ExteriorFacetIntegral.gen_members") 70 return code
71
72 - def gen_constructor(self):
73 sfc_debug("entering ExteriorFacetIntegral.gen_constructor") 74 code = "" 75 76 # FIXME: Update this for multiple integral data, and for precomputation without constructor code 77 for integral in self.itgrep.quadrature_integrals: 78 sfc_assert(len(self.itgrep.quadrature_integrals) <= 1, "Not implemented!") 79 data = self.itgrep.integral_data[integral.measure()] 80 81 member_quad_tokens = self._iter_member_quad_tokens2 # itgrep.iter_member_quad_tokens() 82 member_quad_tokens_assignments = gen_token_assignments(member_quad_tokens) 83 if member_quad_tokens_assignments: 84 # Construct loop over facets with quad loop inside 85 facet_loop_begin = "for(int facet=0; facet<%d; facet++)\n{\n" % self.itgrep.formrep.cell.num_facets 86 quad_loop_start_code = self.gen_quadrature_begin_block(data) 87 88 # Stitch code pieces together 89 code += facet_loop_begin 90 code += quad_loop_start_code 91 code += indent(member_quad_tokens_assignments) 92 code += "\n}\n" 93 code += "\n}\n" 94 95 sfc_debug("leaving ExteriorFacetIntegral.gen_constructor") 96 return code
97
98 - def gen_destructor(self):
99 return ""
100