Noncoprime Subsequences
Source: aizu

Noncoprime Subsequences
Problem Statement

Given a sequence $A = (A_1, A_2, ..., A_N)$, a good subsequence of $A$ is defined as a subsequence, that is not necessarily contiguous, where adjacent elements in the subsequence are not coprime.

Find the maximum length $L$ of a good subsequence of $A$. Also, determine the number of good subsequences of length $L$, modulo $998,244,353$.

Input

The input is given in the following format:

$N$
$A_1$ $A_2$ ... $A_N$

$1 \leq N \leq 200,000$
$1 \leq A_i \leq 10^6$
All input values are integers.
Output

Output 2 lines. On the first line, output $L$. On the second line, output the number of good subsequences of length $L$ of $A$, modulo $998,244,353$.

Sample Input 1
3
2 3 6

Sample Output 1
2
2

Sample Input 2
5
1 1 1 1 1

Sample Output 2
1
5

Sample Input 3
10
631932 735902 895728 78537 723857 330739 286918 329211 539679 238506

Sample Output 3
7
2


C++ solution:
// AOJ #3636 Noncoprime Subsequences
// 2025.4.20

#include <bits/stdc++.h>
using namespace std;
using ll = long long;

#define gc() getchar_unlocked()
#define pc(c) putchar_unlocked(c)

int Cin() {
    int n = 0; int c = gc();
    do n = 10*n + (c & 0xf), c = gc(); while (c >= '0');
    return n;
}

void Couts(const string &s) {
    for (char c : s) pc(c);
    pc('\n');
}

const int MAXA = 1000000;
const int MOD  = 998244353;

int spf[MAXA+1];
int dp_len[MAXA+1];
int dp_cnt[MAXA+1];

int main(){
    int N = Cin();
    vector<int> A(N);
    for(int i = 0; i < N; i++) A[i] = Cin();

    for(int i = 2; i <= MAXA; i++){
        if(spf[i] == 0) for(int j = i; j <= MAXA; j += i) if(spf[j] == 0) spf[j] = i;
    }

    int  L = 0;
    ll ans = 0;

    static int ds[1<<7], dsH[1<<7];
    vector<int> ps; ps.reserve(7);
    vector<int> H;  H.reserve(7);

    for(int x : A){
        ps.clear();
        int t = x;
        while(t > 1){
            int p = spf[t];
            ps.push_back(p);
            while(t % p == 0) t /= p;
        }

        int best = 0;
        for(int p : ps) best = max(best, dp_len[p]);

        ll cnt;
        if(best == 0) cnt = 1;
        else {
            H.clear();
            for(int p : ps) if(dp_len[p] == best) H.push_back(p);

            int k = (int)H.size();
            int M = 1 << k;
            dsH[0] = 1;
            for(int mask = 1; mask < M; mask++){
                int lsb = mask & -mask;
                int i   = __builtin_ctz(lsb);
                dsH[mask] = dsH[mask ^ lsb] * H[i];
            }
            cnt = 0;
            for(int mask = 1; mask < M; mask++){
                int d = dsH[mask];
                if(dp_len[d] == best){
                    if(__builtin_popcount(mask) & 1) cnt = (cnt + dp_cnt[d]) % MOD;
                    else                             cnt = (cnt - dp_cnt[d] + MOD) % MOD;
                }
            }
        }

        int curL = best + 1;
        if(curL > L) L = curL, ans = cnt;
        else if(curL == L) ans = (ans + cnt) % MOD;

        int m = (int)ps.size();
        int M = 1 << m;
        ds[0] = 1;
        for(int mask = 1; mask < M; mask++){
            int lsb = mask & -mask;
            int i   = __builtin_ctz(lsb);
            ds[mask] = ds[mask ^ lsb] * ps[i];
        }
        for(int mask = 1; mask < M; mask++){
            int d = ds[mask];
            if(dp_len[d] < curL){
                dp_len[d] = curL;
                dp_cnt[d] = (int)(cnt % MOD);
            } else if(dp_len[d] == curL) dp_cnt[d] = (dp_cnt[d] + (int)(cnt % MOD)) % MOD;
        }
    }
    cout << L << endl;
    cout << (ans % MOD + MOD) % MOD << endl;
    return 0;
}

Go solution:
package main

import (
	"bufio"
	"fmt"
	"math"
	"os"
	"strconv"
	"strings"
)

const (
	MAXA = 1000000
	MOD  = 998244353
)

var spf [MAXA + 1]int
var dpLen [MAXA + 1]int
var dpCnt [MAXA + 1]int

func main() {
	scanner := bufio.NewScanner(os.Stdin)
	scanner.Split(bufio.ScanWords)

	// Precompute smallest prime factors
	for i := 2; i <= MAXA; i++ {
		if spf[i] == 0 {
			spf[i] = i
			for j := i * 2; j <= MAXA; j += i {
				if spf[j] == 0 {
					spf[j] = i
				}
			}
		}
	}

	scanner.Scan()
	N, _ := strconv.Atoi(scanner.Text())
	A := make([]int, N)
	for i := 0; i < N; i++ {
		scanner.Scan()
		A[i], _ = strconv.Atoi(scanner.Text())
	}

	L := 0
	ans := 0

	ps := make([]int, 0, 7)
	H := make([]int, 0, 7)
	ds := make([]int, 1<<7)
	dsH := make([]int, 1<<7)

	for _, x := range A {
		// Factorize x into its prime factors
		ps = ps[:0]
		t := x
		for t > 1 {
			p := spf[t]
			ps = append(ps, p)
			for t%p == 0 {
				t /= p
			}
		}

		// Find the best length among the prime factors
		best := 0
		for _, p := range ps {
			if dpLen[p] > best {
				best = dpLen[p]
			}
		}

		var cnt int64
		if best == 0 {
			cnt = 1
		} else {
			// Collect primes with dpLen[p] == best
			H = H[:0]
			for _, p := range ps {
				if dpLen[p] == best {
					H = append(H, p)
				}
			}

			k := len(H)
			M := 1 << k
			dsH[0] = 1
			for mask := 1; mask < M; mask++ {
				lsb := mask & -mask
				i := 0
				for (lsb >> i) > 0 {
					i++
				}
				i--
				dsH[mask] = dsH[mask^lsb] * H[i]
			}

			cnt = 0
			for mask := 1; mask < M; mask++ {
				d := dsH[mask]
				if dpLen[d] == best {
					bits := 0
					tmp := mask
					for tmp > 0 {
						bits += tmp & 1
						tmp >>= 1
					}
					if bits%2 == 1 {
						cnt = (cnt + int64(dpCnt[d])) % MOD
					} else {
						cnt = (cnt - int64(dpCnt[d]) + MOD) % MOD
					}
				}
			}
		}

		curL := best + 1
		if curL > L {
			L = curL
			ans = int(cnt % MOD)
		} else if curL == L {
			ans = (ans + int(cnt)) % MOD
		}

		// Update dpLen and dpCnt for all divisors of x (combinations of its primes)
		m := len(ps)
		M := 1 << m
		ds[0] = 1
		for mask := 1; mask < M; mask++ {
			lsb := mask & -mask
			i := 0
			for (lsb >> i) > 0 {
				i++
			}
			i--
			ds[mask] = ds[mask^lsb] * ps[i]
		}

		for mask := 1; mask < M; mask++ {
			d := ds[mask]
			if dpLen[d] < curL {
				dpLen[d] = curL
				dpCnt[d] = int(cnt % MOD)
			} else if dpLen[d] == curL {
				dpCnt[d] = (dpCnt[d] + int(cnt)) % MOD
			}
		}
	}

	fmt.Println(L)
	fmt.Println(ans)
}

Java solution:
import java.io.*;
import java.util.*;

public class Main {
    static final int MAXA = 1000000;
    static final int MOD = 998244353;

    public static void main(String[] args) throws IOException {
        BufferedReader br = new BufferedReader(new InputStreamReader(System.in));
        int N = Integer.parseInt(br.readLine());
        int[] A = new int[N];
        StringTokenizer st = new StringTokenizer(br.readLine());
        for (int i = 0; i < N; i++) {
            A[i] = Integer.parseInt(st.nextToken());
        }

        // Precompute smallest prime factors
        int[] spf = new int[MAXA + 1];
        for (int i = 2; i <= MAXA; i++) {
            if (spf[i] == 0) {
                for (int j = i; j <= MAXA; j += i) {
                    if (spf[j] == 0) {
                        spf[j] = i;
                    }
                }
            }
        }

        int[] dpLen = new int[MAXA + 1];
        int[] dpCnt = new int[MAXA + 1];
        int L = 0;
        long ans = 0;

        for (int x : A) {
            // Factorize x into distinct primes
            List<Integer> ps = new ArrayList<>();
            int t = x;
            while (t > 1) {
                int p = spf[t];
                ps.add(p);
                while (t % p == 0) {
                    t /= p;
                }
            }

            // Find the best length among the primes
            int best = 0;
            for (int p : ps) {
                if (dpLen[p] > best) {
                    best = dpLen[p];
                }
            }

            long cnt;
            if (best == 0) {
                cnt = 1;
            } else {
                List<Integer> H = new ArrayList<>();
                for (int p : ps) {
                    if (dpLen[p] == best) {
                        H.add(p);
                    }
                }
                int k = H.size();
                int M = 1 << k;
                int[] dsH = new int[M];
                dsH[0] = 1;
                for (int mask = 1; mask < M; mask++) {
                    int lsb = mask & -mask;
                    int i = Integer.numberOfTrailingZeros(lsb);
                    dsH[mask] = dsH[mask ^ lsb] * H.get(i);
                }
                cnt = 0;
                for (int mask = 1; mask < M; mask++) {
                    int d = dsH[mask];
                    if (dpLen[d] == best) {
                        if (Integer.bitCount(mask) % 2 == 1) {
                            cnt += dpCnt[d];
                        } else {
                            cnt -= dpCnt[d];
                        }
                    }
                }
                cnt %= MOD;
                if (cnt < 0) {
                    cnt += MOD;
                }
            }

            int curL = best + 1;
            if (curL > L) {
                L = curL;
                ans = cnt;
            } else if (curL == L) {
                ans = (ans + cnt) % MOD;
            }

            // Update dpLen and dpCnt for all divisors formed by the primes of x
            int m = ps.size();
            int M = 1 << m;
            int[] ds = new int[M];
            ds[0] = 1;
            for (int mask = 1; mask < M; mask++) {
                int lsb = mask & -mask;
                int i = Integer.numberOfTrailingZeros(lsb);
                ds[mask] = ds[mask ^ lsb] * ps.get(i);
            }
            for (int mask = 1; mask < M; mask++) {
                int d = ds[mask];
                if (dpLen[d] < curL) {
                    dpLen[d] = curL;
                    dpCnt[d] = (int) cnt;
                } else if (dpLen[d] == curL) {
                    dpCnt[d] = (dpCnt[d] + (int) cnt) % MOD;
                }
            }
        }
        System.out.println(L);
        System.out.println((ans % MOD + MOD) % MOD);
    }
}

JavaScript solution:
const MOD = 998244353;
const MAXA = 1000000;

function main() {
    const input = require('fs').readFileSync(0, 'utf-8').split(/\s+/).filter(Boolean);
    let ptr = 0;
    const N = parseInt(input[ptr++], 10);
    const A = [];
    for (let i = 0; i < N; i++) {
        A.push(parseInt(input[ptr++], 10));
    }

    // Precompute smallest prime factors (SPF)
    const spf = new Array(MAXA + 1).fill(0);
    for (let i = 2; i <= MAXA; i++) {
        if (spf[i] === 0) {
            for (let j = i; j <= MAXA; j += i) {
                if (spf[j] === 0) {
                    spf[j] = i;
                }
            }
        }
    }

    const dp_len = new Array(MAXA + 1).fill(0);
    const dp_cnt = new Array(MAXA + 1).fill(0);

    let L = 0;
    let ans = 0;

    for (const x of A) {
        const ps = [];
        let t = x;
        while (t > 1) {
            const p = spf[t];
            ps.push(p);
            while (t % p === 0) {
                t = Math.floor(t / p);
            }
        }

        let best = 0;
        for (const p of ps) {
            best = Math.max(best, dp_len[p]);
        }

        let cnt;
        if (best === 0) {
            cnt = 1;
        } else {
            const H = [];
            for (const p of ps) {
                if (dp_len[p] === best) {
                    H.push(p);
                }
            }

            const k = H.length;
            const M = 1 << k;
            const dsH = new Array(M).fill(0);
            dsH[0] = 1;
            for (let mask = 1; mask < M; mask++) {
                const lsb = mask & -mask;
                const i = Math.log2(lsb);
                dsH[mask] = dsH[mask ^ lsb] * H[i];
            }

            cnt = 0;
            for (let mask = 1; mask < M; mask++) {
                const d = dsH[mask];
                if (dp_len[d] === best) {
                    if ((mask.toString(2).split('1').length - 1) % 2 === 1) {
                        cnt = (cnt + dp_cnt[d]) % MOD;
                    } else {
                        cnt = (cnt - dp_cnt[d] + MOD) % MOD;
                    }
                }
            }
        }

        const curL = best + 1;
        if (curL > L) {
            L = curL;
            ans = cnt;
        } else if (curL === L) {
            ans = (ans + cnt) % MOD;
        }

        const m = ps.length;
        const M = 1 << m;
        const ds = new Array(M).fill(0);
        ds[0] = 1;
        for (let mask = 1; mask < M; mask++) {
            const lsb = mask & -mask;
            const i = Math.log2(lsb);
            ds[mask] = ds[mask ^ lsb] * ps[i];
        }

        for (let mask = 1; mask < M; mask++) {
            const d = ds[mask];
            if (dp_len[d] < curL) {
                dp_len[d] = curL;
                dp_cnt[d] = cnt % MOD;
            } else if (dp_len[d] === curL) {
                dp_cnt[d] = (dp_cnt[d] + cnt) % MOD;
            }
        }
    }

    console.log(L);
    console.log(ans % MOD);
}

main();

Python solution:
import sys
import math
from math import gcd
from collections import defaultdict

MOD = 998244353

def main():
    input = sys.stdin.read().split()
    ptr = 0
    N = int(input[ptr])
    ptr += 1
    A = list(map(int, input[ptr:ptr+N]))
    
    MAXA = 10**6
    spf = [0] * (MAXA + 1)
    for i in range(2, MAXA + 1):
        if spf[i] == 0:
            spf[i] = i
            for j in range(i * 2, MAXA + 1, i):
                if spf[j] == 0:
                    spf[j] = i
    
    dp_len = defaultdict(int)
    dp_cnt = defaultdict(int)
    
    L = 0
    ans = 0
    
    for x in A:
        # Factorize x into its prime factors
        t = x
        primes = set()
        while t > 1:
            p = spf[t]
            primes.add(p)
            while t % p == 0:
                t //= p
        primes = sorted(primes)
        m = len(primes)
        
        # Find the best previous length
        best = 0
        for p in primes:
            if dp_len[p] > best:
                best = dp_len[p]
        
        # Calculate the count for the new length
        cnt = 0
        if best == 0:
            cnt = 1
        else:
            # Collect primes with dp_len[p] == best
            H = [p for p in primes if dp_len[p] == best]
            k = len(H)
            M = 1 << k
            # Precompute dsH (products of subsets of H)
            dsH = [1] * M
            for mask in range(1, M):
                lsb = mask & -mask
                i = (lsb - 1).bit_length()
                dsH[mask] = dsH[mask ^ (1 << i)] * H[i]
            # Inclusion-exclusion over subsets
            for mask in range(1, M):
                d = dsH[mask]
                if dp_len[d] == best:
                    bits = bin(mask).count('1')
                    if bits % 2 == 1:
                        cnt += dp_cnt[d]
                    else:
                        cnt -= dp_cnt[d]
            cnt %= MOD
        
        curL = best + 1
        if curL > L:
            L = curL
            ans = cnt
        elif curL == L:
            ans = (ans + cnt) % MOD
        
        # Update dp_len and dp_cnt for all divisors formed by primes of x
        M = 1 << m
        ds = [1] * M
        for mask in range(1, M):
            lsb = mask & -mask
            i = (lsb - 1).bit_length()
            ds[mask] = ds[mask ^ (1 << i)] * primes[i]
        
        for mask in range(1, M):
            d = ds[mask]
            if dp_len[d] < curL:
                dp_len[d] = curL
                dp_cnt[d] = cnt
            elif dp_len[d] == curL:
                dp_cnt[d] = (dp_cnt[d] + cnt) % MOD
    
    print(L)
    print(ans % MOD)

if __name__ == '__main__':
    main()

Ruby solution:
require 'set'

MOD = 998244353
MAXA = 1_000_000

def main
  input = $stdin.read.split
  ptr = 0
  n = input[ptr].to_i
  ptr += 1
  a = Array.new(n) { |i| input[ptr + i].to_i }
  
  # Sieve of smallest prime factors
  spf = Array.new(MAXA + 1, 0)
  (2..MAXA).each do |i|
    if spf[i] == 0
      spf[i] = i
      (i * 2).step(MAXA, i) do |j|
        spf[j] = i if spf[j] == 0
      end
    end
  end

  dp_len = Array.new(MAXA + 1, 0)
  dp_cnt = Array.new(MAXA + 1, 0)

  max_len = 0
  total_count = 0

  a.each do |x|
    # Factorize x into distinct primes
    t = x
    primes = Set.new
    while t > 1
      p = spf[t]
      primes.add(p)
      while t % p == 0
        t /= p
      end
    end
    primes = primes.to_a

    # Find the best length among the primes
    best = 0
    primes.each { |p| best = [best, dp_len[p]].max }

    count = 0
    if best == 0
      count = 1
    else
      # Collect primes with dp_len[p] == best
      h = primes.select { |p| dp_len[p] == best }
      k = h.size
      # Inclusion-exclusion over subsets of h
      total = 0
      (1..k).each do |m|
        h.combination(m).each do |combo|
          d = combo.inject(1, :*)
          if dp_len[d] == best
            if m.odd?
              total = (total + dp_cnt[d]) % MOD
            else
              total = (total - dp_cnt[d]) % MOD
            end
          end
        end
      end
      count = total % MOD
    end

    current_len = best + 1
    if current_len > max_len
      max_len = current_len
      total_count = count
    elsif current_len == max_len
      total_count = (total_count + count) % MOD
    end

    # Generate all non-empty subsets of primes
    m = primes.size
    (1..m).each do |s|
      primes.combination(s).each do |combo|
        d = combo.inject(1, :*)
        if dp_len[d] < current_len
          dp_len[d] = current_len
          dp_cnt[d] = count % MOD
        elsif dp_len[d] == current_len
          dp_cnt[d] = (dp_cnt[d] + count) % MOD
        end
      end
    end
  end

  puts max_len
  puts total_count % MOD
end

main