Zależność (10.3) z warunkami brzegowymi (10.2) poprawnie wyznacza podciąg MerLCS dla ciągów Z, A, B

Dowód. Dowód jest identyczny jak dowód twierdzenia 1 z [116] z jednym wyjątkiem. Dla do-pasowania (sytuacji, w której zi= aj lub zi= bk) nie można zapomnieć o możliwości, że skró-cenie jednego z ciągów o jeden symbol może doprowadzić do lepszego wyniku niż podążenie

za dopasowaniem i skrócenie dwu ciągów. 

Podstawową ideą algorytmu Penga i in. [173] jest ograniczenie obliczeń tylko do komó-rek macierzy M reprezentujących dopasowania, których zwykle jest znacznie mniej niż nmr.

10.2. Algorytm równoległości bitowej 195 Zwłaszcza w przypadku dużych alfabetów macierz programowania dynamicznego jest rzadka.

Czas potrzebny na obliczenie pojedynczej komórki macierzy może być wprawdzie większy niż O(1), ale sumarycznie dla całego algorytmu złożoność czasowa jest O(`mr), gdzie ` jest dłu-gością ciągu wynikowego. Zwykle, dla dużych alfabetów ` jest znacznie mniejsze niż n, dzięki czemu algorytm ten może być istotnie szybszy.

10.2. Algorytm równoległości bitowej

10.2.1. Podstawowe idee

Niniejszy podrozdział zawiera opis algorytmu równoległości bitowej, który został zapro-ponowany przez autora w [70]. Zanim jednak zostanie przedstawiony ten algorytm, konieczne jest dowiedzenie kilku lematów.

Lemat 10.2. Wartość komórki M(i, j, k) jest równa bądź większa o 1 niż wartość dowolnego z następujących jej sąsiadów: M(i − 1, j,k), M(i, j − 1,k), M(i, j,k − 1).

Dowód. M(i, j, k) jest długością podciągu MerLCS dla Zi, Aj, Bk. Wymienione komórki są-siednie zawierają długości podciągu MerLCS dla tych samych ciągów, przy czym jeden z nich jest okrojony o ostatni symbol. Niemożliwe jest, aby okrojenie jednego ciągu o ostatni sym-bol spowodowało zmniejszenie długości podciągu MerLCS o więcej niż 1, ponieważ usunięcie ostatniego symbolu z jednego ciągu może co najwyżej spowodować, że ten symbol nie znaj-dzie się w ciągu wynikowym. Podobnie, nie jest możliwe, żeby okrojenie jednego z ciągów wejściowych spowodowało, że długość podciągu MerLCS wzrośnie. 

Na podstawie lematu 10.2 możliwa jest reprezentacja trójwymiarowej macierzy M przez dwuwymiarową macierz G, zawierającą wektory zmian będące wektorami liczb całkowitych z zakresu [1, r]. Wektory te zdefiniowane są następująco:

i∈ G( j,k) ⇐⇒ M(i, j,k) − M(i − 1, j,k) = 1 dla 1 ≤ i ≤ r. (10.4) Równoważność pomiędzy macierzami M i G wynika z faktu, że:

M(i, j, k) = |{x : x ∈ G( j,k) ∧ x ≤ i}|. (10.5)

Łatwo można zauważyć, że definicja wektorów zmian jest analogiczna do sposobu reprezentacji zmian w macierzy programowania dynamicznego dla problemu LCS przedstawionej w [120].

Zanim zostanie sformułowany bezpośredni algorytm wyznaczania G( j, k), dogodne jest za-pisanie zależności (10.3) w następujący sposób:

M⁰(i, j − 1,k) = max









M⁰(i − 1, j − 1,k),

M(i − 1, j − 1,k) + 1, jeśli zi= aj, M(i, j − 1,k),

(10.6)

M⁰⁰(i, j, k − 1) = max









M⁰⁰(i − 1, j,k − 1),

M(i − 1, j,k − 1) + 1, jeśli zi= bk, M(i, j, k − 1),

(10.7)

M(i, j, k) = max

(M⁰(i, j − 1,k),

M⁰⁰(i, j, k − 1). (10.8)

Można zauważyć, że pierwsze składniki funkcji maksimum w (10.6)–(10.7) oznaczają tylko tyle, że wartości M⁰, M⁰⁰, M nie maleją przy wzroście i (por. trzeci składnik funkcji maksimum w (10.3)). Wymaganie to jest inherentne dla definicji G.

Wektor zmian G⁰( j −1,k) dla M⁰(i, j −1,k) może być wyznaczony na podstawie G( j −1,k) następująco:

1. G⁰( j − 1,k) ← G( j − 1,k).

2. Dla każdego i od r do 1 takiego, że zi= ajwykonaj:

a) jeśli i ∈ G⁰( j − 1,k) nie rób nic,

b) w przeciwnym przypadku wstaw i do G⁰( j − 1,k) i usuń z G⁰( j − 1,k) najmniejszą war-tość całkowitą (jeśli taka istnieje) większą niż i.

W podobny sposób obliczany jest wektor G⁰⁰( j, k − 1).

Lemat 10.3. Powyższy algorytm poprawnie oblicza wektory zmian G⁰( j − 1,k) i G⁰⁰( j, k − 1) odpowiadające M⁰(i, j − 1,k) i M⁰⁰(i, j, k − 1) dla wszystkich poprawnych wartości i.

Dowód. Rangą h(x) elementu x w wektorze zmian będzie nazywany indeks tego elementu w posortowanym zbiorze wartości znajdujących się w wektorze zmian. Można zauważyć, że dla każdego i w G⁰( j − 1,k) o randze h(i) zachodzi:

• M(i, j − 1,k) = h(i), a więc i ∈ G( j − 1,k) lub

• M(i − 1, j − 1,k) = h(i) − 1 i zⁱ= aj.

Wobec tego i należy do G⁰( j − 1,k) wtedy i tylko wtedy, gdy:

• i ∈ G( j − 1,k) i nie istnieje i⁰< i takie, że M(i⁰− 1, j − 1,k) = h(i⁰) − 1 i zi⁰ = aj lub

• i jest najmniejszym indeksem takim, że M(i − 1, j − 1,k) = h(i) − 1 i zⁱ= aj.

Dowód dla G⁰⁰( j, k − 1) jest analogiczny. 

10.2. Algorytm równoległości bitowej 197 Łatwo można zauważyć, że powyższy sposób wyznaczania wektorów zmian jest dokładnie taki sam jak pokazany w [120, Observation 4] dla algorytmu LCS-BP-Length rozwiązującego problem LCS.

Dysponując wektorami G⁰( j − 1,k), G⁰⁰( j, k − 1), można wyznaczyć G( j,k). Na podstawie (10.5) oraz (10.8) wiadomo, że:

M(i, j, k) = max

(|{x : x ∈ G⁰( j − 1,k) ∧ x ≤ i}|

|{x : x ∈ G⁰⁰( j, k − 1) ∧ x ≤ i}| . (10.9) W związku z tym, aby obliczyć G( j, k), należy postępować jak poniżej:

1. Dla każdej pary liczb całkowitych z G⁰( j − 1,k) i G⁰⁰( j, k − 1) jednakowej rangi weź mniej-szą z nich i wstaw do G( j, k).

2. Jeśli wśród elementów w G⁰( j − 1,k) i G⁰⁰( j, k − 1) występuje element o unikalnej randze, to wstaw go do G( j, k).

10.2.2. Algorytm

Do efektywnej reprezentacji wektorów zmian mogą być zastosowane wektory bitowe. Po-nieważ wszystkie wektory zmian w proponowanym algorytmie zawierają liczby całkowite z za-kresu [1, r], to długość każdego z używanych wektorów bitowych wynosi r. (Jedynie dla wy-gody prezentacji bity w wektorach bitowych będą numerowane od 1). Dla każdego wektora zmian G( j, k), G⁰( j −1,k), G⁰⁰( j, k−1) istnieje odpowiedni wektor bitowyW ( j,k),W⁰( j −1,k), W⁰⁰( j, k − 1) zdefiniowany następująco: wszystkie bity mają wartość 1, z wyjątkiem tych, któ-rych indeksy występują w odpowiednim wektorze zmian.

Używana będzie także tablica wektorów bitowych Yx, dla każdego x ∈Σ. Każdy wektor Yx

zdefiniowany jest następująco: wszystkie bity mają wartość0, z wyjątkiem tych, które odzwier-ciedlają pozycję symbolu x w T . Te wektory bitowe reprezentują pozycje, na których występują dopasowania przy porównywaniu ciągu T z konkretnym symbolem z A lub B. Jest to dokładnie ta sama reprezentacja jak stosowana w algorytmie LCS-BP-Length (por. podrozdz. 7.3.4).

W celu obliczenia W⁰( j − 1,k) na podstawie W ( j − 1,k) używa się wektora Y^aj, jako że reprezentuje on pozycje, na których w Z występuje aj. Ponieważ operacje wymagane w tym miejscu są identyczne do wykonywanych w algorytmie LCS-BP-Length, więc zostanie użyta następująca sekwencja operacji bitowych z [120]:

W⁰( j − 1,k) ← W ( j − 1,k) & Y^aj,

W⁰( j − 1,k) ← (W ( j − 1,k) +W⁰( j − 1,k)) | (W ( j − 1,k) −W⁰( j − 1,k)).

Analogiczną sekwencję operacji stosuje się do wyznaczenia W⁰⁰( j, k − 1).

Obliczenie W ( j, k) na podstawie W⁰( j − 1,k) i W⁰⁰( j, k − 1) jest nieco bardziej skompliko-wane. Przydatne okażą się tu poniższe lematy.

Lemat 10.4. Dla dowolnego x ∈ G⁰( j − 1,k) o randze h(x) i y ∈ G⁰⁰( j, k − 1) o randze h(y) zachodzi:

a) jeśli x> y, to h(x) ≥ h(y), b) jeśli x< y, to h(x) ≤ h(y), c) jeśli x= y, to |h(x) − h(y)| ≤ 1.

Dowód. Dowód zostanie przeprowadzony przez zaprzeczenie.

Przypadek a: Niech x> y i h(x) < h(y). Z (10.6) i (10.7) wiadomo, że

h(x) = M⁰(x, j − 1,k) = M(x − 1, j − 1,k) + 1, (10.10) h(y) = M⁰⁰(y, j, k − 1) = M(y − 1, j,k − 1) + 1. (10.11) Z założenia h(x) < h(y) wynika, że

M⁰(x, j − 1,k) < M⁰⁰(y, j, k − 1). (10.12)

Na podstawie (10.8) wiadomo

M(y, j, k) ≥ M⁰⁰(y, j, k − 1), (10.13)

a z x > y wynika, że

M(x − 1, j − 1,k) ≥ M(y, j − 1,k). (10.14)

Łącząc wyniki z (10.10)–(10.14), otrzymuje się:

M(y, j, k) ≥ M⁰⁰(y, j, k − 1) > M⁰(x, j − 1,k) > M(x − 1, j − 1,k) ≥ M(y, j − 1,k), (10.15) a więc

M(y, j, k) − M(y, j − 1,k) > 1, (10.16)

co jest niemożliwe z uwagi na lemat 10.2.

Przypadek b: Dowód jest analogiczny jak dla przypadku a.

Przypadek c: Niech|h(x) − h(y)| > 1. Z (10.10) i (10.11) wynika, że

|M(x − 1, j − 1,k) − M(x − 1, j,k − 1)| > 1. (10.17) Ponieważ zarówno M(x − 1, j − 1,k), jak i M(x − 1, j − 1,k) są sąsiadami M(x − 1, j,k), więc na podstawie lematu 10.2 wartość powyższej różnicy nie może być większa niż 1. 

10.2. Algorytm równoległości bitowej 199 Lemat 10.5. W celu obliczenia G( j, k) na podstawie G⁰( j − 1,k) i G⁰⁰( j, k − 1) wystarczy:

1. Znaleźć sumę zbiorów G⁰( j − 1,k) i G⁰⁰( j, k − 1) otrzymując wielozbiór G^∗( j, k), w którym niektóre liczby całkowite mogą występować dwukrotnie.

2. Usunąć z G^∗( j, k) liczby całkowite o parzystych rangach otrzymując zbiór G( j, k).

Dowód. Dowód indukcyjny względem rangi x z G( j, k). Dla x = 1 najmniejsza liczba całkowita w wielozbiorze G^∗( j, k) zostaje przepisana do G( j, k) (jest to element o randze 1 w G⁰( j − 1,k) lub G⁰⁰( j, k − 1)). Kolejny najmniejszy element z G^∗( j, k) (rangi 2) musi być elementem o ran-dze 1 odpowiednio w G⁰⁰( j, k − 1) lub G⁰( j − 1,k) (lemat 10.4), a więc w G( j,k) znajdzie się mniejszy z elementów o randze 1 spośród wejściowych zbiorów.

Niech teraz x > 1 i dla elementów o rangach 2x − 3 i 2x − 2 z G^∗( j, k) jeden z nich jest elementem o randze x − 1 w G⁰( j − 1,k), a drugi elementem o randze x − 1 w G⁰⁰( j, k − 1).

Dla elementów o rangach 2x − 1 i 2x w G^∗( j, k), jeden z nich musi być elementem o randze x w G⁰( j − 1,k), a drugi elementem o randze x w G⁰⁰( j, k − 1) (lemat 10.4). Mniejszy z tych ele-mentów jest wstawiany do G( j, k). Oczywiście, jeśli w G^∗( j, k) znajdują się dwa identyczne ele-menty, ich rangi różnią się parzystością, więc dokładnie jeden z nich jest wstawiany do G( j, k),

dzięki czemu G( j, k) jest zbiorem. 

Ponieważ liczby całkowite występujące w obu wektorach, G⁰( j − 1,k) i G⁰⁰( j, k − 1), muszą należeć do G( j, k), więc można usunąć je z G^∗( j, k) po pierwszym kroku i dodać na samym końcu (tylko po jednej kopii), po usunięciu elementów o parzystych rangach, do G( j, k). Można to wykonać za pomocą operacji różnicy symetrycznej:

W( j, k) ← W⁰( j − 1,k) ⊕ W⁰⁰( j, k − 1).

Kluczowym etapem algorytmu wyznaczania W ( j, k) na podstawie W⁰( j − 1,k) i W⁰⁰( j, k − 1) jest wpisanie wartości 0 do wszystkich bitów wartości 1 o parzystych rangach. W [215, stro-ny 74–77] przedstawiostro-ny jest następujący algorytm wyznaczający „parzystość” w 32-bitowym słowie komputerowym x (jego adaptacja dla innych wartości w jest natychmiastowa):

y← x ⊕ (x  1), y← y ⊕ (y  2), y← y ⊕ (y  4), y← y ⊕ (y  8), y← y ⊕ (y  16).

Bit o indeksie i w y ma wartość1 wtedy i tylko wtedy, gdy liczba bitów o wartościach 1 na pozycjach od 0 do i w x jest parzysta. Dzięki temu, aby usunąć bity 1 o parzystych rangach z wektora bitowego W ( j, k), należy wyznaczyć parzystość tego wektora i wykonać maskowanie bitów1 o parzystych rangach.

Kompletny pseudokod algorytmu rozwiązującego problem MerLCS przedstawiony jest na rys. 10.2. W pierwszej części tego algorytmu (wiersze1–7) rozwiązywany jest problem LCS dla

MerLCS-BP(Z, A, B)

Wejście: Z, A, B – podciągi, dla których wyznaczany jest podciąg MerLCS Wyjście: długość podciągu MerLCS

{Wyznaczanie warunków brzegowych}

1 W(0, 0) ← 0^r 2 for k ← 1 to m do

3 U← W (0,k − 1) & Y^bk

4 W(0, k) ← (W (0,k − 1) +U) | (W (0,k − 1) −U) 5 for j ← 1 to n do

6 U← W ( j − 1,0) & Y^aj

7 W( j, 0) ← (W ( j − 1,0) +U) | (W ( j − 1,0) −U) {Obliczenia właściwe}

8 for j ← 1 to n do 9 for k ← 1 to m do

10 U⁰← W ( j − 1,k) & Yaj

11 W⁰← (W ( j − 1,k) +U⁰) | (W ( j − 1,k) −U⁰) 12 U⁰⁰← W ( j,k − 1) & Y^bj

13 W⁰⁰← (W ( j,k − 1) +U⁰⁰) | (W ( j, k − 1) −U⁰⁰) 14 U← W⁰& W⁰⁰

15 W ← W⁰ ⊕ W⁰⁰

16 V ← W

17 for i ← 0 to dlog2re − 1 do

18 V ← V ⊕ (V  (1 i))

19 W( j, k) ← (W & V ) | U {Wyznaczanie wyniku}

20 ` ← 0; V ← ∼W (n,m) 21 while V 6= 0^rdo

22 V ← V & (V − 1); ` ← ` + 1 23 return `

Rys. 10.2. Algorytm równoległości bitowej wyznaczający długość podciągu MerLCS Fig. 10.2. Bit-parallel algorithm computing MerLCS length

par ciągów Z i A oraz Z i B, w celu wyznaczenia warunków brzegowych (10.2). W wierszach 8–19 wykonywane są obliczenia właściwe, tj. dla każdej pary ( j, k) wyznaczany jest wektor bitowy W ( j, k). Końcowa część algorytmu (wiersze 20–23) to zliczenie bitów o wartości 0 w wektorze W (n, m), dzięki czemu otrzymuje się długość podciągu MerLCS.

10.2.3. Szczegóły implementacyjne

W algorytmie zaprezentowanym powyżej zakłada się, że r ≤ w, tzn. każdy wektor bitowy mieści się w pojedynczym słowie komputerowym, co rzadko jest spełnione. Emulacja wektorów bitowych jako tablic rozmiaru dr/we słów komputerowych nie jest trudna. Uważnej implemen-tacji wymagają jednak niektóre operacje bitowe i arytmetyczne. Wszystkie operacje odejmo-wania (−), bitowe (| , & , ⊕) w wierszach 3–4, 6–7, 10–15 mogą być łatwo zaimplementowane na tablicach słów komputerowych, ponieważ nie występują tu żadne przeniesienia pomiędzy

W dokumencie Sekwencyjne i równoległe algorytmy znajdowania podciągów (Stron 194-200)