diff --git a/.coveragerc b/.coveragerc
deleted file mode 100644
index 753322bd..00000000
--- a/.coveragerc
+++ /dev/null
@@ -1,23 +0,0 @@
-# Configuration for coverage.py (https://pypi.python.org/pypi/coverage)
-
-[run]
-# Enable branch coverage
-branch = True
-
-[report]
-
-# Regexes for lines to exclude from consideration
-exclude_lines =
-
-    pragma: no cover
-
-    # Ignore non-runnable code
-    if __name__ == .__main__.:
-
-    pass
-
-    ...
-
-# Only check coverage for source files
-include =
-    automata/*/*.py
diff --git a/.github/workflows/lint.yml b/.github/workflows/lint.yml
index 1cb65971..6ab7346b 100644
--- a/.github/workflows/lint.yml
+++ b/.github/workflows/lint.yml
@@ -23,8 +23,8 @@ jobs:
         with:
           python-version: "3.11"
 
-      - name: Install graphviz
-        run: sudo apt-get install graphviz
+      - name: Setup Graphviz
+        uses: ts-graphviz/setup-graphviz@v1
 
       - name: Install dependencies
         run: |
diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml
index 62b05cc9..c7b766c1 100644
--- a/.github/workflows/tests.yml
+++ b/.github/workflows/tests.yml
@@ -27,8 +27,8 @@ jobs:
       with:
         python-version: ${{ matrix.python-version }}
 
-    - name: Install graphviz
-      run: sudo apt-get install graphviz
+    - name: Setup Graphviz
+      uses: ts-graphviz/setup-graphviz@v1
 
     - name: Install dependencies
       run: |
diff --git a/automata/fa/dfa.py b/automata/fa/dfa.py
index 89433a26..d7ee9d85 100644
--- a/automata/fa/dfa.py
+++ b/automata/fa/dfa.py
@@ -17,18 +17,15 @@
     Iterable,
     Iterator,
     List,
-    Literal,
     Mapping,
     Optional,
     Set,
     Tuple,
     Type,
-    TypeVar,
     cast,
 )
 
 import networkx as nx
-from pydot import Dot, Edge, Node
 from typing_extensions import Self
 
 import automata.base.exceptions as exceptions
@@ -50,7 +47,6 @@
 class DFA(fa.FA):
     """A deterministic finite automaton."""
 
-    # TODO allow
     __slots__ = (
         "states",
         "input_symbols",
@@ -332,7 +328,7 @@ def minify(self, retain_names: bool = False) -> Self:
         Create a minimal DFA which accepts the same inputs as this DFA.
 
         First, non-reachable states are removed.
-        Then, similiar states are merged using Hopcroft's Algorithm.
+        Then, similar states are merged using Hopcroft's Algorithm.
             retain_names: If True, merged states retain names.
                           If False, new states will be named 0, ..., n-1.
         """
@@ -807,7 +803,7 @@ def predecessors(
         *,
         strict: bool = True,
         key: Optional[Callable[[Any], Any]] = None,
-    ) -> Iterable[str]:
+    ) -> Generator[str, None, None]:
         """
         Generates all strings that come before the input string
         in lexicographical order.
@@ -845,7 +841,7 @@ def successors(
         strict: bool = True,
         key: Optional[Callable[[Any], Any]] = None,
         reverse: bool = False,
-    ) -> Iterable[str]:
+    ) -> Generator[str, None, None]:
         """
         Generates all strings that come after the input string in
         lexicographical order. Passing in None will generate all words. If
@@ -906,7 +902,7 @@ def successors(
             )
             # Traverse to child if candidate is viable
             if candidate_state in coaccessible_nodes:
-                state_stack.append(cast(str, candidate_state))
+                state_stack.append(candidate_state)
                 char_stack.append(cast(str, candidate))
                 candidate = first_symbol
             else:
@@ -1553,38 +1549,35 @@ def get_name_original(states: FrozenSet[DFAStateT]) -> DFAStateT:
             final_states=dfa_final_states,
         )
 
-    def show_diagram(self, path=None):
+    def iter_transitions(
+        self,
+    ) -> Generator[Tuple[DFAStateT, DFAStateT, str], None, None]:
+        return (
+            (from_, to_, symbol)
+            for from_, lookup in self.transitions.items()
+            for symbol, to_ in lookup.items()
+        )
+
+    def _get_input_path(
+        self, input_str
+    ) -> Tuple[List[Tuple[DFAStateT, DFAStateT, DFASymbolT]], bool]:
         """
-        Creates the graph associated with this DFA
+        Calculate the path taken by input.
+
+        Args:
+            input_str (str): The input string to run on the DFA.
+
+        Returns:
+            tuple[list[tuple[DFAStateT, DFAStateT, DFASymbolT], bool]]: A list
+            of all transitions taken in each step and a boolean indicating
+            whether the DFA accepted the input.
+
         """
-        # Nodes are set of states
 
-        graph = Dot(graph_type="digraph", rankdir="LR")
-        nodes = {}
-        for state in self.states:
-            if state == self.initial_state:
-                # color start state with green
-                if state in self.final_states:
-                    initial_state_node = Node(
-                        state, style="filled", peripheries=2, fillcolor="#66cc33"
-                    )
-                else:
-                    initial_state_node = Node(
-                        state, style="filled", fillcolor="#66cc33"
-                    )
-                nodes[state] = initial_state_node
-                graph.add_node(initial_state_node)
-            else:
-                if state in self.final_states:
-                    state_node = Node(state, peripheries=2)
-                else:
-                    state_node = Node(state)
-                nodes[state] = state_node
-                graph.add_node(state_node)
-        # adding edges
-        for from_state, lookup in self.transitions.items():
-            for to_label, to_state in lookup.items():
-                graph.add_edge(Edge(nodes[from_state], nodes[to_state], label=to_label))
-        if path:
-            graph.write_png(path)
-        return graph
+        state_history = list(self.read_input_stepwise(input_str, ignore_rejection=True))
+        path = list(zip(state_history, state_history[1:], input_str))
+
+        last_state = state_history[-1] if state_history else self.initial_state
+        accepted = last_state in self.final_states
+
+        return path, accepted
diff --git a/automata/fa/fa.py b/automata/fa/fa.py
index 2336d609..69cd595e 100644
--- a/automata/fa/fa.py
+++ b/automata/fa/fa.py
@@ -1,12 +1,29 @@
 #!/usr/bin/env python3
 """Classes and methods for working with all finite automata."""
+from __future__ import annotations
 
 import abc
-from typing import Set
+import os
+import pathlib
+import random
+import uuid
+from collections import defaultdict
+from typing import Dict, Generator, List, Literal, Optional, Set, Tuple, Union
 
 from automata.base.automaton import Automaton, AutomatonStateT
 
+# Optional imports for use with visual functionality
+try:
+    import coloraide
+    import pygraphviz as pgv
+except ImportError:
+    _visual_imports = False
+else:
+    _visual_imports = True
+
+
 FAStateT = AutomatonStateT
+LayoutMethod = Literal["neato", "dot", "twopi", "circo", "fdp", "nop"]
 
 
 class FA(Automaton, metaclass=abc.ABCMeta):
@@ -14,6 +31,203 @@ class FA(Automaton, metaclass=abc.ABCMeta):
 
     __slots__ = tuple()
 
+    @staticmethod
+    def get_state_name(state_data: FAStateT) -> str:
+        """
+        Get an string representation of a state. This is used for displaying and
+        uses `str` for any unsupported python data types.
+        """
+        if isinstance(state_data, str):
+            if state_data == "":
+                return "λ"
+
+            return state_data
+
+        elif isinstance(state_data, (frozenset, tuple)):
+            inner = ", ".join(FA.get_state_name(sub_data) for sub_data in state_data)
+            if isinstance(state_data, frozenset):
+                if state_data:
+                    return "{" + inner + "}"
+                else:
+                    return "∅"
+
+            elif isinstance(state_data, tuple):
+                return "(" + inner + ")"
+
+        return str(state_data)
+
+    @staticmethod
+    def get_edge_name(symbol: str) -> str:
+        return "ε" if symbol == "" else str(symbol)
+
+    @abc.abstractmethod
+    def iter_transitions(self) -> Generator[Tuple[FAStateT, FAStateT, str], None, None]:
+        """
+        Iterate over all transitions in the automaton. Each transition is a tuple
+        of the form (from_state, to_state, symbol)
+        """
+
+    def show_diagram(
+        self,
+        input_str: Optional[str] = None,
+        path: Union[str, os.PathLike, None] = None,
+        *,
+        layout_method: LayoutMethod = "dot",
+        horizontal: bool = True,
+        reverse_orientation: bool = False,
+        fig_size: Union[Tuple[float, float], Tuple[float], None] = None,
+        font_size: float = 14.0,
+        arrow_size: float = 0.85,
+        state_separation: float = 0.5,
+    ) -> pgv.AGraph:
+        """
+        Generates the graph associated with the given DFA.
+        Args:
+            input_str (str, optional): String list of input symbols. Defaults to None.
+            - path (str or os.PathLike, optional): Path to output file. If
+              None, the output will not be saved.
+            - horizontal (bool, optional): Direction of node layout. Defaults
+              to True.
+            - reverse_orientation (bool, optional): Reverse direction of node
+              layout. Defaults to False.
+            - fig_size (tuple, optional): Figure size. Defaults to None.
+            - font_size (float, optional): Font size. Defaults to 14.0.
+            - arrow_size (float, optional): Arrow head size. Defaults to 0.85.
+            - state_separation (float, optional): Node distance. Defaults to 0.5.
+        Returns:
+            AGraph corresponding to the given automaton.
+        """
+
+        if not _visual_imports:
+            raise ImportError(
+                "Missing visualization packages; "
+                "please install coloraide and pygraphviz."
+            )
+
+        # Defining the graph.
+        graph = pgv.AGraph(strict=False, directed=True)
+
+        if fig_size is not None:
+            graph.graph_attr.update(size=", ".join(map(str, fig_size)))
+
+        graph.graph_attr.update(ranksep=str(state_separation))
+        font_size_str = str(font_size)
+        arrow_size_str = str(arrow_size)
+
+        if horizontal:
+            rankdir = "RL" if reverse_orientation else "LR"
+        else:
+            rankdir = "BT" if reverse_orientation else "TB"
+
+        graph.graph_attr.update(rankdir=rankdir)
+
+        # we use a random uuid to make sure that the null node has a
+        # unique id to avoid colliding with other states.
+        # To be able to set the random seed, took code from:
+        # https://nathanielknight.ca/articles/consistent_random_uuids_in_python.html
+        null_node = str(
+            uuid.UUID(bytes=bytes(random.getrandbits(8) for _ in range(16)), version=4)
+        )
+        graph.add_node(
+            null_node,
+            label="",
+            tooltip=".",
+            shape="point",
+            fontsize=font_size_str,
+        )
+        initial_node = self.get_state_name(self.initial_state)
+        graph.add_edge(
+            null_node,
+            initial_node,
+            tooltip="->" + initial_node,
+            arrowsize=arrow_size_str,
+        )
+
+        nonfinal_states = map(self.get_state_name, self.states - self.final_states)
+        final_states = map(self.get_state_name, self.final_states)
+        graph.add_nodes_from(nonfinal_states, shape="circle", fontsize=font_size_str)
+        graph.add_nodes_from(final_states, shape="doublecircle", fontsize=font_size_str)
+
+        is_edge_drawn = defaultdict(lambda: False)
+        if input_str is not None:
+            input_path, is_accepted = self._get_input_path(input_str=input_str)
+
+            start_color = coloraide.Color("#ff0")
+            end_color = (
+                coloraide.Color("#0f0") if is_accepted else coloraide.Color("#f00")
+            )
+            interpolation = coloraide.Color.interpolate(
+                [start_color, end_color], space="srgb"
+            )
+
+            # find all transitions in the finite state machine with traversal.
+            for transition_index, (from_state, to_state, symbol) in enumerate(
+                input_path, start=1
+            ):
+                color = interpolation(transition_index / len(input_path))
+                label = self.get_edge_name(symbol)
+
+                is_edge_drawn[from_state, to_state, symbol] = True
+                graph.add_edge(
+                    self.get_state_name(from_state),
+                    self.get_state_name(to_state),
+                    label=f"<{label} <b>[<i>#{transition_index}</i>]</b>>",
+                    arrowsize=arrow_size_str,
+                    fontsize=font_size_str,
+                    color=color.to_string(hex=True),
+                    penwidth="2.5",
+                )
+
+        edge_labels = defaultdict(list)
+        for from_state, to_state, symbol in self.iter_transitions():
+            if is_edge_drawn[from_state, to_state, symbol]:
+                continue
+
+            from_node = self.get_state_name(from_state)
+            to_node = self.get_state_name(to_state)
+            label = self.get_edge_name(symbol)
+            edge_labels[from_node, to_node].append(label)
+
+        for (from_node, to_node), labels in edge_labels.items():
+            graph.add_edge(
+                from_node,
+                to_node,
+                label=",".join(sorted(labels)),
+                arrowsize=arrow_size_str,
+                fontsize=font_size_str,
+            )
+
+        # Set layout
+        graph.layout(prog=layout_method)
+
+        # Write diagram to file. PNG, SVG, etc.
+        if path is not None:
+            save_path_final: pathlib.Path = pathlib.Path(path)
+
+            format = (
+                save_path_final.suffix.split(".")[1] if save_path_final.suffix else None
+            )
+
+            graph.draw(
+                path=save_path_final,
+                format=format,
+            )
+
+        return graph
+
+    @abc.abstractmethod
+    def _get_input_path(
+        self, input_str: str
+    ) -> Tuple[List[Tuple[FAStateT, FAStateT, str]], bool]:
+        """Calculate the path taken by input."""
+
+        raise NotImplementedError(
+            f"_get_input_path is not implemented for {self.__class__}"
+        )
+
+    def _repr_mimebundle_(self, *args, **kwargs) -> Dict[str, Union[bytes, str]]:
+        return self.show_diagram()._repr_mimebundle_(*args, **kwargs)
+
     @staticmethod
     def _add_new_state(state_set: Set[FAStateT], start: int = 0) -> int:
         """Adds new state to the state set and returns it"""
diff --git a/automata/fa/gnfa.py b/automata/fa/gnfa.py
index 551250b3..eafba8f1 100644
--- a/automata/fa/gnfa.py
+++ b/automata/fa/gnfa.py
@@ -4,10 +4,18 @@
 from __future__ import annotations
 
 from itertools import product
-from typing import AbstractSet, Dict, Mapping, Optional, Set, Type, cast
+from typing import (
+    AbstractSet,
+    Dict,
+    Generator,
+    Mapping,
+    Optional,
+    Set,
+    Tuple,
+    Type,
+    cast,
+)
 
-from frozendict import frozendict
-from pydot import Dot, Edge, Node
 from typing_extensions import NoReturn, Self
 
 import automata.base.exceptions as exceptions
@@ -25,17 +33,14 @@
 class GNFA(fa.FA):
     """A generalized nondeterministic finite automaton."""
 
-    # The conventions of using __slots__ state that subclasses automatically
-    # inherit __slots__ from parent classes, so there's no need to redeclare
-    # slotted attributes for each subclass; however, because NFA has a
-    # 'final_states' attribute but GNFA has a 'final_state' attribute, we must
-    # redeclare them below to exclude 'final_states'
+    # Add __dict__ to deal with inheritance issue and the final_states attribute.
     __slots__ = (
         "states",
         "input_symbols",
         "transitions",
         "initial_state",
         "final_state",
+        "__dict__",
     )
 
     final_state: GNFAStateT
@@ -49,15 +54,14 @@ def __init__(
         initial_state: GNFAStateT,
         final_state: GNFAStateT,
     ) -> None:
-        """Initialize a complete NFA."""
+        """Initialize a complete GNFA."""
         super(fa.FA, self).__init__(
-            states=frozenset(states),
-            input_symbols=frozenset(input_symbols),
-            transitions=frozendict(
-                {state: frozendict(paths) for state, paths in transitions.items()}
-            ),
+            states=states,
+            input_symbols=input_symbols,
+            transitions=transitions,
             initial_state=initial_state,
             final_state=final_state,
+            final_states={final_state},
         )
 
     # GNFA should NOT create the lambda closures via NFA.__post_init__()
@@ -327,39 +331,17 @@ def to_regex(self) -> str:
     def read_input_stepwise(self, input_str: str) -> NoReturn:
         raise NotImplementedError
 
-    def show_diagram(self, path=None, show_None=True):
-        """
-        Creates the graph associated with this DFA
-        """
-        # Nodes are set of states
-
-        graph = Dot(graph_type="digraph", rankdir="LR")
-        nodes = {}
-        for state in self.states:
-            if state == self.initial_state:
-                # color start state with green
-                initial_state_node = Node(state, style="filled", fillcolor="#66cc33")
-                nodes[state] = initial_state_node
-                graph.add_node(initial_state_node)
-            else:
-                if state == self.final_state:
-                    state_node = Node(state, peripheries=2)
-                else:
-                    state_node = Node(state)
-                nodes[state] = state_node
-                graph.add_node(state_node)
-        # adding edges
-        for from_state, lookup in self.transitions.items():
-            for to_state, to_label in lookup.items():  # pragma: no branch
-                if to_label is None and show_None:
-                    to_label = "ø"
-                    graph.add_edge(
-                        Edge(nodes[from_state], nodes[to_state], label=to_label)
-                    )
-                elif to_label is not None:
-                    graph.add_edge(
-                        Edge(nodes[from_state], nodes[to_state], label=to_label)
-                    )
-        if path:
-            graph.write_png(path)
-        return graph
+    def iter_transitions(
+        self,
+    ) -> Generator[Tuple[GNFAStateT, GNFAStateT, str], None, None]:
+        return (
+            (from_, to_, symbol)
+            for from_, lookup in self.transitions.items()
+            for to_, symbol in lookup.items()
+            if symbol is not None
+        )
+
+    def _get_input_path(self, input_str: str) -> NoReturn:
+        raise NotImplementedError(
+            f"_get_input_path is not implemented for {self.__class__}"
+        )
diff --git a/automata/fa/nfa.py b/automata/fa/nfa.py
index ac0228c9..4689386b 100644
--- a/automata/fa/nfa.py
+++ b/automata/fa/nfa.py
@@ -23,7 +23,6 @@
 
 import networkx as nx
 from frozendict import frozendict
-from pydot import Dot, Edge, Node
 from typing_extensions import Self
 
 import automata.base.exceptions as exceptions
@@ -35,7 +34,7 @@
 
 NFAPathT = Mapping[str, AbstractSet[NFAStateT]]
 NFATransitionsT = Mapping[NFAStateT, NFAPathT]
-
+InputPathListT = List[Tuple[NFAStateT, NFAStateT, str]]
 DEFAULT_REGEX_SYMBOLS = frozenset(chain(string.ascii_letters, string.digits))
 
 
@@ -808,45 +807,6 @@ def left_quotient(self, other: NFA) -> Self:
             final_states=new_final_states,
         )
 
-    def show_diagram(self, path=None):
-        """
-        Creates the graph associated with this DFA
-        """
-        # Nodes are set of states
-
-        graph = Dot(graph_type="digraph", rankdir="LR")
-        nodes = {}
-        for state in self.states:
-            if state == self.initial_state:
-                # color start state with green
-                if state in self.final_states:
-                    initial_state_node = Node(
-                        state, style="filled", peripheries=2, fillcolor="#66cc33"
-                    )
-                else:
-                    initial_state_node = Node(
-                        state, style="filled", fillcolor="#66cc33"
-                    )
-                nodes[state] = initial_state_node
-                graph.add_node(initial_state_node)
-            else:
-                if state in self.final_states:
-                    state_node = Node(state, peripheries=2)
-                else:
-                    state_node = Node(state)
-                nodes[state] = state_node
-                graph.add_node(state_node)
-        # adding edges
-        for from_state, lookup in self.transitions.items():
-            for to_label, to_states in lookup.items():
-                for to_state in to_states:
-                    graph.add_edge(
-                        Edge(nodes[from_state], nodes[to_state], label=to_label)
-                    )
-        if path:
-            graph.write_png(path)
-        return graph
-
     @staticmethod
     def _load_new_transition_dict(
         state_map_dict: Mapping[NFAStateT, NFAStateT],
@@ -1009,3 +969,79 @@ def add_any_transition(start_state_dict, end_state):
             initial_state=(0, 0),
             final_states=final_states,
         )
+
+    def iter_transitions(
+        self,
+    ) -> Generator[Tuple[NFAStateT, NFAStateT, str], None, None]:
+        return (
+            (from_, to_, symbol)
+            for from_, lookup in self.transitions.items()
+            for symbol, to_lookup in lookup.items()
+            for to_ in to_lookup
+        )
+
+    def _get_input_path(self, input_str: str) -> Tuple[InputPathListT, bool]:
+        """
+        Get best input path. A path is better if (with priority):
+
+        1. It is an accepting path (ends in a final state)
+        2. It reads more of the input (if the input is not accepted we
+            select the path such that we can stay on the nfa the longest)
+        3. It has the fewest jumps (uses less lambda symbols)
+
+        Returns a tuple of:
+        1. the path taken
+        2. whether the path was accepting
+        """
+
+        visited = set()
+        work_queue: Deque[Tuple[InputPathListT, NFAStateT, str]] = deque(
+            [([], self.initial_state, input_str)]
+        )
+
+        last_non_accepting_input: InputPathListT = []
+        least_input_remaining = len(input_str)
+
+        while work_queue:
+            visited_states, curr_state, remaining_input = work_queue.popleft()
+
+            # First final state we hit is the best according to desired criteria
+            if curr_state in self.final_states and not remaining_input:
+                return visited_states, True
+
+            # Otherwise, update longest non-accepting input
+            if len(remaining_input) < least_input_remaining:
+                least_input_remaining = len(remaining_input)
+                last_non_accepting_input = visited_states
+
+            # First, get next states that result from reading from input
+            if remaining_input:
+                next_symbol = remaining_input[0]
+                rest = remaining_input[1:] if remaining_input else ""
+
+                next_states_from_symbol = self.transitions[curr_state].get(
+                    next_symbol, set()
+                )
+
+                for next_state in next_states_from_symbol:
+                    if (next_state, rest) not in visited:
+                        next_visited_states = visited_states.copy()
+                        next_visited_states.append(
+                            (curr_state, next_state, next_symbol)
+                        )
+                        visited.add((next_state, rest))
+                        work_queue.append((next_visited_states, next_state, rest))
+
+            # Next, get next states resulting from lambda transition
+            next_states_from_lambda = self.transitions[curr_state].get("", set())
+
+            for next_state in next_states_from_lambda:
+                if (next_state, remaining_input) not in visited:
+                    next_visited_states = visited_states.copy()
+                    next_visited_states.append((curr_state, next_state, ""))
+                    visited.add((next_state, remaining_input))
+                    work_queue.append(
+                        (next_visited_states, next_state, remaining_input)
+                    )
+
+        return last_non_accepting_input, False
diff --git a/automata/pda/stack.py b/automata/pda/stack.py
index 6c537073..7f534263 100644
--- a/automata/pda/stack.py
+++ b/automata/pda/stack.py
@@ -49,7 +49,7 @@ def __iter__(self) -> Iterator[str]:
         """Return an interator for the stack."""
         return iter(self.stack)
 
-    def __getitem__(self, key):
+    def __getitem__(self, key: int) -> str:
         """Return the stack element at the given index"""
         return self.stack[key]
 
diff --git a/pyproject.toml b/pyproject.toml
index cb52f407..a705eea6 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -13,11 +13,18 @@ maintainers = [
     {name = 'Caleb Evans', email = 'caleb@calebevans.me'}
 ]
 dependencies = [
-    "pydot>=1.4.2",
     "networkx>=2.6.2",
-    "frozendict>=2.3.4", #TODO I think that typing extensions needs to be in here
+    "frozendict>=2.3.4",
+    "typing-extensions>=4.5.0"
 ]
 
+# Per https://setuptools.pypa.io/en/latest/userguide/dependency_management.html#optional-dependencies
+[project.optional-dependencies]
+visual = ["coloraide>=1.8.2", "pygraphviz>=1.10"]
+
+[tool.setuptools.package-data]
+"automata" = ["py.typed"]
+
 [project.urls]
 homepage = "https://github.com/caleb531/automata"
 documentation = "https://github.com/caleb531/automata#readme"
@@ -34,7 +41,7 @@ exclude = ["build"]
 module = [
   "setuptools.*",
   "networkx.*",
-  "pydot.*" # TODO pydot dependency is going to get removed, can delete
+  "pygraphviz.*"
 ]
 ignore_missing_imports = true
 
@@ -46,3 +53,22 @@ extend-ignore = ["E203", "W503"]
 # Per <https://black.readthedocs.io/en/stable/guides/using_black_with_other_tools.html#profile>
 [tool.isort]
 profile = "black"
+
+
+# Configuration for coverage.py (https://pypi.python.org/pypi/coverage)
+
+[tool.coverage.run]
+# Enable branch coverage
+branch = true
+
+[tool.coverage.report]
+# Regexes for lines to exclude from consideration
+exclude_lines = [
+    "pragma: no cover",
+    # Ignore non-runnable code
+    "if __name__ == .__main__.:",
+    "pass",
+]
+
+# Only check coverage for source files
+include =["automata/*/*.py"]
diff --git a/requirements.txt b/requirements.txt
index 411b6c03..9883044c 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -1,11 +1,13 @@
 black==23.3.0
 click==8.1.3
-coverage==6.4.4
+coloraide==1.8.2
+coverage[toml]==6.4.4
 flake8==6.0.0
 flake8-black==0.3.6
 flake8-isort==6.0.0
 Flake8-pyproject==1.2.3
 frozendict==2.3.4
+pygraphviz==1.10
 isort==5.10.1
 mccabe==0.7.0
 mypy==1.1.1
diff --git a/setup.py b/setup.py
index 60684932..7f1a1763 100644
--- a/setup.py
+++ b/setup.py
@@ -1,3 +1,4 @@
 from setuptools import setup
 
-setup()
+if __name__ == "__main__":
+    setup()
diff --git a/tests/test_dfa.py b/tests/test_dfa.py
index 75e268a6..764049fa 100644
--- a/tests/test_dfa.py
+++ b/tests/test_dfa.py
@@ -3,6 +3,7 @@
 
 import os
 import os.path
+import random
 import tempfile
 import types
 from itertools import product
@@ -1519,26 +1520,55 @@ def test_show_diagram_initial_final_different(self) -> None:
         is not a final state.
         """
         graph = self.dfa.show_diagram()
-        self.assertEqual(
-            {node.get_name() for node in graph.get_nodes()}, {"q0", "q1", "q2"}
-        )
-        self.assertEqual(graph.get_node("q0")[0].get_style(), "filled")
-        self.assertEqual(graph.get_node("q1")[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q2")[0].get_peripheries(), None)
-        self.assertEqual(
-            {
-                (edge.get_source(), edge.get_label(), edge.get_destination())
-                for edge in graph.get_edges()
-            },
-            {
-                ("q0", "0", "q0"),
-                ("q0", "1", "q1"),
-                ("q1", "0", "q0"),
-                ("q1", "1", "q2"),
-                ("q2", "0", "q2"),
-                ("q2", "1", "q1"),
-            },
-        )
+        node_names = {node.get_name() for node in graph.nodes()}
+        self.assertTrue(set(self.dfa.states).issubset(node_names))
+        self.assertEqual(len(self.dfa.states) + 1, len(node_names))
+
+        for state in self.dfa.states:
+            node = graph.get_node(state)
+            expected_shape = (
+                "doublecircle" if state in self.dfa.final_states else "circle"
+            )
+            self.assertEqual(node.attr["shape"], expected_shape)
+
+        expected_transitions = {
+            ("q0", "0", "q0"),
+            ("q0", "1", "q1"),
+            ("q1", "0", "q0"),
+            ("q1", "1", "q2"),
+            ("q2", "0", "q2"),
+            ("q2", "1", "q1"),
+        }
+        seen_transitions = {
+            (edge[0], edge.attr["label"], edge[1]) for edge in graph.edges()
+        }
+        self.assertTrue(expected_transitions.issubset(seen_transitions))
+        self.assertEqual(len(expected_transitions) + 1, len(seen_transitions))
+
+        source, symbol, dest = list(seen_transitions - expected_transitions)[0]
+        self.assertEqual(symbol, "")
+        self.assertEqual(dest, self.dfa.initial_state)
+        self.assertTrue(source not in self.dfa.states)
+
+    def test_show_diagram_read_input(self) -> None:
+        """
+        Should construct the diagram for a DFA reading input.
+        """
+        input_strings = ["0111", "001", "01110011", "001011001", "1100", ""]
+
+        for input_string in input_strings:
+            graph = self.dfa.show_diagram(input_str=input_string)
+
+            # Get edges corresponding to input path
+            colored_edges = [
+                edge for edge in graph.edges() if "color" in dict(edge.attr)
+            ]
+            colored_edges.sort(key=lambda edge: edge.attr["label"][2:])
+
+            edge_pairs = [
+                edge[0:2] for edge in self.dfa._get_input_path(input_string)[0]
+            ]
+            self.assertEqual(edge_pairs, colored_edges)
 
     def test_show_diagram_initial_final_same(self) -> None:
         """
@@ -1558,28 +1588,34 @@ def test_show_diagram_initial_final_same(self) -> None:
             initial_state="q0",
             final_states={"q0", "q1"},
         )
+
         graph = dfa.show_diagram()
-        self.assertEqual(
-            {node.get_name() for node in graph.get_nodes()}, {"q0", "q1", "q2"}
-        )
-        self.assertEqual(graph.get_node("q0")[0].get_style(), "filled")
-        self.assertEqual(graph.get_node("q0")[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q1")[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q2")[0].get_peripheries(), None)
-        self.assertEqual(
-            {
-                (edge.get_source(), edge.get_label(), edge.get_destination())
-                for edge in graph.get_edges()
-            },
-            {
-                ("q0", "0", "q0"),
-                ("q0", "1", "q1"),
-                ("q1", "0", "q0"),
-                ("q1", "1", "q2"),
-                ("q2", "0", "q2"),
-                ("q2", "1", "q2"),
-            },
-        )
+        node_names = {node.get_name() for node in graph.nodes()}
+        self.assertTrue(set(dfa.states).issubset(node_names))
+        self.assertEqual(len(dfa.states) + 1, len(node_names))
+
+        for state in self.dfa.states:
+            node = graph.get_node(state)
+            expected_shape = "doublecircle" if state in dfa.final_states else "circle"
+            self.assertEqual(node.attr["shape"], expected_shape)
+
+        expected_transitions = {
+            ("q0", "0", "q0"),
+            ("q0", "1", "q1"),
+            ("q1", "0", "q0"),
+            ("q1", "1", "q2"),
+            ("q2", "0,1", "q2"),
+        }
+        seen_transitions = {
+            (edge[0], edge.attr["label"], edge[1]) for edge in graph.edges()
+        }
+        self.assertTrue(expected_transitions.issubset(seen_transitions))
+        self.assertEqual(len(expected_transitions) + 1, len(seen_transitions))
+
+        source, symbol, dest = list(seen_transitions - expected_transitions)[0]
+        self.assertEqual(symbol, "")
+        self.assertEqual(dest, dfa.initial_state)
+        self.assertTrue(source not in dfa.states)
 
     def test_show_diagram_write_file(self) -> None:
         """
@@ -1596,6 +1632,38 @@ def test_show_diagram_write_file(self) -> None:
         self.assertTrue(os.path.exists(diagram_path))
         os.remove(diagram_path)
 
+    def test_repr_mimebundle_same(self) -> None:
+        """
+        Check that the mimebundle is the same.
+        """
+
+        random.seed(42)
+        first_repr = self.dfa._repr_mimebundle_()
+        random.seed(42)
+        second_repr = self.dfa.show_diagram()._repr_mimebundle_()
+        self.assertEqual(first_repr, second_repr)
+
+    def test_show_diagram_orientations(self) -> None:
+        graph = self.dfa.show_diagram()
+        self.assertEqual(graph.graph_attr["rankdir"], "LR")
+        graph = self.dfa.show_diagram(horizontal=False)
+        self.assertEqual(graph.graph_attr["rankdir"], "TB")
+        graph = self.dfa.show_diagram(reverse_orientation=True)
+        self.assertEqual(graph.graph_attr["rankdir"], "RL")
+        graph = self.dfa.show_diagram(horizontal=False, reverse_orientation=True)
+        self.assertEqual(graph.graph_attr["rankdir"], "BT")
+
+    def test_show_diagram_fig_size(self) -> None:
+        """
+        Testing figure size. Just need to make sure it matches the input
+        (the library handles the rendering).
+        """
+        graph = self.dfa.show_diagram(fig_size=(1.1, 2))
+        self.assertEqual(graph.graph_attr["size"], "1.1, 2")
+
+        graph = self.dfa.show_diagram(fig_size=(3.3,))
+        self.assertEqual(graph.graph_attr["size"], "3.3")
+
     def test_minimal_finite_language(self):
         """Should compute the minimal DFA accepting the given finite language"""
 
diff --git a/tests/test_fa.py b/tests/test_fa.py
index 7b7f3af5..8871389f 100644
--- a/tests/test_fa.py
+++ b/tests/test_fa.py
@@ -4,6 +4,7 @@
 import unittest
 
 from automata.fa.dfa import DFA
+from automata.fa.fa import FA
 from automata.fa.gnfa import GNFA
 from automata.fa.nfa import NFA
 
@@ -56,3 +57,45 @@ def setUp(self) -> None:
             initial_state="q_in",
             final_state="q_f",
         )
+
+
+class TestFAAbstract(unittest.TestCase):
+    def test_abstract_methods_not_implemented(self) -> None:
+        """Should raise NotImplementedError when calling abstract methods."""
+
+        with self.assertRaises(NotImplementedError):
+            getattr(FA, "_get_input_path")(FA, "")
+
+    def test_get_state_name(self) -> None:
+        """Tests get_state_name function."""
+
+        self.assertEqual(FA.get_state_name(""), "λ")
+        self.assertEqual(FA.get_state_name("abc"), "abc")
+        self.assertEqual(FA.get_state_name(123), "123")
+        self.assertEqual(FA.get_state_name(frozenset()), "∅")
+
+        original_sets = [frozenset("abc"), frozenset("dcgf")]
+
+        for original_set in original_sets:
+            set_output = frozenset(
+                FA.get_state_name(original_set)
+                .replace("{", "")
+                .replace("}", "")
+                .replace("'", "")
+                .replace(" ", "")
+                .split(",")
+            )
+            self.assertEqual(set_output, original_set)
+
+        original_tuple = tuple("abc")
+
+        tuple_output = tuple(
+            FA.get_state_name(original_tuple)
+            .replace("(", "")
+            .replace(")", "")
+            .replace("'", "")
+            .replace(" ", "")
+            .split(",")
+        )
+
+        self.assertEqual(tuple_output, original_tuple)
diff --git a/tests/test_gnfa.py b/tests/test_gnfa.py
index 54d1ac81..45d269f4 100644
--- a/tests/test_gnfa.py
+++ b/tests/test_gnfa.py
@@ -16,6 +16,16 @@ class TestGNFA(test_fa.TestFA):
 
     temp_dir_path = tempfile.gettempdir()
 
+    def test_methods_not_implemented(self) -> None:
+        """Should raise NotImplementedError when calling non-implemented methods."""
+        abstract_methods = {
+            "_get_input_path": (GNFA, ""),
+            "read_input_stepwise": (GNFA, ""),
+        }
+        for method_name, method_args in abstract_methods.items():
+            with self.assertRaises(NotImplementedError):
+                getattr(GNFA, method_name)(*method_args)
+
     def test_init_gnfa(self) -> None:
         """Should copy GNFA if passed into NFA constructor."""
         new_gnfa = self.gnfa.copy()
@@ -372,69 +382,43 @@ def test_to_regex(self) -> None:
             # Test equality through DFA regex conversion
             self.assertEqual(dfa_1, dfa_2)
 
-    def test_show_diagram_showNone(self) -> None:
-        """
-        Should construct the diagram for a GNFA when show_None = False
-        """
-
-        gnfa = self.gnfa
-
-        graph = gnfa.show_diagram(show_None=False)
-        self.assertEqual({node.get_name() for node in graph.get_nodes()}, gnfa.states)
-        self.assertEqual(graph.get_node(gnfa.initial_state)[0].get_style(), "filled")
-        self.assertEqual(graph.get_node(gnfa.final_state)[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q2")[0].get_peripheries(), None)
-        self.assertEqual(
-            {
-                (edge.get_source(), edge.get_label(), edge.get_destination())
-                for edge in graph.get_edges()
-            },
-            {
-                ("q0", "a", "q1"),
-                ("q1", "a", "q1"),
-                ("q1", "", "q2"),
-                ("q1", "", "q_f"),
-                ("q2", "b", "q0"),
-                ("q_in", "", "q0"),
-            },
-        )
-
     def test_show_diagram(self) -> None:
         """
-        Should construct the diagram for a GNFA when show_None = True
+        Should construct the diagram for a GNFA.
         """
 
-        gnfa = self.gnfa
+        graph = self.gnfa.show_diagram()
 
-        graph = gnfa.show_diagram()
-        self.assertEqual({node.get_name() for node in graph.get_nodes()}, gnfa.states)
-        self.assertEqual(graph.get_node(gnfa.initial_state)[0].get_style(), "filled")
-        self.assertEqual(graph.get_node(gnfa.final_state)[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q2")[0].get_peripheries(), None)
-        self.assertEqual(
-            {
-                (edge.get_source(), edge.get_label(), edge.get_destination())
-                for edge in graph.get_edges()
-            },
-            {
-                ("q_in", "", "q0"),
-                ("q0", "ø", "q2"),
-                ("q1", "", "q2"),
-                ("q0", "ø", "q_f"),
-                ("q1", "", "q_f"),
-                ("q_in", "ø", "q2"),
-                ("q_in", "ø", "q1"),
-                ("q1", "a", "q1"),
-                ("q2", "b", "q0"),
-                ("q2", "ø", "q2"),
-                ("q_in", "ø", "q_f"),
-                ("q2", "ø", "q1"),
-                ("q0", "ø", "q0"),
-                ("q2", "ø", "q_f"),
-                ("q0", "a", "q1"),
-                ("q1", "ø", "q0"),
-            },
-        )
+        node_names = {node.get_name() for node in graph.nodes()}
+        self.assertTrue(set(self.gnfa.states).issubset(node_names))
+        self.assertEqual(len(self.gnfa.states) + 1, len(node_names))
+
+        for state in self.dfa.states:
+            node = graph.get_node(state)
+            expected_shape = (
+                "doublecircle" if state in self.gnfa.final_states else "circle"
+            )
+            self.assertEqual(node.attr["shape"], expected_shape)
+
+        expected_transitions = {
+            ("q_in", "ε", "q0"),
+            ("q1", "ε", "q2"),
+            ("q1", "ε", "q_f"),
+            ("q1", "a", "q1"),
+            ("q2", "b", "q0"),
+            ("q0", "a", "q1"),
+        }
+        seen_transitions = {
+            (edge[0], edge.attr["label"], edge[1]) for edge in graph.edges()
+        }
+
+        self.assertTrue(expected_transitions.issubset(seen_transitions))
+        self.assertEqual(len(expected_transitions) + 1, len(seen_transitions))
+
+        source, symbol, dest = list(seen_transitions - expected_transitions)[0]
+        self.assertEqual(symbol, "")
+        self.assertEqual(dest, self.gnfa.initial_state)
+        self.assertTrue(source not in self.gnfa.states)
 
     def test_show_diagram_write_file(self) -> None:
         """
diff --git a/tests/test_nfa.py b/tests/test_nfa.py
index 53699c1a..ff3b4ad0 100644
--- a/tests/test_nfa.py
+++ b/tests/test_nfa.py
@@ -584,32 +584,56 @@ def test_show_diagram_initial_final_same(self) -> None:
         is also a final state.
         """
 
-        nfa = NFA(
-            states={"q0", "q1", "q2"},
-            input_symbols={"a", "b"},
-            transitions={
-                "q0": {"a": {"q1"}},
-                "q1": {"a": {"q1"}, "": {"q2"}},
-                "q2": {"b": {"q0"}},
-            },
-            initial_state="q0",
-            final_states={"q0", "q1"},
-        )
-        graph = nfa.show_diagram()
-        self.assertEqual(
-            {node.get_name() for node in graph.get_nodes()}, {"q0", "q1", "q2"}
-        )
-        self.assertEqual(graph.get_node("q0")[0].get_style(), "filled")
-        self.assertEqual(graph.get_node("q0")[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q1")[0].get_peripheries(), 2)
-        self.assertEqual(graph.get_node("q2")[0].get_peripheries(), None)
-        self.assertEqual(
-            {
-                (edge.get_source(), edge.get_label(), edge.get_destination())
-                for edge in graph.get_edges()
-            },
-            {("q0", "a", "q1"), ("q1", "a", "q1"), ("q1", "", "q2"), ("q2", "b", "q0")},
-        )
+        graph = self.nfa.show_diagram()
+        node_names = {node.get_name() for node in graph.nodes()}
+        self.assertTrue(set(self.nfa.states).issubset(node_names))
+        self.assertEqual(len(self.nfa.states) + 1, len(node_names))
+
+        for state in self.dfa.states:
+            node = graph.get_node(state)
+            expected_shape = (
+                "doublecircle" if state in self.nfa.final_states else "circle"
+            )
+            self.assertEqual(node.attr["shape"], expected_shape)
+
+        expected_transitions = {
+            ("q0", "a", "q1"),
+            ("q1", "a", "q1"),
+            ("q1", "ε", "q2"),
+            ("q2", "b", "q0"),
+        }
+        seen_transitions = {
+            (edge[0], edge.attr["label"], edge[1]) for edge in graph.edges()
+        }
+
+        self.assertTrue(expected_transitions.issubset(seen_transitions))
+        self.assertEqual(len(expected_transitions) + 1, len(seen_transitions))
+
+        source, symbol, dest = list(seen_transitions - expected_transitions)[0]
+        self.assertEqual(symbol, "")
+        self.assertEqual(dest, self.nfa.initial_state)
+        self.assertTrue(source not in self.nfa.states)
+
+    def test_show_diagram_read_input(self) -> None:
+        """
+        Should construct the diagram for a NFA reading input.
+        """
+        input_strings = ["ababa", "bba", "aabba", "baaab", "bbaab", ""]
+
+        for input_string in input_strings:
+            graph = self.nfa.show_diagram(input_str=input_string)
+
+            # Get edges corresponding to input path
+            colored_edges = [
+                edge for edge in graph.edges() if "color" in dict(edge.attr)
+            ]
+            colored_edges.sort(key=lambda edge: edge.attr["label"][2:])
+
+            edge_pairs = [
+                edge[0:2] for edge in self.nfa._get_input_path(input_string)[0]
+            ]
+
+            self.assertEqual(edge_pairs, colored_edges)
 
     def test_show_diagram_write_file(self) -> None:
         """
@@ -626,6 +650,71 @@ def test_show_diagram_write_file(self) -> None:
         self.assertTrue(os.path.exists(diagram_path))
         os.remove(diagram_path)
 
+    def test_get_input_path(self) -> None:
+        nfa2 = NFA(
+            states={"q0", "q1", "q2"},
+            input_symbols={"a", "b"},
+            transitions={
+                "q0": {"": {"q1"}, "b": {"q0"}},
+                "q1": {"a": {"q2"}, "": {"q0"}},
+                "q2": {"b": {"q2"}, "a": {"q1"}},
+            },
+            initial_state="q0",
+            final_states={"q1"},
+        )
+        nfa3 = NFA.from_regex("(abab|aabba)*bba*bb")
+
+        input_strings = [
+            "ababa",
+            "bba",
+            "aabba",
+            "baaab",
+            "bbaab",
+            "",
+            "bbbb",
+            "bbaabb",
+        ]
+        nfas = [self.nfa, nfa2, nfa3]
+
+        for input_str, nfa in product(input_strings, nfas):
+            input_path, was_accepted = nfa._get_input_path(input_str)
+            self.assertEqual(was_accepted, nfa.accepts_input(input_str))
+
+            last_vtx = None
+
+            for start_vtx, end_vtx, symbol in input_path:
+                last_vtx = end_vtx
+                self.assertIn(end_vtx, nfa.transitions[start_vtx][symbol])
+
+            self.assertEqual(last_vtx in nfa.final_states, was_accepted)
+
+    def test_input_path_optimality(self) -> None:
+        """
+        A test case for optimality of path found.
+        Checks path length doesn't use the extra epsilon transition.
+        """
+
+        nfa = NFA(
+            states=set(range(6)),
+            input_symbols=set("01"),
+            transitions={
+                0: {"0": {1, 2}},
+                1: {"": {3}},
+                2: {"1": {4}},
+                3: {"1": {4}},
+                4: {"": {5}},
+                5: {},
+            },
+            initial_state=0,
+            final_states={5},
+        )
+
+        input_str = "01"
+
+        input_path, was_accepted = nfa._get_input_path(input_str)
+        self.assertEqual(was_accepted, nfa.accepts_input(input_str))
+        self.assertEqual(len(input_path), 3)
+
     def test_add_new_state_type_integrity(self) -> None:
         """
         Should properly add new state of different type than original states;
diff --git a/tests/test_pdastack.py b/tests/test_pdastack.py
index a1eea107..fce672a3 100644
--- a/tests/test_pdastack.py
+++ b/tests/test_pdastack.py
@@ -8,14 +8,21 @@
 class TestPDAStack(test_pda.TestPDA):
     """A test class for testing stacks of pushdown automata."""
 
-    def test_stack_hashability(self):
-        self.assertEqual(hash(PDAStack(["a", "b"])), hash(PDAStack(["a", "b"])))
+    def setUp(self) -> None:
+        self.stack = PDAStack(["a", "b"])
 
-    def test_stack_iter(self):
+    def test_stack_hashability(self) -> None:
+        self.assertEqual(hash(self.stack), hash(PDAStack(["a", "b"])))
+
+    def test_stack_iter(self) -> None:
         """Should loop through the PDA stack in some manner."""
-        self.assertEqual(list(PDAStack(["a", "b"])), ["a", "b"])
+        self.assertEqual(list(self.stack), ["a", "b"])
+
+    def test_stack_get(self) -> None:
+        """Should retrieve indices in the PDA stack in some manner."""
+        self.assertEqual(self.stack[0], "a")
+        self.assertEqual(self.stack[1], "b")
 
-    def test_stack_repr(self):
+    def test_stack_repr(self) -> None:
         """Should create proper string representation of PDA stack."""
-        stack = PDAStack(["a", "b"])
-        self.assertEqual(repr(stack), "PDAStack(('a', 'b'))")
+        self.assertEqual(repr(self.stack), "PDAStack(('a', 'b'))")