Lab Num1: Uppgiftsbeskrivning
I denna labb skall ni lära er att använda Python för att utföra enklare beräkningar och plotta ersultaten. För detta kommer ni att behöva modulerna Numpy och Matplotlib.
Nyckelkoncept: Vektorer Links to an external site., Matriser Links to an external site., Plottning Links to an external site., Ekvationslösning
Uppgift 0 (förberedelser)
- För att köra Numpy behöver ni installera några extra paket i Pycharm. Det som skall installeras är
- Numpy
- Matplotlib
- Scipy
- Instruktioner för hur man installerar paket finns här
video-instruktion Links to an external site.- PyCharm-hemsidan Links to an external site.
- Notera att man även kan göra på andra sätt om man vill Links to an external site., men vi föreslår Pycharms inbyggda.
- En introduktion till NumPy finns här. Links to an external site.
- Flera av uppgifterna är hämtade från denna tentamen Download tentamen i envariabel, från 2015.
Uppgift 1: Vektorer
Skapa följande variabler. Notera att i uppgift d och e skall ni använda er av vektorn a (och inte skriva in värdena 1 2 3 4 5 för hand). Nyttiga metoder: tolist, reshape, vectorize, arange, array
a: en vektor av längd 5 med värdena 1 till 5
b: en vektor med alla värden från 0 till 2*pi i steg om 0.1
c: en 3x2-matris med alla värden från 1 till 6
d: vektorn a utökad med elementen 6, 7 till en vektor av längd 7
e: vektorn a utökad med elementen -1 till -5 till en 2x5-matris
f: en vektor med sinus-funktionen utvärderad på alla värden i vektorn b
Uppgift 2: Funktioner och Matrisoperationer
I denna uppgift övar vi på att hantera skalärer, vektorer och matriser.
Skapa följande funktioner:
- function_a: en funktion med input x och output x*x (x i kvadrat) (både input och output är skalärer, inte vektorer)
- function_b1: en funktion med input x och output x*x (både input och output är vektorer och x*x skall beräknas för varje element i vektorn
- function_b2: en funktion med input x och output x*x (input är en kolumnvektor och output är en skalär och x*x skall beräknas som en matrismultiplikation (kallas även skalärprodukt, inre produkt eller "dot product")
- function_c1: en funktion med input x och output x*x (både input och output är matriser av samma storlek, och x*x är matrisen där varje element är kvadraten av motsvarande element i x
- function_c2: en funktion med input x och output x*x (både input och output är matriser och x*x skall beräknas med hjälp av matrismultiplikation)
Uppgift 3
Titta på uppgift 1 i Tentamen 2015-10-27 i Envaribelanalys (SF1625) (tentamen som pdf Download tentamen som pdf).
- Plotta Links to an external site. funktionen, begränsa plottfönstret till +/- 10 i både x och y-led.
- Plotta de två asymptoterna (notera att bara den ena går att skriva som en funktion av x, den andra får ni lösa på annat sätt)
Uppgift 4
Skriv kod för att derivera funktionen f numeriskt, och jämför med den analytiska derivatan i en plott.
Funktionen f(x) = 1+sin(x)+0.5*cos(4*x);
Den analytiska derivatan av f blir: f'(x) = cos(x)-2*sin(4*x);
Uppgift 5
Betrakta uppgift 2 i tentan från 2015-10-27. Ni skall nu beräkna de två integralerna numeriskt med hjälp av en Riemansumma Links to an external site.. Hur stort blir felet? Pröva med lite olika steglängder (dx).
Uppgift 6
Betrakta uppgift 6 i tentamen Download tentamen från 2015-10-27. I lösningen till uppgiften beskrivs den allmänna lösningen, och dess beteende efter lång tid. Här skall ni plotta 81 st specialfall av den allmänna lösningen, och se hur de skiljer sig över tiden. Låt t gå från 0 till 10.
- Plotta lösningarna för alla kombinationer av A och B ( A = -4, -3, -2, -1, 0, 1, 2, 3, 4 och B = -4, -3, -2, -1, 0, 1, 2, 3, 4)
- Stämmer påståendet i lösningen till 6B?
Uppgift 7
I denna uppgift skall ni beräkna och plotta Taylorpolynomet av grad 13 kring punkten t=0 till funktionen Sin(x). (minsta termen skall alltså vara x och största vara x^13). Har ni glömt bort hur formeln för Taylorutvecklingen av sin(x) ser ut kan ni titta här http://www.wolframalpha.com/input/?i=taylor+series+sin+x Links to an external site..
- Plotta sin(x) för x från -10 till 10
- Plotta ert taylorpolynom
- För vilka x är approximationen bra respektive dålig?
Uppgift 8
I denna uppgift ska ni plotta och lösa en ekvation numeriskt.
a) Ekvationen x=cos(x) saknar analytisk lösning, men går att lösa numeriskt. Man kan se ungefär vad svaret blir genom att plotta
y1=x och
y2=cos(x) i samma plot, alternativt kan vi skriva om ekvationen till
x−cos(x)=0 och plotta
y=x−cos(x). Rita upp dessa plottar!
b) Lös ekvationen x−cos(x)=0 med intervallhalvering
Links to an external site.. Utgå från plotten ovan, och välj ut två (små) heltal a och b, sådana att y(a) < 0 och y(b) > 0. Upprepa nu följande:
- Välj c till värdet precis i mitten mellan a och b.
- Beräkna y(c)
- Om y(c) < 0, ersätt värdet på a med c, annars ersätter ni värdet på b med c, så att ni har ett nytt par a och b som fortfarande uppfyller y(a) < 0 och y(b) > 0, men ligger hälften så lång ifrån varandra som innan
- Om
|a−b|<10−12 bryt, gå annars tillbaka till 1
Vad är lösningen på x−cos(x)=0, och hur stort är ert fel ?
c) Lös ekvationen x−cos(x)=0 med Newton-Raphsons
Links to an external site. metod. Utgå ifrån samma heltal som ni valde till startgissning för a i uppgift 8b ovan, dvs låt
x0=a. Låt
f(x)=x−cos(x),
f′(x) är dess derivata, och n=0.
Upprepa nu:
xn+1=xn−f(xn)f′(xn)
- om
|xn+1−xn|<10−12 bryt, annars låt n = n+1 och gå till 1.
Vad blir lösningen på x−cos(x)=0, och hur stort är ert fel?
d) Vilken av metoderna i b) och c) kräver flest iterationer innan den kommer fram till svaret? Fungerar båda metoderna lika bra oavsett valet av startvärde? Vad händer om vi låter a börja som ett värde i närheten av −π2, tex
0.000001−π2 ?
e) Testa gärna att lösa ekvationen med en inbyggd ekvationslösare (fsolve Links to an external site.) också!
Uppgift 9
Utgå ifrån följande tutorial: https://matplotlib.org/mpl_toolkits/mplot3d/tutorial.html Links to an external site.
Plotta följande ytor i 3D:
a) En kon Links to an external site.
b) En pyramid Links to an external site.
c) En halvsfär Links to an external site.
d) Två spiraler som går runt varandra Links to an external site.
Notera att era plottar inte behöver vara centrerade, skalade eller riktade på samma sätt som exemplen, så länge det är samma form.
Uppgift 10
Här tar vi en titt på olika sätt att arbeta med ekvationssystem.
a) Läs här om hur man löser linjära ekvationssystem med NumPy:
https://numpy.org/doc/stable/reference/generated/numpy.linalg.solve.html
Links to an external site.
Lös sedan följande ekvationssystem:
4x1−x2−9x3−4x4−6x5=−59x1+x2−x3+4x4−5x5=−21−3x2+4x3+7x4=203x1−5x2−5x3−3x4+7x5=169x1−x2+4x3−8x4−9x5=−11
b) Om man har fler ekvationer än obekanta så säger man att ekvationssystemet är överbestämt. Detta är vanligt förekommande då man genomfört upprepade experiment eller observationer. För att hitta den lösning som ger en så bra anpassning som möjligt till observerade data kan man använda minstakvadratanpassning. Läs här om hur det fungerar i NumPy:
https://docs.scipy.org/doc/stable/reference/generated/numpy.linalg.lstsq.html
Links to an external site.
Stambanan för japanska snabbtåg, Shinkansen, använder ett ganska linjärt system för att beräkna priser baserat på avstånd. Utgå ifrån filen shinkansen.text Download shinkansen.text, och beräkna detta linjära samband. Antag att priset P kan beräknas som P = a + b*x, där x är antalet färdade kilometer, b är priset per kilometer, och a är en grundavgift som inte beror på avståndet. Vad borde det kosta att åka till Sendai, 325.4 km från Tokyo? Givet dagens yenkurs, Links to an external site. vad skulle det kosta att åka mellan Stockholm och Göteborg (455 km), med motsvarande snabbtåg? (Nu antar vi att Japan och Sverige har samma kostnader, kundunderlag, mm. vilket förstås inte stämmer i verkligheten). Tänk på att ni kan testa om ni fått rätt värden på a och b genom att testa vad ni får för priser till alla kända stationer i listan.
c) Här ska ni lösa ett optimeringsproblem. Vi antar att studentkollektivet Snåljåpen har bestämt sig för att laga mat gemensamt, med kraven att den ska vara så billig som möjligt, men samtidigt uppfylla grundläggande krav på näringsriktighet. Ganska förenklat kan vi lista hur mycket av olika näringsämnen den typiske studenten (kvinna, drygt 25 år gammal) behöver varje dag:
Näringsbehov.text Download Näringsbehov.text
Vi antar dessutom att vi har en lista på hur mycket av de olika näringsämnena som olika råvaror innehåller, tillsammans med priset de kostar i närmaste snabbköp:
nutrients.text Download nutrients.text
Vi antar slutgiltligen att studenterna har ett genomsnittlig energibehov på 8710 kJ per dag.
Uppgiften är nu: beräkna vilken som är den billigaste kombination av ingredienser som uppfyller alla kraven på näring.
För att lösa detta kan ni använda linjärprogrammering, läs här:
Texten om hur linprog fungerar kan vara lite kryptisk, så här följer några tips:
- Kostnaden c ska vara en array med alla kostnader på de olika ingredienserna.
- Notera att ordningen på ingredienser och näringsämnen inte får ändras mellan de olika ekvationerna och olikheterna.
- Villkoren för näringsämnen formuleras som olikheter, där A_ub * x <= b_ub. x är det svar som linprog returnerar. Ni måste ange A_ub resp b_ub. Varje kolumn i matrisen A_ub motsvarar ett livsmedel, och varje rad är innehållet av ett näringsämne. Varje element i b_ub motsvarar minimibehovet. Eftersom linprog bara kan räkna med gränser av typen "<=", och ni egentligen vill att varje näringsämne ska vara över eller lika med ett minimikrav, så kan ni utnyttja att
A⋅x≥b⟺−A⋅x≤−b.
- Energiinnehållet ska vara ett specifikt värde - det räcker inte med en miniminivå - så här måste ni ange en likhet, dvs A_eq * x = b_eq. Här är A_eq en matrix med en enda kolumn (det fungerar inte med array) som innehåller energiinnehållet per livsmedel, och b_eq är kravet 8710.
- Eftersom det är väldigt stor skillnad på hur många gram som behövs av t.ex fett resp D-vitamin, så kan linprog få svårt att beräkna ett exakt svar. Lös detta genom att tolerera ett mindre noggrant svar, genom att ange options={"tol": 1e-10}, nu tillåter ni fel i 10:e värdesiffran, vilket fortfarande är ganska nära exakt svar.
- Ni behöver inte ange något för bounds, method eller callback. Defaultvärdena fungerar bra. Bounds utgår från att alla värden är tillåtna så länge de inte är negativa. Ni kan ändra detta om ni t.ex vill tillåta högst ett hekto av varje ingrediens (bounds=(0,1)).
- Hur påverkas svaret om man ansätter en vegansk meny, resp. har med kosttillskott ("supplements")? Eller om ni plockar bort någon annan ingrediens ni inte gillar?
Vad är den optimala ingredienslistan? Hur mycket kostar denna meny per dag?
Eftersom vi tidigare satt tolerance till 1e-10, så kan ingredienslistan i ert svar innehålla fel upp till den storleken. För att få en mer lättläst lista kan ni dels plocka bort alla ingredienser som anges i en mängd under 1e-10 (de är ju ändå inom felmarginalen om de sätts till noll), och sortera listan som är kvar, så att de ingredienser som det ska vara mest av kommer först.
Notera: denna uppgift förutsätter en mycket förenklad syn på näringslära, och har en begränsad mängd ingredienser. I verkligheten finns det förstås många fler faktorer att ta hänsyn till.