visualize solution on boundary

Question

Hello
I have a scalar solution "p" and let

n = unit normal to the boundary

I want to visualize the vector

p*n

on the boundary. Any ideas how I can do this ?

Comments

For unit square I can do this

# If x is on boundary return unit normal, else return zero vector
class UnitNormal(Expression):
   def eval(self, value, x):
       if near(x[0],0.0):
           value[0] = -1.0
           value[1] =  0.0
       elif near(x[0],1.0):
           value[0] =  1.0
           value[1] =  0.0
       elif near(x[1],0.0):
           value[0] =  0.0
           value[1] = -1.0
       elif near(x[1],1.0):
           value[0] =  0.0
           value[1] =  1.0
       else:
           value[0] =  0.0
           value[1] =  0.0
       return value
   def value_shape(self):
      return (2,)

n = UnitNormal()
sigma = project(-p*n, V)
File("ctrb_u.pvd") << sigma

How can I write a function like UnitNormal for a general domain ?