From 33106985e131283df7024a6287701eab86e35904 Mon Sep 17 00:00:00 2001 From: Stephen Mardson McQuay Date: Thu, 3 Feb 2011 09:30:32 -0700 Subject: [PATCH] working test, moving into a real test --- bin/baker_order_test.py | 84 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 84 insertions(+) diff --git a/bin/baker_order_test.py b/bin/baker_order_test.py index 98e7aac..b78253f 100755 --- a/bin/baker_order_test.py +++ b/bin/baker_order_test.py @@ -1,5 +1,89 @@ #!/usr/bin/env python from interp import baker +from interp.grid import grid +import numpy as np + +def exact_func(point): + return 0.5 + point[0] * point[1] + +# a,b,c,d,e,f = (0,1, 1,2, 2,0) +# a,b,c,d,e,f = (2,1, 1,0, 2,0) +# a,b,c,d,e,f = (0,1, 1,2, 2,0) +a,b,c,d,e,f = (2,0, 1,2, 1,0) + +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 + ] +q = [exact_func(p) for p in verts] + +g = grid(verts,q) +R = grid(verts[:3], q[:3]) +S = grid(verts[3:], q[3:]) + +print 'g', g +print 'R', R +print 'S', S + + +X = [ 4,5] +print 'verts', verts +print 'q', q +print 'X', X + +phis,qlin = baker.qlinear(X, R) + +p1, ql1 = baker.qlinear(verts[3], R) +p2, ql2 = baker.qlinear(verts[4], R) +p3, ql3 = baker.qlinear(verts[5], R) +p4, ql4 = baker.qlinear(verts[6], R) + +B = np.array([ + p1[a] * p1[b], p1[c] * p1[d], p1[e] * p1[f], + p2[a] * p2[b], p2[c] * p2[d], p2[e] * p2[f], + p3[a] * p3[b], p3[c] * p3[d], p3[e] * p3[f], + p4[a] * p4[b], p4[c] * p4[d], p4[e] * p4[f], + ]) +B.shape = (4,3) + +q1 = exact_func(verts[3]) +q2 = exact_func(verts[4]) +q3 = exact_func(verts[5]) +q4 = exact_func(verts[6]) + +w = np.array([ + q1 - ql1, + q2 - ql2, + q3 - ql3, + q4 - ql4, + ]) + +print 'WTF ----> q1', q1 +print 'WTF ----> ql1', ql1 +print 'WTF ----> p1', p1 +print 'WTF ----> w', w + +A = np.dot(B.T, B) +rhs = np.dot(B.T, w) +abc = np.linalg.solve(A,rhs) + +print A, b, abc +err = \ + abc[0] * phis[a] * phis[b] + \ + abc[1] * phis[c] * phis[d] + \ + abc[2] * phis[e] * phis[f] + + +answer = baker.run_baker(X,R,S) +print answer['error'] +print answer['final'], exact_func(X) + +np.testing.assert_approx_equal(err, answer['error'])