#hw23.txt Cole Franks 20/04/2014 Help23:=proc() print(`Pw(Lf,x,Lx), Compare(p,q,v,x,x1),BSpline(d,Lx, i),splii(Lx,x,i), TestF3(c,Lx,x), EvalFunct3(p,q,f,x,Lu), RR3(p,q,f,x,Lx)`): end: ##### Problem 1 Pw:=proc(Lf,x,Lx) local f, i: f:=piecewise(Lx[1] nops(Lx) - 2 then return(FAIL): fi: Lu:=[seq(add(a[i1, j1]*x^j1, j1=0..d), i1=0..d)]: var:={seq(seq(a[i1, j1], j1=0..d), i1=0..d)}: eq:={seq(subs(x = Lx[1], diff(Lu[1], [x\$i1])), i1=0..d-1)}: eq:=eq union {seq(subs(x = Lx[d+2], diff(Lu[d+1], [x\$i1])), i1=0..d-1)}: #loop for the spots for j from 2 to d+1 do eq:=eq union {seq(subs(x = Lx[j], diff(Lu[j-1], [x\$i1])) = subs(x = Lx[j], diff(Lu[j], [x\$i1])), i1=0..d-1)}: od: eq:=eq union {subs(x = Lx[i+1], Lu[i]) = 1}: var:=solve(eq, var): if var=NULL then FAIL: else return subs(var, Lu): fi: end: ##### Problem 4 splii:=proc(Lx,x,i) local f, c, Lshift: if 2 < i and i <= nops(Lx) - 2 then Pw(BSpline(3, [seq(Lx[j], j =i-2..i+2)], 2), x, [seq(Lx[j], j =i-2..i +2)]): elif i = 1 then #Lshift is the new list of points for this spline Lshift:=[-Lx[3], -Lx[2], Lx[1], Lx[2], Lx[3]]: f:=Pw(BSpline(3, Lshift, 2), x, Lshift): #then subtract off the right amount to kill it at 0 c:=eval(f, x=0)/eval(f, x=Lx[2]): convert(f - c*subs(x = x+Lx[2], f), piecewise, x): elif i = 2 then Lshift:=[-Lx[2], Lx[1], Lx[2], Lx[3], Lx[4]]: f:=Pw(BSpline(3, Lshift, 2), x, Lshift): c:=eval(f, x=0)/eval(f, x=Lx[3]): convert(f - c*subs(x = x+Lx[3], f), piecewise, x): elif i = nops(Lx) then Lshift:=[Lx[-3], Lx[-2], Lx[-1], 2*Lx[-1] - Lx[-2], 2*Lx[-1] - Lx[-3]]: #print(Lshift): f:=Pw(BSpline(3, Lshift, 2), x, Lshift): c:=eval(f, x=Lx[-1])/eval(f, x=Lx[-2]): convert(f - c*subs(x = x-(Lx[-1]-Lx[-2]), f), piecewise, x): elif i = nops(Lx)-1 then Lshift:=[Lx[-4], Lx[-3], Lx[-2], Lx[-1], 2*Lx[-1] - Lx[-2]]: f:=Pw(BSpline(3, Lshift, 2), x, Lshift): c:=eval(f, x=Lx[-1])/eval(f, x=Lx[-3]): convert(f - c*subs(x = x-(Lx[-1]-Lx[-3]), f), piecewise, x): fi: end: TestF3:=proc(c,Lx,x) local i, f: f:=convert(add(c[i]*splii(Lx,x,i),i=1..nops(Lx)), piecewise, x): convert(piecewise(x<0, 0, x>=0 and x < Lx[-1], f, x >= Lx[-1] ,0), piecewise, x), { seq( c[i], i=1..nops(Lx))}: end: EvalFunct3:=proc(p,q,f,x,Lu) local F,Lu1,i: #now Lu is a piecewise function Lu1:=diff(Lu,x): F:=expand(p*Lu1*Lu1+ q*Lu*Lu-2*f*Lu): int(F, x=0..1): end: RR3:=proc(p,q,f,x,Lx) local c, Lu, F, eq,var, y: Lu:=TestF3(c,Lx,x): var:=Lu[2]: Lu:=Lu[1]: F:=evalf(EvalFunct3(p,q,f,x,Lu)): eq:=evalf({seq(diff(F,y),y in var)}): var:=solve(eq,var): if var=NULL then FAIL: else subs(var,Lu): fi: end: