diff --git a/packages/squiggle-lang/src/rescript/FunctionRegistry/Library/FR_Danger.res b/packages/squiggle-lang/src/rescript/FunctionRegistry/Library/FR_Danger.res index f49f2197..730e1b9f 100644 --- a/packages/squiggle-lang/src/rescript/FunctionRegistry/Library/FR_Danger.res +++ b/packages/squiggle-lang/src/rescript/FunctionRegistry/Library/FR_Danger.res @@ -57,6 +57,187 @@ module Combinatorics = { } } +module Integration = { + let Helpers = { + let integrateFunctionBetweenWithNumIntegrationPoints = ( + aLambda, + min: float, + max: float, + numIntegrationPoints: float, // cast as int? + environment, + reducer, + ) => { + let applyFunctionAtFloatToFloatOption = (point: float) => { + // Defined here so that it has access to environment, reducer + let pointAsInternalExpression = FunctionRegistry_Helpers.Wrappers.evNumber(point) + let resultAsInternalExpression = Reducer_Expression_Lambda.doLambdaCall( + aLambda, + list{pointAsInternalExpression}, + environment, + reducer, + ) + let result = switch resultAsInternalExpression { + | Ok(IEvNumber(x)) => Ok(x) + | Error(_) => + Error( + "Error 1 in Danger.integrate. It's possible that your function doesn't return a number, try definining auxiliaryFunction(x) = mean(yourFunction(x)) and integrate auxiliaryFunction instead", + ) + | _ => Error("Error 2 in Danger.integrate") + } + result + } + // worked example in comments below, assuming min=0, max = 10 + let numTotalPoints = Belt.Float.toInt(numIntegrationPoints) // superflous declaration, but useful to keep track that we are interpreting "numIntegrationPoints" as the total number on which we evaluate the function, not e.g., as the inner integration points. + let numInnerPoints = numTotalPoints - 2 + let numOuterPoints = 2 + let totalWeight = max -. min + let weightForAnInnerPoint = totalWeight /. E.I.toFloat(numTotalPoints - 1) + let weightForAnOuterPoint = totalWeight /. E.I.toFloat(numTotalPoints - 1) /. 2.0 + let innerPointIncrement = (max -. min) /. E.I.toFloat(numTotalPoints - 1) + let innerXs = Belt.Array.makeBy(numInnerPoints, i => + min +. Belt_Float.fromInt(i + 1) *. innerPointIncrement + ) + // Gotcha: makeBy goes from 0 to (n-1): + let ysOptions = Belt.Array.map(innerXs, x => applyFunctionAtFloatToFloatOption(x)) + let okYs = E.A.R.filterOk(ysOptions) + + /* Logging, with a worked example. */ + // Useful for understanding what is happening. + // assuming min = 0, max = 10, numTotalPoints=10, results below: + let verbose = false + if verbose { + Js.Console.log2("numTotalPoints", numTotalPoints) // 5 + Js.Console.log2("numInnerPoints", numInnerPoints) // 3 + Js.Console.log2("numOuterPoints", numOuterPoints) // always 2 + Js.Console.log2("totalWeight", totalWeight) // 10 - 0 = 10 + Js.Console.log2("weightForAnInnerPoint", weightForAnInnerPoint) // 10/4 = 2.5 + Js.Console.log2("weightForAnOuterPoint", weightForAnOuterPoint) // 10/4/2 = 1.25 + Js.Console.log2( + "weightForAnInnerPoint * numInnerPoints + weightForAnOuterPoint * numOuterPoints", + weightForAnInnerPoint *. E.I.toFloat(numInnerPoints) +. + weightForAnOuterPoint *. E.I.toFloat(numOuterPoints), + ) // should be 10 + Js.Console.log2( + "sum of weights == totalWeight", + weightForAnInnerPoint *. E.I.toFloat(numInnerPoints) +. + weightForAnOuterPoint *. E.I.toFloat(numOuterPoints) == totalWeight, + ) // true + Js.Console.log2("innerPointIncrement", innerPointIncrement) // (10-0)/4 = 2.5 + Js.Console.log2("innerXs", innerXs) // 2.5, 5, 7.5 + Js.Console.log2("ysOptions", ysOptions) + Js.Console.log2("okYs", okYs) + } + + //This is pretty hacky. It should use a result type instead of checking that length matches. + let result = if E.A.length(ysOptions) == E.A.length(okYs) { + let innerPointsSum = okYs->E.A.reduce(0.0, (a, b) => a +. b) + let resultWithOuterPoints = switch ( + applyFunctionAtFloatToFloatOption(min), + applyFunctionAtFloatToFloatOption(max), + ) { + | (Ok(yMin), Ok(yMax)) => { + let result = + (yMin +. yMax) *. weightForAnOuterPoint +. innerPointsSum *. weightForAnInnerPoint + let wrappedResult = result->ReducerInterface_InternalExpressionValue.IEvNumber->Ok + wrappedResult + } + | (Error(b), _) => Error(b) + | (_, Error(b)) => Error(b) + } + resultWithOuterPoints + } else { + Error( + "Integration error 3 in Danger.integrate. It's possible that your function doesn't return a number, try definining auxiliaryFunction(x) = mean(yourFunction(x)) and integrate auxiliaryFunction instead", + ) + } + result + } + } + module Lib = { +// Integral in terms of function, min, max, num points + // Note that execution time will be more predictable, because it + // will only depend on num points and the complexity of the function + let integrateFunctionBetweenWithNumIntegrationPoints = Function.make( + ~name="integrateFunctionBetweenWithNumIntegrationPoints", + ~nameSpace, + ~output=EvtNumber, + ~requiresNamespace=false, + ~examples=[`Danger.integrateFunctionBetweenWithNumIntegrationPoints({|x| x+1}, 1, 10, 10)`], + // should be [x^2/2 + x]1_10 = (100/2 + 10) - (1/2 + 1) = 60 - 1.5 = 58.5 + // https://www.wolframalpha.com/input?i=integrate+x%2B1+from+1+to+10 + ~definitions=[ + FnDefinition.make( + ~name="integrateFunctionBetweenWithNumIntegrationPoints", + ~inputs=[FRTypeLambda, FRTypeNumber, FRTypeNumber, FRTypeNumber], + ~run=(inputs, _, env, reducer) => { + let result = switch inputs { + | [_, _, _, IEvNumber(0.0)] => + Error("Integration error 4 in Danger.integrate: Increment can't be 0.") + | [IEvLambda(aLambda), IEvNumber(min), IEvNumber(max), IEvNumber(numIntegrationPoints)] => + Helpers.integrateFunctionBetweenWithNumIntegrationPoints( + aLambda, + min, + max, + numIntegrationPoints, + env, + reducer, + ) + | _ => + Error( + "Integration error 5 in Danger.integrate. Remember that inputs are (function, number (min), number (max), number(increment))", + ) + } + result + }, + (), + ), + ], + (), + ) + // Integral in terms of function, min, max, epsilon (distance between points) + // Execution time will be less predictable, because it + // will depend on min, max and epsilon together, + // as well and the complexity of the function + let integrateFunctionBetweenWithEpsilon = Function.make( + ~name="integrateFunctionBetweenWithEpsilon", + ~nameSpace, + ~output=EvtNumber, + ~requiresNamespace=false, + ~examples=[`Danger.integrateFunctionBetweenWithEpsilon({|x| x+1}, 1, 10, 0.1)`], + ~definitions=[ + FnDefinition.make( + ~name="integrateFunctionBetweenWithEpsilon", + ~inputs=[FRTypeLambda, FRTypeNumber, FRTypeNumber, FRTypeNumber], + ~run=(inputs, _, env, reducer) => { + let result = switch inputs { + | [_, _, _, IEvNumber(0.0)] => + Error("Integration error in Danger.integrate: Increment can't be 0.") + | [IEvLambda(aLambda), IEvNumber(min), IEvNumber(max), IEvNumber(epsilon)] => + Helpers.integrateFunctionBetweenWithNumIntegrationPoints( + aLambda, + min, + max, + (max -. min) /. epsilon, + env, + reducer, + )->E.R2.errMap(_ => + "Integration error 7 in Danger.integrate. Something went wrong along the way" + ) + | _ => + Error( + "Integration error 8 in Danger.integrate. Remember that inputs are (function, number (min), number (max), number(increment))", + ) + } + result + }, + (), + ), + ], + (), + ) + } +} + module Internals = { // Probability functions @@ -65,100 +246,7 @@ module Internals = { ->Belt.Array.map(FunctionRegistry_Helpers.Wrappers.evNumber) ->FunctionRegistry_Helpers.Wrappers.evArray - let integrateFunctionBetweenWithNumIntegrationPoints = ( - aLambda, - min: float, - max: float, - numIntegrationPoints: float, // cast as int? - environment, - reducer, - ) => { - let applyFunctionAtFloatToFloatOption = (point: float) => { - // Defined here so that it has access to environment, reducer - let pointAsInternalExpression = FunctionRegistry_Helpers.Wrappers.evNumber(point) - let resultAsInternalExpression = Reducer_Expression_Lambda.doLambdaCall( - aLambda, - list{pointAsInternalExpression}, - environment, - reducer, - ) - let result = switch resultAsInternalExpression { - | Ok(IEvNumber(x)) => Ok(x) - | Error(_) => - Error( - "Error 1 in Danger.integrate. It's possible that your function doesn't return a number, try definining auxiliaryFunction(x) = mean(yourFunction(x)) and integrate auxiliaryFunction instead", - ) - | _ => Error("Error 2 in Danger.integrate") - } - result - } - // worked example in comments below, assuming min=0, max = 10 - let numTotalPoints = Belt.Float.toInt(numIntegrationPoints) // superflous declaration, but useful to keep track that we are interpreting "numIntegrationPoints" as the total number on which we evaluate the function, not e.g., as the inner integration points. - let numInnerPoints = numTotalPoints - 2 - let numOuterPoints = 2 - let totalWeight = max -. min - let weightForAnInnerPoint = totalWeight /. E.I.toFloat(numTotalPoints - 1) - let weightForAnOuterPoint = totalWeight /. E.I.toFloat(numTotalPoints - 1) /. 2.0 - let innerPointIncrement = (max -. min) /. E.I.toFloat(numTotalPoints - 1) - let innerXs = Belt.Array.makeBy(numInnerPoints, i => - min +. Belt_Float.fromInt(i + 1) *. innerPointIncrement - ) - // Gotcha: makeBy goes from 0 to (n-1): - let ysOptions = Belt.Array.map(innerXs, x => applyFunctionAtFloatToFloatOption(x)) - let okYs = E.A.R.filterOk(ysOptions) - - /* Logging, with a worked example. */ - // Useful for understanding what is happening. - // assuming min = 0, max = 10, numTotalPoints=10, results below: - let verbose = false - if verbose { - Js.Console.log2("numTotalPoints", numTotalPoints) // 5 - Js.Console.log2("numInnerPoints", numInnerPoints) // 3 - Js.Console.log2("numOuterPoints", numOuterPoints) // always 2 - Js.Console.log2("totalWeight", totalWeight) // 10 - 0 = 10 - Js.Console.log2("weightForAnInnerPoint", weightForAnInnerPoint) // 10/4 = 2.5 - Js.Console.log2("weightForAnOuterPoint", weightForAnOuterPoint) // 10/4/2 = 1.25 - Js.Console.log2( - "weightForAnInnerPoint * numInnerPoints + weightForAnOuterPoint * numOuterPoints", - weightForAnInnerPoint *. E.I.toFloat(numInnerPoints) +. - weightForAnOuterPoint *. E.I.toFloat(numOuterPoints), - ) // should be 10 - Js.Console.log2( - "sum of weights == totalWeight", - weightForAnInnerPoint *. E.I.toFloat(numInnerPoints) +. - weightForAnOuterPoint *. E.I.toFloat(numOuterPoints) == totalWeight, - ) // true - Js.Console.log2("innerPointIncrement", innerPointIncrement) // (10-0)/4 = 2.5 - Js.Console.log2("innerXs", innerXs) // 2.5, 5, 7.5 - Js.Console.log2("ysOptions", ysOptions) - Js.Console.log2("okYs", okYs) - } - - //This is pretty hacky. It should use a result type instead of checking that length matches. - let result = if E.A.length(ysOptions) == E.A.length(okYs) { - let innerPointsSum = okYs->E.A.reduce(0.0, (a, b) => a +. b) - let resultWithOuterPoints = switch ( - applyFunctionAtFloatToFloatOption(min), - applyFunctionAtFloatToFloatOption(max), - ) { - | (Ok(yMin), Ok(yMax)) => { - let result = - (yMin +. yMax) *. weightForAnOuterPoint +. innerPointsSum *. weightForAnInnerPoint - let wrappedResult = result->ReducerInterface_InternalExpressionValue.IEvNumber->Ok - wrappedResult - } - | (Error(b), _) => Error(b) - | (_, Error(b)) => Error(b) - } - resultWithOuterPoints - } else { - Error( - "Integration error 3 in Danger.integrate. It's possible that your function doesn't return a number, try definining auxiliaryFunction(x) = mean(yourFunction(x)) and integrate auxiliaryFunction instead", - ) - } - result - } - + // Diminishing returns // Helpers type diminishingReturnsAccumulatorInner = { @@ -272,92 +360,13 @@ module Internals = { } + let library = [ Combinatorics.Lib.laplace, Combinatorics.Lib.factorial, Combinatorics.Lib.choose, Combinatorics.Lib.binomial, - // Integral in terms of function, min, max, num points - // Note that execution time will be more predictable, because it - // will only depend on num points and the complexity of the function - Function.make( - ~name="integrateFunctionBetweenWithNumIntegrationPoints", - ~nameSpace, - ~output=EvtNumber, - ~requiresNamespace=false, - ~examples=[`Danger.integrateFunctionBetweenWithNumIntegrationPoints({|x| x+1}, 1, 10, 10)`], - // should be [x^2/2 + x]1_10 = (100/2 + 10) - (1/2 + 1) = 60 - 1.5 = 58.5 - // https://www.wolframalpha.com/input?i=integrate+x%2B1+from+1+to+10 - ~definitions=[ - FnDefinition.make( - ~name="integrateFunctionBetweenWithNumIntegrationPoints", - ~inputs=[FRTypeLambda, FRTypeNumber, FRTypeNumber, FRTypeNumber], - ~run=(inputs, _, env, reducer) => { - let result = switch inputs { - | [_, _, _, IEvNumber(0.0)] => - Error("Integration error 4 in Danger.integrate: Increment can't be 0.") - | [IEvLambda(aLambda), IEvNumber(min), IEvNumber(max), IEvNumber(numIntegrationPoints)] => - Internals.integrateFunctionBetweenWithNumIntegrationPoints( - aLambda, - min, - max, - numIntegrationPoints, - env, - reducer, - ) - | _ => - Error( - "Integration error 5 in Danger.integrate. Remember that inputs are (function, number (min), number (max), number(increment))", - ) - } - result - }, - (), - ), - ], - (), - ), - // Integral in terms of function, min, max, epsilon (distance between points) - // Execution time will be less predictable, because it - // will depend on min, max and epsilon together, - // as well and the complexity of the function - Function.make( - ~name="integrateFunctionBetweenWithEpsilon", - ~nameSpace, - ~output=EvtNumber, - ~requiresNamespace=false, - ~examples=[`Danger.integrateFunctionBetweenWithEpsilon({|x| x+1}, 1, 10, 0.1)`], - ~definitions=[ - FnDefinition.make( - ~name="integrateFunctionBetweenWithEpsilon", - ~inputs=[FRTypeLambda, FRTypeNumber, FRTypeNumber, FRTypeNumber], - ~run=(inputs, _, env, reducer) => { - let result = switch inputs { - | [_, _, _, IEvNumber(0.0)] => - Error("Integration error in Danger.integrate: Increment can't be 0.") - | [IEvLambda(aLambda), IEvNumber(min), IEvNumber(max), IEvNumber(epsilon)] => - Internals.integrateFunctionBetweenWithNumIntegrationPoints( - aLambda, - min, - max, - (max -. min) /. epsilon, - env, - reducer, - )->E.R2.errMap(_ => - "Integration error 7 in Danger.integrate. Something went wrong along the way" - ) - | _ => - Error( - "Integration error 8 in Danger.integrate. Remember that inputs are (function, number (min), number (max), number(increment))", - ) - } - result - }, - (), - ), - ], - (), - ), + // Diminishing marginal return functions // There are functions diminishingMarginalReturnsForFunctions2 through diminishingMarginalReturnsForFunctions7 // Because of this bug: