3. Domácí úkol: Definitnost matice

Úvod

V třetím, posledním, úkolu vytvoříme balíček implementující metody pro výpočet typů definitností kvadratických forem. Na tuto úlohu jste jistě narazili při vyšetřování typů extrémů funkcí více proměnných.

Teoretický background

Pro čtvercovou symetrickou a reálnou matici $\mathbf{M} \in \mathbb{R}^{n,n}$ definujeme jí příslušející kvadratickou formu jakožto funkci $q\colon \mathbb{R}^n \to \mathbb{R}$ předpisem

q(x) = x^T \cdot \mathbf{M} \cdot x, \quad x \in \mathbb{R}^n,

kde T označuje transpozici a tečky maticové násobení.

NOTE: Matice je symetrická, právě když je rovna své transpozici.

Různé typy definitnosti matice $\mathbf{M}$ vystihují chování znaménka, případně ne/nulovosti, funkčních hodnot jí příslušející kvadratické formy:

$M$ je pozitivně definitní (PD), právě když pro každé nenulové $x \in \mathbb{R}^n$ platí $q(x) > 0$ ,
$M$ je negativně definitní (ND), právě když pro každé nenulové $x \in \mathbb{R}^n$ platí $q(x) < 0$ ,
$M$ je pozitivně semidefinitní (PSD), právě když pro každé $x \in \mathbb{R}^n$ platí $q(x) \geq 0$ ,
$M$ je negativně semidefinitní (NSD), právě když pro každé $x \in \mathbb{R}^n$ platí $q(x) \leq 0$ ,
$M$ je indefinitní (ID), právě když existují $x \in \mathbb{R}^n$ a $y \in \mathbb{R}^n$ splňující $q(x) > 0$ a $q(x) < 0$ .

K určení typu definitnosti matice existuje několik nástrojů:

základní úprava na čtverce,
využití souvislosti s vlastními čísly matice,
obecné Sylvestrovo kritérium využívající determinantů,

existují i další metody založené na faktorizaci matic (tak lze chápat i přístup pomocí vlastních čísel).

Implementační pokyny

Balíček `Definiteness`

Ve svém repozitáři na větvi solution/03-definiteness vytvořte Julia balíček, jehož instalaci bude možné provést v balíčkovacím módu příkazem (v kódu níže USERNAME představuje vaše uživatelské jméno, název vašeho osobního repozitáře)

(@v1.12) pkg> add "git@gitlab.fit.cvut.cz:BI-JUL/B251/USERNAME"#solution/03-definiteness

Pro lokální vývoj samozřejmě doporučuji přidat si balíček ve vývojovém módu pomocí (v kořenovém adresáři repozitáře, případně upravte cestu):

(@v1.12) pkg> devel .

V repozitáři se zadáním jsem nechal pouze tento popis zadání a adresář s testy (viz níže). Vše ostatní nechávám na vaší kreativitě (komentáře, spouštění testů v CI, generátor dokumentace,…).

Metody pro výpočet typu definitnosti

Balíček bude exportovat boolean metody isPD, isPSD, isND, isNSD a isID, které na čtvercové symetrické matici vrací true nebo false podle toho, zda uvedená matice má nebo nemá příslušnou vlastnost.

Metody akceptují matice typu Matrix{T}, kde T je reálného číselného typu (tj. T <: Real). V případě, že matice není čtvercová, nebo není symetrická, pak metody vyvolají výjimku ArgumentError. K porovnání symetrie matic, resp. jejích prvků, použijte metodu isapprox.

Prvky matice mohou být zejména následujících typů:

strojová čísla (např. Float64),
racionální čísla (např. Rational{Int64}),
celá čísla (např. Int64).

Během výpočtu (ať už zvolíte jakoukoliv metodu) budete muset porovnávat čísla s nulou. V exaktní aritmetice racionálních čísel toto není problém. V případě strojových čísel typu T <: AbstractFloat k tomuto účelo použijte porovnání pomocí isapprox s absolutní tolerancí sqrt(eps(T)), tj. strojové x::T budeme považovat za nulové, pokud isapprox(x, zero(T), atol=sqrt(eps(T))) je true. Podobně budeme považovat strojové číslo za kladné, pokud není nulové v uvedeném smyslu a x > zero(T) je true. Analogicky pro zápornost.

Pokud je vstupní matice s celočíselnými prvky (T <: Integer), tak její prvky přetypujeme na Rational{T} a dále počítáme jako s racionálními čísly.

Závěrečné poznámky

Volbu výpočetního algoritmu nechávám zcela na vás (viz výše). Jen dejte pozor, aby algoritmus zvládl pracovat i v exaktní racionální aritmetice.

Metody na výpočet vlastních čísel nebo determinantu samozřejmě nemusíte sami implementovat. Můžete použít metody dostupné v Julia a v jejích balíčcích (jako např. LinearAlgebra).

Pokud to uznáte za vhodné a proveditelné můžete u tohoto úkolu použít i vláknový paralelizmus. Závěrečné porovnání budu spouštět v prostředí s více vlákny.

Spouštění testů

V repozitáři jsem nechal adresář s testy (test), které můžete u správně sestaveného balíčku spustit z balíčkovacího módu takto:

(@v1.12) pkg> test Definiteness
     Testing Definiteness
     ...
     Testing Running tests...
Test Summary: | Pass  Total  Time
isPD          |   13     13  2.2s
Test Summary: | Pass  Total  Time
isPSD         |   14     14  0.2s
Test Summary: | Pass  Total  Time
isND          |   12     12  0.7s
Test Summary: | Pass  Total  Time
isNSD         |   12     12  0.0s
Test Summary: | Pass  Total  Time
isID          |   13     13  0.4s
Test Summary: | Pass  Total  Time
Random        |   15     15  4.0s
     Testing Definiteness tests passed

Součástí vyhodnocení řešení bude porovnání rychlosti i robustnosti vaší implementace vzhledem k ostatními řešením.

Řešení a odevzdání

Opět vytvořte větev odvozenou z větve assignment/03-definiteness a nezvěte ji solution/03-definiteness. Do solution/03-definiteness vložte své řešení editací všeho potřebného. Až budete se svým řešením spokojeni, vytvořte MR (to můžete i dříve, aspoň uvidíte výsledek testů, pokud je nespouštíte lokálně) s cílovou větví assignment/03-definiteness a přiřaďte mě k němu jako assignee. Tímto aktem úkol odevzdáte.