Source Code for Module hedge.models.nd_calculus

  1  # -*- coding: utf8 -*- 
  2  """Canned operators for multivariable calculus.""" 
  3   
  4  from __future__ import division 
  5   
  6  __copyright__ = "Copyright (C) 2009 Andreas Kloeckner" 
  7   
  8  __license__ = """ 
  9  This program is free software: you can redistribute it and/or modify 
 10  it under the terms of the GNU General Public License as published by 
 11  the Free Software Foundation, either version 3 of the License, or 
 12  (at your option) any later version. 
 13   
 14  This program is distributed in the hope that it will be useful, 
 15  but WITHOUT ANY WARRANTY; without even the implied warranty of 
 16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 
 17  GNU General Public License for more details. 
 18   
 19  You should have received a copy of the GNU General Public License 
 20  along with this program.  If not, see U{http://www.gnu.org/licenses/}. 
 21  """ 
 22   
 23   
 24   
 25   
 26  from hedge.models import Operator 
 27   
 28   
 29   
 30   
 31 -class GradientOperator(Operator): 
 32 -    def __init__(self, dimensions): 
 33          self.dimensions = dimensions 
 34   
 35 -    def flux(self): 
 36          from hedge.flux import make_normal, FluxScalarPlaceholder 
 37          u = FluxScalarPlaceholder() 
 38   
 39          normal = make_normal(self.dimensions) 
 40          return u.int*normal - u.avg*normal 
 41   
 42 -    def op_template(self): 
 43          from hedge.mesh import TAG_ALL 
 44          from hedge.optemplate import Field, pair_with_boundary, \ 
 45                  make_nabla, InverseMassOperator, get_flux_operator 
 46   
 47          u = Field("u") 
 48          bc = Field("bc") 
 49   
 50          nabla = make_nabla(self.dimensions) 
 51          flux_op = get_flux_operator(self.flux()) 
 52   
 53          return nabla*u - InverseMassOperator()*( 
 54                  flux_op * u + 
 55                  flux_op * pair_with_boundary(u, bc, TAG_ALL) 
 56                  ) 
 57   
 58 -    def bind(self, discr): 
 59          compiled_op_template = discr.compile(self.op_template()) 
 60   
 61          def op(u): 
 62              from hedge.mesh import TAG_ALL 
 63   
 64              return compiled_op_template(u=u, 
 65                      bc=discr.boundarize_volume_field(u, TAG_ALL)) 
 66   
 67          return op 
 68   
 69   
 70   
 71   
 72 -class DivergenceOperator(Operator): 
 73 -    def __init__(self, dimensions, subset=None): 
 74          self.dimensions = dimensions 
 75   
 76          if subset is None: 
 77              self.subset = dimensions * [True,] 
 78          else: 
 79              # chop off any extra dimensions 
 80              self.subset = subset[:dimensions] 
 81   
 82          from hedge.tools import count_subset 
 83          self.arg_count = count_subset(self.subset) 
 84   
 85 -    def flux(self): 
 86          from hedge.flux import make_normal, FluxVectorPlaceholder 
 87   
 88          v = FluxVectorPlaceholder(self.arg_count) 
 89   
 90          normal = make_normal(self.dimensions) 
 91   
 92          flux = 0 
 93          idx = 0 
 94   
 95          for i, i_enabled in enumerate(self.subset): 
 96              if i_enabled and i < self.dimensions: 
 97                  flux += (v.int-v.avg)[idx]*normal[i] 
 98                  idx += 1 
 99   
100          return flux 
101   
102 -    def op_template(self): 
103          from hedge.mesh import TAG_ALL 
104          from hedge.optemplate import make_vector_field, pair_with_boundary, \ 
105                  get_flux_operator, make_nabla, InverseMassOperator 
106   
107          nabla = make_nabla(self.dimensions) 
108          m_inv = InverseMassOperator() 
109   
110          v = make_vector_field("v", self.arg_count) 
111          bc = make_vector_field("bc", self.arg_count) 
112   
113          local_op_result = 0 
114          idx = 0 
115          for i, i_enabled in enumerate(self.subset): 
116              if i_enabled and i < self.dimensions: 
117                  local_op_result += nabla[i]*v[idx] 
118                  idx += 1 
119   
120          flux_op = get_flux_operator(self.flux()) 
121   
122          return local_op_result - m_inv*( 
123                  flux_op * v + 
124                  flux_op * pair_with_boundary(v, bc, TAG_ALL)) 
125   
126 -    def bind(self, discr): 
127          compiled_op_template = discr.compile(self.op_template()) 
128   
129          def op(v): 
130              from hedge.mesh import TAG_ALL 
131              return compiled_op_template(v=v, 
132                      bc=discr.boundarize_volume_field(v, TAG_ALL)) 
133   
134          return op 
135