Skip to content

Commit

Permalink
Add the streamline-demo example model
Browse files Browse the repository at this point in the history
Add an executable example model used to test and demonstrate the framework.

The example model comprises three main components:
 * ErrorTestingModel
 * DataModel
 * ApproximationModel

The ErrorTestingModel demonstrates error-handling behavior for resources.
Using the `CauseError` activity, we can trigger one of several kinds of errors from the plan,
including effects which fail directly, concurrent effects which conflict, and resources which violate a constraint.
We can also test out different combinations of naming and error handling behaviors, to see what information is propagated back to the user for debugging.
Finally, we can observe that, when logging errors, portions of the model can fail while independent portions of the model continue to function normally.

The DataModel demonstrates a complicated resource modeling problem, leveraging the LinearBoundaryConsistencySolver.
This problem models a data system with bin space shared across multiple buckets, where some buckets have priority over others,
and each bucket independently sets desired write or delete rates.
The problem of allocating bin space to each bucket, allowing higher-priority buckets to overwrite lower-priority ones only if needed,
is phrased as linear inequalities on polynomial resources and given to the solver.
The `ChangeDesiredRate` activity can be used to set up different scenarios for this model.

The ApproximationModel shows three resources, a polynomial, a rational function, and a complicated trig function, approximated various ways.
The `ChangeApproximationInput` activity can be used to change the polynomials feeding this model,
and the simulation configuration contains a setting for the approximation accuracy, affecting those resources which approximate by bounding their maximum error.
In particular, the `approximation/trig/default` resource displays the most complicated resource function, comprising trig and exponentials,
defined through the Java Math library, approximated by the default secant approximation method.
  • Loading branch information
David Legg committed Dec 19, 2023
1 parent 8f98784 commit 2963c74
Show file tree
Hide file tree
Showing 14 changed files with 472 additions and 0 deletions.
42 changes: 42 additions & 0 deletions examples/streamline-demo/build.gradle
Original file line number Diff line number Diff line change
@@ -0,0 +1,42 @@
plugins {
id 'java-library'
id 'jacoco'
}

java {
toolchain {
languageVersion = JavaLanguageVersion.of(19)
}
}

jar {
from {
configurations.runtimeClasspath.filter{ it.exists() }.collect{ it.isDirectory() ? it : zipTree(it) }
} {
exclude 'META-INF/LICENSE.txt', 'META-INF/NOTICE.txt'
}
}

test {
useJUnitPlatform()
}

jacocoTestReport {
dependsOn test
reports {
xml.required = true
}
}

dependencies {
annotationProcessor project(':merlin-framework-processor')

implementation project(':merlin-framework')
implementation project(':contrib')

testImplementation project(':merlin-framework-junit')
testImplementation 'org.assertj:assertj-core:3.23.1'

testImplementation 'org.junit.jupiter:junit-jupiter-engine:5.8.2'
}

Original file line number Diff line number Diff line change
@@ -0,0 +1,102 @@
package gov.nasa.jpl.aerie.streamline_demo;

import gov.nasa.jpl.aerie.contrib.streamline.core.MutableResource;
import gov.nasa.jpl.aerie.contrib.streamline.core.Resource;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.Registrar;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.*;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.linear.Linear;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.Polynomial;

import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.Approximation.approximate;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.Approximation.relative;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.DifferentiableResources.asDifferentiable;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.IntervalFunctions.byBoundingError;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.IntervalFunctions.byUniformSampling;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.SecantApproximation.ErrorEstimates.errorByOptimization;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.SecantApproximation.ErrorEstimates.errorByQuadraticApproximation;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.SecantApproximation.secantApproximation;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.UnstructuredResources.approximateAsLinear;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.UnstructuredResources.approximateAsLinear;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.UnstructuredResources.asUnstructured;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.black_box.monads.UnstructuredResourceApplicative.map;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.PolynomialResources.*;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.PolynomialResources.approximateAsLinear;
import static gov.nasa.jpl.aerie.merlin.protocol.types.Duration.*;
import static gov.nasa.jpl.aerie.merlin.protocol.types.Duration.MINUTE;

public class ApproximationModel {
private static final double EPSILON = 1e-10;

public MutableResource<Polynomial> polynomial;
public MutableResource<Polynomial> divisor;

public Resource<Linear> assumedLinear;
public Resource<Linear> uniformApproximation;
public Resource<Linear> differentiableApproximation;
public Resource<Linear> directApproximation;
public Resource<Linear> defaultApproximation;

public Resource<Unstructured<Double>> rationalFunction;
public Resource<Linear> uniformApproximation2;
public Resource<Linear> directApproximation2;
public Resource<Linear> defaultApproximation2;

public Resource<Unstructured<Double>> trigFunction;
public Resource<Linear> trigDefaultApproximation;

public ApproximationModel(final Registrar registrar, final Configuration config) {
final double tolerance = config.approximationTolerance;

polynomial = polynomialResource(1);
divisor = polynomialResource(1);

assumedLinear = assumeLinear(polynomial);
defaultApproximation = approximateAsLinear(polynomial, tolerance);
uniformApproximation = approximate(
polynomial,
SecantApproximation.<Polynomial>secantApproximation(byUniformSampling(MINUTE)));
differentiableApproximation = approximate(
asDifferentiable(polynomial),
secantApproximation(byBoundingError(
tolerance,
SECOND,
HOUR.times(24),
relative(errorByQuadraticApproximation(), EPSILON))));
directApproximation = approximate(
polynomial,
secantApproximation(byBoundingError(
tolerance,
SECOND,
HOUR.times(24),
Approximation.<Polynomial>relative(errorByOptimization(), EPSILON))));

rationalFunction = map(
asUnstructured(polynomial), asUnstructured(divisor), (p, q) -> p / q);

defaultApproximation2 = approximateAsLinear(rationalFunction, tolerance);
uniformApproximation2 = approximate(rationalFunction,
SecantApproximation.<Unstructured<Double>>secantApproximation(byUniformSampling(MINUTE)));
directApproximation2 = approximate(rationalFunction,
secantApproximation(byBoundingError(
tolerance,
SECOND,
HOUR.times(24),
Approximation.<Unstructured<Double>>relative(errorByOptimization(), EPSILON))));

trigFunction = map(asUnstructured(polynomial), asUnstructured(divisor),
(p, q) -> Math.sin(p * Math.exp(-q / Math.PI)));
trigDefaultApproximation = approximateAsLinear(trigFunction, tolerance);

registrar.real("approximation/polynomial/assumedLinear", assumedLinear);
registrar.real("approximation/polynomial/default", defaultApproximation);
registrar.real("approximation/polynomial/uniform", uniformApproximation);
registrar.real("approximation/polynomial/differentiable", differentiableApproximation);
registrar.real("approximation/polynomial/direct", directApproximation);

registrar.real("approximation/rational/default", defaultApproximation2);
registrar.real("approximation/rational/uniform", uniformApproximation2);
registrar.real("approximation/rational/direct", directApproximation2);

registrar.real("approximation/trig/default", trigDefaultApproximation);
}
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,73 @@
package gov.nasa.jpl.aerie.streamline_demo;

import gov.nasa.jpl.aerie.contrib.streamline.core.MutableResource;
import gov.nasa.jpl.aerie.contrib.streamline.core.Dynamics;
import gov.nasa.jpl.aerie.contrib.streamline.core.DynamicsEffect;
import gov.nasa.jpl.aerie.merlin.framework.annotations.ActivityType;
import gov.nasa.jpl.aerie.merlin.framework.annotations.ActivityType.EffectModel;
import gov.nasa.jpl.aerie.merlin.framework.annotations.Export.Parameter;

import static gov.nasa.jpl.aerie.contrib.streamline.core.monads.DynamicsMonad.effect;
import static gov.nasa.jpl.aerie.contrib.streamline.debugging.Context.contextualized;
import static gov.nasa.jpl.aerie.contrib.streamline.debugging.Context.inContext;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.discrete.DiscreteEffects.set;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.Polynomial.polynomial;
import static gov.nasa.jpl.aerie.merlin.framework.ModelActions.delay;
import static gov.nasa.jpl.aerie.merlin.framework.ModelActions.spawn;
import static gov.nasa.jpl.aerie.merlin.protocol.types.Duration.HOUR;

@ActivityType("CauseError")
public class CauseError {
@Parameter
public ResourceSelection selection;

@Parameter
public String effectName = "";

@Parameter
public String activityName = "";

@EffectModel
public void run(Mission mission) {
inContext("CauseError" + (activityName.isEmpty() ? "" : " " + activityName), () -> {
delay(HOUR);
switch (selection) {
case Bool -> causeError(mission.errorTestingModel.bool);
case Counter -> causeError(mission.errorTestingModel.counter);
case Continuous -> causeError(mission.errorTestingModel.continuous);
case NonCommuting -> {
spawn(contextualized("Branch 1", () -> {
set(mission.errorTestingModel.counter, 5);
}));
spawn(contextualized("Branch 2", () -> {
set(mission.errorTestingModel.counter, 6);
}));
}
case DoomedClamp -> {
MutableResource.set(mission.errorTestingModel.lowerBound, polynomial(-20, 0.001));
MutableResource.set(mission.errorTestingModel.upperBound, polynomial(20, -0.001));
}
}
delay(HOUR);
});
}

private <D extends Dynamics<?, D>> void causeError(MutableResource<D> resource) {
DynamicsEffect<D> effect = effect($ -> {
throw new IllegalStateException("Pretend this is a more informative error message.");
});
if (effectName.isEmpty()) {
resource.emit(effect);
} else {
resource.emit(effectName, effect);
}
}

public enum ResourceSelection {
Bool,
Counter,
Continuous,
NonCommuting,
DoomedClamp
}
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,22 @@
package gov.nasa.jpl.aerie.streamline_demo;

import gov.nasa.jpl.aerie.merlin.framework.annotations.ActivityType;
import gov.nasa.jpl.aerie.merlin.framework.annotations.Export;

import static gov.nasa.jpl.aerie.contrib.streamline.core.MutableResource.set;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.Polynomial.polynomial;

@ActivityType("ChangeApproximationInput")
public class ChangeApproximationInput {
@Export.Parameter
public double[] numeratorCoefficients;

@Export.Parameter
public double[] denominatorCoefficients;

@ActivityType.EffectModel
public void run(Mission mission) {
set(mission.approximationModel.polynomial, polynomial(numeratorCoefficients));
set(mission.approximationModel.divisor, polynomial(denominatorCoefficients));
}
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,29 @@
package gov.nasa.jpl.aerie.streamline_demo;

import gov.nasa.jpl.aerie.merlin.framework.annotations.ActivityType;
import gov.nasa.jpl.aerie.merlin.framework.annotations.ActivityType.EffectModel;
import gov.nasa.jpl.aerie.merlin.framework.annotations.Export.Parameter;

import static gov.nasa.jpl.aerie.contrib.streamline.core.MutableResource.set;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.Polynomial.polynomial;

@ActivityType("ChangeDesiredRate")
public class ChangeDesiredRate {
@Parameter
public Bucket bucket;

@Parameter
public double rate;

@EffectModel
public void run(Mission mission) {
var rateToChange = switch (bucket) {
case A -> mission.dataModel.desiredRateA;
case B -> mission.dataModel.desiredRateB;
case C -> mission.dataModel.desiredRateC;
};
set(rateToChange, polynomial(rate));
}

public enum Bucket { A, B, C }
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
package gov.nasa.jpl.aerie.streamline_demo;

import gov.nasa.jpl.aerie.merlin.framework.annotations.Export.Parameter;
import gov.nasa.jpl.aerie.merlin.protocol.types.Duration;

public final class Configuration {
@Parameter
public boolean traceResources = false;

@Parameter
public double approximationTolerance = 1e-2;

@Parameter
public Duration profilingDumpTime = Duration.ZERO;

}
Original file line number Diff line number Diff line change
@@ -0,0 +1,106 @@
package gov.nasa.jpl.aerie.streamline_demo;

import gov.nasa.jpl.aerie.contrib.streamline.core.MutableResource;
import gov.nasa.jpl.aerie.contrib.streamline.core.Resource;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.Registrar;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.LinearBoundaryConsistencySolver;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.LinearBoundaryConsistencySolver.Domain;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.Polynomial;
import gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.PolynomialResources;

import static gov.nasa.jpl.aerie.contrib.streamline.core.Resources.eraseExpiry;
import static gov.nasa.jpl.aerie.contrib.streamline.core.Resources.forward;
import static gov.nasa.jpl.aerie.contrib.streamline.core.monads.ResourceMonad.*;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.discrete.DiscreteResources.choose;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.LinearBoundaryConsistencySolver.Comparison.*;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.LinearBoundaryConsistencySolver.LinearExpression.lx;
import static gov.nasa.jpl.aerie.contrib.streamline.modeling.polynomial.PolynomialResources.*;

public class DataModel {
public MutableResource<Polynomial> desiredRateA = PolynomialResources.polynomialResource(0);
public MutableResource<Polynomial> desiredRateB = PolynomialResources.polynomialResource(0);
public MutableResource<Polynomial> desiredRateC = PolynomialResources.polynomialResource(0);
public MutableResource<Polynomial> upperBoundOnTotalVolume = PolynomialResources.polynomialResource(10);

public Resource<Polynomial> actualRateA, actualRateB, actualRateC;
public MutableResource<Polynomial> volumeA, volumeB, volumeC;
public Resource<Polynomial> totalVolume;

public DataModel(final Registrar registrar, final Configuration config) {
// Adding a "derivative" operation makes the constraint problem no longer solvable without backtracking.
// Instead of solving for it directly, we'll solve in two steps, one for rate, one for volume.

// Set up the rate solver
var rateSolver = new LinearBoundaryConsistencySolver("DataModel Rate Solver");
var rateA = rateSolver.variable("rateA", Domain::upperBound);
var rateB = rateSolver.variable("rateB", Domain::upperBound);
var rateC = rateSolver.variable("rateC", Domain::upperBound);
this.actualRateA = rateA.resource();
this.actualRateB = rateB.resource();
this.actualRateC = rateC.resource();

// Use a simple feedback loop on volumes to do the integration and clamping.
this.volumeA = PolynomialResources.polynomialResource(0);
this.volumeB = PolynomialResources.polynomialResource(0);
this.volumeC = PolynomialResources.polynomialResource(0);
var clampedVolumeA = clamp(this.volumeA, constant(0), upperBoundOnTotalVolume);
var volumeB_ub = subtract(upperBoundOnTotalVolume, clampedVolumeA);
var clampedVolumeB = clamp(this.volumeB, constant(0), volumeB_ub);
var volumeC_ub = subtract(volumeB_ub, clampedVolumeB);
var clampedVolumeC = clamp(this.volumeC, constant(0), volumeC_ub);
var correctedVolumeA = map(clampedVolumeA, actualRateA, (v, r) -> r.integral(v.extract()));
var correctedVolumeB = map(clampedVolumeB, actualRateB, (v, r) -> r.integral(v.extract()));
var correctedVolumeC = map(clampedVolumeC, actualRateC, (v, r) -> r.integral(v.extract()));
// Use the corrected integral values to set volumes, but erase expiry information in the process to avoid loops:
forward(eraseExpiry(correctedVolumeA), this.volumeA);
forward(eraseExpiry(correctedVolumeB), this.volumeB);
forward(eraseExpiry(correctedVolumeC), this.volumeC);

// Integrate the actual rates.
totalVolume = add(this.volumeA, this.volumeB, this.volumeC);

// Then use the solver to adjust the actual rates

// When full, we never write more than the upper bound will tolerate, in total
var isFull = greaterThanOrEquals(totalVolume, upperBoundOnTotalVolume);
var totalRate_ub = choose(isFull, differentiate(upperBoundOnTotalVolume), constant(Double.POSITIVE_INFINITY));
rateSolver.declare(lx(rateA).add(lx(rateB)).add(lx(rateC)), LessThanOrEquals, lx(totalRate_ub));

// We only exceed the desired rate when we try to delete from an empty bucket.
var isEmptyA = lessThanOrEquals(this.volumeA, 0);
var isEmptyB = lessThanOrEquals(this.volumeB, 0);
var isEmptyC = lessThanOrEquals(this.volumeC, 0);
var rateA_ub = choose(isEmptyA, max(desiredRateA, constant(0)), desiredRateA);
var rateB_ub = choose(isEmptyB, max(desiredRateB, constant(0)), desiredRateB);
var rateC_ub = choose(isEmptyC, max(desiredRateC, constant(0)), desiredRateC);
rateSolver.declare(lx(rateA), LessThanOrEquals, lx(rateA_ub));
rateSolver.declare(lx(rateB), LessThanOrEquals, lx(rateB_ub));
rateSolver.declare(lx(rateC), LessThanOrEquals, lx(rateC_ub));

// We cannot delete from an empty bucket
var rateA_lb = choose(isEmptyA, constant(0), constant(Double.NEGATIVE_INFINITY));
var rateB_lb = choose(isEmptyB, constant(0), constant(Double.NEGATIVE_INFINITY));
var rateC_lb = choose(isEmptyC, constant(0), constant(Double.NEGATIVE_INFINITY));
rateSolver.declare(lx(rateA), GreaterThanOrEquals, lx(rateA_lb));
rateSolver.declare(lx(rateB), GreaterThanOrEquals, lx(rateB_lb));
rateSolver.declare(lx(rateC), GreaterThanOrEquals, lx(rateC_lb));

registerStates(registrar, config);
}

private void registerStates(Registrar registrar, Configuration config) {
registrar.real("desiredRateA", assumeLinear(desiredRateA));
registrar.real("desiredRateB", assumeLinear(desiredRateB));
registrar.real("desiredRateC", assumeLinear(desiredRateC));

registrar.real("actualRateA", assumeLinear(actualRateA));
registrar.real("actualRateB", assumeLinear(actualRateB));
registrar.real("actualRateC", assumeLinear(actualRateC));

registrar.real("volumeA", assumeLinear(volumeA));
registrar.real("volumeB", assumeLinear(volumeB));
registrar.real("volumeC", assumeLinear(volumeC));
registrar.real("totalVolume", assumeLinear(totalVolume));
registrar.real("maxVolume", assumeLinear(upperBoundOnTotalVolume));
}
}
Loading

0 comments on commit 2963c74

Please sign in to comment.