# TMA4320 Introduksjon til vitenskapelige beregninger

## Beskjeder

• 21.08.2018: LF til kont finnes her
• 28.05.2018: Here are possible solutions to today's exam. Please let me know if you see any errors or inconsistencies.
• 14.05.2018:
• I have booked EL6 on 25.05 from 13:15 on for Q&A.
• Here are the results of the evaluation of the last project as well as more general questions about the course. Thanks to everyone who has filled it out and especially for the comments.
• 03.05.2018: There will be a Q&A possibility on 25.05 tentatively 13-15. I will book an auditorium and provide more details later.
• 24.04.2018:
• Curriculum/pensum will be the same as last year (see below)
• The exam is a closed book one, and there are too many formulas to remember (and this is not the point of the course any ways). However, there are some very simple methods which can be derived immediately as long as you understand the idea behind them, and I expect you to be able to apply those at the exam. You can get a very good idea about how the questions are formulated by looking at the previous years' exams, but here is the list to help you:
• Equations solving: bisection; fixed point (and its convergence rate); Newtons (the latter also for systems of equations)
• Interpolation (either Lagrange or Newtons). Do not memorize the error estimate, but be able to utilize it to estimate interpolation errors
• Quadratures: lowest order Newton-Cotes (trapezoid, mid-point) and the idea behind adaptive quadratures. You do not need to memorize the error estimates, those will be given if needed
• Systems of equations: LU with/without pivoting; Cholesky (you do not have to memorize the algorithm, you can always compute the coefficients one-by-one from $A=L L^T$); Jacobi/Gauss-Seidel; Newtons method for systems of non-linear equations
• ODEs: lowest order explicit/implicit RK methods (i.e. explicit/implicit Euler). Have a good idea about how 2nd order RK methods can be derived (i.e. explicit/implicit trapezoid or mid-point). I cannot emphasise enough: be able to rewrite higher order equations/systems as first order systems!
• 15.04.2018: Here is a comment about project #3 reports from Jaakko: "You should provide a simple report where you address the assignments with few illustrative images in each case and answer the posed questions. Please attach your code as well."
• 11.04.2018:
• Please fill out the following form/poor man's doodle ASAP, in order to determine the date for Q&A before the exam.
• Here are the results of your evaluation of project #2.
• 05.04.2018: Project 2 has been graded. Some of you have not marked the project with your student number, and as a result you may not be able to find yourself in the list. If you have questions related to grading please contact Haakon Thømt Simensen
• 04.04.2018: Project 3 is online. Please pay attention to timeplan - it is nontrivial.
• 21.03.2018: Please fill out the evaluation form for Project 2!
• 19.03.2018:
• We are in S5/Tue, EL6/Wed
• Project #3 will be carried out during weeks 15-16. During week 14 we have a lecture on 04.04. I will talk about fast Fourier transform (FFT)
• 14.03.2018: Just wanted to apologize again for deriving a different numerical method at the lecture from what I had in mind, and to clarify the derivations.

The method(s) I have derived during the lecture are midpoint methods. Here we begin with the central difference approximation $(w_{i+1}-w_i)/h \approx y'(t_{i+1/2}) = f(t_{i+1/2}, y(t_{i+1/2}))$. If we then approximate $y(t_{i+1/2}) \approx (w_{i+1}+w_i)/2$ we obtain an implicit method, where we need to solve a system of equations at each time step to find $w_{i+1}$. We could also use forward Euler to approximate $y(t_{i+1/2}) \approx w_i + h/2 f(t_i,w_i)$, which results in the explicit method $w_{i+1} = w_{i} + h f(t_i+h/2,w_i + h/2 f(t_i,w_i))$.

The trapezoid methods begin with the same central difference approximaiton, but then use a slightly different approximation of $y'(t_{i+1/2}) = f(t_{i+1/2}, y(t_{i+1/2}))$. Namely we write $f(t_{i+1/2}, y(t_{i+1/2})) \approx [f(t_{i}, y(t_{i}))+f(t_{i+1}, y(t_{i+1}))]/2 \approx [f(t_i,w_i)+f(t_{i+1},y(t_{i+1}))]/2$. If we use $y(t_{i+1}) \approx w_{i+1}$ in the last formula, we obtain a (semi-)implicit trapezoid method. We could also use forward Euler to approximate the same quantity: $y(t_{i+1}) \approx w_i + h f(t_i,w_i)$, which results in the explicit trapezoid method: $w_{i+1} = w_{i} + h/2 [f(t_i,w_i) + f(w_i + h f(t_i,w_i))]$.

All of these methods have cubic one-step/quadratic global errors. Both explicit methods require only two function evaluations per time step.

• 12.03.2018: Both lectures are in EL6 this week, see https://www.ntnu.no/studier/emner/TMA4320#tab=timeplan
• 08.03.2018:
• The deadline for project #2 is set to the end of week 11, that is to Sunday March 18th 23.59:59.
• There is some confusion in project 2 when it comes to the protein folding model at T=0K. Indeed, the method isn’t properly defined at T=0K as one then tries to divide by 0. The most realistic physical limit is to let T approach 0 from above, but never reach 0. For instance, choose a lower limit at T=0.000001K (or some other very low T). If you have already done the simulations and obtained results using T=0K, I will not reduce your score for this, as it is a consequence of ambiguous formulations in the exercises. However, I encourage all of you to do as I suggest here, as the results will probably be easier to interpret and discuss in this limit. /Haakon
• 05.03.2018:
• The grading for project #1 is done. Amongst those who have submitted the report, the average grade is ca 77%, which is pretty high! I will bring the reports with our comments tomorrow.
• Here is the summary of your feedback for project #1. Thanks to reference group for collecting it, and of course to everyone who has filled out the feedback form. It appears that there are two main issues that many of you raise: the amount of work spent on the project and its imprecise description (many typos etc). With respect to the latter: we work pretty hard to make sure that the description is as accurate as possible and as free of any errors as possible. The text is checked by several "eyes" before it is given to you. Regardless of this, typos make it through (I found and corrected a couple of minor ones on the first day of the project). Personally neither I not my colleagues find the description imprecise, except for the question 8, where it is done on purpose. I would appreciate if you would provide specific examples of imprecisions, and we would try to incorporate your feedback into future project formulations. In any event, I am pretty happy about the fact that the absolute majority of you (87%) have figured the project out, even if with some help! About the amount of work: I do not want to sound dismissive here, but I know from experience with this course that many of you learn quite a bit about programming from project #1, and it goes easier in future projects (including those in year 3 and forward). So hang on, things will feel easier, because you know how to do them now. If you want to discuss project 1 more, please let me know!
• 28.02.2018: Haakon has booked room E3-128 in Realfagbygget 8-10 next Tuesday (March 6th) where he will be available for questions regarding the protein folding project. He encourages you first and foremost to go the exercise hours and ask the student assistants, especially if you have questions concerning e.g. programming technicalities. Additionally, Haakon will come by every exercise session, so you also have the opportunity to meet him there. Rettningslinjer for prosjekt 2 are published now (prosjekter)
• 23.02.2018: Haakon har oppdaget en liten «feil» i oppgaveteksten til prosjekt 2. Alle størrelser bør plottes på temperaturintervallet mellom T=0K og T=1500K, mens det opprinnelig stod opptil T=750K i oppgaveteksten. Dette har ingen store konsekvenser for oppgaven annet enn at det blir lettere å svare på hva som skjer for høye temperaturer i oppgave 2.3. Oppgaveteksten er oppdatert nå.
• 22.02.2018: Just a reminder: weeks 9-10 are dedicated to working on Project #2. We will be in S5/S2 on Tue morning/afternoon, and then in VE1 on Wed.
• 20.02.2018:
• Please take a moment of your time and fill out the evaluation form for project 1. We will discuss the results with the reference group or in general plenum, if necessary. I particularly encourage you to write any comments or suggestions in the evaluation - however constructive suggestions are the most helpful, as it is difficult to know what to change in the course with "this project #\$%^" type comments…
• Project 2 is now published
• These days, a national student survey (Studentenes Helse- og Trivselsundersøkelse 2018) is taking place. In total, 160.000 Norwegian students receive this survey, which focuses on psychosocial issues. Please participate! Link to the survey: https://studenthelse.no/heltaerlig
• 16.02.2018:
• The two remaining lectures on numerical methods for solving linear algebraic systems will be in S5 on Tue and EL6 on Wed. Project 2 will be presented on Tue at 16:15 in S2.
• I was asked to inform you about summer job possibilities, which could be relevant for math/phys students, see here. People from IASTE will tell you more about those during the break on Tue lecture.
• 13.02.2018: Vet ikke om dere har sjekket https://www.ntnu.no/studier/emner/TMA4320#tab=timeplan, men vi er i S5 i dag og de fleste tirsdager.
• 12.02.2018: Several of you have asked what happens if you get "F" on one of the projects. I require that you score a passable grade (at least E - i.e. approximately more than 40% of the questions are answered in a satisfactory fashion) in order to test the learning objectives (3) implement certain numerical algorithms on a computer; and (5) communicate the findings in a form of a report. The rest in my opinion can be assessed at the written exam, to some degree. If you fail all the projects you will not be able to take the exam. If you pass at least one of them, you will be able to take the exam.
• 07.02.2018: Sølve has booked room 734 in SBII from 10 to 14 tomorrow. Just remember there are only two of us, so there is only so much helping we can do.
• 06.02.2018: Just to make it absolutely clear: in all project questions except Q1 and Q8, we use the test example $\rho = \sin(\omega x) \exp(\gamma x)$ and the corresponding measurements $F$, which are implemented in the module test_example.py . This is clearly stated in the last paragraph on page 4 and then pretty much repeated in each question, but since many of you keep wondering, I hereby confirm in writing that this is, indeed, the case. Please see the message from 31.01 about how to load this analytical example from file.
• 05.02.2018: We meet in EL6 tomorrow morning. Because of the long tables it may be difficult for you/us to get in and out, so perhaps we should try to spread out and use the work areas around the auditorium to make the logistics easier? Unfortunately there are no other auditoriums of the right size available this week on Tue 10-12, so please work with what we have. We are in S2 in the afternoon (Tue 16-18) and in VE1 on Wed 10-12.
• 01.02.2018: If you have problems with numpy.geomspace you may need to update your installation of numpy to at least version 1.12 (Kudos to Sindre for sharing this tip). See the instructions on how to update packages using Anaconda Navigator. You could also update directly from command line in Conda terminal: conda update --all --yes
• 31.01.2018:
• Use numpy arrays instead of python lists! They are much faster, and there are plenty of useful functions (inner products, vector norms, matrix-vector and matrix-matrix products, solution of linear systems) which are available. There are plenty of examples on the wiki and also elsewhere on the web on how to use these.
• In order not to perform a high order Taylor series expansion every time you use the analytical solution from the test example, you could write something like this:
import pickle
from test_example import analytical_solution
# define a,b,omega,gamma,Nmax
try:
F = pickle.load( open( "F.pkl", "rb" ) )

except:
F = analytical_solution(a,b,omega,gamma,Nmax)
pickle.dump( F, open( "F.pkl", "wb" ) )
• 30.01.2018:
• There have been some tiny typos relating to the limits of the summation indices in equations (5) and (7) in the project description - these have been corrected. Additionally, some of you found the reference to "y-axis" in question 2 confusing with y as the integration variable in the antiderivative. I replaced the wording with "ordinate axis". Finally, there is one "selv" on line 64 in "test_example.py" which should have been "self". It does not affect the execution of the code (or else I would have seen it), but is a reflection of (1) how long I have lived in Scandinavia and (2) "anything goes" attitude of scripting languages such as Python as opposed to compiled languages such as C++. The updated project description/test_example.py have been uploaded. I emphasize that the changes are really cosmetic, but embarrassing nonetheless.
• I saw that some of you had problems with Anaconda. The whole point of us insisting that you use Anaconda and providing instructions was exactly that you do not run into Python/package installation problems. Still, it is important for us to try to update our instructions for the python installation. I would therefore appreciate if you send me an email briefly describing the issue. I need to know your operating system/version, whether you had any problems with installing anaconda, and what happens when you try to use it. Please put "TMA4320 Anaconda" in the subject line. I appreciate your help with this.
• Regarding the grading: "Prosjekt 1: Den eneste vi skal vurdere er hvorvidt dere har besvaret de spørsmålene som står i teksten på et korrekt og fornuftig måte. Rapporten skal være lesbar, forståelig, og troverdig uten at vi må kikke i koden. Figurene skal inneholde labels og forklaringer til hva er hva. Når det er sagt er det helt mulig at vi skal gi tilbakemeldinger på andre ting (hvis vi har tid) som f.eks. koden, effektiviteten av den, rapportens struktur og presentasjon - men dette skal ikke reflekteres i karakteren. Prosjektene 2 og 3 blir vurdert av FYS. Jeg antar at de skal forklare hvilke vurderingskriterier de bruker."
• Office hours: I will make myself available on Thursdays from ca 12:15 to ca 14:00 at 1038, SBII. See here for more info/updates.
• 24.01.2018:
• Next week there are will be no lectures; instead you will work on project 1 (the description has been published already).
• We are in S5 and not EL6 on 30.01!!!
• I have also booked VE1 on 31.01 and 07.02 - I think it will be better for project work than EL6. The road there is a bit convoluted because of the construction work, so let me know how it works for you.
• I rushed through the derivations today a little, so there have been two mistakes spotted by some of you: the sign mistake in the adaptive estimation of the quadrature error (less dangerous as one is interested in the magnitude and not the sign of the error in the adaptive refinement process), and the weights for a 2-point Gauss quadrature on [-1,1] are of course 1 and not 1/2 (otherwise even the constant polynomial is not integrated correctly). I apologize for these mistakes and thank you for paying attention!
• Some (python) programming tips:
• Do not be afraid of using longer and more descriptive variable and function names - they help a lot when one (including yourself) looks at your code. Ideally, you should also comment your code, but already utilizing descriptive names helps a lot.
• Instead of writing long blocks of code and nested loops, organize chunks of your code into separate functions and modules. This makes the code much more readable, easier to debug and alter. E.g., in this project you could write, debug and test (and then pretty much forget about) helper functions performing integration, assembly of matrices and vectors, evaluation of interpolation polynomials etc. Remember to test the functions as you are writing them, so you could rely upon the result when you call those!
• If you want to use function (or class) foo from a python file bar.py you do not have to copy the code; insted you could simply write from bar import foo
• Use numpy arrays instead of python's lists and tuples for storing vectors/matrices of numbers; they work much faster.
• 09.01.2018:
• There is apparently at least one student who does not speak Norwegian enrolled, so I will try to post in English on this wiki
• I understand that the time-slot on Tue 16:15-18:00 is not convenient for MTFYMA. The unfortunate part is that this is the only possible slot available, which does not collide with other courses for either MTFYMA or BFYS. The only possible solution would be to have different time slots for different programs, but I am not sure we have manpower for this at present. In conclusion, I am aware of the issue, and if it becomes unbearable - please let me know (e.g. through the reference group).
• If there are any medical conditions or travel plans, which are affecting the project deadlines - please send me an email in advance, we will work something out.
• I need at least one representative from BFys in the reference group, and possibly someone from MLREAL (if anyone is taking this course from MLREAL).
• 08.01.2018:
• Velkommen til emnet!
• Vennligst installer Anaconda Python ver. 3 (les instruksjoner her)
• Den første forelesning er i morgen, 09.01 kl 10:15-12:00 i EL6.
• De ikke-obligatoriske øvinger begynner også i morgen, kl 16:15-18:00 i S2
• Husk å stille deg som representant i referansegruppe! Send en email til meg snarest mulig. Jeg trenger representatner fra både MTFYMA og BFys.
• Hvis dere trenger hjelp med å finne en samarbeidspartner(e) for de kommende prosjektene - sen gjerne en email til meg. Jeg forventer at dere arbeider i grupper av 1-3 personer.

## Lærebok, læringsmål, og pensum

• Lærebok: T. Sauer, Numerical analysis: Amazon, Akademika. Dere kan også kikke på H.P. Langtangen, A Primer on Scientific Programming with Python: Springerlink
• Oversettelser av Matlab koder fra læreboken til Python: finnes her
• Læringsmål: finnes her
• Pensum (ver. 24.04.2018):
• Kapittel 0: 0.1, 0.5. Ha en ide om 0.3-0.4.
• Kapittel 1: 1.1-1.5.1
• Kapittel 2: 2.1-2.4; 2.5 (se notatet istedenfor 2.5.3). 2.6.1-2.6.2. 2.7.1
• Kapittel 3: 3.1-3.3
• Kapittel 5: 5.2, 5.4
• Kapittel 6: 6.1-6.6

## Timeplan

• "Teoriuker": uke 2-4, 7-8, 11-12
• "Prosjektuker": uke 5-6, 9-10, 14-16
• Repetisjon: uke 16-17

## Forelesninger (plan med forbehold om endringer)

Dato Tema Avsnitt
Uke 16-17 DFT + spm& svar
24.04 Spm & svar
Uke 15-16 Prosjekt 3
04.04 Definisjon av DFT. FFT algoritmen. jupyternotebooks.zip 10.1
Uke 11-12: Løsning av ordinære diff. likninger (ODEs), Kapittel 6
21.03 Mer om implisitte metoder. Adaptivitet. heat_cond.py, int_flame_rk12.py 6.6, 6.5
20.03 Runge-Kuta metoder. Implisitte metoder. Slides, int_flame.py 6.4, 6.6
14.03 Analyse av numeriske metoder for ODEr, trunkeringsfeil. Systemer av ODEer. 6.2, 6.3
13.03 Initial verdi problemer; eksistens og entydighet of løsninger, osv. Slides 6.1
Uke 9-10 Prosjekt 2
Uke 7-8: Løsning av systemer av likninger, Kapittel 2
21.02 Konvergens av matrix-splitting type metoder: se notatet. Metoder for ikke-lineære systemer 2.5.3, 2.7.1
20.02 Et par ord om Cholesky faktorisering. Intro til sparse matriser. Iterative metoder. matricer.zip, sparsemat_example.py, jacobi_example.py 2.6.1-2.6.2, 2.5.4, 2.5.1-2.5.2
14.02 Feilkilder. LU med pivotering. 2.3-2.4
13.02 Systemer av lineære algebraiske likninger. LU-faktorisering. Python implementasjonen som beregner LU faktoriseringen "radvis" lu1.py og "kolonvis", v.h.a. rad-operasjoner (nettop som i boken og i Gauss eliminasjon) lu.py 2.1-2.2
Uke 5-6 Prosjekt 1
Uke 4: Numerisk integrasjon, Kapittel 5
23.01 Newton-Cotes kvadratur; Trapez, Midtpunkt og Simpsons kvadraturer midpoint.py, trapezoid.py, Jupyter notebooks Simpson.ipynb Example 1.ipynb Middelverdisetningen for integraler 5.2
Uke 3: Interpolasjon, kapittel 3
17.01 Newtons representasjon av interpolasjonspolynomet, interpolasjonsfeil, Runge fenomen. lagrange_poly.py, nest.py, newtdd.py, uniform_interp.py, error_poly.py, chebyshev.py, chebyshev_poly.py 3.2, 3.3
16.01 Interpolasjon: Lagrange polynomet og Newtons representasjon av det 3.1, 3.2 (Teorem 3.5), 0.1
Uke 2: Løsning av likninger, kapittel 1
10.01 Løsning av likninger: konvergens av fikspunktiterasjonen, Newtons metode (og litt om sekant metode), litt om nøyaktighet. fpi.py, newton_div.py 1.2.3, 1.4-1.4.1, 1.3, 1.5.1
09.01 Introduksjon til emnet. Løsning av likninger: halveringsmetode, fikspunktiterasjon. bisect.py, fpi.py, lecture01.pdf 1.1, 1.2

## Prosjekter

Oppgavene skal gjennomføres av grupper på 1-3 studenter.

• Prosjekt 1: **Gravity surveying** (uker 5-6, 29.01.18-11.02.18). Nøkkelord: interpolasjon, numerisk integrasjon, spektrale metoder, inverse problemer, Fredholm integrallikninger, Tikhonov relularisering
• Analytical measurement (used in the project): test_example.py; measurement files + example on how to read those (used in the last question): measurement_files.zip, read_measurement.py
• Rettningslinjer: Oppgaveteksten inneholder konkrete spørsmål. Dere kan skrive rapporten i form av å svare på de ulike oppgavene fortløpende, med etterspurte figurer og korte forklarende tekster til – i en pdf-fil. I tillegg leverer dere inn programfilene (ZIP-fil)
• Innleveringssystemet; frist: 11.02.18, 23:59
• Possible report template: LaTeX source and the resulting PDF.
• Prosjekt 2 Protein folding, Monte-Carlo (weeks 9-10).
• Deadline: end of week 11
• Slides from the project presentation: tma4320_project_2_presentation_1_.pdf
• Rettningslinjer:
• The guidelines for the grading of project 2 will mostly follow the guidelines which are posted on the webpage for project 1. We will assess whether all the questions (marked with bold font) are given a reasonable answer, so be sure that you answer all the questions!
• It is important that you not only present your results, but also discuss them, as is encouraged by the questions. Are your results reasonable, indicating that the implementation of the algorithm is correct? If not, why do you think the results are not reasonable? The background information you should base this discussion on is given in the introduction chapters of the problem text, and you are not expected to have knowledge above this level. We will assess whether your discussion of the results is reasonable, meaning that you may also get points for being critical of your own results if that is reasonable.
• You are allowed to enumerate the answers as we have enumerated the questions (1.1., 1.2., etc.) and give reference directly to each question, and this will in general make it easier for us to recognize the specific answers. The report should otherwise be self-contained, and contain no references to the introduction chapters of the project description, and neither be dependent upon us reading your code. All figures must be labelled, and contain a satisfactory caption which describes the content of the figure. You should also hand in your code, but the quality of this will not affect the grade. You may write the report in either Norwegian or English.
• Prosjekt 3 (uke 15-16): Planetary motion
• Deadline: end of week 16
• Project description: here and here
• Office hours: Thursdays 12-14, Realfagbyg D4-161

## Informasjon om eksamen

• Prosjekter: 30%, skriftlig eksamen: 70%
• Varighet: 4 timer
• Tillatte hjelpemidler (hjelpemiddelkode C):
• Enkel kalkulator
• Rottmann: Matematisk formelsamling

## Gamle eksamensoppgaver og likn.

• Eksamen fra 08.08.2017 med løsningsforslag.
• Eksamen fra 20.05.2017 med løsningsforslag.
• Eksamen fra xx.08.2016 med løsningsforslag.
• Eksamen fra 06.06.2016 med løsningsforslag.
• Her kommer flere treningsoppgaver fra MA2501. Prøv også å løse DFT oppgaver fra læreboken. Løsningene kan finnes her.
• Her ligger et set av de gamle eksamensoppgaver fra TMA4215 som man kan bruke til å forberede til eksamen i TMA4320 (mangler DFT delen; løsning av likninger er "tynn"). De fleste løsninger kan finnes her
• Jeg har lagt sammen noe oppgaver slik at dere kan forberede til den skriftlige eksamen: se her. Løsningsforslaget.