\documentclass[finnish,11pt,twocolumn]{article}
\usepackage{/home/apiola/tex/v3harj}
\include{/home/apiola/tex/defs}
%\usepackage{v23harj}
%\include{defs}

\begin{document}

%\vspace{0.2cm}
%\begin{verbatim}
%\end{verbatim}

\harjoitus{1}{37--38}{13 -- 21.9.2001}

Harjoitussysteemi:

ke 12--14: AV: Liitutaulu + Matlab/Maple \\
to 12--14: KV: Yhdess� ty�skentely�, ohjattu harjoitus.\\
pe 10--12: LV: Liitutaulu+tietokone (Matlab/Maple) (Hiukan enemm�n 
tietokoneharjoitus kuin AV.)

Oppilaat saavat (maahantuojan luvalla) asentaa kurssin ajaksi omiin 
koneisiinsa Matlabin (5.3).
Kotikoneelle kannattaa hakea ainakin \\
{\tt http://www.math.hut.fi/teaching/v/matlab/laode/dfield5.m}

Huomaa, ett� {\tt dfield5}-funktion k�ytt� ei vaadi mit��n muuta
Matlab-taitoa kuin k�ynnistyksen ja polun asettamisen niin, ett�
Matlab l�yt�� dfield5.m:n.

Viikolla 37 emme voi pit�� ylim��r�ist� Matlab-harjoitusta, koska salit
ovat Lapio-kurssien k�yt�ss�.
Opiskelkaapa siis kotona, tarkoitusta varten olen kirjoitellut
teksti� {\tt http://www.math.hut.fi/teaching/v/matlab/opas/}.
(Osa 1 alkaa olla aika valmis ja osa osasta 2). Odotan Matlab-kysmyksi�
erityisesti viikon 37 aikana, teen www-sivulle koosteen Q\&A (tai FAQ).)

Varaamme osan to-pe 20--21.9. harjoituksista neuvottuun Matlab (ja miksei
v�h�n Maplekin)-harjoitteluun.
(Siksi vain 4 teht�v�� pe.)


\subsubsection*{Alkuviikko (AV)}

\begin{enumerate}

\item %1 
Ratkaise differentiaaliyht�l� sijoittamalla ratkaisuehdotus (RE)
annettuun yht�l��n tai esim. integroimalla, arvaamalla tms.:

(a) $y'+y=x^2-2$, RE: $y=Ce^{-x}+x^2-2x$ \quad
(b) $y''+y=0$, RE: $y=a\cos x + b\sin x$ \quad
(c) $y'''=e^x$, \quad
(d) $x+yy'=0$, RE: $x^2+y^2=C$ ($C > 0$, vakio).\\
Vihje (d)-kohtaan: Derivoi implisiittisesti, ts. oleta, ett� on olemassa
derivoituva funktio $x\mapsto y(x)$ s.e. $x^2+y(x)=C$ ja derivoi puolittain.
(T�ss� tapauksessa olemassaolo tiedet��n, onhan $y(x)=\sqrt{C-x^2}$
t�llainen. T�m�n eksplisiittisen lausekkeen k�ytt� ei silti kannata, 
se vain mutkistaa asioita, olkaamme siis implisiittisi�.)


%Sano my�s, mink�tyyppisi� diffyht�l�t ovat ominaisuuksiin
%{\em lineaarinen ja homogeeninen/ep�homog., kertaluku, ep�lineaarinen}
%n�hden.


\item %2
Lis��ntymiskykyinen populaatio, jonka lis��ntymist� rajoittavia tekij�it�
ei ole, noudattaa yleens� likimain {\it Malthus'n} lakia, ts. nykyhetkell�
kasvunopeus on verrannollinen populaation nykykokoon.


Muodosta ilmi�lle differentiaaliyht�l�malli ja ratkaise.

Ratkaisussa esiintyy vakiot $y_0$ = populaation koko alkuhetkell� ja
$k$ = "lis��ntymiskykyvakio". 

T�ydenn�mme vanhaa tuttua USA:n v�kilukutaulukkoa arvoilla, jossa
toinen rivi ilmaisee v�kiluvun miljoonissa.

\begin{verbatim}
1800    1830      1860     1890     1920    1950    1980
5.3      13        31       63      106     150     230
\end{verbatim}

Ota ajan 0-hetkeksi vuosi 1800 ja m��rit� $y_0$ ja $k$ kahden ensimm�isen
sarakkeen perusteella. Miten hyvin malli vastaa todellisuutta?
% Vrt. Heikkil�: Matem. mallinnus s. 11 alh.

%{\bf Otetaan t�h�n liittyv� PNS-sovitus LV-teht�v�ksi.}


%\item %LAODE s. 117
%Talletat $10 000$ mk tilille, josta saat jopa $5.5$ \%:n koron, joka
%lis�t��n p��omaan jatkuvasti. Kirjoita differentiaaliyht�l� ja ratkaise
%saadaksesi selville, paljonko rahaa on tilill�si $\frac{1}{2}$:n vuoden
%kuluttua.
%Paljonko rahaa kertyisi samassa ajassa, jos korko lis�tt�isiin p��omaan
%kuukausittain?


\item %3
T�t� teht�v�� harjoitellaan Matlab-teknisesti loppuviikon harjoituksissa. 
Osattava
esitt�� (ke 19.9.) liitutaululla k�sin piirt�en periaatteessa.

Tarkastellaan diffyht�l��  $y'=y-x$ alueessa $-1 \leq x \leq 1$,

Piirr� $xy$-koordinaatistoon suuntakentt�
ja isokliinej� k�sin ja kokeile my�s Matlabia laskimen roolissa.
Ota hilav�liksi aluksi vaikka $h=0.5$.

Matlab-laskussa kannattaa muodostaa matriisi, sanokaamme K, jonka alkioina
ovat arvot $y_i - x_j$, $i=1\ldots 5, j=1\ldots 5$ t�h�n tapaan:
\begin{verbatim}
  h=0.5; t=-1:h:1;x=0:h:2
  for i=1:5
    for j=1:5
      K(i,j)=y(i)-x(j)
    end
  end
\end{verbatim}

Koska matriisissa rivi-indeksi i juoksee alasp�in, on helpompaa sijoittaa
arvot koordinaatistoon k��nt�m�ll� matriisin sarakkeet yl�salaisin; t�m�
tapahtuu komennolla {\tt flipud}.

Katso siis matriisista {\tt flipud(K)} arvot ja merkitse ne kyn�ll� piirrokseen.
Tarkista k�sin (tai "skalaarilaskimella" (Matlabkin k�y)) muutama alkio 
ainakin. %Mitk� matriisin K alkiot kuuluvat kesken��n samaan isokliiniin?

(Hilan tihent�minen k�y nyt helposti muuttamalla vain yll� $h$:ta.)

My�hemmin opimme, ett� K-matriisin muodostaminen k�y k�tev�mmin,
tehokkaammin ja rutiininomaisemmin n�in:
\begin{verbatim}
x=...;y=...  % Kuten edell�
[X,Y]=meshgrid(x,y);
K=Y-X;    % Jos esiintyy kertolaskua, potenssia ym., 
          % on varustettava pisteell�, siis 
          %     .*, .^  ./  (taulukko-operaatiot)
\end{verbatim}

Nyt meill� on kaikki
data ker�ttyn� suuntakent�n piirt�mist� varten.
Voit katsoa skripti�  
{\tt ...v/3/01/matlab/suuntak1.m}, joka on my�s listattuna teht�v�paperin 
lopussa ja mietti�, mit� siin� tapahtuu.

Lopuksi voit kokeilla LAODE-funktiota {\tt dfield5}, jonka k�ytt� ei
vaadi Matlabin tuntemista lainkaan.


%Seuraava tietokoneharj:

%sitten Matlab-skripti� {\tt suuntak1.m}
%sopivasti editoiden ja Matlabissa ajaen. (Askel $h=0.25$ on aika sopiva,
%voit kokeilla muitakin.)
%Kommentoi j�lkimm�isess� tapauksessa isokliinipiirroskomento 
%pois skriptist�.
%Miksi suuntakentt� ei erotu ennen t�t�? (Perustelu kyn�ll� ja paperilla;
%ei siis riit�: "koska silt� n�ytt��".)


\item %5
Mill� xy-tason k�yrill� on ominaisuus: K�yr�n tangentin kulmakerroin
jokaisessa pisteess� $(x,y)$ on $-\frac{4 x}{y}$ ?

Ratkaise yht�l� muuttujien erottelulla (``separation of variables'').
Piirr� suuntakentt� isokliineja apuna k�ytt�en k�sin vaikkapa alueessa 
$[-2,2] \times [-2,2]$.

\qquad {\tiny Kun p��set t�h�n saakka, voit merkit� rastin.}

Kokeile my�s LAODE-funktiota dfield5. 
T�ss� on k�ytett�v� ahkerasti stop-n�pp�int�, ratkaisu ajautuu aina 
ongelma-alueelle, mik�li x-akseli on mukana.

Voit my�s t�ydent�� kuvaa alussa laskemillasi ratkaisuk�yrill� tyyliin:
\begin{verbatim}
x=linspace(-2,2,30);y=x; [X,Y]=meshgrid(x,y);
Z=... % muista pisteitt�iset laskutoimitukset.
contour(x,y,Z,1:10); shg
\end{verbatim}

Kokeile ja selit�!

\item %5 KRE s. 23, teht. 14
S�ili�ss� on 400 l vett�, johon on liuotettu 100 kg suolaa.
Puhdasta vett� virtaa s�ili��n nopeudella 2 l/min ja t�ysin sekoittunutta
liuosta poistuu s�ili�st� samalla 2 l/min nopeudella.
\begin{verbatim}

         --------|                    |
  vett� 2 l/min ->                    |
         ---------                    |
                 |.................. .|
                 |                    |
                 |    100 l           |
                 |                    -------------
                 |                      2 l/min ->
                  ---------------------------------

\end{verbatim}

Kuinka paljon suolaa on j�ljell� s�ili�ss� 1 tunnin kuluttua ?


\item %6 HUBW s. 50--51 V300, H2AVt2
Mihin yht�l��n liittyy mik�kin suuntakentt�, kuvat k��nt�puolella.


\end{enumerate}


\subsubsection*{Loppuviikko (LV)}

\begin{enumerate}


\item % KRE s. 13 teht. 17
Laskuvarjohypp��j�n yht�l�. Oletetaan, ett� hypp��j�n + varustuksen
massa = $m$ ja ilmanvastus on verrannollinen nopeuden neli��n, olkoon
verrannollisuskerroin $=b$.
T�ll�in Newtonin 2. laki antaa liikeyht�l�n:

$$mv'=mg - b v^2 .$$

Olkoon yksinkertaisuuden vuoksi $m=1, b=1$ ja $g=9.81 m/s^2$. 


Piirr� suuntakentt�.

Oletetaan, ett� laskuvarjo aukeaa, kun $v=10 m/s$, valitaan t�m� 
alkuhetkeksi $t=0$. Piirr� t�m� ratkaisuk�yr� suuntakentt�piirrokseen.
Yrit� n�hd� suuntakent�st�, ett� kaikki ratkaisut n�ytt�v�t l�hestyv�n
rajanopeutta $v\approx 3.13$ ja ett� ratkaisut ovat joko kasvavia tai
pienenevi� (ja mill� alkuarvoilla mit�kin, ja mit� tarkoittaa 
fysikaalisesti)

M��rit� rajanopeus suoraan yht�l�st�.

K�yt� piirroksiin LAODE-funktiota {\tt dfield5}

\item %2
Muodosta edell� olevan laskuvarjoteht�v�n analyyttinen ratkaisu
muuttujien erottelulla. M��rit� edell� mainittu ($v(0)=10$)-ratkaisuk�yr�.
Tarkista ratkaisu Maplella ja kokeile lopuksi Maplen {\tt dsolve}-
komentoa. (Ohje [HAM]-kirjassa.)

Ohje analyyttiseen: Muistathan, ett� osamurtohajoitelma on hy�dyllinen
rationaalilausekkeen integroinnissa 
(Maple: {\tt convert(lauseke,parfrac,muuttuja);} mutta osattava my�s k�sin).


\item
Vaihdamme t�ss� LAODE-tyyliseen notaatioon: $t$ on riippumaton muuttuja,
$x$ on ``riippuva'' muuttuja. Kannattaa totutella eri tyyleihin.

Ratkaise alkuarvoteht�v� $x'=\frac{x}{2}-e^{-t}$, $x(0)=-1$.
Kyseess� on {\em lineaarinen ep�homogeeninen} (EHY). T�m� lasku ei
edellyt� mit��n uutta muuttujien erottelun lis�ksi (ainoastaan uskomista), 
kaikki on t�ss� neuvottu.

Suorita ratkaisu n�in:
\begin{itemize}
\item
Ratkaise ensin vastaava (HY) $x'=\frac{x}{2}$  (yleinen ratkaisu).
\item
Yrit� keksi� jokin (EHY):n erityisratkaisu (siis mik� tahansa (EHY):n 
toteuttava).
Keksiminen on helppoa, kun mietit exp-funktion derivointia. (M��r��m�t�n
kerroin ratkaistaan sijoittamalla yrite (EHY):yyn).
\end{itemize}

Lineaaristen teoria sanoo, ett� (EHY):n yleinen = (HY):n yleinen + 
(EHY):n erikoinen.

Piirr� my�s suuntakentt� ja ratkaisuk�yri� (joko 
{\tt dfield5 tai suuntak}). Miten n�et suuntakent�st�,
ett� yht�l� ei ole autonominen?

\item
Ratkaise edellinen teht�v� integroivan tekij�n menettelyll�.

\end{enumerate}


\subsubsection*{Matlab-avustusta}

Kaikki Matlab-tiedostot sijoitetaan sek� kurssin www-matlab-hakemistoon
{\tt ../matlab} ett� hakemistoon {\tt /p/edu/mat-1.414/matlab/}.
J�lkimm�iset saadaan Matlab-istunnossa k�ytt��n antamalla Matlab-komento
{\tt addpath /p/edu/mat-1.414/matlab/}. T�m� kannattaa sijoittaa omaan
\verb_ $home/matlab/startup.m_ - tiedostoon.

\begin{verbatim}

% Matlab-skripti: ..v/3/01/matlab/suuntak1.m
% Kts. my�s exasuuntak.m
% HA 4.9.2001 
%
% Aihe: Diffyht�l�n  y'=f(x,y)  suuntakentt� ja isokliinit.
%
% Havainnollinen versio, jossa jananp�tk�t piirret��n for-
% silmukassa. Oikeaoppisempi tehoversio: suuntak.m
%
clf
f=inline('y-x','x','y')     % f(x,y)=y-x
%f=inline('x.*y','x','y')   % Toinen esimerkki: f(x,y)=xy

h=0.5;a=-1;b=1;c=0;d=2;   % Muuta tarpeen mukaan
x=a:h:b;y=c:h:d;                          

[X,Y]=meshgrid(x,y);
Z=f(X,Y);

 flipud(Z)   % havainnollistukseen, tosi k�yt�ss� kommentti.

 plot(X(:),Y(:),'x')   % Jos haluat merkit� ristill� pisteet, joissa
                        % derivaatta-arvot on laskettu, poista kommentti
hold on
xx=X(:); yy=Y(:); n=length(xx);
for i=1:n
   v=[1,f(xx(i),yy(i))];v=0.5*h*v/norm(v);
   alkup=[xx(i),yy(i)]-v;loppup=[xx(i),yy(i)]+v;
   plot([alkup(1),loppup(1)],[alkup(2),loppup(2)])
end;

% Jos haluat vain suuntakent�n, lopeta t�h�n.
% 
% Piirret��n edell� olevia korkeusarvoja vastaavat isokliinit.

x=linspace(a,b,30);y=linspace(c,d,30);
[X,Y]=meshgrid(x,y);
ZZ=f(X,Y);
 contour(x,y,ZZ,Z(:))  % Korkeusvektoriksi valitaan Z(:), t�ll�in           
axis equal             %  saadaan kaikki suuntakentt�piirroksen arvot,
axis([a,b,c,d])        %    ovelaa!
grid
title('y''=f(x,y)')     % Jos haluat oikean otsikon, editoi t�h�n
shg











\end{verbatim}
\end{document}



%Tarkastellaan diffyht�l�it� a) $x'=x-t$ alueessa $-1 \leq t \leq 1$,
%$0 \leq x \leq 2$ ja b) $x'=x/t$ alueessa $0 < t \leq 2$, $0 < x \leq 2$.
%(Huomaa, ett� vaihdoimme {\it KRE}-tyylist� {\it LAODE}-tyyliseen notaatioon,
%ts $x\to t$ ja $y\to x$).

%Piirr� kummassakin tapauksessa $tx$-koordinaatistoon suuntakentt�
%ja isokliinej� ensin k�sin ja sitten Matlab-skripti� {\tt suuntak1.m}
%sopivasti editoiden ja Matlabissa ajaen. (Askel $h=0.25$ on aika sopiva,
%voit kokeilla muitakin.)
%Kommentoi j�lkimm�isess� tapauksessa isokliinipiirroskomento 
%pois skriptist�.
%Miksi suuntakentt� ei erotu ennen t�t�? (Perustelu kyn�ll� ja paperilla;
%ei siis riit�: "koska silt� n�ytt��".)

T�m�n voit (soveltuvin osin) esitt�� puhtaasti k�sinlaskuteht�v�n�kin, 
jos Matlab tuntuu viel� vieraalta.
Toisaalta sit�k��n en pane pahakseni, jos kokeilet LAODE-funktiota 
{\tt dfield5}.



%\item % LAODE s. 117, teht. 12: Newton's law of cooling
      % 2000 Harj. 2 teht. 1

%Newtonin j��htymislaki sanoo, ett� kappaleen l�mp�tilan muutosnopeus
%ajanhetkell� $t$ on verrannollinen kappaleen ja ymp�r�iv�n v�liaineen
%l�mp�tilaeroon hetkell� $t$. Differentiaaliyht�l�ll� ilmaistuna:

%$$\frac{dT}{dt}=\alpha (T-T_v),$$

%miss� $T_v$ on ymp�r�iv�n v�liaineen vakiol�mp�tila ja $\alpha$ on
%verrannollisuuskerroin.

%Oletetaan, ett� kappaletta ymp�r�iv�n ilman l�mp�tila on $50^\circ$
%ja kappale j��htyy $100^\circ$:sta $75^\circ$:een $20$ min:ssa.
%Mik� on kappaleen l�mp�tila tunnin kuluttua?





\item %2
Radioaktiivisessa hajoamisessa ainem��r�n muutosnopeus on verrannollinen
nykyhetkell� j�ljell� olevaan ainem��r��n.
Muodosta t�m�n perusteella prosessia kuvaava differentiaaliyht�l�.
 
Radium-isotoopin $_{88}$Ra$^{224}$ puoliintumisaika on $3.6$ vrk. 
Kuinka paljon yhdest� grammasta on j�ljell� a) yhden vrk:n, b)
yhden vuoden kuluttua?

{\bf Ohjetta:}
{\em Puoliintumisaika} tarkoittaa aikaa, jonka kuluttua puolet aineesta on
hajonnut. Kts. my�s KRE$^8$ example 4 ss. 6--7. 




\item %5  HUB-W  

Tarkastellaan diffyht�l�it�

(a) $x'=-\frac{x}{t}$ \quad (b) $x'=\frac{x}{2t}$  \\
(c) $x'=x/t$ \quad %(d)  $x'=\frac{1-t}{1+x}$  \\

Hahmottele k�sin isokliinej� ja muodosta niiden avulla suuntakent�t sek�
hahmottele ratkaisuk�yri�.

Selvit�, miss� yht�l�t on m��ritelty.

Ratkaise lopuksi yht�l�t muuttujien erottelulla.

Vast: (a) $x=\frac{C}{t}$ \quad (b) $x=C\sqrt{|t|}$, 
\quad (c) $x=Ct$ \quad % (d) $(t-1)^2+(x+1)^2=C^2$



\item
Mieti (hetki), miten piirt�isit suuntakent�n, kun k�yt�ss�si on 
{\tt suuntak1.m} skriptiss� saatava data. Katso sitten skriptist�
ja selvit� itsellesi (ja taululla muille), miten for-silmukassa
tapahtuva piirto tapahtuu.



\item %2 KRE 18, Verhulstin logistinen populaatiomalli.
Sopiva my�s t�h�n kohtaan.
Harkitse !


\item KRE t. 20 CAS-project, sopiva my�s, mutta ...

\item
Viime syksyn harj2 LV-teht�viss� on varteenotettavia my�s, valmiiksi
ratkaistuja. Kysell��n lineaarisuutta, autonomisuutta,...

