> # Aufgabe 14
> #Differenzenverfahren
> restart;
> with(linalg): with(plots):
> DiffDirichlet:=proc(a,b,g1,g2,f,c1,c2,N)
> # Löse -y''(x)+g1(x)*y'(x)+g2(x)*y(x) auf Intervall (a,b)
> # Randbedingungen y(a)=c1,y(b)=c2,
> # Diskretisierung in N Teilintervalle
> # Verwendung zentraler Differenzenquotienten
> local A,j,h,x,r,y,z;
> global p;
> # option trace;
> #Aufstellen der Matrix
> h:=(b-a)/N;
> A:=matrix(N-1,N-1,0);
> for j from 1 to N-1 do 
>  x:=a+j*h; A[j,j]:=2/h^2+g2(x)
> od;
> for j from 1 to N-2 do 
>  x:=a+j*h; A[j,j+1]:=-1/h^2+g1(x)/(2*h); A[j+1,j]:=-1/h^2-g1(x+h)/(2*h)
> od;
> # Aufstellen der rechten Seite
> r:=vector(N-1); x:=a+h;
> r[1]:=f(x)+c1/h^2+c1*g1(x)/(2*h);
> for j from 2 to N-2 do x:=x+h; r[j]:=f(x) od;
> x:=b-h; r[N-1]:=f(x)+c2/h^2-c2*b1(x)/(2*h);
> # Berechnung des LGS
> y:=vector(N-1);
> y:=linsolve(A,r);
> #graphische Ausgabe
> p:=vector(N+1); z:=vector(N+1);
> p[1]:=c1; p[N+1]:=c2;
> for j from 2 to N do p[j]:=y[j-1] od;
> for j from 1 to N+1 do z[j]:=a+(j-1)*h od;
> listplot([seq([z[j],p[j]],j=1..N+1)],thickness=0);
> end:
> #Beispiel 1
> # -y''(x)+y(x)=1 in (0,10), y(0)=y(10)=0
> g1:=proc(x)  RETURN(0.0) end:
> g2:=proc(x)  RETURN(1.0) end:
> f:=proc(x)  RETURN(1.0) end:
> N:=20:
> DiffDirichlet(0.0,10.0,g1,g2,f,0.0,0.0,N);
> #exakte Lösung: 
> y:=proc(x) RETURN(evalf((exp(10)-exp(x))*(1-exp(-x))/(exp(10)+1))) end:
> plot(y(x), x=0..10);
> #Fehlerbetrachtung:
> #print(p);
> fehler:=vector(N+1):
> for j from 1 to N+1 do fehler[j]:=evalf(abs(p[j]-y(10*(j-1)/N))) od;
> 
> 



                           fehler[1] := 0.


                       fehler[2] := .0030856095


                       fehler[3] := .0037569159


                       fehler[4] := .0034348140


                       fehler[5] := .0027985362


                       fehler[6] := .0021493924


                       fehler[7] := .0016036052


                       fehler[8] := .0011924655


                       fehler[9] := .0009129585


                      fehler[10] := .0007524851


                      fehler[11] := .0007003295


                      fehler[12] := .0007524849


                      fehler[13] := .0009129584


                      fehler[14] := .0011924654


                      fehler[15] := .0016036052


                      fehler[16] := .0021493924


                      fehler[17] := .0027985362


                      fehler[18] := .0034348140


                      fehler[19] := .0037569159


                      fehler[20] := .0030856093


                           fehler[21] := 0.

>