Internal

This chapter documents internal functions that are typically only for use of the developers of a model library or of Modia.

Variables of built-in Components

The following functions are provided to define and access new variables in built-in components (seee for example model InsulatedRod2 in Modia/models/HeatTransfer.jl).

FunctionsDescription
new_x_segmented_variable!Generate new state variable (x_segmented and der_x_segmented variables)
new_w_segmented_variable!Generate new local variable (w_segmented variable)
new_alias_segmented_variable!Generate new alias variable
new_z_segmented_variable!Generate new zero crossing variables (z_segmented variables)
get_x_startIndex_from_x_segmented_startIndexReturn start index of x_segmented variable with respect to state vector x
copy_scalar_x_segmented_value_from_stateReturn value of scalar x_segmented variable from state vector x
copy_SVector3_x_segmented_value_from_stateReturn value of SVector{3,FloatType} x_segmented variable from state vector x
copy_Vector_x_segmented_value_from_stateReturn value of Vector{FloatType} x_segmented variable from state vector x
copy_der_x_segmented_value_to_stateCopy value of der_x_segmented variable to state derivative vector der(x)
copy_w_segmented_value_to_resultCopy value of local variable (w-segmented) to result
Modia.new_x_segmented_variable!Function
startIndex = new_x_segmented_variable!(
                partiallyInstantiatedModel::InstantiatedModel,
                x_name::String, der_x_name::String, startOrInit, x_unit::String="";
                nominal::Float64 = NaN, unbounded::Bool = false)::Int

Generate new states (x_segmented and der_x_segmented variables) and return the startIndex of the variables in order that actual values can be inquired or copied from the state x and state derivative der(x)vectors via get_x_startIndex_from_x_segmented_startIndex. startOrInit contain the start or init values of the newly generated x_segmented variable.

Actual values of these new variables are stored in:

  • instantiatedModel.x_segmented[startIndex:startIndex+prod(dims(startOrInit))-1]
  • instantiatedModel.der_x_segmented[startIndex:startIndex+prod(dims(startOrInit))-1]

Value startOrInit is the start/init value used during re-initialization of the new segment with initFullRestart!(..).

source
Modia.new_w_segmented_variable!Function
index = new_w_segmented_variable!(
           partiallyInstantiatedModel::InstantiatedModel, name::String,
           w_segmented_default, unit::String="")::Int

Generate new local variable (w_segmented variable) and return the index of the variable in order that actual values can be inquired or copied from the result data structure. New values of w_segmented variables need only to be computed at communication points. Value wsegmenteddefault is stored as default value and defines type and (fixed) size of the variable in this simulation segment.

source
Modia.new_z_segmented_variable!Function
startIndex = new_z_segmented_variable!(instantiatedModel, nz)

Generate nz new zero crossing variables and return the startIndex to of the variables in order that actual values can be copied into the vector of zero crossings.

source
Modia.get_x_startIndex_from_x_segmented_startIndexFunction
x_startIndex = get_x_startIndex_from_x_segmented_startIndex(
                  instantiatedModel::InstantiatedModel, x_segmented_startIndex)

Return the startindex of an x_segmented state with respect to the x-vector, given the startIndex with respect to the x_segmented vector (x_segmented_startIndex is the return value of new_x_segmented_variable!(..)).

source
Modia.copy_scalar_x_segmented_value_from_stateFunction
value = Modia.copy_scalar_x_segmented_value_from_state(instantiatedModel, startIndex)

Return value of scalar x_segmented variable from state vector x by providing its startIndex (returned from new_x_segmented_variable!(..)).

source
Modia.copy_SVector3_x_segmented_value_from_stateFunction
value = Modia.copy_SVector3_x_segmented_value_from_state(instantiatedModel, startIndex)

Return value of SVector{3,FloatType} x_segmented variable from state vector x by providing its startIndex (returned from new_x_segmented_variable!(..)).

source
Modia.copy_Vector_x_segmented_value_from_stateFunction
Modia.copy_Vector_x_segmented_value_from_state(
    instantiatedModel::InstantiatedModel, startIndex, xi::Vector{FloatType})::Nothing

Return value of Vector{FloatType} x_segmented variable from state vector x by providing its startIndex (returned from new_x_segmented_variable!(..)) and copying it into the pre-allocated vector xi.

source
Modia.copy_der_x_segmented_value_to_stateFunction
Modia.copy_der_x_segmented_value_to_state(
   instantiatedModel, startIndex, 
   der_x_segmented_value::[FloatType|Vector{FloatType}])

Copy der_x_segmented_value to state derivative vector der(x) by providing its startIndex (returned from new_x_segmented_variable!(..)) and copying it into the pre-allocated vector der_x_segmented_value.

source
Modia.copy_w_segmented_value_to_resultFunction
Modia.copy_w_segmented_value_to_result(
    instantiatedModel::InstantiatedModel, index::Int, 
    w_segmented_value)::Nothing

Copy value of local variable (w-segmented) to result by providing its index (returned from new_w_segmented_variable!),

source

Inquiries in built-in Components

The following functions are provided to inquire properties in built-in components at the current state of the simulation (see for example model InsulatedRod2 in Modia/models/HeatTransfer.jl).

Modia.isInitialFunction
isInitial(instantiatedModel)

Return true, if initialization phase of simulation (of the current segment of a segmented simulation).

source
Modia.isFirstInitialOfAllSegmentsFunction
isFirstInitialOfAllSegments(instantiatedModel)

Return true, if initialization phase of simulation of the first segment of a segmented simulation.

source
Modia.isTerminalFunction
isTerminal(instantiatedModel)

Return true, if terminal phase of simulation (of the current segment of a segmented simulation).

source
Modia.isTerminalOfAllSegmentsFunction
isTerminalOfAllSegments(instantiatedModel)

Return true, if terminal phase of simulation of the last segment of a segmented simulation.

source
Modia.isEventFunction
isEvent(instantiatedModel)

Return true, if event phase of simulation (including initialization).

source
Modia.isFirstEventIterationFunction
isFirstEventIteration(instantiatedModel)

Return true, if event phase of simulation (including initialization) and during the first iteration of the event iteration.

source
Modia.isFullRestartFunction
isFullRestart(instantiatedModel)

Return true, if FullRestart event of a segmented simulation.

source
Modia.isZeroCrossingFunction
isZeroCrossing(instantiatedModel)

Return true, if event indicators (zero crossings) shall be computed.

source
Modia.getTimeFunction
tCurrent = getTime(instantiatedModel)

Return current simulation time.

source

Code Generation

This section lists internal functions to generate Julia code of the transformed equations.

Modia.InstantiatedModelType
simulationModel = InstantiatedModel{FloatType,TimeType}(
        modelModule, modelName, getDerivatives!, equationInfo, x_startValues,
        parameters, timeName, w_invariant_names;
        vSolvedWithInitValuesAndUnit::OrderedDict{String,Any}(),
        vEliminated::Vector{Int}=Int[],
        vProperty::Vector{Int}=Int[],
        var_name::Function = v->nothing)

Arguments

  • modelModule: Module in which @instantiateModel is invoked (it is used for Core.eval(modelModule, ...)), that is evaluation of expressions in the environment of the user.
  • modelName::String: Name of the model
  • getDerivatives::Function: Function that is used to evaluate the model equations, typically generated with [Modia.generate_getDerivatives!].
  • equationInfo::Modia.EquationInfo: Information about the states and the equations.
  • x_startValues:: Deprecated (is no longer used).
  • parameters: A hierarchical NamedTuple of (key, value) pairs defining the parameter and init/start values.
  • timeName: Name of time (as Symbol)
  • w_invariant_names: A vector of variable names (as vector of symbols or Expr)
source
Modia.generate_getDerivatives!Function
code = generate_getDerivatives!(AST, equationInfo, parameters, timeName, w_invariant_names, functionName;
                                hasUnits=false)

Return the code of the getDerivatives! function as Expr using the Symbol functionName as function name. By eval(code) or fc = @RuntimeGeneratedFunction(code) the function is compiled and can afterwards be called.

Arguments

  • AST::Vector{Expr}: Abstract Syntax Tree of the equations as vector of Expr.

  • equationInfo::Modia.EquationInfo: Data structure returned by `Modia.getSortedAndSolvedAST holding information about the states.

  • parameters: Vector of parameter names (as vector of symbols)

  • timeName: Name of time (as symbol)

  • w_invariant_names: Vector of variable names (as vector of symbols or Expr).

  • functionName::Function: The name of the function that shall be generated.

Optional Arguments

  • pre:Vector{Symbol}: pre-variable names

  • hasUnits::Bool: = true, if variables have units. Note, the units of the state vector are defined in equationinfo.

source
Modia.init!Function
success = init!(simulationModel)

Initialize simulationModel::InstantiatedModel at startTime. In particular:

  • Empty result data structure.

  • Merge parameter and init/start values into simulationModel.

  • Construct x_start.

  • Call simulationModel.getDerivatives! once with isInitial(simulationModel) = true to compute and store all variables in the result data structure at startTime and initialize simulationModel.linearEquations.

  • Check whether explicitly solved variables that have init-values defined, have the required value after initialization (-> otherwise error).

If initialization is successful return true, otherwise false.

source
Modia.outputs!Function
outputs!(x, t, integrator)

DifferentialEquations FunctionCallingCallback function for InstantiatedModel that is used to store results at communication points.

source
Modia.derivatives!Function
derivatives!(derx, x, m, t)

DifferentialEquations callback function to get the derivatives.

source
Modia.DAEresidualsForODE!Function
DAEresidualsForODE!(residuals, derx, x, m, t)

DifferentialEquations callback function for DAE integrator for ODE model

source
Modia.affectEvent!Function
affectEvent!(integrator, stateEvent, eventIndex)

Called when a time event (stateEvent=false) or state event (stateEvent=true) is triggered. In case of stateEvent, eventIndex is the index of the crossing function that triggered the event.

source
Modia.addToResult!Function
addToResult!(instantiatedModel, x, time, w_invariant...)

Add result of current time instant (time, x, der_x, w_invariant, w_segmented) to instantiatedModel.

source
Modia.getFloatTypeFunction
floatType = getFloatType(simulationModel::InstantiatedModel)

Return the floating point type with which simulationModel is parameterized (for example returns: Float64, Float32, DoubleFloat, Measurements.Measurement{Float64}).

source
Modia.measurementToStringFunction
str = measurementToString(v)

Return variable v::Measurements.Measurement{FloatType} or a vector of such variables in form of a string will the full number of significant digits.

source