diff --git a/cmd/radical.go b/cmd/radical.go new file mode 100644 index 0000000..e6621e8 --- /dev/null +++ b/cmd/radical.go @@ -0,0 +1,75 @@ +/* +Copyright © 2025 filifa + +This program is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 3 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with this program. If not, see . +*/ +package cmd + +import ( + "bufio" + "fmt" + "os" + + "github.com/spf13/cobra" + "scm.dairydemon.net/filifa/mathtools/internal/lib" +) + +var radicalN uint + +func radical(cmd *cobra.Command, args []string) { + bufStdout := bufio.NewWriter(os.Stdout) + defer bufStdout.Flush() + + ch := lib.RadicalSieve(radicalN, 1000) + for i := 0; ; i++ { + v, ok := <-ch + if !ok { + break + } + + if i == 0 { + continue + } + + fmt.Fprintln(bufStdout, v) + } +} + +// radicalCmd represents the radical command +var radicalCmd = &cobra.Command{ + Use: "radical -n N", + Short: "Compute the radical of all numbers less than n", + Long: `Compute the radical of all numbers less than n. + +The radical of an integer is the product of all the distinct prime factors of that integer. +For example, 12 = 2^2 * 3, so rad(12) = 2 * 3 = 6.`, + Run: radical, +} + +func init() { + sieveCmd.AddCommand(radicalCmd) + + // Here you will define your flags and configuration settings. + + // Cobra supports Persistent Flags which will work for this command + // and all subcommands, e.g.: + // radicalCmd.PersistentFlags().String("foo", "", "A help for foo") + + // Cobra supports local flags which will only run when this command + // is called directly, e.g.: + // radicalCmd.Flags().BoolP("toggle", "t", false, "Help message for toggle") + + radicalCmd.Flags().UintVarP(&radicalN, "limit", "n", 0, "upper limit") + radicalCmd.MarkFlagRequired("limit") +} diff --git a/internal/lib/radical.go b/internal/lib/radical.go new file mode 100644 index 0000000..1be777f --- /dev/null +++ b/internal/lib/radical.go @@ -0,0 +1,63 @@ +/* +Copyright © 2025 filifa + +This program is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 3 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with this program. If not, see . +*/ +package lib + +func radicalUpdateMultiples(sieve []uint, p uint, n uint) { + for q := p; ; q *= p { + // rad(a*b) = rad(a) * rad(b) if gcd(a,b) = 1 + for i := 2 * q; i < n; i += q { + if i%(p*q) != 0 { + sieve[i] *= sieve[q] + } + } + + if p*q >= n { + break + } + + // rad(p^k) = rad(p) + sieve[p*q] = sieve[q] + } +} + +/* +RadicalSieve computes rad(k) for k=1 to n, where rad(n) is the radical of n. +*/ +func RadicalSieve(n uint, buflen uint) chan uint { + sieve := make([]uint, n) + sieve[0] = 0 + for i := uint(1); i < n; i++ { + sieve[i] = 1 + } + + ch := make(chan uint, buflen) + go func() { + defer close(ch) + for i := uint(0); i < n; i++ { + if i == 0 || i == 1 || sieve[i] != 1 { + ch <- sieve[i] + continue + } + + sieve[i] = i + radicalUpdateMultiples(sieve, i, n) + ch <- sieve[i] + } + }() + + return ch +}