#################################################################
#
#                   Urysohn Lemma
#
#            (Miguel A. Lerma, 7/19/95)
#
# Examples for Urysohn Lemma. Given two sets A and B of closed 
# real intervals, this program computes a function f: R -> [0,1] 
# which takes value 0 on the intervals in A and value 1 on the 
# intervals in B. The function is computed as f(x) = d1/(d1+d2), 
# where d1 = minimum distance from x to the intervals in A, and 
# d2 = minimum distance from x to the intervals in B.
#
#################################################################


# Next function removes every empty interval [a,b], where b < a.

r_empty := (a,b) -> if b < a then NULL else [a,b] fi:


# Next function removes empty intervals from a set.

remove_empty := S -> map(each_interval -> r_empty(op(each_interval)),S):


# Distance from a point x to a closed interval [a,b].

dist := (x,a,b) -> if b < a then infinity  
                   elif x < a then a-x     
                   elif b < x then x-b 
                   else 0 
                   fi:


# Distance from a point to the union of a set of closed intervals.

d := proc(x,S) local T: 
        g_point := x: # a "global" point to be used in the anonymous function
        T := remove_empty(S):
        min(op(map(each_interval -> dist(g_point,op(each_interval)),T))):
     end:


# Distance between two closed sets.

dis := (A,B) -> 
         min(op(map(x->d(x,B),map(op,A))),op(map(x->d(x,A),map(op,B)))):


# Urysohn function.

f := proc(x,A,B) local d1,d2: 
        d1 := d(x,A): 
        d2 := d(x,B): 
        d1/(d1+d2): 
     end:


# Next function computes a suitable range containing all 
# the closed intervals in a set of such intervals.

suitable_range := proc(C) 
                     local lis, a, b, c:
                     lis := map(op,C):
                     a := min(op(lis)):
                     b := max(op(lis)):
                     c := (b-a)/4:
                     a-c..b+c:
                  end:
          

# Plot function f using the range computed above. The intervals in 
# A and B are also plotted as horizontal lines over the X-axis.

plot_urysohn_function := 
proc(A,B) local f1,f2,f3:
   g_A := remove_empty(A):      # Global variables to be used 
   g_B := remove_empty(B):      # in functions f1, f2 and f3.
   if dis(g_A,g_B)=0 then ERROR(`The sets are not disjoint!`) 
   else 
      f1 := x -> f(x,g_A,g_B):
      f2 := x -> if d(x,g_A)=0 then -0.005 else NULL fi:
      f3 := x -> if d(x,g_B)=0 then -0.005 else NULL fi:
      plot({f3,f2,f1},suitable_range([op(g_A),op(g_B)]),-.01..1.01,axes=BOXED):
   fi:
end: