smbinterp/test/baker2dorder.py

#!/usr/bin/env python

import unittest

from interp import baker
from interp.grid import grid
import numpy as np

from interp.grid.simplex import contains

def exact_func(point):
  x = point[0]
  y = point[1]
  return 0.5 + x*x + y

def calculate_error_term(self, a,b,c,d,e,f):
    B = np.array([
          self.p1[a] * self.p1[b], self.p1[c] * self.p1[d], self.p1[e] * self.p1[f],
          self.p2[a] * self.p2[b], self.p2[c] * self.p2[d], self.p2[e] * self.p2[f],
          self.p3[a] * self.p3[b], self.p3[c] * self.p3[d], self.p3[e] * self.p3[f],
          self.p4[a] * self.p4[b], self.p4[c] * self.p4[d], self.p4[e] * self.p4[f],
        ])
    B.shape = (4,3)

    A = np.dot(B.T, B)
    rhs = np.dot(B.T, self.w)
    abc = np.linalg.solve(A,rhs)

    err = \
      abc[0] * self.phis[a] * self.phis[b] + \
      abc[1] * self.phis[c] * self.phis[d] + \
      abc[2] * self.phis[e] * self.phis[f]
    return err

class Test(unittest.TestCase):
  def setUp(self):
    self.verts = [
               [ 2, 3], # 0
               [ 7, 4], # 1
               [ 4, 8], # 2
               [ 0, 7], # 3, 1
               [ 5, 0], # 4, 2
               [10, 5], # 5, 3
               [ 8, 9], # 6, 4
             ]


    self.q = [exact_func(v) for v in self.verts]

    self.g = grid(self.verts,     self.q)
    self.R = grid(self.verts[:3], self.q[:3])
    self.S = grid(self.verts[3:], self.q[3:])

    self.p1, self.ql1 = baker.qlinear(self.verts[3], self.R)
    self.p2, self.ql2 = baker.qlinear(self.verts[4], self.R)
    self.p3, self.ql3 = baker.qlinear(self.verts[5], self.R)
    self.p4, self.ql4 = baker.qlinear(self.verts[6], self.R)

    self.q1 = exact_func(self.verts[3])
    self.q2 = exact_func(self.verts[4])
    self.q3 = exact_func(self.verts[5])
    self.q4 = exact_func(self.verts[6])


    self.w = np.array([
      self.q1 - self.ql1,
      self.q2 - self.ql2,
      self.q3 - self.ql3,
      self.q4 - self.ql4,
    ])

    self.X = [4,5]

    self.g = grid(self.verts, self.q)

    self.phis, self.qlin  = baker.qlinear(self.X, self.R)
    self.exact = exact_func(self.X)
    self.answer = baker.run_baker(self.X,self.R,self.S)


  # def test_R_contains_X(self):
  #   self.assertTrue(contains(self.X, self.R.verts))

  def test_1(self):
    a,b,c,d,e,f = (0,1,   1,2,   2,0)
    err = calculate_error_term(self, a,b,c,d,e,f)
    self.assertAlmostEqual(err, self.answer['error'])
  def test_swap_first_elements(self):
    a,b,c,d,e,f = (1,0,   1,2,   2,0)
    err = calculate_error_term(self, a,b,c,d,e,f)
    self.assertAlmostEqual(err, self.answer['error'])
  def test_swap_two_pairs(self):
    a,b,c,d,e,f = (1,2,   0,1,   2,0)
    err = calculate_error_term(self, a,b,c,d,e,f)
    self.assertAlmostEqual(err, self.answer['error'])
  def test_swap_all_pairs(self):
    a,b,c,d,e,f = (0,2,   0,1,   2,1)
    err = calculate_error_term(self, a,b,c,d,e,f)
    self.assertAlmostEqual(err, self.answer['error'])


if __name__ == '__main__':
  suite = unittest.TestLoader().loadTestsFromTestCase(Test)
  unittest.TextTestRunner(verbosity=3).run(suite)