soedinglab · fnachon · Jul 16, 2026 · Jul 16, 2026 · Jul 16, 2026
diff --git a/CMakeLists.txt b/CMakeLists.txt
@@ -21,13 +21,18 @@ set(USE_SYSTEM_ZSTD 0 CACHE BOOL "Use zstd provided by system instead of bundled
 set(IGNORE_RUST_VERSION 0 CACHE BOOL "Ignore Rust version check")
 set(ENABLE_HIP 0 CACHE BOOL "Enable HIP")
 set(ENABLE_CUDA 0 CACHE BOOL "Enable CUDA")
+set(ENABLE_MPS 0 CACHE BOOL "Enable Metal/MPS GPU backend (Apple Silicon)")
 set(FORCE_STATIC_DEPS 0 CACHE BOOL "Force static linking of deps")
 set(MMSEQS_INT64_IDS 0 CACHE BOOL "Use 64-bit internal sequence/database IDs")
 
 if(ENABLE_HIP)
     set(ENABLE_CUDA 1)
 endif()
 
+if(ENABLE_CUDA AND ENABLE_MPS)
+    message(FATAL_ERROR "ENABLE_CUDA and ENABLE_MPS are mutually exclusive")
+endif()
+
 if(FORCE_STATIC_DEPS)
     if(ENABLE_CUDA)
         set(CMAKE_FIND_LIBRARY_SUFFIXES .a .so CACHE INTERNAL "" FORCE)
@@ -310,6 +315,12 @@ if (ENABLE_CUDA)
     include_directories(lib/libmarv/src)
     add_subdirectory(lib/libmarv/src EXCLUDE_FROM_ALL)
     set_target_properties(marv PROPERTIES POSITION_INDEPENDENT_CODE ON)
+elseif (ENABLE_MPS)
+    # Metal/MPS backend for Apple Silicon; reuses the same marv.h interface
+    # as lib/libmarv so downstream code (GpuUtil, gpuserver, ungappedprefilter)
+    # needs no changes beyond the HAVE_MPS compile definition.
+    include_directories(lib/libmarv/src)
+    add_subdirectory(lib/libmarv-mps EXCLUDE_FROM_ALL)
 endif ()
 
 add_subdirectory(lib/alp)

diff --git a/README.md b/README.md
@@ -97,6 +97,18 @@ To perform searches using GPU acceleration, you can add the `--gpu` flag to `eas
     mmseqs createdb examples/DB.fasta targetDB
     mmseqs makepaddedseqdb targetDB targetDB_padded
     mmseqs easy-search examples/QUERY.fasta targetDB_padded alnRes.m8 tmp --gpu 1
+
+#### Apple Silicon (Metal/MPS) GPU support
+
+> [!NOTE]
+> This is a local, unofficial backend (not part of the upstream project or the precompiled binaries) that ports GPU-accelerated search to Apple Silicon via Metal, as an alternative to the CUDA backend above. It currently only accelerates the default ungapped prefilter (`--prefilter-mode 0`, the same mode used by `--gpu 1` above); the gapped Smith-Waterman GPU kernel (`--prefilter-mode 2`) is CUDA/HIP-only for now.
+
+Build with `-DENABLE_MPS=1` instead of `-DENABLE_CUDA=1` (the two are mutually exclusive):
+
+    cmake -DENABLE_MPS=1 -DHAVE_ARM8=1 -DCMAKE_BUILD_TYPE=Release ..
+    make -j
+
+Usage is identical to the CUDA workflow above (`makepaddedseqdb` + `--gpu 1`); no macOS-specific flags are needed. See [lib/libmarv-mps/marv_mps.mm](lib/libmarv-mps/marv_mps.mm) for the implementation and design notes.
 
 The `databases` workflow provides download and setup procedures for many public reference databases, such as the Uniref, NR, NT, PFAM and many more (see [Downloading databases](https://github.com/soedinglab/mmseqs2/wiki#downloading-databases)). For example, to download and search against a database containing the Swiss-Prot reference proteins run: 
 

diff --git a/lib/libmarv-mps/CMakeLists.txt b/lib/libmarv-mps/CMakeLists.txt
@@ -0,0 +1,17 @@
+cmake_minimum_required(VERSION 3.18)
+
+if (NOT APPLE)
+    message(FATAL_ERROR "libmarv-mps (Metal/MPS GPU backend) requires macOS")
+endif ()
+
+enable_language(OBJCXX)
+
+add_library(marv STATIC marv_mps.mm)
+set_source_files_properties(marv_mps.mm PROPERTIES LANGUAGE OBJCXX)
+target_compile_features(marv PUBLIC cxx_std_17)
+target_compile_options(marv PRIVATE -fobjc-arc)
+set_target_properties(marv PROPERTIES POSITION_INDEPENDENT_CODE ON)
+
+find_library(METAL_FRAMEWORK Metal REQUIRED)
+find_library(FOUNDATION_FRAMEWORK Foundation REQUIRED)
+target_link_libraries(marv PUBLIC ${METAL_FRAMEWORK} ${FOUNDATION_FRAMEWORK})
diff --git a/lib/libmarv-mps/marv_mps.mm b/lib/libmarv-mps/marv_mps.mm
diff --git a/lib/libmarv/src/marv.h b/lib/libmarv/src/marv.h
@@ -47,6 +47,16 @@ class Marv {
     //sequence must be encoded
     Stats scan(const char* sequence,  size_t sequenceLength, int8_t* pssm, Result* results);
 
+    // Optional asynchronous 2-phase form of scan(), letting a caller overlap
+    // building the next query's profile with this query's GPU execution.
+    // submitScan() dispatches GPU work without blocking; it must be followed
+    // by exactly one collectScan() (which blocks until that work completes)
+    // before submitScan() may be called again. Only the Metal/MPS backend
+    // implements these; the CUDA/HIP backend does not define them, so
+    // callers must guard use behind #ifdef HAVE_MPS.
+    void submitScan(const char* sequence, size_t sequenceLength, int8_t* pssm);
+    Stats collectScan(Result* results);
+
 private:
     size_t dbEntries;
     int alphabetSize;

diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
@@ -284,6 +284,11 @@ if (ENABLE_CUDA)
     endif ()
 endif ()
 
+if (ENABLE_MPS)
+    target_compile_definitions(mmseqs-framework PUBLIC -DHAVE_MPS=1)
+    target_link_libraries(mmseqs-framework marv)
+endif ()
+
 if (NOT FRAMEWORK_ONLY)
     include(MMseqsSetupDerivedTarget)
     add_subdirectory(version)

diff --git a/src/commons/GpuUtil.cpp b/src/commons/GpuUtil.cpp
@@ -1,4 +1,4 @@
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
 #include "GpuUtil.h"
 
 #include "Debug.h"

diff --git a/src/commons/Parameters.cpp b/src/commons/Parameters.cpp
@@ -2625,15 +2625,15 @@ void Parameters::setDefaults() {
 
     // gpu
     gpu = 0;
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
     char* gpuEnv = getenv("MMSEQS_FORCE_GPU");
     if (gpuEnv != NULL) {
         gpu = 1;
     }
 #endif
     gpuServer = 0;
     gpuServerWaitTimeout = 10 * 60;
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
     char* gpuServerEnv = getenv("MMSEQS_FORCE_GPUSERVER");
     if (gpuServerEnv != NULL) {
         gpuServer = 1;

diff --git a/src/prefiltering/ungappedprefilter.cpp b/src/prefiltering/ungappedprefilter.cpp
@@ -25,13 +25,13 @@
 #include <omp.h>
 #endif
 // #define HAVE_CUDA 1
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
 #include "GpuUtil.h"
 #include "Alignment.h"
 #include <signal.h>
 #endif
 
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
 
 volatile sig_atomic_t keepRunningClient = 1;
 void intHandlerClient(int) {
@@ -168,6 +168,96 @@ void runFilterOnGpu(Parameters & par, BaseMatrix * subMat,
         sigaction(SIGTERM, &act, NULL);
     }
 
+    // Processes one query's hits (already sitting in `results`) and writes
+    // them out. Factored out so the MPS pipelined loop below can call it for
+    // the *previous* query while the current query's GPU work is in flight.
+    auto processResults = [&](const Marv::Stats& stats, size_t queryKeyForResults, int queryLenForResults) {
+        for(size_t i = 0; i < stats.results; i++){
+            DBKeyType targetKey = tdbr->getDbKey(results[i].id);
+            int score = results[i].score;
+            if(taxonomyHook != NULL){
+                TaxID currTax = taxonomyHook->taxonomyMapping->lookup(targetKey);
+                if (taxonomyHook->expression[0]->isAncestor(currTax) == false) {
+                    continue;
+                }
+            }
+            // check if evalThr != inf
+            // double evalue = 0.0;
+            // if (par.evalThr < std::numeric_limits<double>::max()) {
+            //     evalue = evaluer->computeEvalue(score, queryLenForResults);
+            // }
+            // bool hasEvalue = (evalue <= par.evalThr);
+            bool hasDiagScore = (score > par.minDiagScoreThr);
+
+            const bool isIdentity = (queryKeyForResults == targetKey && (par.includeIdentity || sameDB))? true : false;
+            // --filter-hits
+            if (isIdentity || hasDiagScore) {
+                if(par.prefMode == Parameters::PREF_MODE_UNGAPPED_AND_GAPPED){
+                    Matcher::result_t res;
+                    res.dbKey = targetKey;
+                    res.eval = evaluer->computeEvalue(score, queryLenForResults);
+                    res.dbEndPos = results[i].dbEndPos;
+                    res.dbLen = tdbr->getSeqLen(results[i].id);
+                    res.qEndPos =  results[i].qEndPos;
+                    res.qLen = queryLenForResults;
+                    unsigned int qAlnLen = std::max(static_cast<unsigned int>(res.qEndPos), static_cast<unsigned int>(1));
+                    unsigned int dbAlnLen = std::max(static_cast<unsigned int>(res.dbEndPos), static_cast<unsigned int>(1));
+                    //seqId = (alignment.score1 / static_cast<float>(std::max(dbAlnLen, qAlnLen)))  * 0.1656 + 0.1141;
+                    res.seqId = Matcher::estimateSeqIdByScorePerCol(score, qAlnLen, dbAlnLen);
+                    res.qcov = SmithWaterman::computeCov(0, res.qEndPos, res.qLen );
+                    res.dbcov = SmithWaterman::computeCov(0, res.dbEndPos, res.dbLen );
+                    res.score = evaluer->computeBitScore(score);
+                    if(Alignment::checkCriteria(res, isIdentity, par.evalThr,  par.seqIdThr,  par.alnLenThr,  par.covMode,  par.covThr)){
+                        resultsAln.emplace_back(res);
+                    }
+                } else {
+                    hit_t hit;
+                    hit.seqId = targetKey;
+                    hit.prefScore = score;
+                    hit.diagonal = 0;
+                    shortResults.emplace_back(hit);
+                }
+            }
+        }
+        if(par.prefMode == Parameters::PREF_MODE_UNGAPPED_AND_GAPPED) {
+            SORT_PARALLEL(resultsAln.begin(), resultsAln.end(), Matcher::compareHits);
+            size_t maxSeqs = std::min(par.maxResListLen, resultsAln.size());
+            for (size_t i = 0; i < maxSeqs; ++i) {
+                size_t len = Matcher::resultToBuffer(buffer, resultsAln[i], false);
+                resultBuffer.append(buffer, len);
+            }
+        }else{
+            SORT_PARALLEL(shortResults.begin(), shortResults.end(), hit_t::compareHitsByScoreAndId);
+            size_t maxSeqs = std::min(par.maxResListLen, shortResults.size());
+            for (size_t i = 0; i < maxSeqs; ++i) {
+                size_t len = QueryMatcher::prefilterHitToBuffer(buffer, shortResults[i]);
+                resultBuffer.append(buffer, len);
+            }
+        }
+
+        resultWriter.writeData(resultBuffer.c_str(), resultBuffer.length(), queryKeyForResults, 0);
+        resultBuffer.clear();
+        shortResults.clear();
+        resultsAln.clear();
+        progress.updateProgress();
+    };
+
+    // On the Metal/MPS backend (non-server mode only), pipeline queries:
+    // submit query i's GPU work without blocking, then build query i+1's
+    // profile on the CPU while it runs, and only then collect query i's
+    // results. This overlaps CPU-side profile/composition-bias work with
+    // GPU execution instead of serializing them. The CUDA/HIP path and
+    // gpuserver shared-memory protocol are untouched. These are only
+    // referenced from HAVE_MPS-guarded code below, so they're declared
+    // here under the same guard to avoid an unused-variable error on the
+    // CUDA/HIP build (which compiles with -Werror).
+#ifdef HAVE_MPS
+    const bool pipeline = (serverMode == 0);
+    size_t prevQueryKey = 0;
+    int prevQueryLen = 0;
+    bool havePending = false;
+#endif
+
     // marv.prefetch();
     for (size_t id = 0; id < qdbr->getSize(); id++) {
         if (!keepRunningClient) {
@@ -202,6 +292,24 @@ void runFilterOnGpu(Parameters & par, BaseMatrix * subMat,
                 }
             }
         }
+
+#ifdef HAVE_MPS
+        if (pipeline) {
+            if (havePending) {
+                Marv::Stats stats = marv->collectScan(results.data());
+                if (keepRunningClient == false) {
+                    EXIT(EXIT_FAILURE);
+                }
+                processResults(stats, prevQueryKey, prevQueryLen);
+            }
+            marv->submitScan(reinterpret_cast<const char *>(qSeq.numSequence), qSeq.L, profile);
+            prevQueryKey = queryKey;
+            prevQueryLen = qSeq.L;
+            havePending = true;
+            continue;
+        }
+#endif
+
         Marv::Stats stats;
         if (serverMode == 0) {
             stats = marv->scan(reinterpret_cast<const char *>(qSeq.numSequence), qSeq.L, profile, results.data());
@@ -259,75 +367,17 @@ void runFilterOnGpu(Parameters & par, BaseMatrix * subMat,
             EXIT(EXIT_FAILURE);
         }
 
-        for(size_t i = 0; i < stats.results; i++){
-            DBKeyType targetKey = tdbr->getDbKey(results[i].id);
-            int score = results[i].score;
-            if(taxonomyHook != NULL){
-                TaxID currTax = taxonomyHook->taxonomyMapping->lookup(targetKey);
-                if (taxonomyHook->expression[0]->isAncestor(currTax) == false) {
-                    continue;
-                }
-            }
-            // check if evalThr != inf
-            // double evalue = 0.0;
-            // if (par.evalThr < std::numeric_limits<double>::max()) {
-            //     evalue = evaluer->computeEvalue(score, qSeq.L);
-            // }
-            // bool hasEvalue = (evalue <= par.evalThr);
-            bool hasDiagScore = (score > par.minDiagScoreThr);
-
-            const bool isIdentity = (queryKey == targetKey && (par.includeIdentity || sameDB))? true : false;
-            // --filter-hits
-            if (isIdentity || hasDiagScore) {
-                if(par.prefMode == Parameters::PREF_MODE_UNGAPPED_AND_GAPPED){
-                    Matcher::result_t res;
-                    res.dbKey = targetKey;
-                    res.eval = evaluer->computeEvalue(score, qSeq.L);
-                    res.dbEndPos = results[i].dbEndPos;
-                    res.dbLen = tdbr->getSeqLen(results[i].id);
-                    res.qEndPos =  results[i].qEndPos;
-                    res.qLen = qSeq.L;
-                    unsigned int qAlnLen = std::max(static_cast<unsigned int>(res.qEndPos), static_cast<unsigned int>(1));
-                    unsigned int dbAlnLen = std::max(static_cast<unsigned int>(res.dbEndPos), static_cast<unsigned int>(1));
-                    //seqId = (alignment.score1 / static_cast<float>(std::max(dbAlnLen, qAlnLen)))  * 0.1656 + 0.1141;
-                    res.seqId = Matcher::estimateSeqIdByScorePerCol(score, qAlnLen, dbAlnLen);
-                    res.qcov = SmithWaterman::computeCov(0, res.qEndPos, res.qLen );
-                    res.dbcov = SmithWaterman::computeCov(0, res.dbEndPos, res.dbLen );
-                    res.score = evaluer->computeBitScore(score);
-                    if(Alignment::checkCriteria(res, isIdentity, par.evalThr,  par.seqIdThr,  par.alnLenThr,  par.covMode,  par.covThr)){
-                        resultsAln.emplace_back(res);
-                    }
-                } else {
-                    hit_t hit;
-                    hit.seqId = targetKey;
-                    hit.prefScore = score;
-                    hit.diagonal = 0;
-                    shortResults.emplace_back(hit);
-                }
-            }
-        }
-        if(par.prefMode == Parameters::PREF_MODE_UNGAPPED_AND_GAPPED) {
-            SORT_PARALLEL(resultsAln.begin(), resultsAln.end(), Matcher::compareHits);
-            size_t maxSeqs = std::min(par.maxResListLen, resultsAln.size());
-            for (size_t i = 0; i < maxSeqs; ++i) {
-                size_t len = Matcher::resultToBuffer(buffer, resultsAln[i], false);
-                resultBuffer.append(buffer, len);
-            }
-        }else{
-            SORT_PARALLEL(shortResults.begin(), shortResults.end(), hit_t::compareHitsByScoreAndId);
-            size_t maxSeqs = std::min(par.maxResListLen, shortResults.size());
-            for (size_t i = 0; i < maxSeqs; ++i) {
-                size_t len = QueryMatcher::prefilterHitToBuffer(buffer, shortResults[i]);
-                resultBuffer.append(buffer, len);
-            }
+        processResults(stats, queryKey, qSeq.L);
+    }
+#ifdef HAVE_MPS
+    if (pipeline && havePending) {
+        Marv::Stats stats = marv->collectScan(results.data());
+        if (keepRunningClient == false) {
+            EXIT(EXIT_FAILURE);
         }
-
-        resultWriter.writeData(resultBuffer.c_str(), resultBuffer.length(), queryKey, 0);
-        resultBuffer.clear();
-        shortResults.clear();
-        resultsAln.clear();
-        progress.updateProgress();
+        processResults(stats, prevQueryKey, prevQueryLen);
     }
+#endif
     if (marv != NULL) {
         delete marv;
     } else {
@@ -554,11 +604,11 @@ int prefilterInternal(int argc, const char **argv, const Command &command, int m
         taxonomyHook = new QueryMatcherTaxonomyHook(par.db2, tdbr, par.taxonList, par.threads);
     }
     if(par.gpu){
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
         runFilterOnGpu(par, subMat, qdbr, tdbr, sameDB,
                        resultWriter, evaluer, taxonomyHook);
 #else
-        Debug(Debug::ERROR) << "MMseqs2 was compiled without CUDA support\n";
+        Debug(Debug::ERROR) << "MMseqs2 was compiled without GPU support (CUDA/HIP/MPS)\n";
         EXIT(EXIT_FAILURE);
 #endif
     }else{

diff --git a/src/util/gpuserver.cpp b/src/util/gpuserver.cpp
@@ -7,7 +7,7 @@
 #include "NucleotideMatrix.h"
 
 // #define HAVE_CUDA 1
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
 #include "GpuUtil.h"
 #endif
 
@@ -24,7 +24,7 @@ void intHandler(int) {
 int gpuserver(int argc, const char **argv, const Command& command) {
     Parameters& par = Parameters::getInstance();
     par.parseParameters(argc, argv, command, true, 0, 0);
-#ifdef HAVE_CUDA
+#if defined(HAVE_CUDA) || defined(HAVE_MPS)
     bool touch = (par.preloadMode != Parameters::PRELOAD_MODE_MMAP);
     IndexReader dbrIdx(par.db1, par.threads, IndexReader::SEQUENCES, (touch) ? (IndexReader::PRELOAD_INDEX | IndexReader::PRELOAD_DATA) : 0 );
     DBReader<DBKeyType>* dbr = dbrIdx.sequenceReader;