XMC objects¶
Algorithm¶
-
class
xmc.xmcAlgorithm.
XMCAlgorithm
(**keywordArgs)¶ This top-level class handles the overall algorithm: initialisation as well as everything related to error, tolerance and iterations. It also possesses the necessary methods and attributes to create new types of algorithms. However, the export of results is to be handled outside. This documentation is not yet complete.
Methods
runXMC: run an algorithm with generic structure. It is the method to call to run the algorithm.
runAsynchronousXMC: run an algorithm with generic structure exploiting the asynchronous framework. It is the method to call to run the asynchronous algorithm.
Attributes
estimatorsForHierarchy: List[List]. This is a list of instructions for the indexwise
estimations to be sent to the hierarchy optimiser. estimatorsForHierarchy[0] contains instructions for the first estimation to be sent (estimatorsForHierarchy[1] the second, etc.). estimatorsForHierarchy[0][0] is the index of the estimator concerned, as ordered in MonteCarloIndex.qoiEstimator; estimatorsForHierarchy[0][1] is the list of arguments to be passed to the value method of this estimator. These values are used in XMCAlgorithm.optimalHierarchy and eventually result in a call of the form StatisticalEstimator.value(estimatorsForHierarchy[0], *estimatorsForHierarchy[1]).
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
updateTolerance
(criteriaToUpdate=None)¶
-
indexEstimation
(coordinate, valueMethodArgs)¶
-
indexCostEstimation
(valueMethodArgs)¶
-
predictor
(choice=None)¶
-
costPredictor
()¶
-
tolerances
(*args)¶
-
hierarchy
()¶ Returns current hierarchy of the MC estimator.
-
splitTolerance
(splittingParameter)¶ Method that interfaces with the MultiCriterion class to apply tolerance splitting
-
optimalHierarchy
()¶ Method that interfaces with the HierarchyOptimiser class to compute the optimal hierarchy
-
updateHierarchy
(newHierarchy)¶ Method that interfaces with the monteCarloSample class to execute a given hierarchy
-
estimation
(assemblerCoordinates=None)¶ Method that calls the estimation method of monteCarloSampler
-
errorEstimation
(errorEstimatorCoordinates=None)¶ Method that calls the errorEstimation method of monteCarloSampler
-
updateHierarchySpace
(*args)¶ Method that interfaces with the HierarchyOptimiser class to compute the hierarchy space in which to search for the optimal hierarchy
-
stoppingCriterionFlag
(currentCost=None)¶ Call stoppingCriterion.flag with the proper arguments and return its output (a.k.a flag). Input argument: currentCost is an indication of the cost the algorithm has entailed so far; we usually use the number of iterations. Output argument: criterion flag structure as define in the MultiCriterion class.
-
runXMC
()¶ Run an algorithm with generic structure.
-
asynchronousFinalizeIteration
()¶ Method finalizing an iteration of the asynchornous framework. It synchronizes and calls all relevant methods of one single batch, the first available, before estimating convergence.
-
runAsynchronousXMC
()¶ Run algorithm with asynchronous framework.
-
xmc.methodDefs_xmcAlgorithm.checkInitialisation.
checkInitialisationMC
(XMCAlgorithm)¶ Method checking all attributes of different classes are correctly set to run Monte Carlo algorithm (both standard and asynchronous).
-
xmc.methodDefs_xmcAlgorithm.checkInitialisation.
checkInitialisationCMLMC
()¶
-
xmc.methodDefs_xmcAlgorithm.checkInitialisation.
checkInitialisationAMLMC
()¶
-
xmc.methodDefs_xmcAlgorithm.checkInitialisation.
checkInitialisationMLMC
(XMCAlgorithm)¶ Method checking all attributes of different classes are correctly set to run Multilevel Monte Carlo algorithm (both standard and asynchronous).
-
xmc.methodDefs_xmcAlgorithm.checkInitialisation.
checkInitialisationSolverWrapper
(solverWrapperDictionary)¶ Method checking solver wrapper keys are correctly set. A required key is outputBatchSize, which is needed to organize the quantities of interest into sublists (of future objects). For this reason, outputBatchSize must be smaller than or equal to the total number of quantities of interest. It is required if areSamplesSplit() is true.
-
xmc.methodDefs_xmcAlgorithm.checkInitialisation.
checkMaxNumberIterationsCriterion
(positionMaxNumberIterationsCriterion, tolerances)¶ Method checking maximum number of iterations is set and is the last entry of MonoCriteriaDictionary, which is required for consistency.
-
xmc.methodDefs_xmcAlgorithm.updateTolerance.
updateToleranceConstant
(*args)¶
-
xmc.methodDefs_xmcAlgorithm.updateTolerance.
updateToleranceGeometric
(*args)¶ Return value of tolerance based on iteration number and geometric decrease trend
Stopping criterion¶
-
class
xmc.multiCriterion.
MultiCriterion
(**keywordArgs)¶ This class handles any number of criteria, the associated tolerances (with possible splitting), and the combination of these criteria into an interpred output.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
tolerances
(criteriaReferences=None)¶ Returns the tolerances of the requested criteria.
-
splittingParameter
()¶ Returns the splitting parameter currently applied
-
setTolerance
(criterionReferences, toleranceValues)¶
-
splitTolerance
(splittingParameter)¶
-
updateTolerance
(criteriaToUpdate=None)¶ Update the tolerances of self.criteria entries specified by criteriaToUpdate
-
flag
(values)¶ Return the output of the criterion. It takes the expected values as input, evaluate each elementary criterion on them and combine their boolean outputs into a dictionary flag. This is currently a wrapper for the protected _flagFunction attribute.
-
static
flagStructure
()¶ This function returns the typical output of a MultiCriterion's interpreter, with default values. Such output is meant to be created by first calling this function, then changing specific values according the interpreter. The point is to always have the same data structure.
-
-
xmc.methodDefs_multiCriterion.flag.
plainFlag
(values, criteria, valuesForCriterion, interpreter)¶ Returns the result of the criterion for the input information provided. This result is a dictionary whose entries match possible choices based on that criterion. See the definition of the interpreter for details.
-
xmc.methodDefs_multiCriterion.interpreter.
interpretAsStoppingFlag
(flag)¶ This is the trivial interpreter for an input which is a boolean answering the question `Must I stop the algorithm now'?
-
xmc.methodDefs_multiCriterion.interpreter.
interpretAsConvergenceAndIterationBounds
(flags)¶ This is the common interpreter which expects three booleans as input: 1. Is convergence achieved? 2. Is the number of iterations greater than the lower bound? 3. Is the number of iterations greater than the uppder bound?
-
xmc.methodDefs_multiCriterion.interpreter.
interpretAsMultipleRequiredConvergencesAndIterationBounds
(flags)¶ This is the common interpreter which expects N booleans as input: + 1 to N-2: convergence flags, all necessary (e.g. split tolerance) + N-1: is the number of iteration high enough? + N: is the number of iteration low enough?
-
xmc.methodDefs_multiCriterion.interpreter.
interpretAsMultipleAlternativeConvergencesAndIterationBounds
(flags)¶ This is the common interpreter which expects N booleans as input: + 1 to N-2: alternative convergence flags (at least one necessary) + N-1: is the number of iteration high enough? + N: is the number of iteration low enough?
-
class
xmc.monoCriterion.
MonoCriterion
(criterion_function, tolerance_values=None)¶ A generic, elementary criterion object.
-
__init__
(criterion_function, tolerance_values=None)¶ criterion_function is a necessary input argument. It is a string of the full name of the criterion function, e.g. xmc.methodDefs_monoCriterion.criterionFunctions.isLowerThanOrEqualTo.
tolerance_value is not always necessary; it depends on the chosen criterion_function.
-
flag
(a, b=None)¶ Takes at most two input arguments (a and b) and re- turns a boolean. Calls criterion(a, b, tolerance). Some criterion functions may need only two input arguments, e.g. (a,tolerance) -> a<tolerance.
-
These are functions for elementary tests in criteria. They accept either two or three argument, depending on the operation. The output argument is a boolean.
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isLowerThanOrEqualTo
(a, b)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isStrictlyLowerThan
(a, b)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isGreaterThanOrEqualTo
(a, b)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isStrictlyGreaterThan
(a, b)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isDistanceLowerThanOrEqualTo
(a, b, tolerance)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isDistanceStrictlyLowerThan
(a, b, tolerance)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isDifferenceLowerThanOrEqualTo
(a, b, tolerance)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isDifferenceStrictlyLowerThan
(a, b, tolerance)¶
-
xmc.methodDefs_monoCriterion.criterionFunctions.
isEqualTo
(a, b)¶
Hierarchy optimisation¶
-
class
xmc.hierarchyOptimiser.
HierarchyOptimiser
(**keywordArgs)¶ This class is in charge of deciding the best hierarchy. It is able to: - provide the best index set (according to its implement meaning of
best) within adjustable constraints;
- give the tolerance splitting associated to that choice (be it fixed or
computed)
- provide the sample numbers per Monte Carlo index that minimise the
cost.
It is assumed here that the hierarchy is optimised with respect to a single QoI.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
inputDictionaryTemplate
()¶
-
toleranceSplitting
(inputDict)¶
-
optimalHierarchy
(inputDict)¶ Compute the best hierarchy by individually computing the best index set and best number of samples
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
singleIndexDoubleSampleNumber
(inputDict, newLevels)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
singleIndexConstantSampleNumber
(inputDict, newLevels)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
singleIndexAdaptive
(inputDict, newLevels)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
multiLevelDoubleAllSamples
(inputDict, newLevels)¶ Returns a list of sample numbers of same length as the number of entries in newLevels. Doubles the number of samples from oldHierarchy if an entry of newLevels exists in oldHierarchy. If not, allocate a default newSampleNumber to the entry.
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
multiLevelConstantSampleNumber
(inputDict, newLevels)¶ Returns a list of sample numbers of same length as the number of levels of deault hierarchy. Keeps constant the number of samples from defaultHierarchy if an entry of newLevels exists in deafultHierarchy. If not, allocate a default newSampleNumber to the entry.
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
multiLevelCostMinimisationTErr
(inputDict, newLevels)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalSampleNumbers.
multiLevelCostMinimisationMSE
(inputDict, newLevels)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalIndexSet.
zeroDimensionSamplesOnly
(inputDict)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalIndexSet.
maximalLevelsBiasModel
(inputDict)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalIndexSet.
incrementLevelsByOne
(inputDict)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalIndexSet.
constantNumberLevels
(inputDict)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalIndexSet.
fullIndexSet
(inputDict)¶
-
xmc.methodDefs_hierarchyOptimiser.optimalIndexSet.
optimalIndexSet
(inputDict)¶
-
xmc.methodDefs_hierarchyOptimiser.updateHierarchySpace.
checkIndexSpaceNotTooBig
()¶ Check if index space is of correct size
Monte Carlo hierarchy¶
-
class
xmc.monteCarloSampler.
MonteCarloSampler
(**keywordArgs)¶ This is a generic class for a multi-something Monte Carlo estimator, which also generates its own samples (albeit indirectly). It does not handle errors and tolerances (this is not an algorithm), but it does handle statistical es- timations, and models fitted across Monte Carlo indices.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
hierarchy
()¶
-
indexEstimation
(coordinate, valueMethodArgs)¶
-
indexCostEstimation
(valueMethodArgs)¶
-
estimation
(assemblerCoordinates=None)¶ For an entry of assembler, collect estimations from all indices corresponding to the same entry of estimatorsForAssembler. Then, pass this list to the assembleEstimation method of that entry of assembler.
-
errorEstimation
(errorEstimatorCoordinates=None)¶ For an entry of errorEstimators, collect assembledEstimations from all assemblers corresponding to the same entry of assemblersForError. Then, pass this list to the error method of that entry of errorEstimators.
-
updateIndexSet
(newHierarchy)¶ Remove entries present in self.indices that are no longer present in newHierarchy. Add entries present in newHierarchy that are not there in self.indices
-
_addMonteCarloIndex
(indicesToAdd)¶ Creates and appends MonteCarloIndex instances corresponding to each entry of indicesToAdd.
-
_removeMonteCarloIndex
(indicesToRemove)¶ Removes elements of self.indices corresponding to each entry of indicesToAdd.
-
updatePredictors
(predictorCoordinates=None)¶ For an entry of qoiPredictor, retrieve all index estimations for the corresponding entry of estimatorsForPredictor and pass to the update method of qoiPredictor. If self.isCostUpdated is True, then do the same for costPredictor as well
-
update
(newHierarchy)¶ Add and process new samples. It first runs updateIndexSet, then calls on each indexs update, then calls updatePredictor.
-
asynchronousPrepareBatches
(newHierarchy)¶ Method setting-up batches. If needed, the serialize method is called, to serialize Kratos Model and Kratos Parameters. Otherwise, serialized objects are passed to new batches, to avoid serializing multiple times. For each batch, an index constructor dictionary is built. This way, local estimators may be computed in parallel and then update global estimators.
-
asynchronousSerializeBatchIndices
(batch, index, solver)¶ Method passing serialized Kratos Model and Kratos Parameters and other KratosSolverWrapper members to new batches. Passing them, we avoid serializing for each batch, and we do only once for the first batch.
-
asynchronousUpdateBatches
()¶ Method updating all relevant members when new batches are created.
-
asynchronousInitialize
(newHierarchy)¶ Method initializing new batches.
-
asynchronousLaunchEpoch
(newHierarchy)¶ Method spawning new batches. It first spawn high index batch realizations, since they have a larger computational cost and normally require a larger time-to-solution.
-
asynchronousFinalize
(batch)¶ Method updating all global estimators with the contributions of the new batch.
-
asynchronousUpdate
(newHierarchy)¶ Method called to initialize all new batches to launch, and to launch them.
-
-
class
xmc.monteCarloIndex.
MonteCarloIndex
(**keywordArgs)¶ This is a generic class for an index of an MXMC method. It handles sample generation and index-specific estimators.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
sampleNumber
()¶
-
indexwiseContribution
(qoiEstimatorCoordinate, qoiEstimatorValueArguements)¶
-
costMoment
(costEstimatorValueArguments)¶ Returns the requested moment estimation of the cost. The input arguments are passed to the cost estimator's relevant method.
-
newSample
()¶ Draw a single sample, i.e. set of correlated outputs from the same event.
See the documentation of SampleGenerator.generate for the data structure.
-
indexSetDimension
()¶ Returns the dimension of the Monte Carlo index set.
-
solverIndices
()¶ Compute list of indices of correlated solvers.
-
_addSamples
(newSamples)¶ Add newsamples to stored samples. This is an internal method with no check on whether sample storage is enabled; be careful when using it.
-
update
(newIndexAndSampleNumbers)¶ Update the Monte Carlo index to a new number of samples. First, decide the number of new samples to be generated, then generate each sample and the associated cost with newSample and pass them to the update methods of qoiEstimator and costEstimator.
-
_updateSubsets
(numberOfSamples)¶ Returns the multi-indices of the subsets – of estimators, samples and times – used for update.
Input arguments: - number of samples: int, number of random events = length of list of samples = number of calls to the sampler.
Output arguments: - multi-index of estimator subsets: nested list - multi-index of sample subsets: nested list - multi-index of time (cost) subsets: nested list
We update each estimator by passing it samples. These samples are passed by groups of fixed size (self.eventGroupSize), corresponding to a subset of events. Each estimator needs only some components. Moreover, the outputs from the solvers may be split in sublists. Therefore, the sets of samples and time received from the sample generator, as well as the list of estimator, need to be re-arranged into subsets suitable for this update process. This method generates the multi-indices used to extract the necessary subsets from these multi-dimensional arrays. We use below the notation A_ij := A_(i,j) := A[i][j] and (A_ij) the multi-dimensional array of these components for all acceptable values of i and j.
## Subsets of times (costs)
The array of times T := (T_ij) is bidimensional and its dimensions are (event, solver): T_ij is the computational time returned by the j-th solver for the i-th random event. We want to create from T a tridimensional array gT = (gT_ijk) of dimensions (subset,event,solver) such that gT_i is the i-th subset passed to the cost estimator. This method creates (among others) the multi-index iT such that gT = T_iT; iT is bidimensional. The notation gT = T_iT is to be understood as: iT_ijk is an index of appropriate dimension for T (i.e. 2) and gT_ijk = T_(iT_ijk) for all acceptable values of i, j and k. E.g. if gT_(4,1,3) = T_(2,18), then iT_(4,1,3) = (2,18).
## Subsets of samples and estimators
# Case: solver outputs are not split
If solver outputs are not split, the array of samples S := (S_ijk) is tridimensional. Its dimensions are (event,solver,component): S_ijk is the k-th component (or output, or random variable) of the j-th solver for the i-th random event. We wish to create from the array of estimators E a unidimensional array (gE_i), and from S a tetradimensional array gS = (gS_ijklm), so that gS_ij is the j-th subset to be passed to gE_i. The dimensions of gS are (subset,component,event,solver), and those of gE are (component). In this case, gE = E; not so if solver outputs are split (see below). This method creates the multi-indices iE, iS such that gS = S_iS and gE = E_iE – in addition to the multi-index of the time subset (see above). Since the last dimension (i.e. solver) is unchanged from S to gS, iS is only tridimensional. Obviously, iE is unidimensional.
# Case: solver outputs are split
If solver outputs are split, S := (S_ijkl) is quadridimensional, and its dimensions are (event,solver,split,component). However, in this case the sample subsets are not passed to the update method of a StatisticalEstimator, but to another method (updatePartialQoiEstimators_Task). Unlike StatisticalEstimator.update, who expects an array of dimensions (event,solver), this method expects an array of dimensions (event,split,solver,component). Therefore, we can ignore the component dimension in arranging the subsets of samples: we consider S := (S_ijk) of dimensions (event, solver, split). Unlike the previous case, gE is now a bidimensional array of dimensions (splits,component): gE_ij is the estimator receiving the j-th component of the i-th split. We define iE and iS as previously: gE = E_iE and gS = S_iS.
-
numberOfSolvers
()¶ Returns the numberofsolvers in the sample generator<
-
sampleDimension
(*args)¶ Returns the number of outputs from the sampler. This is the number of outputs from a single solver of the sampler, with no information regarding splitting, if any.
Input arguments: same as SampleGenerator.sampleDimension
Output argument: - sample dimension: integer. See SampleGenerator.sampleDimension
- Return type
int
-
sampleSplitSizes
(*args)¶ Returns the lengths of (sub-)lists into which samples are split, if they are; returns None if they are not.
Input arguments: same as SampleGenerator.samplesSplitSizes.
Output argument: list of integers, or None. If list, it is [len(subSample) for subSample in sample], where sample is the output from a single solver (after processing).
-
areSamplesSplit
(*args)¶ Returns True if samples are split into future (sub-)lists, False otherwise. Accepts the same input arguments as SampleGenerator.areSamplesSplit.
Input arguments: same as SampleGenerator.areSamplesSplit.
Output argument: boolean, True is samples are split and False otherwise.
- Return type
bool
-
Global estimation assemblers¶
-
class
xmc.estimationAssembler.
EstimationAssembler
(**keywordArgs)¶ The class implements what is needed to produce the final estimations desired from the Monte Carlo method, given the statistical estimations provided by each index.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
-
xmc.methodDefs_estimationAssembler.assembleEstimation.
assembleValue
(hierarchy, indexEstimations)¶ Compute a simple sum over all index-wise contributions
-
xmc.methodDefs_estimationAssembler.assembleEstimation.
assembleBias
(hierarchy, indexEstimations)¶ Based on hierarchy, compute the sum over the index estimations that lie on the boundary of the index space as an estimate for the bias. Note - The procedure used in this function to find the boundary of a given hierarchy assumes that the index set is covex
-
xmc.methodDefs_estimationAssembler.assembleEstimation.
assembleStatisticalError
(hierarchy, indexEstimations)¶ Add together the variance over all indices to estimate global variance
-
xmc.methodDefs_estimationAssembler.assembleEstimation.
assembleInterpolationError
(hierarchy, indexEstimations)¶
-
class
xmc.errorEstimator.
ErrorEstimator
(**keywordArgs)¶ This class is responsible for the error estimation of the multi-something MC estimation.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
error
(*args)¶
-
-
xmc.methodDefs_errorEstimator.errorEstimation.
errorEstimationStatError
(cdfValue, globalEstimations)¶ Accept the summation over the variances divided by number of samples over all levels and return its square root as the statistical error
-
xmc.methodDefs_errorEstimator.errorEstimation.
errorEstimationTerr
(cdfValue, globalEstimations)¶ Return total error given moments/variances
-
xmc.methodDefs_errorEstimator.errorEstimation.
errorEstimationMSE
(_, globalEstimations)¶ Return mean squared error given moments/variances
Predictors¶
-
class
xmc.modelEstimator.
ModelEstimator
(**keywordArgs)¶ Generic class for model with a known analytical expression, dependent on parameters to be fitted on data. A data point is an observed (antecedent ; image) couple. Some methods descriptions. - value: return model evaluation at given point. - valueForParameters: Returns the model evaluation at given point with given parameter values. - fitDistance: computes the distance between two sets of data points; this is the distance used to fit the model to new observations.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
value
(antecedentList)¶
-
update
(*data)¶
-
-
xmc.methodDefs_modelEstimator.update.
updatePredictorGeometric
(data)¶ Fit a model of the form y = C*exp(-w.x) where x is the antecedent and y is the corresponding image.
-
xmc.methodDefs_modelEstimator.valueForParameters.
geometricModel
(parameters, point)¶ Compute y = C*exp(-w.x) where C = parameters[0] and w = parameters[1:]
-
class
xmc.bayesianEstimator.
BayesianEstimator
(**keywordArgs)¶ This class is for Bayesian estimation, where an a priori estimation is updated by observations.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
blend
(newLevels, inputDict)¶
-
-
xmc.methodDefs_bayesianEstimator.blend.
noBlending
(parameters, newLevels, inputDict)¶ Return list of variances based purely on the variance model
-
xmc.methodDefs_bayesianEstimator.blend.
bayesianUpdate
(parameters, newLevels, inputDict)¶ Use bayesian estimation to blend the variance estimations and variance model such that the value of the estimation is used if there are a large number of samples and the value of the model is used if there is a small number of samples.
Statistical estimators¶
-
class
xmc.statisticalEstimator.
StatisticalEstimator
(**keywordArgs)¶ Estimators of statistics of a given random variable. The estimations are based on a set of realisations of (and possibly some prior assumptions on) that random variable. The estimator is able to update its estimation when given new observations.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
reset
()¶ Forgets about all observations and return to its initial state. Useful when samples are not recycled (i.e. the update process actually discards the previous observations).
-
_addData
(newData)¶ Add given value(s) to data (if self.isDataStored) and increment self._sampleCounter
-
value
()¶ Returns the current estimation.
-
update
(newData)¶ Update the estimation with new observations.
-
sampleNumber
()¶ Returns the current number of samples known to the estimator.
-
-
class
xmc.momentEstimator.
MomentEstimator
(**keywordArgs)¶ Bases:
xmc.statisticalEstimator.StatisticalEstimator
This class estimates raw and central moments up to a given order (including expectations, obviously). These estimations are computed using h-statistics, so we store (and update) the power sums from which they are computed.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
update
(samples)¶ Method updating the power sums given new samples.
Inputs: - self: an instance of the class - samples: list containing new samples (dimension at least 3)
-
value
(order, isCentral=None, isErrorEstimationRequested=False)¶ Method returning a specific order moment, its error estimate or both.
Inputs: - self: an instance of the class - order: order of the moment - isCentral: setting if the working moment is central or raw - computeErrorEstimation: setting if computing or not the errorEstimation
-
errorEstimation
(order, isCentral)¶ Method returning the error estimation of working moment from stored power sums.
Inputs: - self: an instance of the class - order: order of the moment
-
_computeCentralMoment
(order, isErrorEstimationRequested)¶ Method returning the central moment of working moment from stored power sums.
Inputs: - self: an instance of the class - order: order of the moment
-
_computeRawMoment
(order, isErrorEstimationRequested)¶ Compute raw moment of requested order from stored power sums.
-
powerSumsOrder
()¶ Returns the maximum order to which the power sums are stored to compute estimation and error of moments up to order self.order.
-
reset
()¶ Reset the value of power sums to zero. Reset the number of samples
-
-
class
xmc.momentEstimator.
CombinedMomentEstimator
(**keywordArgs)¶ Bases:
xmc.momentEstimator.MomentEstimator
This class estimates raw and central moments up to a given order (including expectations, obviously). These estimations are computed using h-statistics, so we store (and update) the power sums from which they are computed. The computed statistics are the composition of sampling and time statistical processes. The first difference with respect to MomentEstimator is the update method, which is supposed to receive power sums instead of sample values. Another difference with respect to MomentEstimator is that, for multi-index moment estimators (as Multilevel Monte Carlo), we can compute only power sums of order a of the form (Sa0,S0a), which requires to estimate h-statistics exploiting equations as eq. (4) of [S. Krumscheid et al. / Journal ofComputational Physics 414 (2020) 109466].
Methods: - update: method updating power sums and number of realizations.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
update
(samples)¶ Method updating the power sums and number of realizations, given new samples. The new data is passed through as a list containing power sums, up to a given order, and the number of samples from which it was computed. For example, let us discretise a scalar time signal u into M time steps to obtain a vector U of length M. A power sum of order a is computed as S_a = sum_{i=1}^{M} U[i]^a. Samples will have the following shape: [[[[S1], [S2], M]]] , where S1 and S2 are power sums of order one and two, respectively, and M is the total number of terms in the power sum (equivalent to number of time steps times realizations).
Inputs: - self: an instance of the class - samples: list containing new sample power sums.
-
-
class
xmc.momentEstimator.
MultiMomentEstimator
(**keywordArgs)¶ Bases:
xmc.statisticalEstimator.StatisticalEstimator
This class estimates statistical moments of multi-valued real random variables, from independent realisations. It inherits from StatisticalEstimator; the description below describes only its differences with it.
order : int
Mandatory. All statistical moments up to this order are estimated. It determines the power of the sums of realisations to be computed; see MultiMomentEstimator.powerSumsOrder.
_referenceComponent : int
Default: 0. MultiMomentEstimator.value must return a float. Of all the components of the multi-valued random variable, MultiMomentEstimator.value will return the requested statistic of the one of index MultiMomentEstimator._referenceComponent.
_isEstimationParallel : bool
Default: True. Whether the estimation process is to run in parallel. If not, power sums are synchronised synchronised before computing the estimation.
_isEstimationParallel : bool
Default: True. Whether the update process is to run in parallel. If not, samples are synchronised before the update.
_variableDimension : int
Mandatory if self._isParallel; otherwise, automatically assigned. The dimension of the image space (i.e. number of scalar components) of the multi-valued random variable. If self._isParallel, then it must be supplied to the constructor; if not, it can be guessed from the first samples.
_indexSetDimension : int
Automatically assigned. Local dimension of the Monte Carlo index set, i.e. the number of solvers sending samples to this estimator is 2^self._indexSetDimenion. It is not supplied at construction but deduced from the first samples.
_powerSums : DefaultDict[str, List[float]]
Default: see self._initialisePowerSums. This dictionary stores sums of powers of sample values, and can be supplied at construction. The summation and exponentiation are component-wise, so that the values in each entry is an array of dimension self._variableDimension. The keys are strings of integers, of length 2^self._indexSetDimension; each integer is a power of a multivariate sum or difference of multi-valued random variables across solver levels. Examples: - if self._indexSetDimension is 0: there is one solver level (multi-valued random variable X), so the keys are of the form 'p' where p is the integer power of the sums stored in this entry. self._powerSums['p'] = sum_{e in events} X(e)^p - if self._indexSetDimension is 1: there are two solver levels (multi-valued random variables X and Y), so the keys are of the form 'pq' where p and q are respectively the integer powers of the sum and difference of multi-valued random variables across solver levels. self._powerSums['pq'] = sum_{e in events} (X(e)+Y(e))^p (X(e)-Y(e))^q
_powerSumsUpdater : callable
Default: xmc.methodDefs_momentEstimator.powerSums.updatedPowerSums. Function which updates power sums with additional sample values; called by self.update (see also its documentation). It must abide by the following specifications; see default function for an example. Input arguments: 1. dictionary of power sums: DefaultDict[str, List[float]] 2. array of samples: List[List[List[float]]] Output arguments: 1. dictionary of power sums (updated): DefaultDict[str, List[float]]
_estimation: callable
Default: xmc.methodDefs_momentEstimator.hStatistics.hStatistics. Function which computes the desired statistics from the available power sums; called by self.multiValue (see also its documentation). It must abide by the following specifications; see default function for an example. Input arguments: 1. components of the multi-valued random variable to consider: Union[int, slice, Iterable[int]] 2. dictionary of power sums: DefaultDict[str, List[float]]. 3. number of random events having contributed to the power sums: int. 4. statistical order of the moment estimated: int. 5. whether the error on the estimation is requested, instead of the estimation: bool. 6. whether the moment estimated is central or raw: bool. Output arguments: 1. statistics for every component considered: Union[List[float],float] (float if single component)
-
__init__
(**keywordArgs)¶ Keyword arguments required: - order
Keyword arguments required for parallel update: - variableDimension
Optional keyword arguments: - referenceComponent - isUpdateParallel - isEstimationParallel - isParallel - estimationFunction - powerSumsUpdater - powerSums
Notes: - isParallel: value is assigned to both self._isUpdateParallel and self._isEstimationParallel. If the keyword arguments isUpdateParallel or isEstimationParallel are given, they take priority.
-
_initialisePowerSums
()¶
-
_guessIndexSetDimension
(singleSample=None)¶ Return the dimension of the index set considered by this estimator. If unknown, returns None. Relies on the keys of dictionary self._powerSums being of length dimension+1
Alternatively: if input argument is provided, guess index set dimension from it.
Input argument: - singleSample (optional): set of correlated samples (list) for a single event. Expected to be of length 2^dimension
Output argument: integer, or None if no dimension could be computed.
- Return type
int
-
_guessVariableDimension
(singleSample=None)¶ Return the dimension of the space of values of the multi-valued random variable.
Input argument: - singleSample (optional): set of correlated samples (list) for a single event. Every element of that list is expected to be a list of length equal to the dimension of the multi-valued random variable.
Output argument: integer, or None if no dimension could be computed
An error will be raised if the sample is a Future object. In that case, the variable dimension must be given. Use the variableDimension keyword at instantiation..
- Return type
int
-
update
(samples)¶ Method to update the power sums from new samples.
Input arguments: 1. samples: nested list of depth ≥ 3. (S_ijk) with: i the random event; j the solver level; k the component of the multi-valued random variable.
-
value
(order, isError=False, isCentral=None)¶ Function returning the current estimation of a statistical moment, or its error estimation. This function only returns estimations for the component indicated by self._referenceComponent.
Inputs: same as MultiMomentEstimator.multiValue, except for the 'component', which is not an input argument of MultiMomentEstimator.value.
- Return type
float
-
multiValue
(order, isError=False, component=slice(None, None, None), isCentral=None)¶ Function returning a specific order moment, its error estimate or both.
Inputs: 1. order: order of the moment 2. isError: whether to return the error of the estimation (as opposed to the estimation itself). 3. component: index (or indices) of components of the multi-valued random variable to consider 4. isCentral: whether the moment is central (as opposed to raw)
Outputs: float if component is an integer, list of floats otherwise.
- Return type
Union
[float
,List
[float
]]
-
class
xmc.momentEstimator.
MultiCombinedMomentEstimator
(**keywordArgs)¶ Bases:
xmc.momentEstimator.MultiMomentEstimator
This class estimates contributions to statistical 'combined' moments of multi-valued real random variables, from independent realisations. These contribution take the form of h-statistics, or differences of such. They are computed from monovariate power sums. The computed statistics are the composition of sampling and time statistical processes. This class inherits from MultiMomentEstimator; the description below describes only its differences with it.
_powerSums: Union[DefaultDict[str, PowerSumsDict], PowerSumsDict] with
PowerSumsDict = DefaultDict[str, Union[float, List[float]]] See the documentation of this attribute of the parent class. The only difference is if self._indexSetDimension is 1: since only monovariate power sums are used, self._powerSums then has exactly two entries. Each of these entries is a dictionary of monovariate power sums as described for the parent class. These entries have respective keys 'upper' and 'lower' for the first and second solver levels in the array received by self._update. See self._initialisePowerSums for details.
_powerSumsUpdater : callable
Default: xmc.methodDefs_momentEstimator.powerSums.addPowerSumsAndCounter. Function which updates current power sums with additional power sums, and likewise for the sample count; called by self.update (see also its documentation). It must abide by the following specifications; see default function for an example. Input arguments: 1. dictionary of power sums: DefaultDict[str, List[float]] 2. current counter of samples: int. 3. array of power sums: List[List[Union[List[float], int]]] 4. keys of power sum (in argument 1) in the order of argument 3 : tuple (Default: ('1','2', ...)) Output arguments: 1. dictionary of power sums (updated): DefaultDict[str, List[float]]
_estimation: callable
Default: xmc.methodDefs_momentEstimator.hStatistics.hStatisticsMonoPowerSums. See the description of this attribute of the parent class. The only difference is that the power sums are always monovariate, hence the different default value.
-
__init__
(**keywordArgs)¶ Keyword arguments required: - order
Keyword arguments required for parallel update: - variableDimension
Optional keyword arguments: - referenceComponent - isUpdateParallel - isEstimationParallel - isParallel - estimationFunction - powerSumsUpdater - powerSums
Notes: - isParallel: value is assigned to both self._isUpdateParallel and self._isEstimationParallel. If the keyword arguments isUpdateParallel or isEstimationParallel are given, they take priority.
-
update
(samples)¶ Method updating the power sums and number of realizations, given new samples. The new data is passed through as a list containing power sums, up to a given order, and the number of samples from which it was computed. For example, let us discretise a scalar time signal u into M time steps to obtain a vector U of length M. A power sum of order a is computed as S_a = sum_{i=1}^{M} U[i]^a. Samples will have the following shape: [[[[S1], [S2], M]]] , where S1 and S2 are power sums of order one and two, respectively, and M is the total number of terms in the power sum (equivalent to number of time steps times realizations).
Inputs: - self: an instance of the class - samples: nest list of dimension 4. Axes are: event, MC index, power, component.
-
_guessOrder
(singleSample)¶ Guesses the statistical order (value of MultiCombinedMomentEstimator.order) from the power sums received. This is done under the assumption that the order of power sums must be twice the statistical order.
- Return type
int
-
_initialisePowerSums
()¶
-
xmc.methodDefs_momentEstimator.hStatistics.
hStatistics
(component, powerSums, numberOfSamples, order, isError, isCentral)¶ This is the generic function to compute h-statistics using the hStatistics module.
It extracts the relevant power sums for the selected component(s), and calls the right function for the requested h-statistics. It abides by the specifications laid out in the documentation of xmc.momentEstimator.MultiMomentEstimator.
- Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
hStatisticsMonoPowerSums
(component, powerSums, numberOfSamples, order, isError, isCentral)¶ Compute h-statistics from monovariate power sums. Obviously, this makes a difference only for differences of h-statistics (e.g. for MLMC). In this case, error estimation is not supported. This is useful for combined moment estimators, since they store only monovariate power sums.
The input and output arguments are the same as for hStatistics, which it calls.
- Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
rawMomentOrder1Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
rawMomentErrorOrder1Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentOrder2Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentErrorOrder2Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentOrder3Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentErrorOrder3Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentOrder4Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentErrorOrder4Dimension0
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
rawMomentOrder1Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
rawMomentErrorOrder1Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentOrder2Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentErrorOrder2Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentOrder3Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentErrorOrder3Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentOrder4Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
centralMomentErrorOrder4Dimension1
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
_uglyD0O3E
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
_uglyD0O4E
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
_uglyD1O3E
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
_uglyD1O4V
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.hStatistics.
_uglyD1O4E
(powerSums, numberOfSamples)¶ - Return type
Union
[float
,List
[float
]]
-
xmc.methodDefs_momentEstimator.powerSums.
updatedPowerSums
(powerSums, samples)¶ Increments power sums from an array of samples.
Supports multi-valued random variables, and Monte Carlo index sets of dimension up to 1.
Input arguments: - powerSums: dictionary following the format of MomentEstimator._powerSums - samples: array of samples (see below).
Output arguments: - powerSums: same as input, with updated values.
The array of samples is expected to be an array of dimension at least three, viz. 1. random event; 2. Monte Carlo index (e.g. solver, model or level); 3. (and beyond) random variable component.
Any dimension beyond the third will be collapsed into this latter one.
- Return type
DefaultDict
[str
,Union
[float
,List
[float
]]]
-
xmc.methodDefs_momentEstimator.powerSums.
addPowerSumsAndCounter
(psDict, counter, psArray, keyOrder=None)¶ Increment existing dictionary of power sums and sample counter for an array of such. This function is meant as a defintion for MultiCombinedPowerSums._powerSumsUpdater.
Input argument: - psDict: dictionary of power sums to be incremented; in the format of MultiCombinedPowerSums._powerSums - psArray: multidimensional array of power sums and sample count to add; in the format expected by MultiCombinedMomentEstimator.update. - counter: sample counter to be incremented. - keyOrder: (optional) keys of entries of power sums in psDict, in their order in psArray. Default value: ('1', '2', ..., 'p'), with 'p' the number of distinct power sums in psArray.
Output arguments: - psDict: same format as input (depends on dimension of the MC index set) - counter: same format as input
See the documentation of MultiCombinedPowerSums for more details.
- Return type
(typing.Union[typing.DefaultDict[str, typing.Union[float, typing.List[float]]], typing.DefaultDict[str, typing.DefaultDict[str, typing.Union[float, typing.List[float]]]]], <class 'int'>)
-
xmc.methodDefs_momentEstimator.updatePowerSums.
updatePowerSumsOrder2Dimension0
(samples, power_sum_1, power_sum_2)¶
-
xmc.methodDefs_momentEstimator.updatePowerSums.
updatePowerSumsOrder4Dimension0
(samples, power_sum_1, power_sum_2, power_sum_3, power_sum_4)¶
-
xmc.methodDefs_momentEstimator.updatePowerSums.
updatePowerSumsOrder10Dimension0
(samples, power_sum_1, power_sum_2, power_sum_3, power_sum_4, power_sum_5, power_sum_6, power_sum_7, power_sum_8, power_sum_9, power_sum_10)¶
-
xmc.methodDefs_momentEstimator.updatePowerSums.
updatePowerSumsOrder2Dimension1
(samples, power_sum_10, power_sum_01, power_sum_20, power_sum_11, power_sum_02)¶
-
xmc.methodDefs_momentEstimator.updateCombinedPowerSums.
updatePowerSumsOrder1Dimension0
()¶
-
xmc.methodDefs_momentEstimator.updateCombinedPowerSums.
updatePowerSumsOrder2Dimension0
(old_sample_counter, samples, power_sum_1, power_sum_2)¶
-
xmc.methodDefs_momentEstimator.updateCombinedPowerSums.
updatePowerSumsOrder10Dimension0
(old_sample_counter, samples, power_sum_1, power_sum_2, power_sum_3, power_sum_4, power_sum_5, power_sum_6, power_sum_7, power_sum_8, power_sum_9, power_sum_10)¶
-
xmc.methodDefs_momentEstimator.updateCombinedPowerSums.
updatePowerSumsOrder2Dimension1
(old_sample_counter, samples, power_sum_upper_1, power_sum_lower_1, power_sum_upper_2, power_sum_lower_2)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
centralMomentWrapper
(dimension, order, *args)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
centralCombinedMomentWrapper
(dimension, order, *args)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderOneDimensionZero
(power_sum_1, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderTwoDimensionZero
(power_sum_1, power_sum_2, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderThreeDimensionZero
(power_sum_1, power_sum_2, power_sum_3, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderFourDimensionZero
(power_sum_1, power_sum_2, power_sum_3, power_sum_4, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderTwoDimensionZeroBiased
(power_sum_1, power_sum_2, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderOneDimensionOne
(power_sum_10, power_sum_01, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeCentralMoments.
computeCentralMomentsOrderTwoDimensionOne
(power_sum_10, power_sum_01, power_sum_20, power_sum_11, power_sum_02, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeErrorEstimation.
centralMomentErrorWrapper
(dimension, order, *args)¶
-
xmc.methodDefs_momentEstimator.computeErrorEstimation.
centralCombinedMomentErrorWrapper
(dimension, order, *args)¶
-
xmc.methodDefs_momentEstimator.computeErrorEstimation.
computeCentralMomentsErrorEstimationOrderOneDimensionZero
(power_sum_1, power_sum_2, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeErrorEstimation.
computeCentralMomentsErrorEstimationOrderTwoDimensionZero
(power_sum_1, power_sum_2, power_sum_3, power_sum_4, number_samples)¶
-
xmc.methodDefs_momentEstimator.computeErrorEstimation.
computeCentralMomentsErrorEstimationOrderOneDimensionOne
(power_sum_10, power_sum_01, power_sum_20, power_sum_11, power_sum_02, number_samples)¶
Random sampling¶
-
class
xmc.randomGeneratorWrapper.
RandomGeneratorWrapper
(**keywordArgs)¶ This is a class is responsible for the generation of random inputs: scalars, fields, whatever. It is very similar to SolverWrapper: it handles calls to whatever tool (most likely external) is chosen for this task. Therefore, this is a barebone template from which tool-specific implementation will be fashioned, either by inheritance or pointers.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
realisation
()¶
-
-
xmc.methodDefs_randomGeneratorWrapper.generator.
normal
(mean, deviation, shape=None)¶ Draws samples from the normal distribution defined by arguments.
Input parameters: - mean: float, mean of the normal distribution. - deviation: float, standard deviation of the normal distribution. - shape (optional): tuple, shape of the desired output. If omitted, the output is a float.
Output: float if shape is omitted. Otherwise, nested list of required shape, whose every element is a float. These floats are drawn independantly from the normal distribution specified by the input arguments.
Note: differences with returnStandardNormal. - random number generator: this method uses numpy.random.Generator.normal, whereas returnStandardNormal uses numpy.random.normal. - output type and shape: this method returns either a float or a nested list of floats of arbitrary shape, whereas returnStandardNormal returns a list containing a single float.
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnStandardNormal
(*args)¶ Return a normally-distributed variable with given parameters
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnStandardNormalandUniform
(*args)¶ Return one standard normal random variable and one uniformly distributed random variable in the interval [0,1]
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnBeta
(*args)¶ Return a beta distributed random variables with given parameters
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnUniform
(*args)¶ Return one uniformly distributed random variable in the interval [0,1]
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnRayleigh
(*args)¶ Return one rayleigh random variable with given parameters use numpy.random.rayleigh(scale, size) already returns as a numpy array
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnRayleighAndUniform
(*args)¶ Return one rayleigh and uniformly distributed random variables with given parameters
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnUniformAndNormal
(*args)¶ Return one integer uniformly distributed random variable and one normal random variable
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnUniformAndTwoNormal
(*args)¶ Return one integer uniformly distributed random variable and two normal random variables
-
xmc.methodDefs_randomGeneratorWrapper.generator.
returnIntegerUniform
(*args)¶ Return one integer uniformly distributed random variable
-
class
xmc.sampleGenerator.
SampleGenerator
(**keywordArgs)¶ Class used for sample generation. Interacts with a sample generator and solvers to generate a set of correlated samples.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
qoiFromRaw
(rawSolutions)¶ Function that generates additionally QoI from given QoI. For example, if one wished to study the statisics of the product or sum or difference of two or more QoI
-
randomEvent
()¶
-
rawSolutions
(randomEvent)¶ Return the array of raw solutions and computation time produced by the solvers.
Input argument: random event, as expect by the solve method of the solvers.
Output argument: - solutions: list of solver outputs; solutions[j] is the solution returned by self.solvers[j]. - times: list of computational times; times[j] is the time returned by self.solvers[j].
See the documentation of SampleGenerator.generate for more details on the different possible data structures of solutions.
-
generate
(randomEvent=None)¶ Generate new single sample. Calls randomEvent (if no evaluation point is provided), then rawSolutions, then pass the outputs to qoiFromRaw and return its output.
Input arguments: - randomEvent (optional): random event to be passed to solvers. By default, the solvers receive a random event generated internally by the method SampleGenerator.randomEvent.
Output arguments: - sample: list of length the number of solvers. Each element of that list is either a list of length SampleGenerator.sampleDimension() or a list of future lists of length len(SampleGenerator.sampleSplitSizes()). For the length of these future lists, see SampleGenerator.sampleSplitSizes. - time: list of floats of length the number of solvers. Measurement of the cost of each solver call.
Let us note S = sample, S_j = sample[j], etc.. All such objects are lists. S is a nested list of length len(SampleGenerator.solvers) and depth at least 2. Elements of S_j may be future objects.
If solver outputs are not split, (i.e. not SampleGenerator.areSamplesSplit()), S_j is of length SampleGenerator.sampleDimension(). Example. If len(SampleGenerator.solvers) = 2 and SampleGenerator.sampleDimension() = 3: S_1 = [s1_o1, s1_o2, s1_o3] (s: solver, o: output); idem for S_2.
If solver outputs are split, (i.e. SampleGenerator.areSamplesSplit()), S_j is of length len(SampleGenerator.sampleSplitSizes()) and S_jk is of length SampleGenerator.sampleSplitSizes()[k]. Example. If len(SampleGenerator.solvers) = 2 and SampleGenerator.sampleSplitSizes() = [3, 2]: S_1 = [ S_(1,1), S_(1,2) ] = [ [s1_o1, s1_o2, s1_o3], [s1_o4, s1_o5] ] (s: solver, o: output); idem for solver S_2.
-
sampleDimension
(*args)¶ Returns the number of outputs of a single solver, as integer. The time measurement is not counted. This gives no information on whether the samples are split into future (sub-lists). Use SampleGenerator.areSamplesSplit for this.
- Return type
int
-
sampleSplitSizes
(*args)¶ Returns the size of lists in which solver outputs are split (None if they are not). Denoting s the first output argument of self.solvers[any].solve, splitSizes[j]=len(s[j]).
Input arguments: same as SampleGenerator._solverOutputDimensions.
Output arguments: - splitSizes: list of integers if SampleGenerator.areSamplesSplit(); None otherwise.
-
areSamplesSplit
(*args)¶ Returns True if samples are split into future (sub-)lists, False otherwise. Accepts the same input arguments as _solverOutputDimensions.
- Return type
bool
-
_solverOutputDimension
(solverChoice=0)¶ Returns the dimension of the output of a solver. See SolverWrapper documentation for details.
Input arguments: - solverChoice (optional): index of solver from which to fetch. Default is 0.
-
-
xmc.methodDefs_sampleGenerator.qoiProcessor.
doNothing
(rawSolutions)¶
Solvers interfaces¶
-
class
xmc.solverWrapper.
SolverWrapper
(**keywordArgs)¶ Class used to call the relevant solver. Since this is very problem- and tool- specific, sub-classes (inheriting from that one) should be implement for each solver. Therefore, this is a template of what SampleGenerator objects need.
Attributes: - solverWrapperIndex: list, place in the index set. Useful for the constructor. - outputDimension: integer or list of integer. If integer, equals to len(sample), where sample is the first output argument of self.solve(). If list of integers, then it means that samples is split in future lists, and outputDimension is [len(subSample) for subSample in sample].
Methods: - solve. See its documentation.
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
solve
(randomEvent)¶ Input: random event in the form expect by the solver wrapper.
Outputs: 1. sample: list of structure respecting self.outputDimension 2. cost: float, measure of time to generate the sample.
-
-
class
xmc.classDefs_solverWrapper.singleLevelRNGSolverWrapper.
SingleLevelRNGSolverWrapper
(**kwargs)¶ Bases:
xmc.solverWrapper.SolverWrapper
solveWrapper type whose solve method accepts a random number and returns it in a list of size self.ouptutDimension
Constructor arguments - solverWrapperIndex - Index-space position of the solverWrapper instance
-
__init__
(**kwargs)¶ Initialize self. See help(type(self)) for accurate signature.
-
_drawSample_Task
(**kwargs)¶
-
solve
(randomInput)¶ Input: random event in the form expect by the solver wrapper.
Outputs: 1. sample: list of structure respecting self.outputDimension 2. cost: float, measure of time to generate the sample.
-
-
class
xmc.classDefs_solverWrapper.multiLevelRNGSolverWrapper.
MultiLevelRNGSolverWrapper
(**keywordArgs)¶ Bases:
xmc.solverWrapper.SolverWrapper
solveWrapper type whose solve method accepts a standard normal random number and scales it by an index-dependent mean and variance
Constructor arguements - solverWrapperIndex - Index-space position of the solverWrapper instance rates - list of length 4 whose first and last two entries are respectively the constant and rate required to compute the mean and variance using the equation mean (or) variance = constant * exp (-rate*level)
-
__init__
(**keywordArgs)¶ Initialize self. See help(type(self)) for accurate signature.
-
_drawSample_Task
(**kwargs)¶
-
solve
(randomInput)¶ Input: random event in the form expect by the solver wrapper.
Outputs: 1. sample: list of structure respecting self.outputDimension 2. cost: float, measure of time to generate the sample.
-
Miscellaneous functions¶
-
xmc.tools.
dynamicImport
(fullName)¶ This function returns the requested attribute of a module; it manages the necessary imports. The input argument is the full name as a string, e.g. 'module.package.method'. Alternatively, the input argument may be a list of such full names. The output will then be a list of same size containing the requested attributes in the same order. This only works for modern Python (at least Python > 3).
-
xmc.tools.
instantiateObject
(fullName, **kwargs)¶
-
xmc.tools.
doNothing
(*args, **kwArgs)¶ This is a function that does nothing. Useful to hold a place without taking any action.
-
xmc.tools.
returnInput
(*args)¶ This is a function that returns exactly what is passed to it. Useful to hold a place in the dataflow without modifying the data.
-
xmc.tools.
getUnionAndMap
(listOfLists)¶ Given a list of lists, computes union of all lists, as well as map from list of lists to union lists.
-
xmc.tools.
splitOneListIntoTwo
(inputList)¶ Function that takes a list, each of whose entries is a list containing two entries, and assembles two separate lists, one for each of these entries
-
xmc.tools.
mergeTwoListsIntoOne
(list1, list2)¶ Function that takes two lists and assembles them into one list such that each entry of the final list is a list containing the two corresponding entries of the two input lists
-
xmc.tools.
strictlyPositiveBoundaryBooleans
(indexSet)¶ Accepts a convex set of non-negative multi-indices and returns a list of booleans of the same size that states whether the corresponding index lies on the surface of the set where none of the indices have zero-values.
-
xmc.tools.
summation
(*args)¶
-
xmc.tools.
splitList
(listToSplit, num_sublists=1)¶ Function that splits a list into a given number of sublists Reference: https://stackoverflow.com/questions/752308/split-list-into-smaller-lists-split-in-half
-
xmc.tools.
normalInverseCDF
(y0)¶ Computes inverse of normal distribution function (percent point function) Sources: https://github.com/scipy/scipy/blob/59347ae8b86bcc92c339efe213128f64ab6df98c/scipy/special/cephes/ndtri.c, https://stackoverflow.com/questions/41338539/how-to-calculate-a-normal-distribution-percent-point-function-in-python References: Pauli Virtanen, Ralf Gommers, Travis E. Oliphant, Matt Haberland, Tyler Reddy, David Cournapeau, Evgeni Burovski, Pearu Peterson, Warren Weckesser, Jonathan Bright, Stéfan J. van der Walt, Matthew Brett, Joshua Wilson, K. Jarrod Millman, Nikolay Mayorov, Andrew R. J. Nelson, Eric Jones, Robert Kern, Eric Larson, CJ Carey, İlhan Polat, Yu Feng, Eric W. Moore, Jake VanderPlas, Denis Laxalde, Josef Perktold, Robert Cimrman, Ian Henriksen, E.A. Quintero, Charles R Harris, Anne M. Archibald, Antônio H. Ribeiro, Fabian Pedregosa, Paul van Mulbregt, and SciPy 1.0 Contributors. (2020) SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python. Nature Methods, in press.
-
xmc.tools.
doubleFactorial
(n)¶ Returns double factorial of an integer.
-
xmc.tools.
binomial
(n, k)¶ Returns binomial coefficient (n,k) or 'n choose k' with n ≥ k.
-
xmc.tools.
flatten
(sequence)¶ Returns a generator of all non-iterable elements in the input sequence. This can be used to flatten nested heterogeneous arrays. It may crash on deeply-nested array.
Framework for distributed computing¶
Defines the framework to compute in a distributed environment.
This module provides the definition of the API used throughout the library. There are several implementation of this API. The default one is standalone, without dependencies, and provides no parallelisation. The other ones allow distributed computing, but depend on external libraries. Edit this module to choose a different implementation.
-
class
xmc.distributedEnvironmentFramework.
ExaquteTask
(*args, **kwargs)¶ -
__init__
(*args, **kwargs)¶ Initialize self. See help(type(self)) for accurate signature.
-
-
xmc.distributedEnvironmentFramework.
barrier
()¶
-
xmc.distributedEnvironmentFramework.
compute
(obj)¶
-
class
xmc.distributedEnvironmentFramework.
constraint
(*args, **kwargs)¶ -
__init__
(*args, **kwargs)¶ Initialize self. See help(type(self)) for accurate signature.
-
-
xmc.distributedEnvironmentFramework.
delete_file
(file_path)¶
-
xmc.distributedEnvironmentFramework.
delete_object
(*objs)¶
-
xmc.distributedEnvironmentFramework.
get_value_from_remote
(obj)¶