docs: def get_mutations(self, pdb_id_1, pdb_id_2, alignment_tool = 'clustalw', gap_opening_penalty = 0.2):
          '''Returns a mapping chain_of_pdb_1 -> List[PDBMutationPair] representing the mutations needed to transform each chain of pdb_1 into the respective chains of pdb_2.
             This function also sets self.seqres_sequence_maps, a mapping Tuple(pdb_id_1, chain ID in pdb_id_1) -> Tuple(pdb_id_2, chain ID in pdb_id_2) -> a SequenceMap representing the mapping of residues between the chains based on the alignment.
  
             Warning: This function does not consider what happens if a chain in pdb_id_1 matches two chains in pdb_id_2
                      which have differing sequences. In this case, an self-inconsistent set of mutations is returned.
                      One solution would be to extend the mapping to:
  
                          chain_m_of_pdb_1 -> common_mutations -> List[PDBMutationPair]
                                           -> chain_x_of_pdb_2 -> List[PDBMutationPair]
                                           -> chain_y_of_pdb_2 -> List[PDBMutationPair]
                                           ...
  
                      where common_mutations contains the set of mutations common to the mapping m->x and m->y etc. whereas the
                      other two mappings contain the set of mutations from m->x or m->y etc. respectively which do not occur in
                      common_mutations. In general, both chain_x_of_pdb_2 and chain_y_of_pdb_2 will be empty excepting the
                      considered case where x and y differ in sequence.
          '''
  
          # Set up the objects
          p1 = self.add_pdb(pdb_id_1)
          p2 = self.add_pdb(pdb_id_2)
          self.chain_map[pdb_id_1] = self.chain_map.get(pdb_id_1, {})
  
          # Determine which chains map to which
          alignment = self.get_alignment(pdb_id_1, pdb_id_2, alignment_tool = alignment_tool, gap_opening_penalty = gap_opening_penalty)
          best_matches = alignment['best_matches']
  
          # Create the list of mutations
          mutations = {}
          for from_chain, mtches in sorted(best_matches.iteritems()):
  
              from_pdb_id, from_chain_id = from_chain.split('_')
  
              # Only consider matches from pdb_id_1 to pdb_id_2
              if from_pdb_id == pdb_id_1:
  
                  self.seqres_sequence_maps[(from_pdb_id, from_chain_id)] = {}
  
                  self.chain_map[from_pdb_id][from_chain_id] = self.chain_map[from_pdb_id].get(from_chain_id, {})
  
                  # Do not consider matches from pdb_id_1 to itself or matches with poor sequence similarity
                  restricted_mtchs = {}
                  for to_chain, similarity in sorted(mtches.iteritems()):
                      if to_chain.split('_')[0] == pdb_id_2 and similarity >= self.acceptable_sequence_similarity:
                          restricted_mtchs[to_chain] = similarity
  
                  # Take the best matching chains and create a list of mutations needed to transform from_chain to those chains
                  # Warning: This does NOT take into account whether the sequences of the best matches differ.
                  if restricted_mtchs:
                      top_similarity = max(restricted_mtchs.values())
  
                      #todo: if the sequences of the best matches differ, raise an Exception. Use 2ZNW and 1DQJ as an example (2ZNW chain A matches with 48% to both 1DQJ chain A and chain B)
                      #top_matches = [to_chain for to_chain, similarity in sorted(restricted_mtchs.iteritems()) if similarity == top_similarity]
                      #pprint.pprint(restricted_mtchs)
                      #print(from_pdb_id, from_chain, 'top_matches', top_matches)
                      #sys.exit(0)
  
  
                      for to_chain, similarity in sorted(restricted_mtchs.iteritems()):
                          to_pdb_id, to_chain_id = to_chain.split('_')
                          if similarity == top_similarity:
                              #print(from_pdb_id, from_chain_id)
                              #print(restricted_mtchs)
                              #print(to_pdb_id, to_chain, similarity)
                              self.chain_map[from_pdb_id][from_chain_id][to_pdb_id] = self.chain_map[from_pdb_id][from_chain_id].get(to_pdb_id, set())
                              self.chain_map[from_pdb_id][from_chain_id][to_pdb_id].add(to_chain)
                              mutations[from_chain_id] = mutations.get(from_chain_id, [])
                              chain_mutations = self.get_chain_mutations(from_pdb_id, from_chain_id, to_pdb_id, to_chain_id)
                              mutations[from_chain_id].extend(chain_mutations)
  
          # mutations can contain duplicates so we remove those
          for chain_id, mlist in mutations.iteritems():
              mutations[chain_id] = sorted(set(mlist))
          return mutations