NunoSempere/fermi
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Compare commits

base: NunoSempere:9dc6e1443d98c0e7a32277c7cd351872ecfa0c49
Branches Tags
NunoSempere:master
..
compare: NunoSempere:4c6fb35ac2656d74b43bc8f0c337e2185691e898
Branches Tags
NunoSempere:master

No commits in common. "9dc6e1443d98c0e7a32277c7cd351872ecfa0c49" and "4c6fb35ac2656d74b43bc8f0c337e2185691e898" have entirely different histories.
9dc6e1443d ... 4c6fb35ac2

6 changed files with 172 additions and 716 deletions
Show Stats Expand all files Collapse all files

68
README.md
View File

@ -1,6 +1,6 @@
# A minimalist calculator for f estimation
# A minimalist calculator for fermi estimation
This project is a minimalist, stack-based DSL for f estimation. It can multiply and divide scalars, lognormals and beta distributions.
This project is a minimalist, stack-based DSL for Fermi estimation. It can multiply and divide scalars, lognormals and beta distributions.
## Motivation
@ -11,7 +11,7 @@ Sometimes, [Squiggle](https://github.com/quantified-uncertainty/squiggle), [simp
Here is an example
```
$ go run f.go
$ go run fermi.go
5000000 12000000
=> 5.0M 12.0M
* beta 1 200
@ -31,7 +31,7 @@ $ go run f.go
Perhaps this example is more understandable with comments and better units:
```
$ sed -u "s|#.*||" | sed -u 's|M|000000|g' | go run f.go
$ sed -u "s|#.*||" | sed -u 's|M|000000|g' | go run fermi.go
5M 12M # number of people living in Chicago
=> 5.0M 12.0M
* beta 1 200 # fraction of people that have a piano
@ -74,40 +74,12 @@ x
The difference between `=: x` and `=. y` is that `=.` clears the stack after the assignment.
If you type "help", you can see a small grammar:
```
help
Operation | Variable assignment | Special
Operation: operator operand
operator: (empty) | * | / | + | -
operand: scalar | lognormal | beta | variable
lognormal: low high
beta: beta alpha beta
Variable assignment: =: variable_name
Clear stack: . | c | clear
Variable assignment and clear stack: =. variable_name
Other special operations: help | debug | exit
Examples:
+ 2
/ 2.5
* 1 10 (interpreted as lognormal)
+ 1 10
* beta 1 10
1 10 (multiplication taken as default operation)
=: x
.
1 100
+ x
exit
```
## Installation
```
make build
sudo make install
f # rather than the previous go run f.go
f # rather than the previous go run fermi.go
```
Why use make instead of the built-in go commands? Because the point of make is to be able to share command-line recipes.
@ -121,8 +93,8 @@ sed -u "s|#.*||" | sed -u 's|M|000000|g' | f
cat more/piano-tuners.f | f
cat more/piano-tuners-commented.f | sed -u "s|#.*||" | sed -u 's|M|000000|g' | f
tee -a input.log | go run f.go | tee -a output.log
tee -a io.log | go run f.go | tee -a io.log
tee -a input.log | go run fermi.go | tee -a output.log
tee -a io.log | go run fermi.go | tee -a io.log
function f(){
sed -u "s|#.*||" |
@ -142,17 +114,8 @@ Note that these sed commands are just hacks, and won't parse e.g., `3.5K` correc
- Sums and divisions now also supported
- For things between 0 and 1, consider using a beta distribution
## Different levels of complexity
The top level f.go file (400 lines) has a bunch of complexity: variables, parenthesis, samples, beta distributions. In the simple/ folder:
- f_simple.go (370 lines) strips variables and parenthesis, but keeps beta distributions, samples, and addition and substraction
- f_minimal.go (140 lines) strips everything that isn't lognormal and scalar multiplication and addition, plus a few debug options.
## Roadmap
Done:
- [x] Write README
- [x] Add division?
- [x] Read from file?
@ -161,10 +124,13 @@ Done:
- [x] Use a sed filter?
- [x] Add show more info version
- [x] Scalar multiplication and division
- [ ] Program into a small device, like a calculator?
- [-] Think of some way of calling bc
- [x] Think how to integrate with squiggle.c to draw samples
- [x] Copy the time to botec go code
- [x] Define samplers
- [x] Call those samplers when operating on distributions that can't be operted on algebraically
- [ ] Think about how to draw a histogram from samples
- [x] Display output more nicely, with K/M/B/T
- [x] Consider the following: make this into a stack-based DSL, with:
- [x] Variables that can be saved to and then displayed
@ -174,20 +140,8 @@ Done:
- Joint types
- Enums
- [x] Fix correlation problem, by spinning up a new randomness thing every time some serial computation is done.
- [x] Clean up error code. Right now only needed for division
- [x] Maintain *both* a more complex thing that's more featureful *and* the more simple multiplication of lognormals thing.
To (possibly) do:
- [ ] Document parenthesis syntax
- [ ] Allow input with K/M/T
- [ ] Add functions. Now easier to do with an explicit representation of the stakc
- [ ] Think about how to draw a histogram from samples
- [ ] Dump samples to file
- [ ] Represent samples/statistics in some other way
- [ ] Perhaps use qsort rather than full sorting
- [ ] Program into a small device, like a calculator?
Discarded:
- [ ] ~~Think of some way of calling bc~~
Some possible syntax for a more expressive stack-based DSL (now implemented)

299
f.go → fermi.go
View File

@ -4,8 +4,8 @@ import (
"bufio"
"errors"
"fmt"
"git.nunosempere.com/NunoSempere/fermi/pretty"
"git.nunosempere.com/NunoSempere/fermi/sample"
"git.nunosempere.com/NunoSempere/fermi/pretty"
"math"
"os"
"sort"
@ -14,6 +14,7 @@ import (
)
/* Types and interfaces */
type Stack struct {
old_dist Dist
vars map[string]Dist
@ -51,12 +52,14 @@ func (p Scalar) Samples() []float64 {
func (ln Lognormal) Samples() []float64 {
sampler := func(r sample.Src) float64 { return sample.Sample_to(ln.low, ln.high, r) }
// return sample.Sample_parallel(sampler, N_SAMPLES)
// Can't do parallel because then I'd have to await throughout the code
return sample.Sample_serially(sampler, N_SAMPLES)
}
func (beta Beta) Samples() []float64 {
sampler := func(r sample.Src) float64 { return sample.Sample_beta(beta.a, beta.b, r) }
// return sample.Sample_parallel(sampler, N_SAMPLES)
return sample.Sample_serially(sampler, N_SAMPLES)
}
@ -65,7 +68,7 @@ func (fs FilledSamples) Samples() []float64 {
}
/* Constants */
const HELP_MSG = " Operation | Variable assignment | Special\n" +
const GENERAL_ERR_MSG = " Operation | Variable assignment | Special\n" +
" Operation: operator operand\n" +
" operator: (empty) | * | / | + | -\n" +
" operand: scalar | lognormal | beta | variable\n" +
@ -91,55 +94,131 @@ const NORMAL90CONFIDENCE = 1.6448536269514727
const INIT_DIST Scalar = Scalar(1)
const N_SAMPLES = 100_000
/* Printers */
/* Pretty print for distributions */
// Needs types
func prettyPrintDist(dist Dist) {
switch v := dist.(type) {
case Lognormal:
fmt.Printf("=> ")
pretty.PrettyPrint2Floats(v.low, v.high)
fmt.Println()
case FilledSamples:
tmp_xs := make([]float64, N_SAMPLES)
copy(tmp_xs, v.xs)
sort.Slice(tmp_xs, func(i, j int) bool {
return tmp_xs[i] < tmp_xs[j]
})
low_int := N_SAMPLES / 20
low := tmp_xs[low_int]
high_int := N_SAMPLES * 19 / 20
high := tmp_xs[high_int]
fmt.Printf("=> ")
pretty.PrettyPrintFloat(low)
fmt.Printf(" ")
pretty.PrettyPrintFloat(high)
fmt.Printf(" (")
pretty.PrettyPrintInt(N_SAMPLES)
fmt.Printf(" samples)\n")
case Beta:
fmt.Printf("=> beta ")
pretty.PrettyPrint2Floats(v.a, v.b)
fmt.Println()
case Scalar:
fmt.Printf("=> scalar ")
w := float64(v)
pretty.PrettyPrintFloat(w)
fmt.Println()
case FilledSamples:
sorted_xs := make([]float64, N_SAMPLES)
copy(sorted_xs, v.xs)
sort.Slice(sorted_xs, func(i, j int) bool {
return sorted_xs[i] < sorted_xs[j]
})
low := sorted_xs[N_SAMPLES/20]
high := sorted_xs[N_SAMPLES*19/20]
fmt.Printf("=> ")
pretty.PrettyPrint2Floats(low, high)
fmt.Printf(" (")
pretty.PrettyPrintInt(N_SAMPLES)
fmt.Printf(" samples)")
fmt.Println()
default:
fmt.Printf("%v\n", v)
fmt.Printf("%v", v)
}
}
func printAndReturnErr(err_msg string) error {
// Parse line into Distribution
func parseLineErr(err_msg string) (string, Dist, error) {
fmt.Println(GENERAL_ERR_MSG)
fmt.Println(err_msg)
fmt.Println(HELP_MSG)
return errors.New(err_msg)
var errorDist Dist
return "", errorDist, errors.New(err_msg)
}
func parseLineIntoOpAndDist(line string, vars map[string]Dist) (string, Dist, error) {
words := strings.Split(strings.TrimSpace(line), " ")
op := ""
var dist Dist
switch words[0] {
case "*", "/", "+", "-":
op = words[0]
words = words[1:]
default:
op = "*" // later, change the below to
}
switch len(words) {
case 0:
return parseLineErr("Operator must have operand; can't operate on nothing")
case 1:
var_word, var_word_exists := vars[words[0]]
single_float, err1 := strconv.ParseFloat(words[0], 64) // abstract this away to search for K/M/B/T/etc.
switch {
case var_word_exists:
dist = var_word
case err1 == nil:
dist = Scalar(single_float)
case err1 != nil && !var_word_exists:
return parseLineErr("Trying to operate on a scalar, but scalar is neither a float nor an assigned variable")
}
case 2:
new_low, err1 := strconv.ParseFloat(words[0], 64)
new_high, err2 := strconv.ParseFloat(words[1], 64)
if err1 != nil || err2 != nil {
return parseLineErr("Trying to operate by a distribution, but distribution is not specified as two floats")
}
dist = Lognormal{low: new_low, high: new_high}
case 3:
if words[0] == "beta" || words[0] == "b" {
a, err1 := strconv.ParseFloat(words[1], 64)
b, err2 := strconv.ParseFloat(words[2], 64)
if err1 != nil || err2 != nil {
return parseLineErr("Trying to specify a beta distribution? Try beta 1 2")
}
dist = Beta{a: a, b: b}
} else {
return parseLineErr("Input not understood or not implemented yet")
}
default:
return parseLineErr("Input not understood or not implemented yet")
}
return op, dist, nil
}
func multiplyLogDists(l1 Lognormal, l2 Lognormal) Lognormal {
logmean1 := (math.Log(l1.high) + math.Log(l1.low)) / 2.0
logstd1 := (math.Log(l1.high) - math.Log(l1.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean2 := (math.Log(l2.high) + math.Log(l2.low)) / 2.0
logstd2 := (math.Log(l2.high) - math.Log(l2.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean_product := logmean1 + logmean2
logstd_product := math.Sqrt(logstd1*logstd1 + logstd2*logstd2)
h := logstd_product * NORMAL90CONFIDENCE
loglow := logmean_product - h
loghigh := logmean_product + h
return Lognormal{low: math.Exp(loglow), high: math.Exp(loghigh)}
}
func multiplyBetaDists(beta1 Beta, beta2 Beta) Beta {
return Beta{a: beta1.a + beta2.a, b: beta1.b + beta2.b}
}
/* Operations */
// Generic operations with samples
func operateDistsAsSamples(dist1 Dist, dist2 Dist, op string) (Dist, error) {
xs := dist1.Samples()
ys := dist2.Samples()
// fmt.Printf("xs: %v\n", xs)
// fmt.Printf("ys: %v\n", ys)
zs := make([]float64, N_SAMPLES)
for i := 0; i < N_SAMPLES; i++ {
@ -160,31 +239,10 @@ func operateDistsAsSamples(dist1 Dist, dist2 Dist, op string) (Dist, error) {
}
}
// fmt.Printf("%v\n", zs)
return FilledSamples{xs: zs}, nil
}
// Multiplication
func multiplyLogDists(l1 Lognormal, l2 Lognormal) Lognormal {
logmean1 := (math.Log(l1.high) + math.Log(l1.low)) / 2.0
logstd1 := (math.Log(l1.high) - math.Log(l1.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean2 := (math.Log(l2.high) + math.Log(l2.low)) / 2.0
logstd2 := (math.Log(l2.high) - math.Log(l2.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean_product := logmean1 + logmean2
logstd_product := math.Sqrt(logstd1*logstd1 + logstd2*logstd2)
h := logstd_product * NORMAL90CONFIDENCE
loglow := logmean_product - h
loghigh := logmean_product + h
return Lognormal{low: math.Exp(loglow), high: math.Exp(loghigh)}
}
func multiplyBetaDists(beta1 Beta, beta2 Beta) Beta {
return Beta{a: beta1.a + beta2.a, b: beta1.b + beta2.b}
}
func multiplyDists(old_dist Dist, new_dist Dist) (Dist, error) {
switch o := old_dist.(type) {
@ -232,16 +290,10 @@ func divideDists(old_dist Dist, new_dist Dist) (Dist, error) {
{
switch n := new_dist.(type) {
case Lognormal:
if n.high == 0 || n.low == 0 {
fmt.Println("Error: Can't divide by 0.0")
return nil, errors.New("Error: division by zero")
}
// to do: check division by zero
return multiplyLogDists(o, Lognormal{low: 1.0 / n.high, high: 1.0 / n.low}), nil
case Scalar:
if n == 0.0 {
fmt.Println("Error: Can't divide by 0.0")
return nil, errors.New("Error: division by zero scalar")
}
// to do: check division by zero
return multiplyLogDists(o, Lognormal{low: 1.0 / float64(n), high: 1.0 / float64(n)}), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "/")
@ -253,10 +305,7 @@ func divideDists(old_dist Dist, new_dist Dist) (Dist, error) {
case Lognormal:
return multiplyLogDists(Lognormal{low: float64(o), high: float64(o)}, Lognormal{low: 1.0 / n.high, high: 1.0 / n.low}), nil
case Scalar:
if n == 0.0 {
fmt.Println("Error: Can't divide by 0.0")
return nil, errors.New("Error: division by zero scalar")
}
// to do: check division by zero
return Scalar(float64(o) / float64(n)), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "/")
@ -267,133 +316,95 @@ func divideDists(old_dist Dist, new_dist Dist) (Dist, error) {
}
}
// Generic distribution operations
func operateDists(old_dist Dist, new_dist Dist, op string) (Dist, error) {
/* Combine old dist and new line */
// We want this as a function to be able to have parenthesis/recusion, possibly functions
func operateStackWithDist(stack Stack, new_dist Dist, op string) Stack {
var combined_dist Dist
var err error
switch op {
case "*":
return multiplyDists(old_dist, new_dist)
if combined_dist, err = multiplyDists(stack.old_dist, new_dist); err == nil {
stack.old_dist = combined_dist
}
case "/":
return divideDists(old_dist, new_dist)
if combined_dist, err = divideDists(stack.old_dist, new_dist); err == nil {
stack.old_dist = combined_dist
}
case "+":
return operateDistsAsSamples(old_dist, new_dist, "+")
if combined_dist, err = operateDistsAsSamples(stack.old_dist, new_dist, "+"); err == nil {
stack.old_dist = combined_dist
}
case "-":
return operateDistsAsSamples(old_dist, new_dist, "-")
default:
return nil, printAndReturnErr("Can't combine distributions in this way")
}
}
/* Parser and repl */
func parseWordsErr(err_msg string) (string, Dist, error) {
return "", nil, printAndReturnErr(err_msg)
}
func parseWordsIntoOpAndDist(words []string, vars map[string]Dist) (string, Dist, error) {
op := ""
var dist Dist
switch words[0] {
case "*", "/", "+", "-":
op = words[0]
words = words[1:]
default:
op = "*" // later, change the below to
}
switch len(words) {
case 0:
return parseWordsErr("Operator must have operand; can't operate on nothing")
case 1:
var_word, var_word_exists := vars[words[0]]
single_float, err1 := strconv.ParseFloat(words[0], 64) // abstract this away to search for K/M/B/T/etc.
switch {
case var_word_exists:
dist = var_word
case err1 == nil:
dist = Scalar(single_float)
case err1 != nil && !var_word_exists:
return parseWordsErr("Trying to operate on a scalar, but scalar is neither a float nor an assigned variable")
}
case 2:
new_low, err1 := strconv.ParseFloat(words[0], 64)
new_high, err2 := strconv.ParseFloat(words[1], 64)
if err1 != nil || err2 != nil {
return parseWordsErr("Trying to operate by a distribution, but distribution is not specified as two floats")
}
dist = Lognormal{low: new_low, high: new_high}
case 3:
if words[0] == "beta" || words[0] == "b" {
a, err1 := strconv.ParseFloat(words[1], 64)
b, err2 := strconv.ParseFloat(words[2], 64)
if err1 != nil || err2 != nil {
return parseWordsErr("Trying to specify a beta distribution? Try beta 1 2")
}
dist = Beta{a: a, b: b}
} else {
return parseWordsErr("Input not understood or not implemented yet")
if combined_dist, err = operateDistsAsSamples(stack.old_dist, new_dist, "-"); err == nil {
stack.old_dist = combined_dist
}
default:
return parseWordsErr("Input not understood or not implemented yet")
fmt.Println("Can't combine distributions in this way")
}
return op, dist, nil
}
return stack
/* Combine old dist and new line */
// We want this as a function (rather than just be in main)
// to be able to have parenthesis/recusion, possibly functions
}
func runRepl(stack Stack, reader *bufio.Reader) Stack {
replForLoop:
replForLoop:
for {
new_line, _ := reader.ReadString('\n')
words := strings.Split(strings.TrimSpace(new_line), " ")
switch {
case strings.TrimSpace(new_line) == "": /* Empty line case */
/* Empty line case */
case strings.TrimSpace(new_line) == "":
continue replForLoop
/* Parenthesis */
case len(words) == 2 && (words[0] == "*" || words[0] == "+" || words[0] == "-" || words[0] == "/") && words[1] == "(":
new_stack := runRepl(Stack{old_dist: INIT_DIST, vars: stack.vars}, reader)
combined_dist, err := operateDists(stack.old_dist, new_stack.old_dist, words[0])
if err == nil {
stack.old_dist = combined_dist
}
new_stack := runRepl(Stack{old_dist: INIT_DIST, vars: stack.vars}, reader) // to do: think if I want to shadow variables or not (right now, variables persist, since I'm copying the map(.
stack = operateStackWithDist(stack, new_stack.old_dist, words[0])
prettyPrintDist(stack.old_dist)
case len(words) == 1 && words[0] == ")":
return stack
/* Special operations */
case words[0] == "exit" || words[0] == "e":
os.Exit(0)
case words[0] == "help" || words[0] == "h":
fmt.Println(HELP_MSG)
fmt.Println(GENERAL_ERR_MSG)
continue replForLoop
case words[0] == "debug" || words[0] == "d":
fmt.Printf("%v", stack)
fmt.Printf("Old dist: %v\n", stack.old_dist)
fmt.Printf("Vars: %v\n", stack.vars)
continue replForLoop
case words[0] == "clear" || words[0] == "c" || words[0] == ".":
stack.old_dist = INIT_DIST
fmt.Println()
/* Variable assignment */
continue replForLoop
/* Variable assignment */
case words[0] == "=:" && len(words) == 2:
stack.vars[words[1]] = stack.old_dist
fmt.Printf("%s ", words[1])
prettyPrintDist(stack.old_dist)
continue replForLoop
case words[0] == "=." && len(words) == 2:
stack.vars[words[1]] = stack.old_dist
fmt.Printf("%s ", words[1])
prettyPrintDist(stack.old_dist)
stack.old_dist = INIT_DIST
// fmt.Println()
// continue replForLoop
fmt.Println()
continue replForLoop
default:
op, new_dist, err := parseWordsIntoOpAndDist(words, stack.vars)
op, new_dist, err := parseLineIntoOpAndDist(new_line, stack.vars)
if err != nil {
continue replForLoop
}
combined_dist, err := operateDists(stack.old_dist, new_dist, op)
if err == nil {
stack.old_dist = combined_dist
}
stack = operateStackWithDist(stack, new_dist, op)
prettyPrintDist(stack.old_dist)
}
prettyPrintDist(stack.old_dist)
}
}
/* Main event loop */
func main() {
reader := bufio.NewReader(os.Stdin)
stack := Stack{old_dist: INIT_DIST, vars: make(map[string]Dist)}
runRepl(stack, reader)

10
makefile
View File

@ -1,9 +1,9 @@
run:
go run f.go
go run fermi.go
build:
go build f.go
go build fermi.go
install: f
sudo mv f /usr/bin/f
sudo ln -s /usr/bin/f /usr/bin/fermi
install: fermi
sudo mv fermi /usr/bin/fermi
sudo ln -s /usr/bin/fermi /usr/bin/f

1
pretty/pretty.go
View File

@ -48,4 +48,5 @@ func PrettyPrint2Floats(low float64, high float64) {
PrettyPrintFloat(low)
fmt.Printf(" ")
PrettyPrintFloat(high)
fmt.Printf("\n")
}

140
simple/f_minimal.go
View File

@ -1,140 +0,0 @@
package main
import (
"bufio"
"errors"
"fmt"
"git.nunosempere.com/NunoSempere/fermi/pretty"
"math"
"os"
"strconv"
"strings"
)
type Lognormal struct {
low float64
high float64
}
const HELP_MSG = " Operation | Variable assignment | Special\n" +
" Operation: operator operand\n" +
" operator: (empty) | * | / | + | -\n" +
" operand: scalar | lognormal\n" +
" lognormal: low high\n" +
" Clear stack: . | c | clear\n" +
" Other special operations: help | debug | exit\n" +
" Examples: \n" +
" / 2.5\n" +
" * 1 10 (interpreted as lognormal)\n" +
" / 1 10\n" +
" 1 10 (multiplication taken as default operation)\n" +
" .\n" +
" exit\n"
const NORMAL90CONFIDENCE = 1.6448536269514727
const N_SAMPLES = 100_000
func prettyPrintLognormal(l Lognormal) {
fmt.Printf("=> ")
pretty.PrettyPrint2Floats(l.low, l.high)
fmt.Println()
}
func printAndReturnErr(err_msg string) error {
fmt.Println(err_msg)
fmt.Println(HELP_MSG)
return errors.New(err_msg)
}
func multiplyLogDists(l1 Lognormal, l2 Lognormal) Lognormal {
logmean1 := (math.Log(l1.high) + math.Log(l1.low)) / 2.0
logstd1 := (math.Log(l1.high) - math.Log(l1.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean2 := (math.Log(l2.high) + math.Log(l2.low)) / 2.0
logstd2 := (math.Log(l2.high) - math.Log(l2.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean_product := logmean1 + logmean2
logstd_product := math.Sqrt(logstd1*logstd1 + logstd2*logstd2)
h := logstd_product * NORMAL90CONFIDENCE
loglow := logmean_product - h
loghigh := logmean_product + h
return Lognormal{low: math.Exp(loglow), high: math.Exp(loghigh)}
}
func divideLogDists(l1 Lognormal, l2 Lognormal) (Lognormal, error) {
if l2.high == 0 || l2.low == 0 {
fmt.Println("Error: Can't divide by 0.0")
return Lognormal{}, errors.New("Error: division by zero")
}
return multiplyLogDists(l1, Lognormal{low: 1.0 / l2.high, high: 1.0 / l2.low}), nil
}
func parseWordsErr(err_msg string) (string, Lognormal, error) {
return "", Lognormal{}, printAndReturnErr(err_msg)
}
func parseWordsIntoOpAndDist(words []string) (string, Lognormal, error) {
op := ""
var dist Lognormal
switch words[0] {
case "*", "/":
op = words[0]
words = words[1:]
default:
op = "*"
}
switch len(words) {
case 0:
return parseWordsErr("Operator must have operand; can't operate on nothing")
case 1:
single_float, err := strconv.ParseFloat(words[0], 64) // abstract this away to search for K/M/B/T/etc.
if err != nil {
return parseWordsErr("Trying to operate on a scalar, but scalar is neither a float nor an assigned variable")
}
dist = Lognormal{low: single_float, high: single_float}
case 2:
new_low, err1 := strconv.ParseFloat(words[0], 64)
new_high, err2 := strconv.ParseFloat(words[1], 64)
if err1 != nil || err2 != nil {
return parseWordsErr("Trying to operate by a distribution, but distribution is not specified as two floats")
}
dist = Lognormal{low: new_low, high: new_high}
default:
return parseWordsErr("Input not understood or not implemented yet")
}
return op, dist, nil
}
func main() {
reader := bufio.NewReader(os.Stdin)
old_dist := Lognormal{low: 1, high: 1}
replForLoop:
for {
new_line, _ := reader.ReadString('\n')
words := strings.Split(strings.TrimSpace(new_line), " ")
if strings.TrimSpace(new_line) == "" {
continue replForLoop
}
op, new_dist, err := parseWordsIntoOpAndDist(words)
if err != nil {
continue replForLoop
}
switch op {
case "*":
old_dist = multiplyLogDists(old_dist, new_dist)
case "/":
result_dist, err := divideLogDists(old_dist, new_dist)
if err != nil {
continue replForLoop
}
old_dist = result_dist
}
prettyPrintLognormal(old_dist)
}
}

370
simple/f_simple.go
View File

@ -1,370 +0,0 @@
package main
import (
"bufio"
"errors"
"fmt"
"git.nunosempere.com/NunoSempere/fermi/pretty"
"git.nunosempere.com/NunoSempere/fermi/sample"
"math"
"os"
"sort"
"strconv"
"strings"
)
/* Types and interfaces */
type Dist interface {
Samples() []float64
}
type Scalar float64
type Lognormal struct {
low float64
high float64
}
type Beta struct {
a float64
b float64
}
type FilledSamples struct {
xs []float64
}
/* Dist interface functions */
// https://go.dev/tour/methods/9
func (p Scalar) Samples() []float64 {
xs := make([]float64, N_SAMPLES)
for i := 0; i < N_SAMPLES; i++ {
xs[i] = float64(p)
}
return xs
}
func (ln Lognormal) Samples() []float64 {
sampler := func(r sample.Src) float64 { return sample.Sample_to(ln.low, ln.high, r) }
return sample.Sample_serially(sampler, N_SAMPLES)
}
func (beta Beta) Samples() []float64 {
sampler := func(r sample.Src) float64 { return sample.Sample_beta(beta.a, beta.b, r) }
return sample.Sample_serially(sampler, N_SAMPLES)
}
func (fs FilledSamples) Samples() []float64 {
return fs.xs
}
/* Constants */
const HELP_MSG = " Operation | Variable assignment | Special\n" +
" Operation: operator operand\n" +
" operator: (empty) | * | / | + | -\n" +
" operand: scalar | lognormal | beta\n" +
" lognormal: low high\n" +
" beta: beta alpha beta\n" +
" Clear stack: . | c | clear\n" +
" Other special operations: help | debug | exit\n" +
" Examples: \n" +
" + 2\n" +
" / 2.5\n" +
" * 1 10 (interpreted as lognormal)\n" +
" + 1 10\n" +
" * beta 1 10\n" +
" 1 10 (multiplication taken as default operation)\n" +
" .\n" +
" 1 100\n" +
" exit\n"
const NORMAL90CONFIDENCE = 1.6448536269514727
const INIT_DIST Scalar = Scalar(1)
const N_SAMPLES = 100_000
/* Printers */
func prettyPrintDist(dist Dist) {
switch v := dist.(type) {
case Lognormal:
fmt.Printf("=> ")
pretty.PrettyPrint2Floats(v.low, v.high)
fmt.Println()
case Beta:
fmt.Printf("=> beta ")
pretty.PrettyPrint2Floats(v.a, v.b)
fmt.Println()
case Scalar:
fmt.Printf("=> scalar ")
w := float64(v)
pretty.PrettyPrintFloat(w)
fmt.Println()
case FilledSamples:
sorted_xs := make([]float64, N_SAMPLES)
copy(sorted_xs, v.xs)
sort.Slice(sorted_xs, func(i, j int) bool {
return sorted_xs[i] < sorted_xs[j]
})
low := sorted_xs[N_SAMPLES/20]
high := sorted_xs[N_SAMPLES*19/20]
fmt.Printf("=> ")
pretty.PrettyPrint2Floats(low, high)
fmt.Printf(" (")
pretty.PrettyPrintInt(N_SAMPLES)
fmt.Printf(" samples)")
fmt.Println()
default:
fmt.Printf("%v\n", v)
}
}
func printAndReturnErr(err_msg string) error {
fmt.Println(err_msg)
fmt.Println(HELP_MSG)
return errors.New(err_msg)
}
/* Operations */
// Generic operations with samples
func operateDistsAsSamples(dist1 Dist, dist2 Dist, op string) (Dist, error) {
xs := dist1.Samples()
ys := dist2.Samples()
zs := make([]float64, N_SAMPLES)
for i := 0; i < N_SAMPLES; i++ {
switch op {
case "*":
zs[i] = xs[i] * ys[i]
case "/":
if ys[0] != 0 {
zs[i] = xs[i] / ys[i]
} else {
fmt.Println("Error: When dividing as samples, division by zero")
return nil, errors.New("Division by zero")
}
case "+":
zs[i] = xs[i] + ys[i]
case "-":
zs[i] = xs[i] - ys[i]
}
}
// fmt.Printf("%v\n", zs)
return FilledSamples{xs: zs}, nil
}
// Multiplication
func multiplyLogDists(l1 Lognormal, l2 Lognormal) Lognormal {
logmean1 := (math.Log(l1.high) + math.Log(l1.low)) / 2.0
logstd1 := (math.Log(l1.high) - math.Log(l1.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean2 := (math.Log(l2.high) + math.Log(l2.low)) / 2.0
logstd2 := (math.Log(l2.high) - math.Log(l2.low)) / (2.0 * NORMAL90CONFIDENCE)
logmean_product := logmean1 + logmean2
logstd_product := math.Sqrt(logstd1*logstd1 + logstd2*logstd2)
h := logstd_product * NORMAL90CONFIDENCE
loglow := logmean_product - h
loghigh := logmean_product + h
return Lognormal{low: math.Exp(loglow), high: math.Exp(loghigh)}
}
func multiplyBetaDists(beta1 Beta, beta2 Beta) Beta {
return Beta{a: beta1.a + beta2.a, b: beta1.b + beta2.b}
}
func multiplyDists(old_dist Dist, new_dist Dist) (Dist, error) {
switch o := old_dist.(type) {
case Lognormal:
{
switch n := new_dist.(type) {
case Lognormal:
return multiplyLogDists(o, n), nil
case Scalar:
return multiplyLogDists(o, Lognormal{low: float64(n), high: float64(n)}), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "*")
}
}
case Scalar:
{
if o == 1 {
return new_dist, nil
}
switch n := new_dist.(type) {
case Lognormal:
return multiplyLogDists(Lognormal{low: float64(o), high: float64(o)}, n), nil
case Scalar:
return Scalar(float64(o) * float64(n)), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "*")
}
}
case Beta:
switch n := new_dist.(type) {
case Beta:
return multiplyBetaDists(o, n), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "*")
}
default:
return operateDistsAsSamples(old_dist, new_dist, "*")
}
}
func divideDists(old_dist Dist, new_dist Dist) (Dist, error) {
switch o := old_dist.(type) {
case Lognormal:
{
switch n := new_dist.(type) {
case Lognormal:
// to do: check division by zero
if n.high == 0 || n.low == 0 {
fmt.Println("Error: Can't divide by 0.0")
return nil, errors.New("Error: division by zero")
}
return multiplyLogDists(o, Lognormal{low: 1.0 / n.high, high: 1.0 / n.low}), nil
case Scalar:
// to do: check division by zero
if n == 0.0 {
fmt.Println("Error: Can't divide by 0.0")
return nil, errors.New("Error: division by zero scalar")
}
return multiplyLogDists(o, Lognormal{low: 1.0 / float64(n), high: 1.0 / float64(n)}), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "/")
}
}
case Scalar:
{
switch n := new_dist.(type) {
case Lognormal:
return multiplyLogDists(Lognormal{low: float64(o), high: float64(o)}, Lognormal{low: 1.0 / n.high, high: 1.0 / n.low}), nil
case Scalar:
// to do: check division by zero
if n == 0.0 {
fmt.Println("Error: Can't divide by 0.0")
return nil, errors.New("Error: division by zero scalar")
}
return Scalar(float64(o) / float64(n)), nil
default:
return operateDistsAsSamples(old_dist, new_dist, "/")
}
}
default:
return operateDistsAsSamples(old_dist, new_dist, "/")
}
}
// Generic distribution operations
func operateDists(old_dist Dist, new_dist Dist, op string) (Dist, error) {
switch op {
case "*":
return multiplyDists(old_dist, new_dist)
case "/":
return divideDists(old_dist, new_dist)
case "+":
return operateDistsAsSamples(old_dist, new_dist, "+")
case "-":
return operateDistsAsSamples(old_dist, new_dist, "-")
default:
return nil, printAndReturnErr("Can't combine distributions in this way")
}
}
/* Parser and repl */
func parseWordsErr(err_msg string) (string, Dist, error) {
return "", nil, printAndReturnErr(err_msg)
}
func parseWordsIntoOpAndDist(words []string) (string, Dist, error) {
op := ""
var dist Dist
switch words[0] {
case "*", "/", "+", "-":
op = words[0]
words = words[1:]
default:
op = "*" // later, change the below to
}
switch len(words) {
case 0:
return parseWordsErr("Operator must have operand; can't operate on nothing")
case 1:
single_float, err := strconv.ParseFloat(words[0], 64) // abstract this away to search for K/M/B/T/etc.
if err != nil {
return parseWordsErr("Trying to operate on a scalar, but scalar is neither a float nor an assigned variable")
}
dist = Scalar(single_float)
case 2:
new_low, err1 := strconv.ParseFloat(words[0], 64)
new_high, err2 := strconv.ParseFloat(words[1], 64)
if err1 != nil || err2 != nil {
return parseWordsErr("Trying to operate by a distribution, but distribution is not specified as two floats")
}
dist = Lognormal{low: new_low, high: new_high}
case 3:
if words[0] == "beta" || words[0] == "b" {
a, err1 := strconv.ParseFloat(words[1], 64)
b, err2 := strconv.ParseFloat(words[2], 64)
if err1 != nil || err2 != nil {
return parseWordsErr("Trying to specify a beta distribution? Try beta 1 2")
}
dist = Beta{a: a, b: b}
} else {
return parseWordsErr("Input not understood or not implemented yet")
}
default:
return parseWordsErr("Input not understood or not implemented yet")
}
return op, dist, nil
}
func main() {
reader := bufio.NewReader(os.Stdin)
var old_dist Dist = INIT_DIST
replForLoop:
for {
new_line, _ := reader.ReadString('\n')
words := strings.Split(strings.TrimSpace(new_line), " ")
switch {
case strings.TrimSpace(new_line) == "": /* Empty line case */
continue replForLoop
case words[0] == "exit" || words[0] == "e":
os.Exit(0)
case words[0] == "help" || words[0] == "h":
fmt.Println(HELP_MSG)
continue replForLoop
case words[0] == "debug" || words[0] == "d":
fmt.Printf("%v", old_dist)
continue replForLoop
case words[0] == "clear" || words[0] == "c" || words[0] == ".":
old_dist = INIT_DIST
fmt.Println()
continue replForLoop
}
op, new_dist, err := parseWordsIntoOpAndDist(words)
if err != nil {
continue replForLoop
}
combined_dist, err := operateDists(old_dist, new_dist, op)
if err != nil {
continue replForLoop
}
old_dist = combined_dist
prettyPrintDist(old_dist)
}
}