\334bung 7 - 24.05.04

Length := proc(ll) return nops(ll); end proc: Last := proc(ll) return op(-1, ll); end proc: Append := proc(ll, x) return [op(ll), x]; end proc: Prepend := proc(ll, x) return [x, op(ll)]; end proc: Member := proc(ll, p) return op(p, ll); end proc: Take := proc(ll, n) if n > 0 then return [ op(1..n, ll) ] elif n < 0 then return [ op(n..-1, ll) ] end if; end proc: Total := proc(ll) local i; return sum(Member(ll,i),i=1..Length(ll)); end proc: Count := proc(ll, x) local i,r; r := 0; for i from 1 by 1 to Length(ll) do if Member(ll,i) = x then r := r+1 end if; end do; return r; end proc: Times := proc(l1, l2) local i, r; r := []; if Length(l1) = Length(l2) then for i from 1 to Length(l1) do r := Append(r, Member(l1, i)*Member(l2, i)); end do: end if; return r; end proc:with(plots):

# trapez: integriert fnc im Intervall [a,b] mit Schrittweite n nach der # Trapezformel #trapez := proc(fnc,n,a,b) local h,i,res; h := (b-a)/n; res := (fnc(a)+fnc(b))/2; res := res+sum(fnc(a+i*h), i=1..n-1); return evalf(res*h); end proc:# mitte: integriert fnc im Intervall [a,b] mit Schrittweite 2n nach der # Mittelpunktsformel #mitte := proc(fnc,n,a,b) local h,i,res; h := (b-a)/(2*n); res := 0; for i from 1 by 2 to 2*n-1 do res := res+fnc(a+i*h); end do; return evalf(res*2*h); end proc:# simpson: integriert fnc im Intervall [a,b] mit Schrittweite 2n nach der # Simpson'schen Regel #simpson := proc(fnc,n,a,b) local h,i,res,s; h := (b-a)/(2*n); res := fnc(a)+fnc(b); s := 0; for i from 1 by 2 to 2*n-1 do s := s+fnc(a+i*h); end do; res := res+4*s; s := 0; for i from 2 by 2 to 2*n-2 do s := s+fnc(a+i*h); end do; res := res+2*s; return evalf(res*h/3); end proc:

# rnl, cnl: Nullstellen und Koeffizienten der Lagrange-Ploynome f\374r n=2..5 #rnl := [ [-1/sqrt(3), 1/sqrt(3)], [-sqrt(3/5), 0, sqrt(3/5)], [-sqrt((3 + 2*sqrt(6/5))/7), -sqrt((15 - 2*sqrt(30))/35), sqrt((15 - 2*sqrt(30))/35), sqrt((3 + 2*sqrt(6/5))/7)], [-sqrt((5 + 2*sqrt(10/7))/9), -sqrt((35 - 2*sqrt(70))/63), 0, sqrt((35 - 2*sqrt(70))/63), sqrt((5 + 2*sqrt(10/7))/9)] ]:cnl := [ [1, 1], [0.5555555556, 0.8888888889, 0.5555555556], [0.3478548451, 0.6521451549, 0.6521451549, 0.3478548451], [0.2369268850, 0.4786286705, 0.5688888889, 0.4786286705, 0.2369268850] ]:# gauss: integriert fnc im Intervall [a,b] mit Guass'scher Quadratur der # Ordnung n #gauss := proc(fnc,n,a,b) local res,i,rl,cl,s; if n<2 or n>5 then return 0; end if; res := (b-a)/2; rl := Member(rnl,n-1); cl := Member(cnl,n-1); s := sum(Member(cl,i)*fnc(((b-a)*Member(rl,i)+a+b)/2), i=1..n); return evalf(res*s); end proc:

# romb: Integriert nach Romberg-Schema # romb := proc(fnc,a,b,k,j) local r,h,res; h := u -> (b-a)/2^(u-1); r := proc(u,v) option remember; local i; if u=1 and v=1 then return( (fnc(a)+fnc(b))*(b-a)/2 ); elif v=1 then return( (1/2)*(r(u-1,1)+h(u-1)*sum(fnc(a+(2*i-1)*h(u)),i=1..2^(u-2))) ); else return( r(u,v-1)+(r(u,v-1)-r(u-1,v-1))/(4^(v-1)-1) ); end if; end proc; res :=r(k,j); return evalf(res); end proc:

AWP: NiMvLCgmJSJ5RzYjJSN0dEciIiIqJiYlImNHNiNGKUYpJkYmNiMlInRHRilGKSokRiYiIiRGKSomJkYsNiMiIiNGKS0lJHNpbkdGL0Yp , NiMvLSUieUc2IyIiISIiIg== , NiMvLSYlInlHNiMlInRHNiMiIiFGKg== NiMvJiUiY0c2IyIiIiomIiImRicpIiM1LCQiIiMhIiJGJw== , NiMvJiUiY0c2IyIiIy0lJkZsb2F0RzYkIiN2ISIi , 0 <= t <= 100DGL-System: NiMvJiUjdTFHNiMlInRHJSN1Mkc= , NiMvJiUjdTJHNiMlInRHLCgqJiYlImNHNiMiIiMiIiItJSRzaW5HRiZGLkYuKiYmRis2I0YuRi5GJUYuISIiKiQlI3UxRyIiJEY0 , NiMvLSUjdTFHNiMiIiEiIiI= , NiMvLSUjdTJHNiMiIiFGJw== # Approximation #c1 := 0.05: c2 := 7.5: f1 := (t,u1,u2) -> u2: f2 := (t,u1,u2) -> c2*sin(t)-c1*u2-u1^3:a := 0: b := 100: n := 3000: h := (b-a)/n:au1 := [1]: au2 := [0]:for i from 0 to n-1 do k11 := evalf( h*f1(a+i*h, Last(au1), Last(au2))); k12 := evalf( h*f2(a+i*h, Last(au1), Last(au2))); k21 := evalf( h*f1( a+(i+1/2)*h, Last(au1)+k11/2, Last(au2)+k12/2)); k22 := evalf( h*f2( a+(i+1/2)*h, Last(au1)+k11/2, Last(au2)+k12/2)); k31 := evalf( h*f1( a+(i+1/2)*h, Last(au1)+k21/2, Last(au2)+k22/2)); k32 := evalf( h*f2( a+(i+1/2)*h, Last(au1)+k21/2, Last(au2)+k22/2)); k41 := evalf( h*f1( a+(i+1)*h, Last(au1)+k31, Last(au2)+k32)); k42 := evalf( h*f2( a+(i+1)*h, Last(au1)+k31, Last(au2)+k32)); au1 := Append(au1, Last(au1)+(1/6)*(k11+2*k21+2*k31+k41)); au2 := Append(au2, Last(au2)+(1/6)*(k12+2*k22+2*k32+k42)); end do:plist := []: for i from 1 to n do plist := Append( plist, [ Member(au1,i), Member(au2,i) ] ); end do:listplot(plist);

DGL's: Beute: NiMvJiUieEc2IyUidEcqJkYlIiIiLCYqJiYlImNHNiNGKUYpJSJ5R0YpRikmRi02IyIiIyEiIkYp => u1 = x J\344ger: NiMvJiUieUc2IyUidEcqJkYlIiIiLCYmJSJjRzYjIiIkRikqJiZGLDYjIiIlRiklInhHRikhIiJGKQ== => u2 = yStabile L\366sungen: NiQvLSUieEc2IyIiIUYnLy0lInlHRiZGJw== und NiQvLSUieEc2IyIiISIlKzovLSUieUdGJiIkKyY= # Approximation: u1(t) = x(t), u2(t) = y(t) #c1 := 0.0008: c2 := 0.4: c3 := 3.0: c4 := 0.002: f1 := (t,u1,u2) -> u1*(c1*u2-c2): f2 := (t,u1,u2) -> u2*(c3-c4*u1):a := 0: b := 100: n := 1000: h := (b-a)/n:au1 := [1500]: au2 := [500]:for i from 0 to n-1 do k11 := evalf( h*f1(a+i*h, Last(au1), Last(au2))); k12 := evalf( h*f2(a+i*h, Last(au1), Last(au2))); k21 := evalf( h*f1( a+(i+1/2)*h, Last(au1)+k11/2, Last(au2)+k12/2)); k22 := evalf( h*f2( a+(i+1/2)*h, Last(au1)+k11/2, Last(au2)+k12/2)); k31 := evalf( h*f1( a+(i+1/2)*h, Last(au1)+k21/2, Last(au2)+k22/2)); k32 := evalf( h*f2( a+(i+1/2)*h, Last(au1)+k21/2, Last(au2)+k22/2)); k41 := evalf( h*f1( a+(i+1)*h, Last(au1)+k31, Last(au2)+k32)); k42 := evalf( h*f2( a+(i+1)*h, Last(au1)+k31, Last(au2)+k32)); au1 := Append(au1, Last(au1)+(1/6)*(k11+2*k21+2*k31+k41)); au2 := Append(au2, Last(au2)+(1/6)*(k12+2*k22+2*k32+k42)); end do:plist := []: for i from 1 to n do plist := Append( plist, [ Member(au1,i), Member(au2,i) ] ); end do:listplot(plist);# Stabile L\366sung #solve( {x*(c1*y-c2) = 0, y*(c3-c4*x) = 0}, {x,y} );

Zu Aufgabe (b):Lineare Interpolation: St\374tzwerte (NiQmJSJ4RzYjIiIiJiUieUdGJQ==) und (NiQmJSJ4RzYjIiIjJiUieUdGJQ==) Interpolationsgerade im Intervall [NiQmJSJ4RzYjIiIiJkYkNiMiIiM=] : NiMvLSUieUc2IyUieEcsJi0lImFHNiMsJkYnIiIiJkYnNiNGLSEiIkYtJkYlRi9GLQ== , wobei NiMvJSJhRyomLCYmJSJ5RzYjIiIjIiIiJkYoNiNGKyEiIkYrLCYmJSJ4R0YpRismRjFGLUYuRi4= # (a) #f := x -> sqrt(1+cos(x)^2):g := x -> simpson(f,32,0,x):plot(g(x),x=1..50);# Lineare Interpolation zwischen x1=31 und x2=32 f\374r x=31.5 #y1 := g(31); y2 := g(32); lp := (y2-y1)*0.5+y1; g(31.5);(g(32)+g(31))/2;# Vergleich Maple-intern #l := seq( evalf( int( f(x),x=0..k)), k=1..50):listplot([l]);Member([l],32);