lampoyhtalo1.html

Lämpöyhtälö muuttujat erotellen

KP3-II 9.12.2008

Alustukset

> restart:

> with(plots): setoptions3d(orientation=[-120,50],axes=box):

> #read("e:ns05.mpl");

Esim

KRE s. 603 Exa 1

> f:=x->100*sin(Pi*x/80); plot(f,0..80);

Plot_2d

> u:=(x,t)->100*sin(Pi*x/80)*exp(-0.001785*t);

(2.1)

> plot3d(u(x,t),x=0..80,t=0..400,axes=box);

Plot

> plot([seq(u(x,t),t=[seq(i*40,i=0..10)])],x=0..80);

Plot_2d

> with(plots):animate(u(x,t),x=0..80,t=0..400,frames=30);

Plot_2d

Exa 2. Muuten sama, mutta nyt

> setoptions3d(orientation=[-120,50],axes=box):

> f:=x->100*sin(3*Pi*x/80); plot(f,0..80);

Plot_2d

> lambda[3]:=c*3*Pi/L;

(2.2)

> c:=sqrt(1.158):L:=80: lambda[3];

(2.3)

> u:=(x,t)->100*sin(3*Pi*x/L)*exp(-lambda[3]^2*t);

(2.4)

> plot3d(u(x,t),x=0..80,t=0..200,axes=box);

Plot

"Aikasiivuparvi" :

> #[seq(u(x,t),t=[seq(i*0.1,i=0..5)])];

Tällainen lausekkeiden lista voidaan antaa suoraan plot:lle, kuten tiedämme: Otetaan samantien enemmän aikaviipaleita.

> plot([seq(u(x,t),t=[seq(i*20,i=0..10)])],x=0..80);

Plot_2d

> animate(u(x,t),x=0..80,t=0..200,frames=30);

Plot_2d

Esim. 3

> summaf:=(x,t,N)->400/Pi*add(1/(2*k-1)*sin((2*k-1)*Pi*x/80)*exp(-(lambda(2*k-1))^2*t),k=1..N);

proc (x, t, N) options operator, arrow; `+`(`/`(`*`(400, `*`(add(`/`(`*`(sin(`+`(`*`(`/`(1, 80), `*`(`+`(`*`(2, `*`(k)), `-`(1)), `*`(Pi, `*`(x)))))), `*`(exp(`+`(`-`(`*`(`^`(lambda(`+`(`*`(2, `*`(k))... (2.5)

> lambda:=n->c*n*Pi/L;

(2.6)

> c:=1.158;L:=80;

(2.7)

> summaf(x,t,5);

(2.8)

> plot3d(summaf(x,t,10),x=0..L,t=0.1..300,axes=box);

Plot

Summafunktion "aikasiivuparvi" :

> plot([seq(summaf(x,t,10),t=[seq(10*i,i=0..40)])],x=0..L);

Plot_2d

> animate(summaf(x,t,10),x=0..L,t=0..400,frames=60,numpoints=200);

Plot_2d

Muuttujien erottelu

> restart:
with(LinearAlgebra): with(plots):

> setoptions3d(orientation=[-120,50],axes=box):

> #read("c:\\usr\\heikki\\numsym05\\maple\\ns05.mpl");

> read("e:\\ns05.mpl");

> # read("/p/edu/mat-1.192/ns05.mpl");

> osdy:=c^2*diff(u(x,t),x$2)=diff(u(x,t),t);eval(subs(u(x,t)=F(x)*G(t), osdy));

> separoitu:=simplify(%/(c^2*F(x)*G(t)));

> xyht:=lhs(separoitu)=-p^2;
tyht:=rhs(separoitu)=-p^2;

> FDY:=lhs(xyht)*F(x)=rhs(xyht)*F(x);

> GDY:=c^2*lhs(tyht)*G(t)=c^2*rhs(tyht)*G(t);

> ff:=rhs(dsolve(FDY,F(x)));
dsolve(GDY,G(t)): gg:=rhs(%);

> ff:=subs(_C2=A,_C1=B,ff); gg:=subs(_C1=1,gg);

RE1 & RE2 ==>

> ff:=subs(A=0,B=1,p=n*Pi/L,ff);

> gg:=subs(p=n*Pi/L,gg);

> uu:=ff*gg;

> u:=sum(B[n]*ff*gg,n=1..infinity);

Sinisarjaksi:

Alkuehto saadaan toteutumaan valitsemalla B[n]:t AE-funktion f (2L-jaksoisen) sinisarjan kertoimiksi.

Määritellään sinisarjan kertoimet ja annetaan alkuarvofunktion f sekä jakson puolikkaan L olla

vapaita, myöhemmin annettavia symboleja (vrt. fourier.mws-työarkki).

> f:='f':L:='L':n:='n':Bn:=(2/L)*Int(f(x)*sin(n*Pi*x/L),x=0..L);

Otetaan esimerkiksi KRE s. 648, Exa 3

> f:=x->piecewise(x<L/2,x,x>L/2,L-x);L:=Pi;

> plot(f,0..Pi);

:n pituinen sauva sopiii matemaatikoille.

> bn:=simplify(value(Bn));

> bn:=trigsiev(%,n);

> uu;

> un:=bn*uu;

> summaf:=(x,t,N)->add(un,n=1..N);

> summaf(x,t,4);

> summaf(x,t,10);

> c:=1:

> plot3d(summaf(x,t,5),x=0..Pi,t=0..2,axes=box);

Summafunktion "aikasiivuparvi" :

> [seq(summaf(x,t,5),t=[seq(i*0.1,i=0..5)])];

Tällainen lausekkeiden lista voidaan antaa suoraan plot:lle, kuten tiedämme: Otetaan samantien enemmän aikaviipaleita.

> plot([seq(summaf(x,t,40),t=[seq(i*0.1,i=0..10)])],x=0..Pi);

> animate(summaf(x,t,40),x=0..Pi,t=0..10,frames=30);

Eristetyt päät, adiabaattiset reunaehdot

> restart:
with(plots):setoptions3d(orientation=[-120,50],axes=box):

> #read("c:\\usr\\heikki\\numsym05\\maple\\ns05.mpl");

> read("e:ns05.mpl");

> #uu:=(x,t,n)->cos(n*Pi*x/L)*exp(-(c*n*Pi/L)^2*t);

> #usum:=(x,t,N)->add(a(n)*uu(x,t,n),n=0..N);

> #usum(x,t,3);

> #f:='f':L:='L':n:='n':an:=(2/L)*Int(f(x)*cos(n*Pi*x/L),x=0..L);

Esimerkiksi KRE s. 64x, Exa 4

> f:=x->piecewise(x<L/2,x,x>L/2,L-x);L:=Pi;

> plot(f,0..Pi);

> L:='L':

> a:=n->2/L*int(x*cos(n*Pi*x/L),x=0..L/2)+2/L*int((L-x)*cos(n*Pi*x/L),x=L/2..L);

> a(0):=1/L*(int(x,x=0..L/2)+int(L-x,x=L/2..L));

> a(n);

> trigsiev(a(n),n);

> a(0);

> seq(a(n),n=0..10);

> u:=(x,t,N)->L/4-8/L/Pi^2*sum(1/(4*k-2)^2*cos((4*k-2)*Pi/L*x)*exp(-(((4*k-2)*c*Pi/L)^2)*t),k=1..N);

> u(x,t,5);

> c:=1: L:=Pi: u(x,t,5);

> plot3d(u(x,t,10),x=0..Pi,t=0..2,axes=box);

Summafunktion "aikasiivuparvi" :

> #[seq(usum(x,t,20),t=[seq(i*0.1,i=0..5)])];

Tällainen lausekkeiden lista voidaan antaa suoraan plot:lle, kuten tiedämme: Otetaan samantien enemmän aikaviipaleita.

> plot([seq(u(x,t,10),t=[seq(i*0.1,i=0..10)])],x=0..Pi);

> animate(u(x,t,20),x=0..Pi,t=0..2,frames=30);