Aaltoyhtälöpohja
29.4.2005
Alustukset
> | restart:with(plots):setoptions3d(orientation=[-30,50],axes=box): |
Warning, the name changecoords has been redefined
> | #with(LinearAlgebra): |
> | read("c:/usr/heikki/numsym05/maple/ns05.mpl"); |
> | # read("/p/edu/mat-1.192/ns05.mpl"); #read("/home/apiola/opetus/numsym/05/maple/ns05.mpl"); |
D'Alembertin ratkaisu: Aaltoyhtälö ilman R-ehtoja
> | u:=(x,t)->1/2*(f(x+c*t)+f(x-c*t))+1/2/c*int(g(y),y=x-c*t..x+c*t); |
> | u(x,t); |
Coo Exe 5.3 s. 175
> | c:=1: f:=x->exp(-x^2); g:=x->0: |
> | animate(u(x,t),x=-10..10,t=0..9,frames=10,numpoints=300); |
> | plot3d(u(x,t),x=-10..10,t=0..9,grid=[50,50]); pinta:=%: |
> | display(pinta,seq(spacecurve([x,t,u(x,t),x=-10..10],thickness=3,color=red),t=[0,1,2,3,5,9])); |
Karakteristikat:
> | T1:=a->-1/c*x+a: T2:=a->1/c*x+a: |
> | karkuva:=a->plot([T1(a),T2(a)],x=-10..10); |
> | karkuva(1); |
> | display(seq(karkuva(k),k=[-2,-1,0,3,4,5])); |
Aaltoyhtälö reunaehdoin
Jos on L:n pituinen sauva, ja Dirichlet'n RE't, niin D'Alembertin kaava pätee, kun jatketaan funktiot
f ja g parittomina,
- jaksoisina.
KRE-kaava (17) s. 592, f:n paritonta 2L-jaksoista jatketta merkitty f* (tässä ftahti ). KRE-kaavassa alkunopeus =0, joten integraaliosa puuttuu.
> | u:=(x,t)->1/2*(ftahti(x+c*t)+ftahti(x-c*t))+1/2/c*int(gtahti(y),y=x-c*t..x+c*t); |
> |
Esimerkkejä, joissa alkunopeus g=0
Esim: g=0
Näissä voi suuresti nautiskellen hyödyntää suoraan ns05.mpl-funktioitamme ParitonJ ja JJ.
Esim 1)
Kitaran kieltä näpsäytetään venyttämällä keskeltä ilman alkunopeutta.
> | restart:with(plots):setoptions3d(orientation=[-30,50],axes=box): #with(LinearAlgebra): read("c:/usr/heikki/numsym05/maple/ns05.mpl"); # read("/p/edu/mat-1.192/ns05.mpl"); #read("/home/apiola/opetus/numsym/05/maple/ns05.mpl"); |
Warning, the name changecoords has been redefined
> | f:=x->piecewise(x<L/2,x,x>L/2,L-x): g:=x->0:L:=10: |
> | paritonf:=ParitonJ(f): |
> | plot(paritonf(x),x=-20..20); |
> | ftahti:=JJ(paritonf,-10..10); |
> | plot(ftahti(x),x=-20..20); |
> | u:=(x,t)->1/2*(ftahti(x+c*t)+ftahti(x-c*t)); |
> | c:=1: |
> | animate(u(x,t),x=0..10,t=0..20,frames=50,numpoints=300); |
> | plot3d(u(x,t),x=0..L,t=0..20,grid=[30,30]); |
> | pinta:=%: |
Kuva tulee vielä havainnollisemmaksi, kun lisätään sopivia aika-arvoja vastaavat aaltomuodot, jolloin animaatiota on mukava verrata 3d-kuvaan.
> | display(pinta,seq(spacecurve([x,t,u(x,t),x=0..L],thickness=2,color=red),t=[0,1,2,3,5,10,20])); |
> |
Tehtäviä: Tee tämän mallin mukaisesti KRE8-kirjan s. 594 Problem set 11.3 teht. 3 ja 9.
(Kannattaa laskea toinen myös käsin (käyttäen tarvittaessa Maplea "laskimena".)
Alkunopeus
Tällöin D'Alembertin kaavassa tulee myös integraalitermi mukaan.
> | u:=(x,t)->1/2*(f(x+c*t)+f(x-c*t))+1/2/c*int(g(y),y=x-c*t..x+c*t); |
Sekä f että g ovat nyt parillisen jatkon 2L-jaksoisia jatkeita. Ongelmaksi tulee, että Maple ei mitenkään voi osata integroida JJ-operaattorin palauttamaa funktiota (proseduuria) symbolisesti. Numeerinen integrointi toki onnistuu. Jos halutaan symbolista, jaksollinen jatko pitää rakentaa sellaiseksi, että Maplen
int osaa sitä käsitellä.
Ajatus: Välillä [-L,L] määritelty pariton laajennus kerrotaan välin karakterisitisella funktiolla (funktiolla, joka on 1, kun ollaan välillä ja 0 välin ulkopuolella). Tällaisia funktioita siirretään jakson monikerran verran vasemmalle ja oikealle ja lasketaan yhteen riittävän suureen n:n arvoon saakka. Näin saadaan jaksollinen jatke, jota Maple taatusti osaa integroida.
> | restart:with(plots):setoptions3d(orientation=[-30,50],axes=box): #with(LinearAlgebra): read("c:/usr/heikki/numsym05/maple/ns05.mpl"); # read("/p/edu/mat-1.192/ns05.mpl"); #read("/home/apiola/opetus/numsym/05/maple/ns05.mpl"); |
Warning, the name changecoords has been redefined
> | alias(H=Heaviside): |
> | CHI:=(a,b,x)->H(x-a)-H(x-b); # Välin [a,b] karakteristinen funktio. |
> | JaksJ:=(f,L,N,x)->sum(f(x+nnn*2*L)*CHI(-L,L,x+nnn*2*L),nnn=-N..N); |
Esim:
> | f:=x->x*(1-x);plot(f,-1..1): |
> | fPton:=ParitonJ(f); |
> | plot(fPton,-1..1); |
> | ftahti:=x->JaksJ(fPton,1,10,x); |
> | plot(ftahti,-5..5); |
Vrt:
> | plot(JJ(fPton,-1..1),-5..5); |
Kumpikin toimii, mutta JJ:llä tuotettua ei niin vain voida integroida symbolisesti. Sensijaan JaksJ:llä tuotettu voidaan ilman muuta integroida millä välillä vain
(kunhan huolehditaan siitä, että JaksJ:n kolmas argumentti N on ao. väliä ajatellen riittävän iso (piirrä, niin näet, mitä tarkoitan!)).
> | int(JaksJ(fPton,1,5,x),x=-2..0); |
> | int(f(x),x=-1..1); |
> | Jf:=JJ(fPton,-1..1); |
> | Jf(x); |
Numeerinen integrointi käy:
> | evalf(Int(Jf(x),x=-2..0)); |
> | int(Jf(x),x=-2..0);evalf(%); |
No nyt ei ole ongelmaa D' Alembert'n kaavan käytölle yleisessä muodossa.
Mallitehtävä: Vaihdetaan edellisessä f:n ja g:n roolit. Ts. Alkuhetkellä kitaran kieli on lepotilassa, mutta alkunopeus on "kolmiomainen".
> | restart:with(plots):setoptions3d(orientation=[-30,50],axes=box): #with(LinearAlgebra): read("c:/usr/heikki/numsym05/maple/ns05.mpl"); # read("/p/edu/mat-1.192/ns05.mpl"); #read("/home/apiola/opetus/numsym/05/maple/ns05.mpl"); |
Warning, the name changecoords has been redefined
> | f:=x->0: g:=x->piecewise(x<L/2,x,x>L/2,L-x):L:=10: |
> | paritong:=ParitonJ(g): |
> | plot(paritong(x),x=-20..20): |
> | gtahti:=x->JaksJ(paritong,10,2,x): |
> | plot(gtahti,-200..200); # Kun otettiin N=2 yllä, voidaan ajassa mennä yli 50:n, jotta gtahti(x-t) ja gtahti(x+t) laskevat oikein. Jos halutaan mennä sen yli on Jaksj:n 3:tta argumenttia N kasvatettava. |
Koska nyt f=0, D'Alembert'n kaavassa on vain integraalitermi:
> | u:=(x,t)->1/2/c*int(gtahti(y),y=x-c*t..x+c*t); |
> | c:=1: |
> | animate(u(x,t),x=0..10,t=0..20,frames=50,numpoints=300); |
> | plot3d(u(x,t),x=0..L,t=0..20,grid=[30,30]); |
> | pinta:=%: |
> |
Kuva tulee vielä havainnollisemmaksi, kun lisätään sopivia aika-arvoja vastaavat aaltomuodot, jolloin animaatiota on mukava verrata 3d-kuvaan.
> | display(pinta,seq(spacecurve([x,t,u(x,t),x=0..L],thickness=2,color=red),t=[0,1,10,20])); |
Jospa vilkaistaan u(x,t):n lauseketta. Aikamoinen taidonnäyte ohjelmalta, tällaisten hallitseminen käsin ei kyllä taitaisi onnistua.
> | u(x,t): |
Entä se JJ-tyyli, onnistuuko sittenkään?
> | gtahti2:=JJ(paritong,-10..10); |
> | u2:=(x,t)->1/2/c*(evalf@Int)(gtahti2(y),y=x-c*t..x+c*t): |
> | u2(1,2);u(1,2); |
> | plot3d(u2(x,t),x=0..L,t=0..1); |
> | plot3d(u(x,t),x=0..L,t=0..1); |
> | u(2,3);u2(2,3); |
> | u2(4.9,1); |
Error, (in evalf/int) unable to convert to pwlist
> |
Kun ollaan sauvan puolenvälin lähellä, niin JJ-tyylillä tehty ei osaa laskea, kuten kuvastakin näkyy.
Siten tuo uusi summatyyli jaksolliseen jatkamiseen onkin ainoa toimiva näissä tehtävissä.
Huom! Mikään ei estä laskemasta sellaista tehtävää, jossa on sekä alkupoikkeutusta (f<>0) ja alkunopeutta (g<>0). Oppikirjoista nämä yleensä puuttuvat johtuen käsinlaskun vaivalloisuudesta, mutta sitä vaivaa meillä ei ole.
Ehdotus: Kokeilkaa jotain esimerkkiä!