Home My Page Projects Code Snippets Project Openings diderot
Summary Activity Tracker Tasks SCM

SCM Repository

[diderot] Diff of /branches/ein16/synth/d2/test_eval.py
ViewVC logotype

Diff of /branches/ein16/synth/d2/test_eval.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 3946, Sat Jun 11 00:39:19 2016 UTC revision 3998, Sun Jun 19 17:12:03 2016 UTC
# Line 19  Line 19 
19  def fn_subtract(exp1,exp2):  def fn_subtract(exp1,exp2):
20      return exp1-exp2      return exp1-exp2
21  # scaling operator  # scaling operator
22  def fn_scaling(fld1, fld2):  def fn_multiplication(exp_s, t):
     def get_sca():  
         if(field.is_Scalar(fld1)):  
             return (fld1, fld2)  
         else:  
             return (fld2, fld1)  
     (s, t) = get_sca()  
     exp_s = field.get_data(s)  
23      exp_t = field.get_data(t)      exp_t = field.get_data(t)
24      ityp1 = field.get_ty(t)      ityp1 = field.get_ty(t)
25      shape1 = fty.get_shape(ityp1)      shape1 = fty.get_shape(ityp1)
# Line 63  Line 56 
56      else:      else:
57          raise "unsupported scaling"          raise "unsupported scaling"
58    
59    #scaling of a field
60    def fn_scaling(fld1, fld2):
61        def get_sca():
62            if(field.is_Scalar(fld1)):
63                return (fld1, fld2)
64            else:
65                return (fld2, fld1)
66        (s, t) = get_sca()
67        exp_s = field.get_data(s)
68        return fn_multiplication(exp_s, t)
69    
70    #division of a field
71    def fn_division(t, s):
72        if(field.is_Scalar(s)):
73            exp_s = (1/field.get_data(s))
74            return fn_multiplication(exp_s, t)
75        else:
76            raise Exception ("err second arg in division should be a scalar")
77    
78  # negation  of field  # negation  of field
79  def fn_negation(exp):  def fn_negation(exp):
80    
# Line 151  Line 163 
163                   [diff(exp[1],x), diff(exp[1],y), diff(exp[1],z)],                   [diff(exp[1],x), diff(exp[1],y), diff(exp[1],z)],
164                   [diff(exp[2],x), diff(exp[2],y), diff(exp[2],z)]]                   [diff(exp[2],x), diff(exp[2],y), diff(exp[2],z)]]
165      else:      else:
166          raise "unsupported type for divergence"          raise "unsupported type for jacob"
167    
168    #evaluate norm
169    def fn_norm(fld, dim):
170        exp = field.get_data(fld)
171        ityp = field.get_ty(fld)
172        dim = field.get_dim(fld)
173        #print " exp1: ",exp1," exp2: ",exp2
174        # vectors
175        def iter (es):
176            sum = 0
177            for i in es:
178                sum+=i*i
179            return sqrt(sum)
180        if(field.is_Scalar(fld)):
181            [] = fty.get_shape(ityp)
182            return exp
183        elif(field.is_Vector(fld)):
184            [n] = fty.get_shape(ityp)
185            rtn = []
186            for i in range(n):
187                rtn.append(exp[i])
188            return iter(rtn)
189        elif(field.is_Matrix(fld)):
190            [n, m] = fty.get_shape(ityp)
191            rtn = []
192            for i in range(n):
193                for j in range(m):
194                    rtn.append(exp[i][j])
195            return iter(rtn)
196        elif(field.is_Ten3(fld)):
197            [n, m, p] = fty.get_shape(ityp)
198            rtn = []
199            for i in range(n):
200                for j in range(m):
201                    for k in range(p):
202                        rtn.append(exp[i][j][k])
203            return iter(rtn)
204        else:
205            raise "unsupported type for norm"
206    
207    #evaluate norm
208    def fn_normalize(fld, dim):
209        exp = field.get_data(fld)
210        ityp = field.get_ty(fld)
211        dim = field.get_dim(fld)
212        #print " exp1: ",exp1," exp2: ",exp2
213        norm = fn_norm(fld, dim)
214        if(field.is_Scalar(fld)):
215            #print "scal",exp
216            return exp
217        elif(field.is_Vector(fld)):
218            [n] = fty.get_shape(ityp)
219            rtn = []
220            for i in range(n):
221                rtn.append(exp[i]/norm)
222            #print "vec",rtn
223            return rtn
224        elif(field.is_Matrix(fld)):
225            [n, m] = fty.get_shape(ityp)
226            rtn = []
227            for i in range(n):
228                rtni = []
229                for j in range(m):
230                    rtni.append(exp[i][j]/norm)
231                rtn.append(rtni)
232                #print "matrix:",rtn
233            return rtn
234        elif(field.is_Ten3(fld)):
235            [n, m, p] = fty.get_shape(ityp)
236            rtn = []
237            for i in range(n):
238                rtni = []
239                for j in range(m):
240                    rtnj = []
241                    for k in range(p):
242                        rtnj.append(exp[i][j][k]/norm)
243                    rtni.append( rtnj)
244                rtn.append( rtni)
245    #print "ten3",rtn
246            return rtn
247        else:
248            raise "unsupported type for norm"
249    
250    #evaluate slice
251    def fn_slice(fld1):
252        exp1 = field.get_data(fld1)
253        ityp1 = field.get_ty(fld1)
254        print exp1
255        rtn=[]
256        if(fty.is_Matrix(ityp1)):
257            [n2,n3] = fty.get_shape(ityp1)
258            for j in range(n3):
259                rtn.append(exp1[j][0])
260            print rtn
261            return rtn
262        else:
263            raise "unsupported type for slice"
264    
265    #evaluate trace
266    def fn_trace(fld):
267        exp = field.get_data(fld)
268        ityp = field.get_ty(fld)
269        rtn=[]
270        if(field.is_Matrix(fld)):
271            [n, m] = fty.get_shape(ityp)
272            if (n!=m):
273                raise Exception("matrix is not identitical")
274            rtn = exp[0][0]+exp[1][1]
275            if(n==2):
276                return rtn
277            elif(n==3):
278                return rtn+exp[2][2]
279            else:
280                raise "unsupported matrix field"
281        else:
282            raise "unsupported trace"
283    
284    #evaluate transpose
285    def fn_transpose(fld):
286        exp = field.get_data(fld)
287        ityp = field.get_ty(fld)
288        if(field.is_Matrix(fld)):
289            [n, m] = fty.get_shape(ityp)
290            rtn = []
291            for i in range(n):
292                rtni = []
293                for j in range(m):
294                    rtni.append(exp[j][i])
295                rtn.append(rtni)
296            return rtn
297        else:
298            raise "unsupported transpose"
299    
300    #evaluate det
301    def fn_det(fld):
302        exp = field.get_data(fld)
303        ityp = field.get_ty(fld)
304        rtn=[]
305        if(field.is_Matrix(fld)):
306            [n, m] = fty.get_shape(ityp)
307            if (n!=m):
308                raise Exception("matrix is not identitical")
309            a = exp[0][0]
310            d = exp[1][1]
311            c = exp[1][0]
312            b = exp[0][1]
313            if(n==2):
314                return a*d-b*c
315            elif(n==3):
316                a = exp[0][0]
317                b = exp[0][1]
318                c = exp[0][2]
319                d = exp[1][0]
320                e = exp[1][1]
321                f = exp[1][2]
322                g = exp[2][0]
323                h = exp[2][1]
324                i = exp[2][2]
325                return a*(e*i-f*h)-b*(d*i-f*g)+c*(d*h-e*g)
326            else:
327                raise "unsupported matrix field"
328        else:
329            raise "unsupported trace"
330    
331    
332  #evaluate outer product  #evaluate outer product
# Line 373  Line 548 
548      if(op_probe==fn_name): #probing      if(op_probe==fn_name): #probing
549          return exp          return exp
550      elif(op_negation==fn_name): #negation      elif(op_negation==fn_name): #negation
   
551          return applyToM(exp, fn_negation)          return applyToM(exp, fn_negation)
552      elif(op_jacob==fn_name): #jacob      elif(op_jacob==fn_name): #jacob
553          return fn_jacob(fld)          return fn_jacob(fld)
554        elif(op_slice==fn_name):
555            return fn_slice(fld)
556        elif(op_trace == fn_name):
557            return fn_trace(fld)
558        elif(op_transpose==fn_name):
559            return fn_transpose(fld)
560        elif(op_det==fn_name):
561            return fn_det(fld)
562      else:      else:
563          raise Exception("unsupported unary operator:"+ fn_name.name)          raise Exception("unsupported unary operator:"+ fn_name.name)
564    
# Line 406  Line 588 
588          return fn_subtract(exp1,exp2)          return fn_subtract(exp1,exp2)
589      elif(op_scale==fn_name): #scaling      elif(op_scale==fn_name): #scaling
590          return fn_scaling(fld1,fld2)          return fn_scaling(fld1,fld2)
591        elif(op_division==fn_name): #division
592            return fn_division(fld1,fld2)
593      else:      else:
594          raise Exception("unsupported binary operator on scalar fields:"+ fn_name.name)          raise Exception("unsupported binary operator on scalar fields:"+ fn_name.name)
595    
# Line 449  Line 633 
633      fn_name=e.opr      fn_name=e.opr
634      exp = field.get_data(fld)      exp = field.get_data(fld)
635      dim = field.get_dim(fld)      dim = field.get_dim(fld)
636      if (field.is_Scalar(fld)): # input is a scalar field      if(op_norm==fn_name):#norm
637            return fn_norm(fld, dim)
638        if(op_normalize==fn_name):#normalize
639            x= fn_normalize(fld, dim)
640            return x
641        elif (field.is_Scalar(fld)): # input is a scalar field
642          return applyToExp_U_S(fn_name, fld)          return applyToExp_U_S(fn_name, fld)
643      elif(field.is_Vector(fld)): # input is a vector field      elif(field.is_Vector(fld)): # input is a vector field
644          return applyToExp_U_V(fn_name, fld)          return applyToExp_U_V(fn_name, fld)
# Line 460  Line 649 
649    
650  def binary(e):  def binary(e):
651      (f, g) =apply.get_binary(e)      (f, g) =apply.get_binary(e)
652      #print "binary- f: ",str(field.get_data(f)), " g: ",str(field.get_data(g))      fn_name = e.opr
653      if (field.is_Scalar(f) and field.is_Scalar(g)): # input is a scalar field      # type is checked elsewhere or does not matter
654        if(op_division==fn_name): #division
655            return fn_division(f, g)
656        elif (field.is_Scalar(f) and field.is_Scalar(g)): # input is a scalar field
657          return applyToExp_B_S(e)          return applyToExp_B_S(e)
658      elif (field.is_Vector(f)):# input is a vector field      elif (field.is_Vector(f)):# input is a vector field
659          if(field.is_Vector(g)):          if(field.is_Vector(g)):
# Line 472  Line 664 
664           return applyToExp_B_uneven(e)           return applyToExp_B_uneven(e)
665    
666  def applyUnaryOnce(oexp_inner,app_inner,app_outer):  def applyUnaryOnce(oexp_inner,app_inner,app_outer):
667      print "applyUnaryOnce"      #print "applyUnaryOnce"
668      #apply.toStr(app_inner)      #apply.toStr(app_inner)
669      oty_inner = apply.get_oty(app_inner)      oty_inner = apply.get_oty(app_inner)
670      oty_outer = apply.get_oty(app_outer)      oty_outer = apply.get_oty(app_outer)
# Line 554  Line 746 
746  # ***************************  evaluate at positions ***************************  # ***************************  evaluate at positions ***************************
747  #evaluate scalar field exp  #evaluate scalar field exp
748  def eval_d2(pos0, pos1, exp):  def eval_d2(pos0, pos1, exp):
749        #print "exp:",exp
750      #evaluate field defined by coefficients at position      #evaluate field defined by coefficients at position
751      exp0 = exp.subs(x,pos0)      exp0 = exp.subs(x,pos0)
752      exp1 = exp0.subs(y,pos1)      exp1 = exp0.subs(y,pos1)

Legend:
Removed from v.3946  
changed lines
  Added in v.3998

root@smlnj-gforge.cs.uchicago.edu
ViewVC Help
Powered by ViewVC 1.0.0