refactor genetic/errors.go out to mu8 root; errors overhaul; modify G…

…enomeGrad interface

README.md

@@ -74,8 +74,9 @@ pop := genetic.NewPopulation(individuals, rand.NewSource(1), func() *mygenome {
 })
 
 const Ngeneration = 100
+ctx := context.Background()
 for i := 0; i < Ngenerations; i++ {
-		err := pop.Advance()
+		err := pop.Advance(ctx)
 		if err != nil {
 			panic(err.Error())
 		}
@@ -88,6 +89,7 @@ fmt.Printf("champ fitness=%.3f\n", pop.ChampionFitness())
 ```
 The final fitness should be close to 1.0 if the algorithm did it's job. For the code see 
 [`mu8_test.go`](./mu8_test.go)
+
 ### Rocket stage optimization example
 
 See [`rocket`](./examples/rocket/main.go) for a demonstration on rocket stage optimization. 
@@ -109,6 +111,44 @@ Stage 0: coast=281.2s, propMass=0.0kg, Δm=99.35kg/s, totalMass=200.0
 Stage 1: coast=0.0s, propMass=1.6kg, Δm=0.01kg/s, totalMass=21.6
 ```
 
+### Gradient "ascent" example
+```go
+src := rand.NewSource(1)
+const (
+	genomelen      = 6
+	gradMultiplier = 10.0
+	epochs         = 6
+)
+// Create new individual and mutate it randomly.
+individual := newGenome(genomelen)
+rng := rand.New(src)
+for i := 0; i < genomelen; i++ {
+	individual.GetGene(i).Mutate(rng)
+}
+// Prepare for gradient descent.
+grads := make([]float64, genomelen)
+ctx := context.Background()
+// Champion will harbor our best individual.
+champion := newGenome(genomelen)
+for epoch := 0; epoch < epochs; epoch++ {
+	// We calculate the gradients of the individual passing a nil
+	// newIndividual callback since the GenomeGrad type we implemented
+	// does not require blank-slate initialization.
+	err := mu8.Gradient(ctx, grads, individual, nil)
+	if err != nil {
+		panic(err)
+	}
+	// Apply gradients.
+	for i := 0; i < individual.Len(); i++ {
+		gene := individual.GetGeneGrad(i)
+		grad := grads[i]
+		gene.SetValue(gene.Value() + grad*gradMultiplier)
+	}
+	mu8.CloneGrad(champion, individual)
+	fmt.Printf("fitness=%f with grads=%f\n", individual.Simulate(ctx), grads)
+}
+```
+
 ## Contributing
 Contributions very welcome! I myself have no idea what I'm doing so I welcome
 issues on any matter :)

errors.go

@@ -0,0 +1,68 @@
+package mu8
+
+import (
+	"context"
+	"errors"
+	"fmt"
+	"math"
+	"math/rand"
+)
+
+var (
+	ErrNegativeFitness = errors.New("negative fitness")
+	ErrInvalidFitness  = errors.New("got infinite or NaN fitness")
+	ErrCodependency    = errors.New("codependency between individuals")
+)
+
+// FindCodependecy returns error if inconsistency detected in newIndividual function
+// for use with mu8.Genome genetic algorithm implementations.
+//
+// Users can check if a codependency is found by checking the error:
+//
+//	if errors.Is(err, ErrCodependency) {
+//		// Handle codependency case.
+//	}
+//
+// The error will have a descriptive text of how the Genome is codependent.
+func FindCodependecy[G Genome](src rand.Source, newIndividual func() G) error {
+	ctx := context.Background()
+	starter1 := newIndividual()
+	starter2 := newIndividual()
+	// These two fitnesses should be equal if there is not codependency.
+	fit1 := starter1.Simulate(ctx)
+	fit2 := starter2.Simulate(ctx)
+	switch {
+	case fit1 != fit2:
+		return fmt.Errorf("%w: during subsequent calls to newIndividual which should return identical fitnesses. check for closure variable capture modification or preserved slice reference?", ErrCodependency)
+	case fit1 == 0:
+		return errors.New("cannot reliably determine codependency with zero fitness simulation results")
+	case fit1 < 0 || fit2 < 0:
+		return ErrNegativeFitness
+	case math.IsNaN(fit1) || math.IsNaN(fit2) || math.IsInf(fit1, 0) || math.IsInf(fit2, 0):
+		return ErrInvalidFitness
+	}
+
+	rng := rand.New(src)
+	var codependents []int
+	for i := 0; i < starter1.Len(); i++ {
+		parent1 := newIndividual()
+		parent2 := newIndividual()
+		fit1 := parent1.Simulate(ctx)
+		g := parent1.GetGene(i)
+		// This line should have no effect on parent2's simulation (should be "initial" fitness).
+		g.Mutate(rng)
+		fit2 := parent2.Simulate(ctx)
+		if math.IsNaN(fit1) || math.IsNaN(fit2) || math.IsInf(fit1, 0) || math.IsInf(fit2, 0) {
+			if codependents != nil {
+				return fmt.Errorf("invalid fitness and %w: detected in genes indices: %v", ErrCodependency, codependents)
+			}
+			return ErrInvalidFitness
+		} else if fit1 != fit2 {
+			codependents = append(codependents, i)
+		}
+	}
+	if codependents == nil {
+		return nil // No codependency detected.
+	}
+	return fmt.Errorf("%w: genes indices: %v", ErrCodependency, codependents)
+}

genes/constrainedfloat.go

@@ -10,6 +10,7 @@ import (
 
 // NewConstrainedFloat returns a mu8.Gene implementation for a number
 // that should be kept within bounds [min,max] during mutation.
+// start is the initial value of the gene.
 func NewConstrainedFloat(start, min, max float64) *ConstrainedFloat {
 	if min > max {
 		panic(errBadConstraints)
@@ -42,23 +43,30 @@ func (c *ConstrainedFloat) SetValue(f float64) {
 	c.gene = c.clamp(f)
 }
 
+// Mutate changes the gene's value by a random amount within constraints.
+// Mutate implements the [mu8.Gene] interface.
 func (c *ConstrainedFloat) Mutate(rng *rand.Rand) {
 	// Uniform mutation distribution.
 	random := rng.Float64()
 	random = c.min + random*c.rangeLength()
 	c.gene = c.clamp(random)
 }
 
+// CloneFrom copies the argument gene into the receiver. CloneFrom implements
+// the [mu8.Gene] interface. If g is not of type *ConstrainedFloat, CloneFrom panics.
 func (c *ConstrainedFloat) CloneFrom(g mu8.Gene) {
 	co := castGene[*ConstrainedFloat](g)
 	c.gene = co.gene
 }
 
+// Copy returns a new ConstrainedFloat with the same value as the receiver.
 func (c *ConstrainedFloat) Copy() *ConstrainedFloat {
 	clone := *c
 	return &clone
 }
 
+// Splice performs a crossover operation between the receiver and g.
+// Splice implements the [mu8.Gene] interface. If g is not of type *ConstrainedFloat, Splice panics.
 func (c *ConstrainedFloat) Splice(rng *rand.Rand, g mu8.Gene) {
 	co := castGene[*ConstrainedFloat](g)
 	random := rng.Float64()
@@ -96,6 +104,7 @@ func (c *ConstrainedFloat) rangeLength() float64 {
 	return c.maxMinus1 + 1 - c.min
 }
 
+// String implements fmt.Stringer interface.
 func (c *ConstrainedFloat) String() string {
 	return fmt.Sprintf("%f", c.gene)
 }

genes/constrainedint.go

@@ -47,21 +47,29 @@ func (c *ConstrainedInt) SetValue(f int) {
 	c.gene = f
 }
 
+// Mutate changes the gene's value by a random amount within constraints.
+// Mutate implements the [mu8.Gene] interface.
 func (c *ConstrainedInt) Mutate(rng *rand.Rand) {
 	// Uniform mutation distribution.
 	c.gene = c.min + rng.Intn(c.rangeMinus1+1)
 }
 
+// CloneFrom copies the argument gene into the receiver. CloneFrom implements
+// the [mu8.Gene] interface. If g is not of type *ConstrainedInt, CloneFrom panics.
 func (c *ConstrainedInt) CloneFrom(g mu8.Gene) {
 	co := castGene[*ConstrainedInt](g)
 	c.gene = co.gene
 }
 
+// Copy returns a copy of the gene.
 func (c *ConstrainedInt) Copy() *ConstrainedInt {
 	clone := *c
 	return &clone
 }
 
+// Splice performs a crossover between the argument and the receiver genes
+// and stores the result in the receiver. It implements the [mu8.Gene] interface.
+// If g is not of type *ConstrainedInt, Splice panics.
 func (c *ConstrainedInt) Splice(rng *rand.Rand, g mu8.Gene) {
 	co := castGene[*ConstrainedInt](g)
 	diff := c.gene - co.gene
@@ -77,6 +85,7 @@ func (c *ConstrainedInt) Splice(rng *rand.Rand, g mu8.Gene) {
 	c.gene = minGene + random
 }
 
+// String returns a string representation of the gene.
 func (c *ConstrainedInt) String() string {
 	return fmt.Sprintf("%d", c.gene)
 }

genes/constrainednormal.go

@@ -7,11 +7,19 @@ import (
 	"github.com/soypat/mu8"
 )
 
+// ConstrainedFloat is a float64 gene that is constrained to a range and whose
+// mutations are drawn from a uniform distribution.
 type ConstrainedNormalDistr struct {
 	NormalDistribution
 	minPlus3sd, maxMinus3sd float64
 }
 
+// NewConstrainedNormalDistr returns a new ConstrainedNormalDistr.
+// The arguments are:
+//   - val: the initial value of the gene.
+//   - stdDeviation: the standard deviation of mutations.
+//   - min: the minimum value the gene can take.
+//   - max: the maximum value the gene can take.
 func NewConstrainedNormalDistr(val, stdDeviation, min, max float64) *ConstrainedNormalDistr {
 	if min > max {
 		panic(errBadConstraints)
@@ -23,34 +31,44 @@ func NewConstrainedNormalDistr(val, stdDeviation, min, max float64) *Constrained
 	}
 }
 
+// Mutate mutates the gene according to a random normal distribution.
+// It implements the [mu8.Gene] interface.
 func (cn *ConstrainedNormalDistr) Mutate(rng *rand.Rand) {
 	cn.NormalDistribution.Mutate(rng)
 	cn.clamp()
 }
 
+// CloneFrom copies the argument gene into the receiver. CloneFrom implements the
+// [mu8.Gene] interface. If g is not of type *ConstrainedNormalDistr, CloneFrom panics.
 func (cn *ConstrainedNormalDistr) CloneFrom(g mu8.Gene) {
 	co := castGene[*ConstrainedNormalDistr](g)
 	cn.gene = co.gene
 }
 
+// Splice performs a crossover between the argument and the receiver genes
+// and stores the result in the receiver. It implements the [mu8.Gene] interface.
+// If g is not of type *ConstrainedNormalDistr, Splice panics.
 func (cn *ConstrainedNormalDistr) Splice(rng *rand.Rand, g mu8.Gene) {
 	co := castGene[*ConstrainedNormalDistr](g)
 	cn.NormalDistribution.Splice(rng, &co.NormalDistribution)
 	cn.clamp()
 }
 
+// Copy returns a copy of the gene.
 func (cn *ConstrainedNormalDistr) Copy() *ConstrainedNormalDistr {
 	clone := *cn
 	return &clone
 }
 
+// clamp clamps the gene value to the constraints.
 func (cn *ConstrainedNormalDistr) clamp() {
 	sd3 := cn.StdDev() * 3
 	min := cn.minPlus3sd - sd3
 	max := cn.maxMinus3sd + sd3
 	cn.gene = math.Max(min, math.Min(max, cn.gene))
 }
 
+// SetValue sets the value of the gene.
 func (cn *ConstrainedNormalDistr) SetValue(f float64) {
 	cn.gene = f
 	cn.clamp()

genes/genes.go

@@ -7,6 +7,7 @@ import (
 	"github.com/soypat/mu8"
 )
 
+// Gene errors.
 var (
 	ErrMismatchedGeneType = errors.New("mu8.Gene argument in Splice or CloneFrom not same type as receiver")
 	errStartOutOfBounds   = errors.New("start value should be contained within bounds [min,max] for Contrained types")

genes/grad.go

@@ -4,34 +4,60 @@ import "github.com/soypat/mu8"
 
 const defaultStep = 5e-7
 
+// compile time check that these types implement the GeneGrad interface.
 var (
 	_ mu8.GeneGrad = (*ConstrainedFloatGrad)(nil)
+	_ mu8.GeneGrad = (*ConstrainedNormalDistrGrad)(nil)
 )
 
+// ConstrainedFloatGrad is a ConstrainedFloat that implements the GeneGrad interface
+// with a programmable step size.
+// It implements the [mu8.GeneGrad] interface.
 type ConstrainedFloatGrad struct {
 	ConstrainedFloat
 	stepMinusDefaultStep float64
 }
 
+// Step returns the step size that should be performed during gradient descent.
+// It implements the [mu8.GeneGrad] interface.
 func (cf *ConstrainedFloatGrad) Step() float64 {
 	return cf.stepMinusDefaultStep + defaultStep
 }
 
+// NewConstrainedFloatGrad returns a new ConstrainedFloatGrad.
+// The arguments are:
+//   - start: the initial value of the gene.
+//   - min: the minimum value the gene can take.
+//   - max: the maximum value the gene can take.
+//   - step: the step size that should be used gradient descent.
 func NewConstrainedFloatGrad(start, min, max, step float64) *ConstrainedFloatGrad {
 	return &ConstrainedFloatGrad{
 		ConstrainedFloat:     *NewConstrainedFloat(start, min, max),
 		stepMinusDefaultStep: step - defaultStep,
 	}
 }
 
+// ConstrainedNormalDistrGrad is a ConstrainedNormalDistr that implements the GeneGrad interface
+// with a programmable step size.
+// It implements the [mu8.GeneGrad] interface.
 type ConstrainedNormalDistrGrad struct {
 	ConstrainedNormalDistr
 	stepMinusDefaultStep float64
 }
 
+// Step returns the step size that should be performed during gradient descent.
+// It implements the [mu8.GeneGrad] interface.
 func (cf *ConstrainedNormalDistrGrad) Step() float64 {
 	return cf.stepMinusDefaultStep + defaultStep
 }
+
+// NewConstrainedNormalGrad returns a new ConstrainedNormalDistrGrad
+// where arguments are:
+//   - start: the initial value of the gene.
+//   - stddev: the standard deviation of mutations.
+//   - min: the minimum value the gene can take.
+//   - max: the maximum value the gene can take.
+//   - step: the step size that should be used during gradient descent.
 func NewConstrainedNormalGrad(start, stddev, min, max, step float64) *ConstrainedNormalDistrGrad {
 	return &ConstrainedNormalDistrGrad{
 		ConstrainedNormalDistr: *NewConstrainedNormalDistr(start, stddev, min, max),