SHARE
TWEET

Base Functions on arbitrary Reference Element

a guest Jul 9th, 2015 163 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. dimension, order = 3, 2
  2. from numpy.linalg import inv
  3.  
  4. degreesOfFreedom = []
  5. monomials = []
  6.  
  7. if dimension == 2:
  8.         for i in range(order+1):
  9.                 for j in range(order+1):
  10.                         if i + j > order:
  11.                                 break
  12.                         t = map(lambda x: x/float(order), [i,j] )
  13.                         degreesOfFreedom.append( tuple(t) )
  14.                         monomials.append( (i,j) )
  15.  
  16. elif dimension == 3:
  17.         for i in range(order+1):
  18.                 for j in range(order+1):
  19.                         for k in range(order+1):
  20.                                 if i + j + k > order:
  21.                                         break
  22.                                 t = map(lambda x: x/float(order), [i,j,k] )
  23.                                 degreesOfFreedom.append( tuple(t) )
  24.                                 monomials.append( (i,j,k) )
  25.  
  26.  
  27. referenceBaseFunctions = []
  28.  
  29. if dimension == 2:
  30.         for i in range( len(degreesOfFreedom) ):
  31.                 matrix = []
  32.                 for (x,y) in degreesOfFreedom:
  33.                         row = []
  34.                         for (ex,ey) in monomials:
  35.                                 row.append( pow(x,ex) * pow(y,ey) )
  36.                         matrix.append(row)
  37.                 inverse = inv(matrix)
  38.                 column = list( inverse[:,i] )       # inverse * vector
  39.                 def baseFunc(x,y):
  40.                         s = 0
  41.                         for i in range( len(monomials) ):
  42.                                 C = column[i]
  43.                                 (ex,ey) = monomials[i]
  44.                                 s += C * (pow(x,ex) * pow(y,ey))
  45.                         return s
  46.                 referenceBaseFunctions.append(baseFunc)
  47.  
  48. elif dimension == 3:
  49.         for i in range( len(degreesOfFreedom) ):
  50.                 matrix = []
  51.                 for (x,y,z) in degreesOfFreedom:
  52.                         row = []
  53.                         for (ex,ey,ez) in monomials:
  54.                                 row.append( pow(x,ex) * pow(y,ey) * pow(z,ez) )
  55.                         matrix.append(row)
  56.                 inverse = inv(matrix)
  57.                 column = list( inverse[:,i] )       # inverse * vector
  58.                 def baseFunc(x,y,z):
  59.                         s = 0
  60.                         for i in range( len(monomials) ):
  61.                                 C = column[i]
  62.                                 (ex,ey,ez) = monomials[i]
  63.                                 s += C * (pow(x,ex) * pow(y,ey) * pow(z,ez))
  64.                         return s
  65.                 referenceBaseFunctions.append(baseFunc)
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top