FIT ČVUT Course Pages
BI-JUL — Programování v jazyku Julia
Jdi na navigaci předmětu

1. Domácí úkol: Interpolace

Úvod

V tomto úkolu budeme implementovat jednoduchý interpolační algoritmus. Na problém interpolace narazíme v mnoha formách v počítačové grafice, počítačově asistovaném designu (CAD) nebo v datové analýze. V rámci tohoto dokumentu nejprve popíšeme konkrétní algoritmus a poté shrneme požadavky na implementaci.

Teoretický základ

Nechť je dána konečná neprázdná množina bodů {(xj,yj)}j=1n⊂R2\{(x_j, y_j)\}_{j=1}^n \subset \mathbb{R}^{2}, kde n≥2n \geq 2, a dále liché s∈Ns \in \mathbb{N}. Bez újmy na obecnosti nyní předpokládejme, že x1<x2<⋯<xnx_1 < x_2 < \cdots < x_n. Naším cílem je nalézt m:=n−1m := n - 1 polynomů stupně nejvýše ss tvaru fi(x)=∑k=0sαikxkf_i(x) = \sum_{k=0}^s \alpha_{ik} x^k, i=1,…,mi=1,\ldots,m, a kde αi0,…,αis\alpha_{i0},\ldots,\alpha_{is} jsou zatím neznámé koeficienty. Každý z těchto polynomů fif_i uvažujeme jako definovaný na intervalu ⟨xi,xi+1⟩\langle x_{i}, x_{i+1}\rangle. Na tyto polynomy postupně naklademe následující podmínky:

  1. pro každé i=1,…,mi=1,\ldots,m platí fi(xi)=yif_i(x_i) = y_i a fi(xi+1)=yi+1f_i(x_{i+1}) = y_{i+1} (tj. 2m2m rovnic; polynomy prochází na svých definičních oborech předepsanými body),
  2. pro každé i=1,…,m−1i=1,\ldots,m-1 a pro každé k=1,2,…,s−1k=1,2,\ldots,s-1 platí fi(k)(xi+1)=fi+1(k)(xi+1)f_i^{(k)}(x_{i+1}) = f_{i+1}^{(k)}(x_{i+1}) (tj. (m−1)×(s−1)(m-1)\times(s-1) rovnic; polynomy na sebe hezky navazují),
  3. v krajních bodech (tj. x1x_1 a xnx_n) explicitně předepíšeme pro f1f_1 a fmf_m derivace stupně 11 až (s−1)/2(s-1)/2.

Neznámých koeficientů (α\alpha) máme m×(s+1)m\times(s+1), musí splňovat výše uvedených

2m+(m−1)×(s−1)+s−1=m×(s+1)2m \mathbf{+} (m - 1) \times (s - 1) \mathbf{+} s - 1 = m\times(s+1)

lineárních rovnic pro neznámé koeficienty. Jde o tzv. splajny (spline) lichého stupně ss s předepsanými okrajovými derivacemi.

Tuto soustavu je možné vyřešit pomocí nástrojů lineární algebry. Za výše uvedených podmínek je i známo, že tato soustava má právě jedno řešení (to zde podrobněji rozebírat nebudeme). Tím získáme neznámé koeficienty polynomů fif_i. Pokud následně dostaneme hodnotu x∈Rx\in\mathbb{R}, tak stačí určit do jakého intervalu ⟨xi,xi+1⟩\langle x_i, x_{i+1}\rangle, i=1,…,n−1i=1,\ldots,n-1, patří a následně vyhodnotit funkční hodnotu fi(x)f_i(x).

Implementační pokyny

Typ Interpolate{T}

V souboru src/interpol.jl je předpřipravený modul Interpolation v rámci něhož je potřeba implementovat parametrický typ Interpolate{T}. Typový parametr T omezme na reálná čísla (abstraktní typ Real). Konstruktor tohoto typu bude popořadě přijímat poziční argumenty

  • x: vektor xx-ových souřadnic interpolovaných bodů, hodnoty typu T,
  • y: vektor yy-ových souřadnic interpolovaných bodů, hodnoty typu T

a klíč-hodnota argumenty

  • order (role ss výše): liché přirozené číslo představující maximální stupeň interpolačních polynomů; výchozí hodnota 1,
  • left: vektor předepsaných derivací v levém krajním bodě (od nejmenší po největší; tedy od 11 do (s−1)/2(s - 1) / 2, pro s=1s = 1 je tato podmínka prázdná); výchozí hodnota T[],
  • right: vektor předepsaných derivací v pravém krajním bodě (od nejmenší po největší; tedy od 11 do (s−1)/2(s - 1) / 2, pro s=1s = 1 je tato podmínka prázdná); výchozí hodnota T[].

Vektory x, y, left a right jsou stejného typu Vector{T}. Pokud dostaneme vektory celočíselného nebo racionálního typu, tak všechny hodnoty převedeme do typu Float64.

Pokud argumenty nesplňují uvedené podmínky, dojde k vyvolání výjimky ArgumentError s vhodným textovým popisem. Viz testy v souboru test/test_constructor.jl, konkrétně jde o:

  • typové požadavky uvedené výše,
  • pole x a y musí být stejné délky, stejně tak pole left a right,
  • pole x nesmí obsahovat duplicitní hodnoty (kontrola stačí vzhledem k ==),
  • order je omezen pouze na liché hodnoty Int64,
  • délka polí left a right závisí správně na hodnotě order.

Hodnoty v poli x nemusí být seřazené podle velikosti, v rámci přípravy dat je nutné si tyto hodnoty seřadit. Nezapomeňte na odpovídající seřazení i pole hodnot y!

Objekty typu Interpolate{T} bude možné vyhodnotit v zadaném bodě z::T. Pokud je hodnota z mimo rozsah daný minimální a maximální hodnotou v poli x (včetně), tak dojde k vyvolání výjimky DomainError s vhodným textovým popisem. V opačném případě vrátíme uživateli hodnotu příslušného interpolačního polynomu v bodě z. Pokud dostaneme na vstupu hodnotu jiného typu než T, tak se ji do tohoto typu nejprve pokusíme konvertovat.

Například:

julia> include("src/interpol.jl");
julia> f = Interpolation.Interpolate([1.,2.],[3.,4.], order=1);
julia> f(1.3)
3.3

Další ukázky lze nalézt v adresáři test.

Další náležitosti

Vedle výše uvedeného typu a metod (konstruktor, volání) můžete v modulu implementovat i další pomocné funkce pro zpřehlednění kódu. Je to velmi vhodné k zpřehlednění vašeho kódu. Modul a implementované metody opatřete vhodnými docstringy, aby případný čtenář dokázal pochopit jejich smysl. Snažte se svůj kód psát konzistentně s některým ze stylů nabízených balíčkem JuliaFormatter (např. BlueStyle).

Pokud použijete nástroje umělé inteligence (LLM; doplňování kódu v editoru, chat, nebo agentní mód), tak do textového popisu Merge requestu v Gitlabu sepište stručný report (reflexi použití) obsahující:

  • Jaký model a v jaké verzi jste při práci použili.
  • V jakém módu jste tyto nástroje využívali (doplňování kódu v editoru, web chat, v případě agentního módu uveďte i použitý postroj).
  • Jakým způsobem jste nástroj používali a ověřovali správnou funkčnost kódu (přímé generování, vlastní implementace a LLM review, atp.).
  • Uveďte problémy, na které jste při práci s nástrojem narazili. Například navrhnuté chybné řešení a jak jste ho následně odhalili.

V opačném případě do textového popisu uveďte, že se jedná o "ruční práci". Tento report bude i vstupem pro výsledné hodnocení, na základě tohoto reportu také bude klást důraz na některé aspekty odevzdaného kódu.

Lokální spouštění testů

Toto zadání je úmyslně vytvořeno lehce "na koleni" a nevyužívá například Julia projekty ke správě závislostí (látka, ke které se dostaneme později během semestru).

K lokálnímu spuštění budete proto pravděpodobně potřebovat doinstalovat balíček Glob.

import Pkg; Pkg.add("Glob")

Poté v kořenovém adresáři tohoto zadání stačí z příkazové řádky spustit

$ julia --color=yes test/runtests.jl

Tímto příkazem spustíte všechny testy v souborech test/test_*.jl. Případně tak do adresáře test můžete snadno přidat i další pomocné vlastní testy vlastních pomocných metod, jen je musíte pojmenovat ve tvaru test_*.jl.

Pokud vše dopadne dobře, měli byste vidět následující standardní výstup:

┌ Info:
│ Running tests in test_constructor.jl...
└ ────────────────────────────────────────────────────────────
Test Summary:                                                | Pass  Total  Time
Interpolation: constructor                                   |   15     15  3.4s
┌ Info:
│ Running tests in test_cubic.jl...
└ ────────────────────────────────────────────────────────────
Test Summary:                                                | Pass  Total  Time
Interpolation: cubic, 2 points                               |   22     22  0.0s
Test Summary:                                                | Pass  Total  Time
Interpolation: cubic, 3 points                               |   21     21  0.1s
Test Summary:                                                | Pass  Total  Time
Interpolation: cubic, 25 points                              |   51     51  0.1s
┌ Info:
│ Running tests in test_linear.jl...
└ ────────────────────────────────────────────────────────────
Test Summary:                                                | Pass  Total  Time
Interpolation: handling types during evaluation              |    7      7  0.2s
Test Summary:                                                | Pass  Total  Time
Interpolation: linear, 2 points                              |   10     10  0.0s
Test Summary:                                                | Pass  Total  Time
Interpolation: linear, 3 points                              |    7      7  0.0s

V opačném případě dostanete od Julie vynadáno!

Řešení a odevzdání

Ve svém repozitáři vytvořte větev odvozenou z větve assignment/01-interpolation a nazvěte ji například solution/01-interpolation. Do solution/01-interpolation vložte své řešení editací souboru src/interpol.jl případně přidáním testů do složky test. Až budete se svým řešením spokojeni, vytvořte MR (to můžete i dříve, alespoň uvidíte výsledek testů, pokud je nespouštíte lokálně) a přiřaďte mě k němu jako assignee. Tímto aktem úkol odevzdáte.