How do I find the most similar sequences from a large set of short sequences?
130k * 20bp is a small data set. At this scale, SSE2 Smith-Waterman may work well:

git clone https://github.com/attractivechaos/klib
cd klib &amp;&amp; gcc -O2 -D_KSW_MAIN ksw.c -o ksw -lz
./ksw -a1 -b1 -q1 -r1 -t14 20bp.fa 20bp.fa &gt; out.tsv


I simulated 130k 20bp reads from E. coli K-12. The last command line takes about ~2 CPU hours based on the partial output.

It is possible to make this much faster by using 5bp seeds. You can also try ssearch and fasta from the fasta aligner package. They are slower, probably because they are not optimized for such input. Computing Hamming distance with SSE2 will be much faster, but is probably of little use as it does not allow gaps. There are also Gene Myers' edit-distance based algorithms, which can be faster than SSE2-SW.

Bowtie2 (as well as most short-read mappers) won't work well at the default setting. Bowtie2 uses long exact seeds. It will miss many 1-mismatch hits, let alone 3-mismatch ones. The number you get from bowtie2 is an underestimate. You might be able to tune bowtie2, but to find 3-mismatch hits, you have to make it a lot slower.