diff --git a/parse_magic.py b/parse_magic.py
index 840739a..8e0bb76 100644
--- a/parse_magic.py
+++ b/parse_magic.py
@@ -4,7 +4,6 @@ Parses the Magic Ugly Data File Format
 Assumes the base logic is R with no extra connectives
 """
 import argparse
-import re
 import sys
 from typing import TextIO, List, Optional, Tuple, Set, Dict
 
@@ -24,31 +23,14 @@ class SourceFile:
         self.current_line = 0
         self.reststr = ""
 
-    def next_line(self):
+    def __next__(self):
         if self.reststr != "":
-            reststr = self.reststr
-            self.reststr = ""
-            return reststr
+            return self.reststr
 
-        contents = next(self.fileobj).strip()
+        contents = next(self.fileobj)
         self.current_line += 1
         return contents
 
-    def __next__(self):
-        """
-        Grabs the next word token from the stream
-        """
-        if self.reststr == "":
-            self.reststr = next(self.fileobj).strip()
-            self.current_line += 1
-
-        tokens = self.reststr.split(" ")
-        next_token = tokens[0]
-        self.reststr = " ".join(tokens[1:])
-
-        return next_token
-
-
     def set_reststr(self, reststr: str):
         self.reststr = reststr
 
@@ -62,95 +44,49 @@ class ModelBuilder:
     def __init__(self):
         self.size : int = 0
         self.carrier_set : Set[ModelValue] = set()
+        self.num_negation: int = 0
         self.mnegation: Optional[ModelFunction] = None
+        self.num_order: int = 0
         self.mconjunction: Optional[ModelFunction] = None
         self.mdisjunction: Optional[ModelFunction] = None
+        self.num_designated: int = 0
         self.designated_values: Set[ModelValue] = set()
+        self.num_implication: int = 0
         self.mimplication: Optional[ModelFunction] = None
+        self.num_necessitation: int = 0
         self.mnecessitation: Optional[ModelFunction] = None
-        self.custom_model_functions: Dict[str,  ModelFunction] = {}
+        self.custom_model_functions: List[Tuple[int, ModelFunction]]
 
 class Stage:
     def __init__(self, name: str):
         self.name = name
         self.next: Optional['Stage'] = None
         self.previous: Optional['Stage'] = None
-        self.num = 0
-
-    def increment(self):
-        self.num += 1
-
-    def reset(self):
-        self.num = 0
-
     def __str__(self):
         return self.name
 
 class Stages:
     def __init__(self):
         self.stages: Dict[str, Stage] = {}
-        self.last_added_stage: Optional[Stage] = None
-        self.first_stage: Optional[Stage] = None
-
+        self.last_stage: Optional[Stage] = None
     def add(self, name: str):
         stage = Stage(name)
         stage.next = stage
-
-        if self.last_added_stage is not None:
-            stage.previous = self.last_added_stage
-            self.last_added_stage.next = stage
+        if self.last_stage is not None:
+            stage.previous = self.last_stage
+            self.last_stage.next = stage
         else:
-            # The previous of the first stage
-            # is the end
-            self.first_stage = stage
             stage.previous = Stage("end")
 
         self.stages[name] = stage
-        self.last_added_stage = stage
-
+        self.last_stage = stage
     def next_stage(self, name):
         return self.stages[name].next
-
     def previous_stage(self, name):
         return self.stages[name].previous
-
-    def reset_after(self, name):
-        stage = self.stages[name]
-        stage.reset()
-        if stage.next.name == "process_model":
-            return
-        next_stage = stage.next
-        while next_stage is not None:
-            next_stage.reset()
-            if next_stage.next.name == "process_model":
-                next_stage = None
-            else:
-                next_stage = next_stage.next
-
     def get(self, name):
         return self.stages[name]
 
-    def name(self):
-        result = ""
-        stage = self.first_stage
-        if stage is None:
-            return ""
-
-        result = f"{stage.num}"
-        if stage.next == "process_model":
-            return result
-
-        stage = stage.next
-        while stage is not None:
-            result += f".{stage.num}"
-            if stage.next.name != "process_model":
-                stage = stage.next
-            else:
-                stage = None
-
-        return result
-
-
 def derive_stages(header: UglyHeader) -> Stages:
     stages = Stages()
     stages.add("size")
@@ -161,8 +97,9 @@ def derive_stages(header: UglyHeader) -> Stages:
     stages.add("implication")
     if header.necessitation:
         stages.add("necessitation")
-    for (adicity, symbol) in header.custom_model_functions:
-        stages.add(f"custom--{adicity}--{symbol}")
+    # TODO: Can't use -- as that can be a custom symbol
+    # for (acidity, symbol) in header.custom_model_functions:
+        # stages.add(f"custom--{acidity}--{symbol}")
     stages.add("process_model")
 
     # After processing the model, go to the previous stage
@@ -170,6 +107,7 @@ def derive_stages(header: UglyHeader) -> Stages:
 
     return stages
 
+
 def parse_matrices(infile: SourceFile) -> List[Tuple[Model, Dict]]:
     solutions = [] # Reset
     header = parse_header(infile)
@@ -178,85 +116,34 @@ def parse_matrices(infile: SourceFile) -> List[Tuple[Model, Dict]]:
     current_model_parts = ModelBuilder()
     stage = stages.get("size")
     while True:
+        print("Current stage:", stage)
         match stage.name:
             case "end":
                 break
             case "process_model":
-                process_model(stages.name(), current_model_parts, solutions)
+                process_model(current_model_parts, solutions)
                 stage = stage.next
             case "size":
                 processed = process_sizes(infile, current_model_parts, first_run)
                 first_run = False
-                if processed:
-                    stage.num = current_model_parts.size + 1
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
+                stage = stage.next if processed else stage.previous
             case "negation":
                 processed = process_negations(infile, current_model_parts)
-                if processed:
-                    stage.increment()
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
+                stage = stage.next if processed else stage.previous
             case "order":
                 processed = process_orders(infile, current_model_parts)
-                if processed:
-                    stage.increment()
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
+                stage = stage.next if processed else stage.previous
             case "designated":
                 processed = process_designateds(infile, current_model_parts)
-                if processed:
-                    stage.increment()
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
+                stage = stage.next if processed else stage.previous
             case "implication":
                 processed = process_implications(infile, current_model_parts)
-                if processed:
-                    stage.increment()
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
+                stage = stage.next if processed else stage.previous
             case "necessitation":
                 processed = process_necessitations(infile, current_model_parts)
-                if processed:
-                    stage.increment()
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
+                stage = stage.next if processed else stage.previous
             case _:
-                custom_stage = re.search(r"custom--(\d+)--(\S+)", stage.name)
-                if custom_stage is None or len(custom_stage.groups()) != 2:
-                    raise NotImplementedError(f"Unrecognized Stage: {stage.name}")
-                adicity, symbol = custom_stage.groups()
-                adicity = int(adicity)
-                processed = True
-                if adicity == 0:
-                    # We don't need to do anything here
-                    stage = stage.next
-                elif adicity == 1:
-                    processed = process_custom_monadic_connective(infile, symbol, current_model_parts)
-                elif adicity == 2:
-                    processed = process_custom_dyadic_connective(infile, symbol, current_model_parts)
-                else:
-                    raise NotImplementedError("Unable to process connectives of adicity greater than c2")
-
-                if processed:
-                    stage.increment()
-                    stage = stage.next
-                else:
-                    stages.reset_after(stage.name)
-                    stage = stage.previous
-
+                raise NotImplementedError("Custom Connectives are not yet supported")
 
     return solutions
 
@@ -269,8 +156,15 @@ def process_sizes(infile: SourceFile, current_model_parts: ModelBuilder, first_r
         return False
 
     carrier_set = carrier_set_from_size(size)
-    current_model_parts.carrier_set = carrier_set
     current_model_parts.size = size
+    current_model_parts.carrier_set = carrier_set
+
+    # Reset counts
+    current_model_parts.num_negation = 0
+    current_model_parts.num_order = 0
+    current_model_parts.num_designated = 0
+    current_model_parts.num_implication = 0
+    current_model_parts.num_necessitation = 0
 
     return True
 
@@ -280,7 +174,15 @@ def process_negations(infile: SourceFile, current_model_parts: ModelBuilder) ->
     if mnegation is None:
         return False
 
+    current_model_parts.num_negation += 1
     current_model_parts.mnegation = mnegation
+
+    # Reset counts
+    current_model_parts.num_order = 0
+    current_model_parts.num_designated = 0
+    current_model_parts.num_implication = 0
+    current_model_parts.num_necessitation = 0
+
     return True
 
 def process_orders(infile: SourceFile, current_model_parts: ModelBuilder) -> bool:
@@ -290,9 +192,15 @@ def process_orders(infile: SourceFile, current_model_parts: ModelBuilder) -> boo
         return False
 
     mconjunction, mdisjunction = result
+    current_model_parts.num_order += 1
     current_model_parts.mconjunction = mconjunction
     current_model_parts.mdisjunction = mdisjunction
 
+    # Reset counts
+    current_model_parts.num_designated = 0
+    current_model_parts.num_implication = 0
+    current_model_parts.num_necessitation = 0
+
     return True
 
 def process_designateds(infile: SourceFile, current_model_parts: ModelBuilder) -> bool:
@@ -301,48 +209,58 @@ def process_designateds(infile: SourceFile, current_model_parts: ModelBuilder) -
     if designated_values is None:
         return False
 
+    current_model_parts.num_designated += 1
     current_model_parts.designated_values = designated_values
+
+    # Reset counts
+    current_model_parts.num_implication = 0
+    current_model_parts.num_necessitation = 0
+
     return True
 
 def process_implications(infile: SourceFile, current_model_parts: ModelBuilder) -> bool:
     """Stage 5"""
-    mimplication = parse_single_implication(infile, current_model_parts.size)
+    if infile.reststr == "":
+        infile.reststr = next(infile).strip()
+
+    mimplication, reststr = parse_single_implication(infile.reststr, infile.current_line, current_model_parts.size)
+    infile.reststr = reststr
     if mimplication is None:
         return False
 
+    current_model_parts.num_implication += 1
     current_model_parts.mimplication = mimplication
+
+     # Reset counts
+    current_model_parts.num_necessitation = 0
+
     return True
 
 def process_necessitations(infile: SourceFile, current_model_parts: ModelBuilder) -> bool:
-    mnecessitation = parse_single_necessitation(infile, current_model_parts.size)
+    # TODO: In progress
+    if infile.reststr != "":
+        infile.reststr = next(infile).strip()
+
+    # TODO: This should match the implication way of doing things when
+    # custom conenctives are involved (returning reststr)
+    mnecessitation, reststr = parse_single_necessitation_from_str(infile.reststr, infile.current_line, current_model_parts.size)
+    infile.reststr = reststr
     if mnecessitation is None:
         return False
 
+    current_model_parts.num_necessitation += 1
     current_model_parts.mnecessitation = mnecessitation
+
     return True
 
-def process_custom_monadic_connective(infile: SourceFile, symbol: str, current_model_parts: ModelBuilder) -> bool:
-    mfunction = parse_single_monadic_connective(infile, symbol, current_model_parts.size)
-    if mfunction is None:
-        return False
-
-    current_model_parts.custom_model_functions[symbol] = mfunction
-    return True
-
-def process_custom_dyadic_connective(infile: SourceFile, symbol: str, current_model_parts: ModelBuilder) -> bool:
-    mfunction = parse_single_dyadic_connective(infile, symbol, current_model_parts.size)
-    if mfunction is None:
-        return False
-
-    current_model_parts.custom_model_functions[symbol] = mfunction
-    return True
-
-def process_model(model_name: str, mp: ModelBuilder,  solutions: List[Tuple[Model, Dict]]):
+def process_model(mp: ModelBuilder, solutions: List[Tuple[Model, Dict]]):
     """Create Model"""
     assert mp.mimplication is not None
+    assert len(mp.carrier_set) > 0
     assert mp.size + 1 == len(mp.carrier_set)
 
     logical_operations = { mp.mimplication }
+    model_name = f"{mp.size + 1}{'.' + str(mp.num_negation) if mp.num_negation != 0 else ''}.{mp.num_order}.{mp.num_designated}.{mp.num_implication}{'.' + str(mp.num_necessitation) if mp.num_necessitation != 0 else ''}"
     model = Model(mp.carrier_set, logical_operations, mp.designated_values, name=model_name)
     interpretation = {
         Implication: mp.mimplication
@@ -360,12 +278,6 @@ def process_model(model_name: str, mp: ModelBuilder,  solutions: List[Tuple[Mode
         logical_operations.add(mp.mnecessitation)
         interpretation[Necessitation] = mp.mnecessitation
 
-    for custom_mf in mp.custom_model_functions.values():
-        if custom_mf is not None:
-            logical_operations.add(custom_mf)
-            # NOTE: No need to assign interpretation
-            # for VSP check
-
     solutions.append((model, interpretation))
     print(f"Parsed Matrix {model.name}")
 
@@ -373,7 +285,7 @@ def parse_header(infile: SourceFile) -> UglyHeader:
     """
     Parse the header line from the ugly data format.
     """
-    header_line = infile.next_line()
+    header_line = next(infile).strip()
     header_tokens = header_line.split(" ")
     assert header_tokens[0] in ["0", "1"]
     assert header_tokens[6] in ["0", "1"]
@@ -402,12 +314,12 @@ def parse_size(infile: SourceFile, first_run: bool) -> Optional[int]:
     Parse the line representing the matrix size.
     """
 
-    size = int(infile.next_line())
+    size = int(next(infile))
     # HACK: When necessitation and custom connectives are enabled
     # MaGIC may produce -1s at the beginning of the file
     if first_run:
         while size == -1:
-            size = int(infile.next_line())
+            size = int(next(infile))
 
     if size == -1:
         return None
@@ -418,7 +330,20 @@ def parse_single_negation(infile: SourceFile, size: int) -> Optional[ModelFuncti
     """
     Parse the line representing the negation table.
     """
-    return parse_single_monadic_connective(infile, "¬", size)
+    line = next(infile).strip()
+    if line == '-1':
+        return None
+
+    row = line.split(" ")
+    assert len(row) == size + 1, f"Negation table doesn't match size at line {infile.current_line}"
+    mapping = {}
+
+    for i, j in zip(range(size + 1), row):
+        x = mvalue_from_index(i)
+        y = parse_mvalue(j)
+        mapping[(x, )] = y
+
+    return ModelFunction(1, mapping, "¬")
 
 def mvalue_from_index(i: int) -> ModelValue:
     """
@@ -485,7 +410,7 @@ def parse_single_order(infile: SourceFile, size: int) -> Optional[Tuple[ModelFun
     """
     Parse the line representing the ordering table
     """
-    line = infile.next_line()
+    line = next(infile).strip()
     if line == '-1':
         return None
 
@@ -536,7 +461,7 @@ def parse_single_designated(infile: SourceFile, size: int) -> Optional[Set[Model
     """
     Parse the line representing which model values are designated.
     """
-    line = infile.next_line()
+    line = next(infile).strip()
     if line == '-1':
         return None
 
@@ -552,49 +477,17 @@ def parse_single_designated(infile: SourceFile, size: int) -> Optional[Set[Model
 
     return designated_values
 
-def parse_single_implication(infile: SourceFile, size: int) -> Tuple[ModelFunction]:
+def parse_single_implication(instr: str, line: int, size: int) -> Tuple[ModelFunction, str]:
     """
-    Take the current string, parse an implication table from it.
+    Take the current string, parse an implication table from it,
+    and return along with it the remainder of the string
     """
-    return parse_single_dyadic_connective(infile, "→", size)
+    if instr == "-1":
+        return None, ""
 
-def parse_single_necessitation(infile: SourceFile, size: int) -> Optional[ModelFunction]:
-    """
-    Parse the line representing the necessitation table.
-    """
-    return parse_single_monadic_connective(infile, "!", size)
+    table = instr.split(" ")
 
-def parse_single_monadic_connective(infile: SourceFile, symbol: str, size: int) -> Optional[ModelFunction]:
-    line = infile.next_line()
-    if line == '-1':
-        return None
-
-    row = line.split(" ")
-    assert len(row) == size + 1, f"{symbol} table doesn't match size at line {infile.current_line}"
-    mapping = {}
-
-    for i, j in zip(range(size + 1), row):
-        x = mvalue_from_index(i)
-        y = parse_mvalue(j)
-        mapping[(x, )] = y
-
-    return ModelFunction(1, mapping, symbol)
-
-def parse_single_dyadic_connective(infile: SourceFile, symbol: str, size: int) -> Optional[ModelFunction]:
-    try:
-        first_token = next(infile)
-        if first_token == "-1":
-            return None
-    except StopIteration:
-        return None
-
-    table = []
-    try:
-        table = [first_token] + [next(infile) for _ in range((size + 1)**2 - 1)]
-    except StopIteration:
-        pass
-
-    assert len(table) == (size + 1)**2, f"{symbol} table does not match expected size at line {infile.current_line}"
+    assert len(table) >= (size + 1)**2, f"Implication table does not match expected size at line {line}"
 
     mapping = {}
     table_i = 0
@@ -608,7 +501,30 @@ def parse_single_dyadic_connective(infile: SourceFile, symbol: str, size: int) -
 
             mapping[(x, y)] = r
 
-    return ModelFunction(2, mapping, symbol)
+    mimplication = ModelFunction(2, mapping, "→")
+    reststr = " ".join(table[(size + 1)**2:])
+    return mimplication, reststr
+
+def parse_single_necessitation_from_str(instr: str, line: int, size: int) -> Optional[ModelFunction]:
+    """
+    Parse the line representing the necessitation table.
+    """
+    if instr == "-1":
+        return None
+
+    row = instr.split(" ")
+    assert len(row) >= size + 1, f"Necessitation table doesn't match size at line {line}"
+    mapping = {}
+
+    for i, j in zip(range(size + 1), row):
+        x = mvalue_from_index(i)
+        y = parse_mvalue(j)
+        mapping[(x, )] = y
+
+    mnecessitation = ModelFunction(1, mapping, "!")
+    reststr = " ".join(row[(size + 1):])
+
+    return mnecessitation, reststr
 
 
 if __name__ == "__main__":