text:
#
# Copyright (c) 2021 Tobias Thummerer, Lars Mikelsons, Josef Kircher
# Licensed under the MIT license. See LICENSE file in the project root for details.
#

# What is included in the file `FMI3_c.jl`?
# - default implementations for the `fmi3CallbackFunctions`
# - enum definition for `fmi3InstanceState` to protocol the current FMU-state
# - julia-implementaions of the functions inside the FMI-specification
# Any c-function `f(c::fmi3Instance, args...)` in the spec is implemented as `f(c::FMU3Instance, args...)`.
# Any c-function `f(args...)` without a leading `fmi3Instance`-arguemnt is implented as `f(c_ptr, args...)` where `c_ptr` is a pointer to the c-function (inside the DLL).

import FMICore: fmi3InstantiateCoSimulation, fmi3InstantiateModelExchange, fmi3InstantiateScheduledExecution, fmi3FreeInstance!, fmi3GetVersion
import FMICore: fmi3SetDebugLogging, fmi3EnterInitializationMode, fmi3ExitInitializationMode, fmi3Terminate, fmi3Reset
import FMICore: fmi3GetFloat32!, fmi3SetFloat32, fmi3GetFloat64!, fmi3SetFloat64
import FMICore: fmi3GetInt8!, fmi3SetInt8, fmi3GetInt16!, fmi3SetInt16,fmi3GetInt32!, fmi3SetInt32, fmi3GetInt64!, fmi3SetInt64
import FMICore: fmi3GetUInt8!, fmi3SetUInt8, fmi3GetUInt16!, fmi3SetUInt16,fmi3GetUInt32!, fmi3SetUInt32, fmi3GetUInt64!, fmi3SetUInt64
import FMICore: fmi3GetBoolean!, fmi3SetBoolean, fmi3GetString!, fmi3SetString, fmi3GetBinary!, fmi3SetBinary, fmi3GetClock!, fmi3SetClock
import FMICore: fmi3GetFMUState!, fmi3SetFMUState, fmi3FreeFMUState!, fmi3SerializedFMUStateSize!, fmi3SerializeFMUState!, fmi3DeSerializeFMUState!
import FMICore: fmi3SetIntervalDecimal, fmi3SetIntervalFraction, fmi3GetIntervalDecimal!, fmi3GetIntervalFraction!, fmi3GetShiftDecimal!, fmi3GetShiftFraction!, fmi3ActivateModelPartition
import FMICore: fmi3GetNumberOfVariableDependencies!, fmi3GetVariableDependencies!
import FMICore: fmi3GetDirectionalDerivative!, fmi3GetAdjointDerivative!, fmi3GetOutputDerivatives!
import FMICore: fmi3EnterConfigurationMode, fmi3ExitConfigurationMode
import FMICore: fmi3GetNumberOfContinuousStates!, fmi3GetNumberOfEventIndicators!
import FMICore: fmi3DoStep!, fmi3EnterStepMode
import FMICore: fmi3SetTime, fmi3SetContinuousStates, fmi3EnterEventMode, fmi3UpdateDiscreteStates, fmi3EnterContinuousTimeMode, fmi3CompletedIntegratorStep!
import FMICore: fmi3GetContinuousStateDerivatives, fmi3GetEventIndicators!, fmi3GetContinuousStates!, fmi3GetNominalsOfContinuousStates!, fmi3EvaluateDiscreteStates

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.1. Super State: FMU State Setable

Function that is called in the FMU, usually if an fmi3XXX function, does not behave as desired. If “logger” is called with “status = fmi3OK”, then the message is a pure information message.
instanceEnvironment - is the instance name of the model that calls this function.
category - is the category of the message. The meaning of “category” is defined by the modeling environment that generated the FMU. Depending on this modeling environment, none, some or all allowed values of “category” for this FMU are defined in the modelDescription.xml file via element “<fmiModelDescription><LogCategories>”, see section 2.4.5. Only messages are provided by function logger that have a category according to a call to fmi3SetDebugLogging (see below).
message - is provided in the same way and with the same format control as in function “printf” from the C standard library. [Typically, this function prints the message and stores it optionally in a log file.]
"""
# TODO error in the specification??
function fmi3CallbackLogger(_instanceEnvironment::Ptr{Cvoid},
    _status::Cuint,
    _category::Ptr{Cchar},
    _message::Ptr{Cchar})

message = unsafe_string(_message)
category = unsafe_string(_category)
status = fmi3StatusToString(_status)
instanceEnvironment = unsafe_load(_instanceEnvironment)

if status == fmi3StatusOK && instanceEnvironment.logStatusOK
@info "[$status][$category][$instanceName]: $message"
elseif (status == fmi3StatusWarning && instanceEnvironment.logStatusWarning)
@warn "[$status][$category][$instanceName]: $message"
elseif (status == fmi3StatusDiscard && instanceEnvironment.logStatusDiscard) ||
   (status == fmi3StatusError   && instanceEnvironment.logStatusError) ||
   (status == fmi3StatusFatal   && instanceEnvironment.logStatusFatal)
@error "[$status][$category][$instanceName]: $message"
end

return nothing
end

"""
Source: FMISpec3.0, Version D5ef1c1: 4.2.2. State: Intermediate Update Mode

When a Co-Simulation FMU provides values for its output variables at intermediate points between two consecutive communication points, and is able to receive new values for input variables at these intermediate points, the Intermediate Update Callback function is called. This is typically required when the FMU uses a numerical solver to integrate the FMU's internal state between communication points in fmi3DoStep.
The callback function switches the FMU from Step Mode (see 4.2.1.) in the Intermediate Update Mode (see 4.2.2.) and returns to Step Mode afterwards. The parameters of this function are:

instanceEnvironment - is the instance name of the model that calls this function.

intermediateUpdateTime - is the internal value of the independent variable [typically simulation time] of the FMU at which the callback has been called for intermediate and final steps. If an event happens or an output Clock ticks, intermediateUpdateTime is the time of event or output Clock tick. In Co-Simulation, intermediateUpdateTime is restricted by the arguments to fmi3DoStep as follows:
currentCommunicationPoint ≤ intermediateUpdateTime ≤ (currentCommunicationPoint + communicationStepSize).
The FMU must not call the callback function fmi3CallbackIntermediateUpdate with an intermediateUpdateTime that is smaller than the intermediateUpdateTime given in a previous call of fmi3CallbackIntermediateUpdate with intermediateStepFinished == fmi3True.

If intermediateVariableSetRequested == fmi3True, the co-simulation algorithm may provide intermediate values for continuous input variables with intermediateUpdate = true by calling fmi3Set{VariableType}. The set of variables for which the co-simulation algorithm will provide intermediate values is declared through the requiredIntermediateVariables argument to fmi3InstantiateXXX. If a co-simulation algorithm does not provide a new value for any of the variables contained in the set it registered, the last value set remains.

If intermediateVariableGetAllowed == fmi3True, the co-simulation algorithm may collect intermediate output variables by calling fmi3Get{VariableType} for variables with intermediateUpdate = true. The set of variables for which the co-simulation algorithm can get values is supplied through the requiredIntermediateVariables argument to fmi3InstantiateXXX.

If intermediateStepFinished == fmi3False, the intermediate outputs of the FMU that the co-simulation algorithm inquires with fmi3Get{VariableType} resulting from tentative internal solver states and may still change for the same intermediateUpdateTime [e.g., if the solver deems the tentative state to cause a too high approximation error, it may go back in time and try to re-estimate the state using smaller internal time steps].
If intermediateStepFinished == fmi3True, intermediate outputs inquired by the co-simulation algorithm with fmi3Get{VariableType} correspond to accepted internal solver step.

When canReturnEarly == fmi3True the FMU signals to the co-simulation algorithm its ability to return early from the current fmi3DoStep.

earlyReturnRequested - If and only if canReturnEarly == fmi3True, the co-simulation algorithm may request the FMU to return early from fmi3DoStep by setting earlyReturnRequested == fmi3True.

earlyReturnTime is used to signal the FMU at which time to return early from the current fmi3DoStep, if the return value of earlyReturnRequested == fmi3True. If the earlyReturnTime is greater than the last signaled intermediateUpdateTime, the FMU may integrate up to the time instant earlyReturnTime.

If the ModelDescription has the "providesIntermediateUpdate" flag, the Intermediate update callback function is called. That flag is ignored in ModelExchange and ScheduledExecution.
"""
function fmi3CallbackIntermediateUpdate(instanceEnvironment::Ptr{Cvoid},
    intermediateUpdateTime::fmi3Float64,
    intermediateVariableSetRequested::fmi3Boolean,
    intermediateVariableGetAllowed::fmi3Boolean,
    intermediateStepFinished::fmi3Boolean,
    canReturnEarly::fmi3Boolean,
    earlyReturnRequested::Ptr{fmi3Boolean},
    earlyReturnTime::Ptr{fmi3Float64})
    @debug "To be implemented!"
end

"""
Source: FMISpec3.0, Version D5ef1c1: 5.2.2. State: Clock Activation Mode

A model partition of a Scheduled Execution FMU calls fmi3CallbackClockUpdate to signal that a triggered output Clock ticked or a new interval for a countdown Clock is available.
fmi3CallbackClockUpdate switches the FMU itself then into the Clock Update Mode (see 5.2.3.). The callback may be called from several model partitions.

instanceEnvironment - is the instance name of the model that calls this function.
"""
function fmi3CallbackClockUpdate(instanceEnvironment::Ptr{Cvoid})
    @debug "to be implemented!"
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.1. Super State: FMU State Setable

Disposes the given instance, unloads the loaded model, and frees all the allocated memory and other resources that have been allocated by the functions of the FMU interface. If a NULL pointer is provided for argument instance, the function call is ignored (does not have an effect).
"""
function fmi3FreeInstance!(c::FMU3Instance)

    ind = findall(x->x==c, c.fmu.instances)
    deleteat!(c.fmu.instances, ind)
    ccall(c.fmu.cFreeInstance, Cvoid, (Ptr{Cvoid},), c.compAddr)

    nothing
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.4. Inquire Version Number of Header Files

This function returns fmi3Version of the fmi3Functions.h header file which was used to compile the functions of the FMU. This function call is allowed always and in all interface types.

The standard header file as documented in this specification has version "3.0-beta.2", so this function returns "3.0-beta.2".
"""
function fmi3GetVersion(fmu::FMU3)

    fmi3Version = fmi3GetVersion(fmu.cGetVersion)

    unsafe_string(fmi3Version)
end
# special case
function fmi3GetVersion(c::FMU3Instance)
    fmi3GetVersion(c.fmu)
end

# helper
function checkStatus(c::FMU3Instance, status::fmi3Status)
    @assert (status != fmi3StatusWarning) || !c.fmu.executionConfig.assertOnWarning "Assert on `fmi3StatusWarning`. See stack for errors."

    if status == fmi3StatusError
        c.state = fmi3InstanceStateError
        @assert !c.fmu.executionConfig.assertOnError "Assert on `fmi3StatusError`. See stack for errors."

    elseif status == fmi3StatusFatal
        c.state = fmi3InstanceStateFatal
        @assert false "Assert on `fmi3StatusFatal`. See stack for errors."
    end
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.1. Super State: FMU State Setable

The function controls debug logging that is output via the logger function callback. If loggingOn = fmi3True, debug logging is enabled, otherwise it is switched off.
"""
function fmi3SetDebugLogging(c::FMU3Instance, logginOn::fmi3Boolean, nCategories::UInt, categories::Ptr{Nothing})
    status = fmi3SetDebugLogging(c.fmu.cSetDebugLogging, c.compAddr, logginOn, nCategories, categories)
    checkStatus(c, status)
    return status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.2. State: Instantiated

Informs the FMU to enter Initialization Mode. Before calling this function, all variables with attribute <Datatype initial = "exact" or "approx"> can be set with the “fmi3SetXXX” functions (the ScalarVariable attributes are defined in the Model Description File, see section 2.4.7). Setting other variables is not allowed.
Also sets the simulation start and stop time.
"""
function fmi3EnterInitializationMode(c::FMU3Instance, toleranceDefined::fmi3Boolean,
    tolerance::fmi3Float64,
    startTime::fmi3Float64,
    stopTimeDefined::fmi3Boolean,
    stopTime::fmi3Float64)
    if c.state != fmi3InstanceStateInstantiated
        @warn "fmi3EnterInitializationMode(...): Needs to be called in state `fmi3IntanceStateInstantiated`."
    end
    status = fmi3EnterInitializationMode(c.fmu.cEnterInitializationMode, c.compAddr, toleranceDefined, tolerance, startTime, stopTimeDefined, stopTime)
    checkStatus(c, status)
    if status == fmi3StatusOK
        c.state = fmi3InstanceStateInitializationMode
    end
    return status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.3. State: Initialization Mode

Informs the FMU to exit Initialization Mode.
"""
function fmi3ExitInitializationMode(c::FMU3Instance)
    if c.state != fmi3InstanceStateInitializationMode
        @warn "fmi3ExitInitializationMode(...): Needs to be called in state `fmi3InstanceStateInitializationMode`."
    end

    status = fmi3ExitInitializationMode(c.fmu.cExitInitializationMode, c.compAddr)
    checkStatus(c, status)
    if status == fmi3StatusOK
        if fmi3IsCoSimulation(c.fmu) && c.fmu.modelDescription.coSimulation.hasEventMode
            c.state = fmi3InstanceStateStepMode
        elseif fmi3IsScheduledExecution(c.fmu)
            c.state = fmi3InstanceStateClockActivationMode
        else
            c.state = fmi3InstanceStateEventMode
        end
    end
    return status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.4. Super State: Initialized

Informs the FMU that the simulation run is terminated.
"""
function fmi3Terminate(c::FMU3Instance; soft::Bool=false)
    if c.state != fmi3InstanceStateContinuousTimeMode && c.state != fmi3InstanceStateEventMode && c.state != fmi3InstanceStateClockActivationMode && c.state != fmi3InstanceStateStepMode
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3Terminate(_): Needs to be called in state `fmi3InstanceStateContinuousTimeMode`, `fmi3InstanceStateEventMode`, `fmi3InstanceStateClockActivationMode` or `fmi3InstanceStateStepMode`."
        end
    end

    status = fmi3Terminate(c.fmu.cTerminate, c.compAddr)
    checkStatus(c, status)
    if status == fmi3StatusOK
        c.state = fmi3InstanceStateTerminated
    end
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.1. Super State: FMU State Setable

Is called by the environment to reset the FMU after a simulation run. The FMU goes into the same state as if fmi3InstantiateXXX would have been called.
"""
function fmi3Reset(c::FMU3Instance)
    if c.state != fmi3InstanceStateTerminated && c.state != fmi3InstanceStateError
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3Reset(_): Needs to be called in state `fmi3InstanceStateTerminated` or `fmi3InstanceStateError`."
        end
    end

    if c.fmu.cReset == C_NULL
        fmi3FreeInstance!(c.fmu.cFreeInstance, c.compAddr)
        if fmi3IsCoSimulation(c.fmu)
            compAddr = fmi3InstantiateCoSimulation!(c.fmu)
        elseif fmi3IsModelExchange(c.fmu)
            compAddr = fmi3InstantiateModelExchange!(c.fmu)
        elseif fmi3IsScheduledExecution(c.fmu)
            compAddr = fmi3InstantiateScheduledExecution!(c.fmu)
        end

        if compAddr == Ptr{Cvoid}(C_NULL)
            @error "fmi3Reset(...): Reinstantiation failed!"
            return fmi3StatusError
        end

        c.compAddr = compAddr
        return fmi3StatusOK
    else
        status = fmi3Reset(c.fmu.cReset, c.compAddr)
        checkStatus(c, status)
        if status == fmi3StatusOK
            c.state = fmi3InstanceStateInstantiated
        end
        return status
    end
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetFloat32!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Float32}, nvalue::Csize_t)
    status = fmi3GetFloat32!(c.fmu.cGetFloat32,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetFloat32(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Float32}, nvalue::Csize_t)
    status = fmi3SetFloat32(c.fmu.cSetFloat32,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetFloat64!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Float64}, nvalue::Csize_t)
    status = fmi3GetFloat64!(c.fmu.cGetFloat64,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetFloat64(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Float64}, nvalue::Csize_t)
    status = fmi3SetFloat64(c.fmu.cSetFloat64,
               c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetInt8!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int8}, nvalue::Csize_t)
    status = fmi3GetInt8!(c.fmu.cGetInt8,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetInt8(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int8}, nvalue::Csize_t)
    status = fmi3SetInt8(c.fmu.cSetInt8,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetUInt8!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt8}, nvalue::Csize_t)
    status = fmi3GetUInt8!(c.fmu.cGetUInt8,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetUInt8(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt8}, nvalue::Csize_t)
    status = fmi3SetUInt8(c.fmu.cSetUInt8,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetInt16!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int16}, nvalue::Csize_t)
    status = fmi3GetInt16!(c.fmu.cGetInt16,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetInt16(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int16}, nvalue::Csize_t)
    status = fmi3SetInt16(c.fmu.cSetInt16,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetUInt16!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt16}, nvalue::Csize_t)
    status = fmi3GetUInt16!(c.fmu.cGetUInt16,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetUInt16(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt16}, nvalue::Csize_t)
    status = fmi3SetUInt16(c.fmu.cSetUInt16,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetInt32!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int32}, nvalue::Csize_t)
    status = fmi3GetInt32!(c.fmu.cGetInt32,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetInt32(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int32}, nvalue::Csize_t)
    status = fmi3SetInt32(c.fmu.cSetInt32,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetUInt32!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt32}, nvalue::Csize_t)
    status = fmi3GetUInt32!(c.fmu.cGetUInt32,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetUInt32(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt32}, nvalue::Csize_t)
    status = fmi3SetUInt32(c.fmu.cSetUInt32,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetInt64!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int64}, nvalue::Csize_t)
    status = fmi3GetInt64!(c.fmu.cGetInt64,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetInt64(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Int64}, nvalue::Csize_t)
    status = fmi3SetInt64(c.fmu.cSetInt64,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO test, no variable in FMUs
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetUInt64!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt64}, nvalue::Csize_t)
    status = fmi3GetUInt64!(c.fmu.cGetUInt64,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetUInt64(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3UInt64}, nvalue::Csize_t)
    status = fmi3SetUInt64(c.fmu.cSetUInt64,
            c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetBoolean!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Boolean}, nvalue::Csize_t)
    status = fmi3GetBoolean!(c.fmu.cGetBoolean,
          c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetBoolean(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Boolean}, nvalue::Csize_t)
    status = fmi3SetBoolean(c.fmu.cSetBoolean,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO change to fmi3String when possible to test
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetString!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::Vector{Ptr{Cchar}}, nvalue::Csize_t)
    status = fmi3GetString!(c.fmu.cGetString,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetString(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::Union{AbstractArray{Ptr{Cchar}}, AbstractArray{Ptr{UInt8}}}, nvalue::Csize_t)
    status = fmi3SetString(c.fmu.cSetString,
                c.compAddr, vr, nvr, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValues - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetBinary!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, valueSizes::AbstractArray{Csize_t}, value::AbstractArray{fmi3Binary}, nvalue::Csize_t)
    status = fmi3GetBinary!(c.fmu.cGetBinary,
                c.compAddr, vr, nvr, valueSizes, value, nvalue)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetBinary(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, valueSizes::AbstractArray{Csize_t}, value::AbstractArray{fmi3Binary}, nvalue::Csize_t)
    status = fmi3SetBinary(c.fmu.cSetBinary,
                c.compAddr, vr, nvr, valueSizes, value, nvalue)
    checkStatus(c, status)
    status
end

# TODO, Clocks not implemented so far thus not tested
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3GetClock!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Clock})
    status = fmi3GetClock!(c.fmu.cGetClock,
                c.compAddr, vr, nvr, value)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.2. Getting and Setting Variable Values

Functions to get and set values of variables idetified by their valueReference.

nValue - is different from nvr if the value reference represents an array and therefore are more values tied to a single value reference.
"""
function fmi3SetClock(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, value::AbstractArray{fmi3Clock})
    status = fmi3SetClock(c.fmu.cSetClock,
                c.compAddr, vr, nvr, value)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.4. Getting and Setting the Complete FMU State

fmi3GetFMUstate makes a copy of the internal FMU state and returns a pointer to this copy
"""
function fmi3GetFMUState!(c::FMU3Instance, FMUstate::Ref{fmi3FMUState})
    status = fmi3GetFMUState!(c.fmu.cGetFMUState,
                c.compAddr, FMUstate)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.4. Getting and Setting the Complete FMU State

fmi3SetFMUstate copies the content of the previously copied FMUstate back and uses it as actual new FMU state.
"""
function fmi3SetFMUState(c::FMU3Instance, FMUstate::fmi3FMUState)
    status = fmi3SetFMUState(c.fmu.cSetFMUState,
                c.compAddr, FMUstate)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.4. Getting and Setting the Complete FMU State

fmi3FreeFMUstate frees all memory and other resources allocated with the fmi3GetFMUstate call for this FMUstate.
"""
function fmi3FreeFMUState!(c::FMU3Instance, FMUstate::Ref{fmi3FMUState})
    status = fmi3FreeFMUState!(c.fmu.cFreeFMUState,
                c.compAddr, FMUstate)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.4. Getting and Setting the Complete FMU State

fmi3SerializedFMUstateSize returns the size of the byte vector which is needed to store FMUstate in it.
"""
function fmi3SerializedFMUStateSize!(c::FMU3Instance, FMUstate::fmi3FMUState, size::Ref{Csize_t})
    status = fmi3SerializedFMUStateSize!(c.fmu.cSerializedFMUStateSize,
                c.compAddr, FMUstate, size)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.4. Getting and Setting the Complete FMU State

fmi3SerializeFMUstate serializes the data which is referenced by pointer FMUstate and copies this data in to the byte vector serializedState of length size
"""
function fmi3SerializeFMUState!(c::FMU3Instance, FMUstate::fmi3FMUState, serialzedState::AbstractArray{fmi3Byte}, size::Csize_t)
    status = fmi3SerializeFMUState!(c.fmu.cSerializeFMUState,
                c.compAddr, FMUstate, serialzedState, size)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.6.4. Getting and Setting the Complete FMU State

fmi3DeSerializeFMUstate deserializes the byte vector serializedState of length size, constructs a copy of the FMU state and returns FMUstate, the pointer to this copy.
"""
function fmi3DeSerializeFMUState!(c::FMU3Instance, serialzedState::AbstractArray{fmi3Byte}, size::Csize_t, FMUstate::Ref{fmi3FMUState})
    status = fmi3DeSerializeFMUState!(c.fmu.cDeSerializeFMUState,
                c.compAddr, serialzedState, size, FMUstate)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.9. Clocks

fmi3SetIntervalDecimal sets the interval until the next clock tick
"""
# TODO Clocks and dependencies functions
function fmi3SetIntervalDecimal(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, intervals::AbstractArray{fmi3Float64})
    status = fmi3SetIntervalDecimal(c.fmu.cSetIntervalDecimal,
                c.compAddr, vr, nvr, intervals)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.9. Clocks

fmi3SetIntervalFraction sets the interval until the next clock tick
Only allowed if the attribute 'supportsFraction' is set.
"""
function fmi3SetIntervalFraction(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, intervalCounters::AbstractArray{fmi3UInt64}, resolutions::AbstractArray{fmi3UInt64})
    status = fmi3SetIntervalFraction(c.fmu.cSetIntervalFraction,
                c.compAddr, vr, nvr, intervalCounters, resolutions)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.9. Clocks

fmi3GetIntervalDecimal retrieves the interval until the next clock tick.

For input Clocks it is allowed to call this function to query the next activation interval.
For changing aperiodic Clock, this function must be called in every Event Mode where this clock was activated.
For countdown aperiodic Clock, this function must be called in every Event Mode.
Clock intervals are computed in fmi3UpdateDiscreteStates (at the latest), therefore, this function should be called after fmi3UpdateDiscreteStates.
For information about fmi3IntervalQualifiers, call ?fmi3IntervalQualifier
"""
function fmi3GetIntervalDecimal!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, intervals::AbstractArray{fmi3Float64}, qualifiers::fmi3IntervalQualifier)
    status = fmi3GetIntervalDecimal!(c.fmu.cGetIntervalDecimal,
                c.compAddr, vr, nvr, intervals, qualifiers)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.9. Clocks

fmi3GetIntervalFraction retrieves the interval until the next clock tick.

For input Clocks it is allowed to call this function to query the next activation interval.
For changing aperiodic Clock, this function must be called in every Event Mode where this clock was activated.
For countdown aperiodic Clock, this function must be called in every Event Mode.
Clock intervals are computed in fmi3UpdateDiscreteStates (at the latest), therefore, this function should be called after fmi3UpdateDiscreteStates.
For information about fmi3IntervalQualifiers, call ?fmi3IntervalQualifier
"""
function fmi3GetIntervalFraction!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, intervalCounters::AbstractArray{fmi3UInt64}, resolutions::AbstractArray{fmi3UInt64}, qualifiers::fmi3IntervalQualifier)
    status = fmi3GetIntervalFraction!(c.fmu.cGetIntervalFraction,
                c.compAddr, vr, nvr, intervalCounters, resolutions, qualifiers)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.9. Clocks

fmi3GetShiftDecimal retrieves the delay to the first Clock tick from the FMU.
"""
function fmi3GetShiftDecimal!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, shifts::AbstractArray{fmi3Float64})
    status = fmi3GetShiftDecimal!(c.fmu.cGetShiftDecimal,
                c.compAddr, vr, nvr, shifts)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.9. Clocks

fmi3GetShiftFraction retrieves the delay to the first Clock tick from the FMU.
"""
function fmi3GetShiftFraction!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nvr::Csize_t, shiftCounters::AbstractArray{fmi3UInt64}, resolutions::AbstractArray{fmi3UInt64})
    status = fmi3GetShiftFraction!(c.fmu.cGetShiftFraction,
                c.compAddr, vr, nvr, shiftCounters, resolutions)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 5.2.2. State: Clock Activation Mode

During Clock Activation Mode (see 5.2.2.) after fmi3ActivateModelPartition has been called for a calculated, tunable or changing Clock the FMU provides the information on when the Clock will tick again, i.e. when the corresponding model partition has to be scheduled the next time.

Each fmi3ActivateModelPartition call is associated with the computation of an exposed model partition of the FMU and therefore to an input Clock.
"""
function fmi3ActivateModelPartition(c::FMU3Instance, vr::fmi3ValueReference, activationTime::AbstractArray{fmi3Float64})
    status = fmi3ActivateModelPartition(c.fmu.cActivateModelPartition,
                c.compAddr, vr, activationTime)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.10. Dependencies of Variables

The number of dependencies of a given variable, which may change if structural parameters are changed, can be retrieved by calling fmi3GetNumberOfVariableDependencies.

This information can only be retrieved if the 'providesPerElementDependencies' tag in the ModelDescription is set.
"""
# TODO not tested
function fmi3GetNumberOfVariableDependencies!(c::FMU3Instance, vr::fmi3ValueReference, nvr::Ref{Csize_t})
    status = fmi3GetNumberOfVariableDependencies!(c.fmu.cGetNumberOfVariableDependencies,
                c.compAddr, vr, nvr)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.10. Dependencies of Variables

The actual dependencies (of type dependenciesKind) can be retrieved by calling the function fmi3GetVariableDependencies:

dependent - specifies the valueReference of the variable for which the dependencies should be returned.

nDependencies - specifies the number of dependencies that the calling environment allocated space for in the result buffers, and should correspond to value obtained by calling fmi3GetNumberOfVariableDependencies.

elementIndicesOfDependent - must point to a buffer of size_t values of size nDependencies allocated by the calling environment.
It is filled in by this function with the element index of the dependent variable that dependency information is provided for. The element indices start with 1. Using the element index 0 means all elements of the variable. (Note: If an array has more than one dimension the indices are serialized in the same order as defined for values in Section 2.2.6.1.)

independents -  must point to a buffer of fmi3ValueReference values of size nDependencies allocated by the calling environment.
It is filled in by this function with the value reference of the independent variable that this dependency entry is dependent upon.

elementIndicesIndependents - must point to a buffer of size_t values of size nDependencies allocated by the calling environment.
It is filled in by this function with the element index of the independent variable that this dependency entry is dependent upon. The element indices start with 1. Using the element index 0 means all elements of the variable. (Note: If an array has more than one dimension the indices are serialized in the same order as defined for values in Section 2.2.6.1.)

dependencyKinds - must point to a buffer of dependenciesKind values of size nDependencies allocated by the calling environment.
It is filled in by this function with the enumeration value describing the dependency of this dependency entry.
For more information about dependenciesKinds, call ?fmi3DependencyKind

If this function is called before the fmi3ExitInitializationMode call, it returns the initial dependencies. If this function is called after the fmi3ExitInitializationMode call, it returns the runtime dependencies.
The retrieved dependency information of one variable becomes invalid as soon as a structural parameter linked to the variable or to any of its depending variables are set. As a consequence, if you change structural parameters affecting B or A, the dependency of B becomes invalid. The dependency information must change only if structural parameters are changed.

This information can only be retrieved if the 'providesPerElementDependencies' tag in the ModelDescription is set.
"""
function fmi3GetVariableDependencies!(c::FMU3Instance, vr::fmi3ValueReference, elementIndiceOfDependents::AbstractArray{Csize_t}, independents::AbstractArray{fmi3ValueReference},
    elementIndiceOfInpendents::AbstractArray{Csize_t}, dependencyKind::AbstractArray{fmi3DependencyKind}, ndependencies::Csize_t)
    status = fmi3GetVariableDependencies!(c.fmu.cGetVariableDependencies,
               c.compAddr, vr, elementIndiceOfDependents, independents, elementIndiceOfInpendents, dependencyKind, ndependencies)
    checkStatus(c, status)
    status
end

# TODO not tested
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.11. Getting Partial Derivatives

This function computes the directional derivatives v_{sensitivity} = J ⋅ v_{seed} of an FMU.

unknowns - contains value references to the unknowns.

nUnknowns - contains the length of argument unknowns.

knowns - contains value references of the knowns.

nKnowns - contains the length of argument knowns.

seed - contains the components of the seed vector.

nSeed - contains the length of seed.

sensitivity - contains the components of the sensitivity vector.

nSensitivity - contains the length of sensitivity.

This function can only be called if the 'ProvidesDirectionalDerivatives' tag in the ModelDescription is set.
"""
function fmi3GetDirectionalDerivative!(c::FMU3Instance,
                                       unknowns::AbstractArray{fmi3ValueReference},
                                       nUnknowns::Csize_t,
                                       knowns::AbstractArray{fmi3ValueReference},
                                       nKnowns::Csize_t,
                                       seed::AbstractArray{fmi3Float64},
                                       nSeed::Csize_t,
                                       sensitivity::AbstractArray{fmi3Float64},
                                       nSensitivity::Csize_t)
    status = fmi3GetDirectionalDerivative!(c.fmu.cGetDirectionalDerivative,
          c.compAddr, unknowns, nUnknowns, knowns, nKnowns, seed, nSeed, sensitivity, nSensitivity)
    checkStatus(c, status)
    status

end

# TODO not tested
"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.11. Getting Partial Derivatives

This function computes the adjoint derivatives v^T_{sensitivity}= v^T_{seed} ⋅ J of an FMU.

unknowns - contains value references to the unknowns.

nUnknowns - contains the length of argument unknowns.

knowns - contains value references of the knowns.

nKnowns - contains the length of argument knowns.

seed - contains the components of the seed vector.

nSeed - contains the length of seed.

sensitivity - contains the components of the sensitivity vector.

nSensitivity - contains the length of sensitivity.

This function can only be called if the 'ProvidesAdjointDerivatives' tag in the ModelDescription is set.
"""
function fmi3GetAdjointDerivative!(c::FMU3Instance,
                                       unknowns::AbstractArray{fmi3ValueReference},
                                       nUnknowns::Csize_t,
                                       knowns::AbstractArray{fmi3ValueReference},
                                       nKnowns::Csize_t,
                                       seed::AbstractArray{fmi3Float64},
                                       nSeed::Csize_t,
                                       sensitivity::AbstractArray{fmi3Float64},
                                       nSensitivity::Csize_t)
    status = fmi3GetAdjointDerivative!(c.fmu.cGetAdjointDerivative,
          c.compAddr, unknowns, nUnknowns, knowns, nKnowns, seed, nSeed, sensitivity, nSensitivity)
    checkStatus(c, status)
    status

end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.2.12. Getting Derivatives of Continuous Outputs

Retrieves the n-th derivative of output values.

valueReferences - is a vector of value references that define the variables whose derivatives shall be retrieved. If multiple derivatives of a variable shall be retrieved, list the value reference multiple times.

nValueReferences - is the dimension of the arguments valueReferences and orders.

orders - contains the orders of the respective derivative (1 means the first derivative, 2 means the second derivative, …, 0 is not allowed).
If multiple derivatives of a variable shall be retrieved, provide a list of them in the orders array, corresponding to a multiply occurring value reference in the valueReferences array.
The highest order of derivatives retrievable can be determined by the 'maxOutputDerivativeOrder' tag in the ModelDescription.

values - is a vector with the values of the derivatives. The order of the values elements is derived from a twofold serialization: the outer level corresponds to the combination of a value reference (e.g., valueReferences[k]) and order (e.g., orders[k]), and the inner level to the serialization of variables as defined in Section 2.2.6.1. The inner level does not exist for scalar variables.

nValues - is the size of the argument values. nValues only equals nValueReferences if all corresponding output variables are scalar variables.
"""
function fmi3GetOutputDerivatives!(c::FMU3Instance, vr::AbstractArray{fmi3ValueReference}, nValueReferences::Csize_t, order::AbstractArray{fmi3Int32}, values::AbstractArray{fmi3Float64}, nValues::Csize_t)
    status = fmi3GetOutputDerivatives!(c.fmu.cGetOutputDerivatives,
               c.compAddr, vr, nValueReferences, order, values, nValues)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.2. State: Instantiated

If the importer needs to change structural parameters, it must move the FMU into Configuration Mode using fmi3EnterConfigurationMode.
"""
function fmi3EnterConfigurationMode(c::FMU3Instance)
    if c.state != fmi3InstanceStateInstantiated && (c.state != fmi3InstanceStateStepMode && fmi3IsCoSimulation(c.fmu)) && (c.state != fmi3InstanceEventMode && fmi3IsModelExchange(c.fmu))
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3EnterConfigurationMode(...): Called at the wrong time."
        end
    end

    status = fmi3EnterConfigurationMode(c.fmu.cEnterConfigurationMode,
    c.compAddr)
    checkStatus(c, status)
    if status == fmi3StatusOK
        if c.state == fmi3InstanceStateInstantiate
            c.state = fmi3InstanceStateConfigurationMode
        else
            c.state = fmi3InstanceStateReconfigurationMode
        end
    end
    return status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.6. State: Configuration Mode

Exits the Configuration Mode and returns to state Instantiated.
"""
function fmi3ExitConfigurationMode(c::FMU3Instance)
    if c.state != fmi3InstanceStateConfigurationMode && c.state != fmi3InstanceStateReconfigurationMode
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3ExitConfigurationMode(...): Called at the wrong time."
        end
    end

    status = fmi3ExitConfigurationMode(c.fmu.cExitConfigurationMode,
         c.compAddr)
    checkStatus(c, status)
    if status == fmi3StatusOK
        if c.state == fmi3InstanceStateConfigurationMode
            c.state = fmi3InstanceStateInstantiate
        elseif fmi3IsCoSimulation(c.fmu)
            c.state = fmi3InstanceStateStepMode
        elseif fmi3IsModelExchange(c.fmu)
            c.state = fmi3InstanceStateEventMode
        elseif fmi3IsScheduledExecution(c.fmu)
            c.state = fmi3InstanceStateClockActivationMode
        end
    end
    return status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.2. State: Instantiated

This function returns the number of continuous states.
This function can only be called in Model Exchange.

fmi3GetNumberOfContinuousStates must be called after a structural parameter is changed. As long as no structural parameters changed, the number of states is given in the modelDescription.xml, alleviating the need to call this function.
"""
function fmi3GetNumberOfContinuousStates!(c::FMU3Instance, nContinuousStates::Ref{Csize_t})
    status = fmi3GetNumberOfContinuousStates!(c.fmu.cGetNumberOfContinuousStates,
           c.compAddr, nContinuousStates)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.2. State: Instantiated

This function returns the number of event indicators.
This function can only be called in Model Exchange.

fmi3GetNumberOfEventIndicators must be called after a structural parameter is changed. As long as no structural parameters changed, the number of states is given in the modelDescription.xml, alleviating the need to call this function.
"""
function fmi3GetNumberOfEventIndicators!(c::FMU3Instance, nEventIndicators::Ref{Csize_t})
    status = fmi3GetNumberOfEventIndicators!(c.fmu.cGetNumberOfEventIndicators,
            c.compAddr, nEventIndicators)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.3. State: Initialization Mode

Return the states at the current time instant.

This function must be called if fmi3UpdateDiscreteStates returned with valuesOfContinuousStatesChanged == fmi3True. Not allowed in Co-Simulation and Scheduled Execution.
"""
function fmi3GetContinuousStates!(c::FMU3Instance, nominals::AbstractArray{fmi3Float64}, nContinuousStates::Csize_t)
    status = fmi3GetContinuousStates!(c.fmu.cGetContinuousStates,
            c.compAddr, nominals, nContinuousStates)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.3. State: Initialization Mode

Return the nominal values of the continuous states.

If fmi3UpdateDiscreteStates returned with nominalsOfContinuousStatesChanged == fmi3True, then at least one nominal value of the states has changed and can be inquired with fmi3GetNominalsOfContinuousStates.
Not allowed in Co-Simulation and Scheduled Execution.
"""
function fmi3GetNominalsOfContinuousStates!(c::FMU3Instance, x_nominal::AbstractArray{fmi3Float64}, nx::Csize_t)
    status = fmi3GetNominalsOfContinuousStates!(c.fmu.cGetNominalsOfContinuousStates,
                    c.compAddr, x_nominal, nx)
    checkStatus(c, status)
    status
end

# TODO not testable not supported by FMU
"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.3. State: Initialization Mode

This function is called to trigger the evaluation of fdisc to compute the current values of discrete states from previous values.
The FMU signals the support of fmi3EvaluateDiscreteStates via the capability flag providesEvaluateDiscreteStates.
"""
function fmi3EvaluateDiscreteStates(c::FMU3Instance)
    status = fmi3EvaluateDiscreteStates(c.fmu.cEvaluateDiscreteStates,
            c.compAddr)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.5. State: Event Mode

This function is called to signal a converged solution at the current super-dense time instant. fmi3UpdateDiscreteStates must be called at least once per super-dense time instant.
"""
function fmi3UpdateDiscreteStates(c::FMU3Instance, discreteStatesNeedUpdate::Ref{fmi3Boolean}, terminateSimulation::Ref{fmi3Boolean},
                                    nominalsOfContinuousStatesChanged::Ref{fmi3Boolean}, valuesOfContinuousStatesChanged::Ref{fmi3Boolean},
                                    nextEventTimeDefined::Ref{fmi3Boolean}, nextEventTime::Ref{fmi3Float64})
    status = fmi3UpdateDiscreteStates(c.fmu.cUpdateDiscreteStates,
            c.compAddr, discreteStatesNeedUpdate, terminateSimulation, nominalsOfContinuousStatesChanged, valuesOfContinuousStatesChanged, nextEventTimeDefined, nextEventTime)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.5. State: Event Mode

The model enters Continuous-Time Mode and all discrete-time equations become inactive and all relations are “frozen”.
This function has to be called when changing from Event Mode (after the global event iteration in Event Mode over all involved FMUs and other models has converged) into Continuous-Time Mode.
"""
function fmi3EnterContinuousTimeMode(c::FMU3Instance; soft::Bool=false)
    if c.state != fmi3InstanceStateEventMode
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3EnterContinuousTimeMode(...): Needs to be called in state `fmi3InstanceStateEventMode`."
        end
    end

    status = fmi3EnterContinuousTimeMode(c.fmu.cEnterContinuousTimeMode,
          c.compAddr)
    checkStatus(c, status)
    if status == fmi3StatusOK
        c.state = fmi3InstanceStateContinuousTimeMode
    end
    return status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 2.3.5. State: Event Mode

This function must be called to change from Event Mode into Step Mode in Co-Simulation (see 4.2.).
"""
function fmi3EnterStepMode(c::FMU3Instance)
    if c.state != fmi3InstanceStateEventMode
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3EnterContinuousTimeMode(...): Needs to be called in state `fmi3InstanceStateEventMode`."
        end
    end

    status = fmi3EnterStepMode(c.fmu.cEnterStepMode,
    c.compAddr)
    checkStatus(c, status)
    if status == fmi3StatusOK
        c.state = fmi3InstanceStateStepMode
    end
    return status

end

"""
Source: FMISpec3.0, Version D5ef1c1: 3.2.1. State: Continuous-Time Mode

Set a new time instant and re-initialize caching of variables that depend on time, provided the newly provided time value is different to the previously set time value (variables that depend solely on constants or parameters need not to be newly computed in the sequel, but the previously computed values can be reused).
"""
function fmi3SetTime(c::FMU3Instance, time::fmi3Float64)
    c.t = time
    status = fmi3SetTime(c.fmu.cSetTime,
          c.compAddr, time)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 3.2.1. State: Continuous-Time Mode

Set a new (continuous) state vector and re-initialize caching of variables that depend on the states. Argument nx is the length of vector x and is provided for checking purposes
"""
function fmi3SetContinuousStates(c::FMU3Instance,
                                 x::AbstractArray{fmi3Float64},
                                 nx::Csize_t)
    status = fmi3SetContinuousStates(c.fmu.cSetContinuousStates,
         c.compAddr, x, nx)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 3.2.1. State: Continuous-Time Mode

Compute first-oder state derivatives at the current time instant and for the current states.
"""
function fmi3GetContinuousStateDerivatives(c::FMU3Instance,
                            derivatives::AbstractArray{fmi3Float64},
                            nx::Csize_t)
    status = fmi3GetContinuousStateDerivatives(c.fmu.cGetContinuousStateDerivatives,
          c.compAddr, derivatives, nx)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 3.2.1. State: Continuous-Time Mode

Compute event indicators at the current time instant and for the current states. EventIndicators signal Events by their sign change.
"""
function fmi3GetEventIndicators!(c::FMU3Instance, eventIndicators::AbstractArray{fmi3Float64}, ni::Csize_t)
    status = fmi3GetEventIndicators!(c.fmu.cGetEventIndicators,
                   c.compAddr, eventIndicators, ni)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 3.2.1. State: Continuous-Time Mode

This function must be called by the environment after every completed step of the integrator provided the capability flag needsCompletedIntegratorStep = true.
If enterEventMode == fmi3True, the event mode must be entered
If terminateSimulation == fmi3True, the simulation shall be terminated
"""
function fmi3CompletedIntegratorStep!(c::FMU3Instance,
                                      noSetFMUStatePriorToCurrentPoint::fmi3Boolean,
                                      enterEventMode::Ref{fmi3Boolean},
                                      terminateSimulation::Ref{fmi3Boolean})
    status = fmi3CompletedIntegratorStep!(c.fmu.cCompletedIntegratorStep,
         c.compAddr, noSetFMUStatePriorToCurrentPoint, enterEventMode, terminateSimulation)
    checkStatus(c, status)
    status
end

"""
Source: FMISpec3.0, Version D5ef1c1: 3.2.1. State: Continuous-Time Mode

The model enters Event Mode from the Continuous-Time Mode in ModelExchange oder Step Mode in CoSimulation and discrete-time equations may become active (and relations are not “frozen”).
"""
function fmi3EnterEventMode(c::FMU3Instance, stepEvent::fmi3Boolean, stateEvent::fmi3Boolean, rootsFound::AbstractArray{fmi3Int32}, nEventIndicators::Csize_t, timeEvent::fmi3Boolean; soft::Bool=false)
    if c.state != fmi3InstanceStateContinuousTimeMode && c.state != fmi3InstanceStateStepMode
        if soft
            return fmi3StatusOK
        else
            @warn "fmi3EnterEventMode(...): Called at the wrong time."
        end
    end

    status =  fmi3EnterEventMode(c.fmu.cEnterEventMode,
    c.compAddr, stepEvent, stateEvent, rootsFound, nEventIndicators, timeEvent)
    checkStatus(c, status)
    if status == fmi3StatusOK
        c.state = fmi3InstanceStateEventMode
    end
    return status

end

"""
Source: FMISpec3.0, Version D5ef1c1: 4.2.1. State: Step Mode

The computation of a time step is started.
"""
function fmi3DoStep!(c::FMU3Instance, currentCommunicationPoint::fmi3Float64, communicationStepSize::fmi3Float64, noSetFMUStatePriorToCurrentPoint::fmi3Boolean,
                    eventEncountered::Ref{fmi3Boolean}, terminateSimulation::Ref{fmi3Boolean}, earlyReturn::Ref{fmi3Boolean}, lastSuccessfulTime::Ref{fmi3Float64})
    status = fmi3DoStep!(c.fmu.cDoStep,
          c.compAddr, currentCommunicationPoint, communicationStepSize, noSetFMUStatePriorToCurrentPoint, eventEncountered, terminateSimulation, earlyReturn, lastSuccessfulTime)
    checkStatus(c, status)
    status
end

ext: jl


text:
#!/usr/bin/julia

# Trizen
# Date: 28 August 2020
# https://github.com/trizen

# A new algorithm for generating Fermat super-pseudoprimes to base 2.

# OEIS:
#   https://oeis.org/A050217 -- Super-Poulet numbers: Poulet numbers whose divisors d all satisfy d|2^d-2.
#   https://oeis.org/A214305 -- Fermat pseudoprimes to base 2 with two prime factors.

# See also:
#   https://en.wikipedia.org/wiki/Fermat_pseudoprime
#   https://trizenx.blogspot.com/2020/08/pseudoprimes-construction-methods-and.html

using Primes
using Combinatorics

function divisors(n)

    d = Int64[1]

    for (p,e) in factor(n)
        t = Int64[]
        r = 1

        for i in 1:e
            r *= p
            for u in d
                push!(t, u*r)
            end
        end

        append!(d, t)
    end

    return d
end

function fermat_pseudoprimes(limit)

    table = Dict{Int64, Array{Int64}}()

    for p in primes(3,limit)
        for d in (divisors(p-1))
            if (powermod(2, d, p) == 1)
                if (haskey(table, d))
                    push!(table[d], p)
                else
                    table[d] = [p]
                end
            end
        end
    end

    for (k,v) in table
        for k in 2:length(v)
            for c in combinations(v,k)
                println(prod(map(x -> BigInt(x), c)))
            end
        end
    end

    return true
end

fermat_pseudoprimes(10^5)

max_stars_repo_name: trizen/julia-scripts

max_issues_repo_name: trizen/julia-scripts

max_forks_repo_name: trizen/julia-scripts


text:
macro c(ex)
    Meta.isexpr(ex, :call) || throw(ArgumentError("not a function call expression."))
    prologue = Expr(:block)
    epilogue = Expr(:block)
    func_expr = Expr(:call, first(ex.args))
    for arg in ex.args[2:end]
        if Meta.isexpr(arg, :&)
            refee = arg.args[]
            if Meta.isexpr(refee, :ref) && length(refee.args) == 2
                # &a[n] => pointer(a) + n * Core.sizeof(eltype(a))
                array_name, n = refee.args
                push!(func_expr.args, :(pointer($array_name) + $n * Core.sizeof(eltype($array_name))))
            else
                ref_sym = gensym("cref")
                push!(prologue.args, Expr(:(=), ref_sym, Expr(:call, :Ref, refee)))
                push!(func_expr.args, ref_sym)
                push!(epilogue.args, Expr(:(=), refee, Expr(:ref, ref_sym)))
            end
        else
            push!(func_expr.args, arg)
        end
    end
    func_ret = gensym("cref_ret")
    push!(prologue.args, Expr(:(=), func_ret, func_expr))
    append!(prologue.args, epilogue.args)
    push!(prologue.args, func_ret)
    return esc(prologue)
end


text:
# organize/sort this later.

"""
Info on run_optim usage: minf, minx, ret, N_evals = run_optim(p_β)
"""
function setupcostnestedλd(Es::Vector{NMRSpectraSimulator.κCompoundFIDType{T, SST}},
    As::Vector{NMRSpectraSimulator.CompoundFIDType{T, SST}},
    fs::T,
    SW::T,
    LS_inds,
    U_rad_cost,
    y_cost::Vector{Complex{T}},
    shift_constants;
    w = ones(T, length(Es)),
    optim_algorithm = :GN_DIRECT_L,
    κ_lb_default = 0.2,
    κ_ub_default = 5.0,
    max_iters = 500,
    xtol_rel = 1e-9,
    ftol_rel = 1e-9,
    maxtime = Inf) where {T <: Real, SST}

    # model.
    f = uu->NMRSpectraSimulator.evalitpproxymixture(uu, Es; w = w)

    ##### update functions.
    N_d = sum( getNd(As[n]) for n = 1:length(As) )
    N_λ = sum( getNλ(As[n]) for n = 1:length(As) )

    st_ind_d = 1
    fin_ind_d = st_ind_d + N_d - 1
    updatedfunc = pp->updatemixtured!(As, pp, st_ind_d, fs, SW, shift_constants)

    #λupdate.
    st_ind_λ = fin_ind_d + 1
    fin_ind_λ = st_ind_λ + N_λ -1
    updateλfunc = pp->updateλ!(As, pp, st_ind_λ)

    N_vars_set = [N_d; N_λ]

    # β, κ update.
    run_optim, obj_func_β, E_BLS, κ_BLS, b_BLS, updateβfunc, updateκfunc,
    q_β = setupβLSsolver(optim_algorithm,
        Es, As, LS_inds, U_rad_cost, y_cost;
        κ_lb_default = κ_lb_default,
        κ_ub_default = κ_ub_default,
        max_iters = max_iters,
        xtol_rel = xtol_rel,
        ftol_rel = ftol_rel,
        maxtime = maxtime)

    #### extract parameters from p.
    getshiftfunc = pp->pp[st_ind_d:fin_ind_d]
    getλfunc = pp->pp[st_ind_λ:fin_ind_λ]

    return f, updatedfunc, updateλfunc, getshiftfunc, getλfunc, N_vars_set,
    run_optim, obj_func_β, E_BLS, κ_BLS, b_BLS, updateβfunc, q_β
end

# function setupcostnestedλdwarp(Es::Vector{NMRSpectraSimulator.κCompoundFIDType{T, SST}},
#     Bs,
#     As,
#     fs::T,
#     SW::T,
#     LS_inds,
#     U_rad_cost,
#     y_cost::Vector{Complex{T}},
#     Δsys_cs::Vector{Vector{T}},
#     itp_a,
#     itp_b;
#     w = ones(T, length(Es)),
#     optim_algorithm = :GN_DIRECT_L,
#     κ_lb_default = 0.2,
#     κ_ub_default = 5.0,
#     max_iters = 500,
#     xtol_rel = 1e-9,
#     ftol_rel = 1e-9,
#     maxtime = Inf) where {T <: Real, SST}
#
#     # model.
#     f = uu->NMRSpectraSimulator.evalitpproxymixture(uu, As, Es; w = w)
#
#     ##### update functions.
#     N_d = sum( getNd(Bs[n]) for n = 1:length(Bs) )
#     N_λ = sum( getNλ(Bs[n]) for n = 1:length(Bs) )
#
#     st_ind_d = 1
#     fin_ind_d = st_ind_d + N_d - 1
#     updatedfunc = pp->updatemixturedwarp!(Bs, pp, st_ind_d, fs, SW,
#         Δsys_cs, itp_a, itp_b)
#
#     #λupdate.
#     st_ind_λ = fin_ind_d + 1
#     fin_ind_λ = st_ind_λ + N_λ -1
#     updateλfunc = pp->updateλ!(Bs, pp, st_ind_λ)
#
#     N_vars_set = [N_d; N_λ]
#
#     # β, κ update.
#     run_optim, obj_func_β, E_BLS, κ_BLS, b_BLS, updateβfunc,
#     q_β = setupβLSsolver(optim_algorithm,
#         Es, Bs, As, LS_inds, U_rad_cost, y_cost,
#         κs_β_DOFs, κs_β_orderings;
#         κ_lb_default = κ_lb_default,
#         κ_ub_default = κ_ub_default,
#         max_iters = max_iters,
#         xtol_rel = xtol_rel,
#         ftol_rel = ftol_rel,
#         maxtime = maxtime)
#
#     #### extract parameters from p.
#     getshiftfunc = pp->pp[st_ind_d:fin_ind_d]
#     getλfunc = pp->pp[st_ind_λ:fin_ind_λ]
#
#     return f, updatedfunc, updateλfunc, getshiftfunc, getλfunc, N_vars_set,
#     run_optim, obj_func_β, E_BLS, κ_BLS, b_BLS, updateβfunc, q_β
# end

function costnestedλd(U,
    S_U::Vector{Complex{T}},
    updatedfunc, updateλfunc,
    p::Vector{T}, Es, As, f,
    run_optim_β_κ::Function,
    E_BLS::Matrix{Complex{T}}, κ_BLS::Vector{T}, b_BLS,
    p_β::Vector{T})::T where T <: Real

    updatedfunc(p)
    updateλfunc(p)

    ### minimize inner problem.
    fill!(p_β, zero(T)) # always start from 0-phase?
    minf, minx, ret, N_evals = run_optim_β_κ(p_β)
    p_β[:] = minx

    # ensure Es/As is updated with the latest β and κ.
    updateβ!(As, minx, 1)
    parseκ!(Es, κ_BLS)

    # evaluate cost.
    cost = norm(reinterpret(T, E_BLS)*κ_BLS - b_BLS)^2

    return cost
end

##### only d. calibrate.

function setupcostnesteddwarp(Es,
    Bs,
    As,
    fs::T,
    SW::T,
    LS_inds,
    U_rad_cost,
    y_cost::Vector{Complex{T}},
    Δsys_cs::Vector{Vector{T}},
    itp_a,
    itp_b,
    κs_β_DOFs, κs_β_orderings;
    w = ones(T, length(Es)),
    optim_algorithm = :GN_DIRECT_L,
    κ_lb_default = 0.2,
    κ_ub_default = 5.0,
    max_iters = 500,
    xtol_rel = 1e-9,
    ftol_rel = 1e-9,
    maxtime = Inf) where T <: Real

    # model.
    f = uu->NMRSpectraSimulator.evalitpproxymixture(uu, As, Es; w = w)

    ##### update functions.
    #N_β = sum( getNβ(κs_β_DOFs[n], Bs[n]) for n = 1:length(Bs) )
    N_d = sum( getNd(Bs[n]) for n = 1:length(Bs) )

    st_ind_d = 1
    fin_ind_d = st_ind_d + N_d - 1
    updatedfunc = pp->updatemixturedwarp!(Bs, pp, st_ind_d, fs, SW,
        Δsys_cs, itp_a, itp_b)

    N_vars_set = [N_d; ]

    # β, κ update.
    run_optim, obj_func_β, E_BLS, κ_BLS, b_BLS, updateβfunc, updateκfunc,
    q_β = setupβLSsolver(optim_algorithm,
        Es, Bs, As, LS_inds, U_rad_cost, y_cost,
        κs_β_DOFs, κs_β_orderings;
        κ_lb_default = κ_lb_default,
        κ_ub_default = κ_ub_default,
        max_iters = max_iters,
        xtol_rel = xtol_rel,
        ftol_rel = ftol_rel,
        maxtime = maxtime)

    #### extract parameters from p.
    getshiftfunc = pp->pp[st_ind_d:fin_ind_d]

    return f, updatedfunc, getshiftfunc, N_vars_set,
    run_optim, obj_func_β, E_BLS, κ_BLS, b_BLS, updateβfunc, updateκfunc, q_β
end

##### cost.
"""
Note that NLopt handle exceptions gracefully. Check the termination status to see if an exception occurs (would return FORCED_STOP)
"""
function costnestedd(U,
    S_U::Vector{Complex{T}},
    updatedfunc, updateκfunc,
    p::Vector{T}, Es, Bs, #κs_β_orderings, κs_β_DOFs, f,
    run_optim_β_κ::Function,
    E_BLS::Matrix{Complex{T}}, κ_BLS::Vector{T}, b_BLS,
    p_β::Vector{T})::T where T <: Real

    updatedfunc(p)
    # try run_optim_β_κ(p_β)
    # catch e
    #     println("error in run_optim_β_κ(p_β)")
    # 	bt = catch_backtrace()
    # 	showerror(stdout, e, bt)
    # 	rethrow(e)
    # end

    ### minimize inner problem.
    fill!(p_β, zero(T)) # always start from 0-phase?
    minf, minx, ret, N_evals = run_optim_β_κ(p_β)
    p_β[:] = minx

    # ensure Es/As is updated with the latest β and κ.
    #updateβ!(Bs, κs_β_orderings, κs_β_DOFs, minx, 1)
    updateβ!(Bs, minx, 1)
    updateκfunc(0.0)
    parseκ!(Es, κ_BLS)

    # evaluate cost.
    cost = norm(reinterpret(T, E_BLS)*κ_BLS - b_BLS)^2

    return cost
end

##### only d. quantify.

function setupcostnesteddwarpw(Es,
    Bs,
    As,
    fs::T,
    SW::T,
    LS_inds,
    U_rad_cost,
    y_cost::Vector{Complex{T}},
    Δsys_cs::Vector{Vector{T}},
    a_setp, b_setp,
    κs_β_DOFs, κs_β_orderings;
    β_optim_algorithm = :GN_DIRECT_L,
    w_lb_default = 0.2,
    w_ub_default = 5.0,
    β_max_iters = 500,
    β_xtol_rel = 1e-9,
    β_ftol_rel = 1e-9,
    β_maxtime = Inf) where T <: Real

    ##### update functions.
    N_β = sum( getNβ(κs_β_DOFs[n], Bs[n]) for n = 1:length(Bs) )
    N_d = sum( getNd(Bs[n]) for n = 1:length(Bs) )

    st_ind_d = 1
    fin_ind_d = st_ind_d + N_d - 1
    updatedfunc = pp->updatemixturedwarp!(Bs, pp, st_ind_d, fs, SW,
        Δsys_cs, a_setp, b_setp)

    N_vars_set = [N_d; ]

    # β, κ update.
    run_optim, obj_func_β, E_BLS, w_BLS, b_BLS, updateβfunc, updatewfunc,
    q_β = setupβLSsolverw(β_optim_algorithm,
        Es, Bs, As, LS_inds, U_rad_cost, y_cost, κs_β_DOFs, κs_β_orderings;
        w_lb_default = w_lb_default,
        w_ub_default = w_ub_default,
        β_max_iters = β_max_iters,
        β_xtol_rel = β_xtol_rel,
        β_ftol_rel = β_ftol_rel,
        β_maxtime = β_maxtime)

    #### extract parameters from p.
    getshiftfunc = pp->pp[st_ind_d:fin_ind_d]

    # model. (optional)
    #f = uu->NMRSpectraSimulator.evalitpproxymixture(uu, As, Es; w = w_BLS)
    f = q_β

    return f, updatedfunc, getshiftfunc, N_vars_set,
    run_optim, obj_func_β, E_BLS, w_BLS, b_BLS, updateβfunc, updatewfunc, q_β
end

##### cost.
"""
Note that NLopt handle exceptions gracefully. Check the termination status to see if an exception occurs (would return FORCED_STOP)
"""
function costnesteddw(U,
    S_U::Vector{Complex{T}},
    updatedfunc,
    updatewfunc,
    p::Vector{T}, Es, Bs, κs_β_orderings, κs_β_DOFs, #f,
    run_optim_β_κ::Function,
    E_BLS::Matrix{Complex{T}}, w_BLS::Vector{T}, b_BLS,
    p_β::Vector{T})::T where T <: Real

    updatedfunc(p)

    ### minimize inner problem.
    fill!(p_β, zero(T)) # always start from 0-phase?
    minf, minx, ret, N_evals = run_optim_β_κ(p_β)
    p_β[:] = minx # take out?

    # ensure Es/As is updated with the latest β.
    #updateβ!(Bs, κs_β_orderings, κs_β_DOFs, minx, 1)
    updateβ!(Bs, minx, 1)
    updatewfunc(1.0)

    # evaluate cost.
    cost = norm(reinterpret(T, E_BLS)*w_BLS - b_BLS)^2

    return cost
end


text:
using Metatheory
using Test
using Metatheory.Library
using Metatheory.EGraphs
using Metatheory.EGraphs.Schedulers

t = @theory begin
    f(a...) |> (println.(a); 42)
end

dump(t)

g = EGraph(:(f(1,2,3)))
saturate!(g, t)


text:
module DryRun

using Mocking

export @dryer

include("dryer.jl")

# Create the initial definition of `activated` which defaults DryRun to be disabled
activated() = false

"""
    DryRun.activate()

Enable `@mock` call sites to allow for calling patches instead of the original function.
"""
function activate()
    # Avoid redefining `activated` when it's already set appropriately
    !activated() && @eval activated() = true
    return nothing
end

"""
    DryRun.deactivate()

Disable `@mock` call sites to only call the original function.
"""
function deactivate()
    # Avoid redefining `activated` when it's already set appropriately
    activated() && @eval activated() = false
    return nothing
end

const NULLIFIED = Ref{Bool}(false)

"""
    DryRun.nullify()

Force any packages loaded after this point to treat the `@mock` macro as a no-op. Doing so
will maximize performance by eliminating any runtime checks taking place at the `@mock` call
sites but will break any tests that require patches to be applied.

Note to ensure that all `@mock` macros are inoperative be sure to call this function before
loading any packages which depend on DryRun.jl.
"""
function nullify()
    global NULLIFIED[] = true
    return nothing
end

end # module


text:
#Algorithm 12.2
#Heat equation with backwards-difference
#Boundary conditions u(0, t) = u(l, t) = 0, 0 < t < T

using Printf

function heatbd(f, L, T, α, m, N)
	w = zeros(Float64, m)
	l = zeros(Float64, m)
	u = zeros(Float64, m)
	z = zeros(Float64, m)
	h = L/m
	k = T/N
	λ = k*α^2/h^2
	for i ∈ 1:(m - 1)
		w[i] = f(i * h)
	end
	l[1] = 1 + 2λ
	u[1] = -λ/l[1]
	for i ∈ 2:(m - 2)
		l[i] = 1 + 2λ + λ*u[i - 1]
		u[i] = -λ/l[i]
	end
	l[m - 1] = 1 + 2λ + λ*u[m - 2]
	for j ∈ 1:N
		t = j*k
		z[1] = w[1] / l[1]
		for i ∈ 2:(m - 1)
			z[i] = (w[i] + λ*z[i - 1])/l[i]
		end
		w[m - 1] = z[m - 1]
		for i ∈ (m - 2):1
			w[i] = z[i][t] - u[i]w[i+1]
		end
		for i ∈ 1:(m - 1)
			x = i*h
			println(x, w[i])
			x = i * h
			@printf("%3d %11.8f %14.8f\n", i, x, w[i])
		end
	end
end

function distrib(x)
	return sin(π*x)
end

heatbd(distrib, 1, 10, 0.5, 100, 100)


text:
# Inspired by:
#
# BW Bader and TG Kolda. Efficient MATLAB Computations with Sparse
# and Factored Tensors, SIAM J Scientific Computing 30:205-231, 2007.
# DOI: 10.1137/060676489. [BibTeX]

export SparseTensor, ttt, sptenrand

import SparseArrays

mutable struct SparseTensor{T,N} <: AbstractArray{T,N}
    dict::Dict{Dims{N},T}
    # dict::Dict{Tuple{Int,Vararg{Int}},T} # Julia doesn't like this
    dims::Dims{N}
    function SparseTensor(v::Array{E,1},subs::Array{Int,2},dims="auto") where E <: Number
        newsubs = unique(subs,dims=1)

        # This feels like it should be slow
        dict = Dict(s => zero(eltype(v)) for s in Tuple.(eachrow(newsubs)))
        for i in 1:length(v)
            dict[Tuple(subs[i,:])] += v[i]
        end

        #dims = CartesianIndex(maximum(newsubs, dims=1)...)
        SparseTensor(dict,dims == "auto" ? maximum(subs,dims=1) |> Tuple : dims)
    end
    SparseTensor(d::Dict,dims) = new{d |> values |> eltype, length(dims)}(d,dims)
    SparseTensor(d::Dict) = SparseTensor(d,maximum(_indarray(d),dims=1) |> Tuple)
    SparseTensor(a::AbstractArray) = SparseTensor(Dict(Tuple(ind) => a[ind] for ind in findall(x->x!==zero(x),a)),size(a))
end

function _indarray(d::Dict)
    reduce(hcat,d |> keys .|> collect) |> permutedims
end
_indarray(t::SparseTensor) = _indarray(t.dict)

Base.size(t::SparseTensor) = t.dims

# TODO: add error checking - out of bounds etc.
function Base.getindex(A::SparseTensor{T,N}, inds::Vararg{Int,N}) where {N,T}
    get(A.dict,inds,zero(eltype(A)))
end
# Bug: Julia asks for (1,1) indices for rank-1 tensors - this is a workaround but not a pretty one
# and this doesn't seem to do anything
# Base.IndexStyle(::SparseTensor{T,1}) where T = IndexLinear()
function Base.getindex(A::SparseTensor{T,1}, inds::Vararg{Int,N}) where {N,T}
    get(A.dict,(inds[1],),zero(eltype(A)))
end

function Base.setindex!(A::SparseTensor, value, inds::Vararg{Int,N}) where N
    A.dict[inds] = value
end

# If you need exactly n non-zeroes, provide subrng: n,dims -> Array{Int,2}
# E.g.
# ```
# subrng=(n,dims)->begin
#     reduce(hcat, StatsBase.sample(
#         Iterators.product([range(1;length=e) for e in dims]...) |> collect,
#         n;
#         replace=false
#     ) .|> collect) |> permutedims
# end
# ``` (NB: very slow for large numbers of dimensions - ideally we want a sampler that can operate on an iterator)
function sptenrand(dims::Array{Int,1},n::Int;rng=rand,subrng=_subrng)
    subs = subrng(n,dims)
    vals = rng(size(subs)[1])
    SparseTensor(vals,subs,dims|>Tuple)
end

# Generates n indices into dims and discards repeats
_subrng(n,dims::Array{Int,1}) = unique(round.((rand(Float64,n,length(dims))) .* (dims' .-1)) .+ 1 .|> Int, dims=1)

sptenrand(dims,d::AbstractFloat;kwargs...) = sptenrand(dims,d*prod(dims)|>round|>Int;kwargs...)

# Broadcasting - inspired by https://github.com/andyferris/Dictionaries.jl/blob/9b22a254683260354fda8e32291727dfad040106/src/broadcast.jl#L94
# See also: https://docs.julialang.org/en/v1/manual/interfaces/index.html
Base.Broadcast.broadcastable(t::SparseTensor) = t

struct SparseTensorBroadcastStyle{N} <: Base.Broadcast.AbstractArrayStyle{N}
end

Base.Broadcast.BroadcastStyle(::Type{<:SparseTensor{T,N}}) where {T,N} = SparseTensorBroadcastStyle{N}()

# Fall back to dense array if one of the broadcast arguments is a dense array
Base.Broadcast.BroadcastStyle(::SparseTensorBroadcastStyle{D}, a::Base.Broadcast.DefaultArrayStyle{N}) where {D,N} = Base.Broadcast.DefaultArrayStyle{max(N,D)}()

# TODO:
# - sparse .* dense = sparse
# - sparse .+ dense = dense, but it can be sped up by being s .+ d = begin c = copy(d); for (k,v) in s.dict; c[k...] .+ v; end
# (need to make own SparseDense broadcast style to replace the DefaultArrayStyle above and then make the function below specialise on f={Base.:*,Base.:+}

function Base.Broadcast.broadcasted(::SparseTensorBroadcastStyle, f, args...)
    # max() currently fails for (2,2,2),(3,1,3) - should return (3,2,3) or maybe just an error
    SparseTensor(merge_default(f,getfield.(args,:dict)...),max(getfield.(args,:dims)...))
end

" Merge two dictionaries such that ans[k] = f((get(l,k,0),get(r,k,0)) for all k in keys(l) ∪ keys(r)"
function merge_default(f,l,r)
    merged = promote_type(typeof(l),typeof(r))()
    for k in keys(l)
        merged[k] = f(l[k],get(r,k,zero(valtype(r))))
    end
    for k in keys(r)
        haskey(merged, k) && continue
        merged[k] = f(zero(valtype(l)),r[k])
    end
    merged
end

# I'm not sure what these do but they're not good
# Base.Broadcast.materialize(t::SparseTensor) = copy(t)
# Base.Broadcast.materialize!(out::SparseTensor, t::SparseTensor) = copyto!(out, t)

# Base.:* - see ttt.m
# Outer product: (l⊗r)_{i_1,...,i_n,j_1,...,j_n} = l_{i_1,...,i_n} r_{j_1,...,j_n} (i.e. not the Kronecker product)
function ttt(l::SparseTensor, r::SparseTensor)
    SparseTensor(Dict((a..., b...) => l[a...] * r[b...] for (a, b) in Iterators.product(l.dict |> keys,r.dict |> keys)),(l.dims...,r.dims...))
end

function Base.Array(t::SparseTensor)
    a = zeros(size(t))
    for (k,v) in t.dict
        a[k...] = v
    end
    a
end

# To implement for ttv:
# - permutedims
# - reshape
# - sparse-vector multiplication
# function ttv2(X::AbstractArray{<:Number,N},V::VectorCell,modes::AbstractVector{<:Integer}) where N
#   remmodes=setdiff(1:N,modes)'
#   if N > 1
#     X=permutedims(X,[remmodes modes'])
#     @show [remmodes modes']
#     @show X
#   end
#   sz=size(X)
#   @show sz
#   if length(modes) < length(V)
#     V=V[modes]
#   end
#   M=N
#   for n=length(modes):-1:1
#     X=reshape(X,prod(sz[1:M-1]),sz[M])
#     @show X
#     X=X*V[n]
#     @show X
#     M-=1
#   end
#   if M>0
#     X=reshape(X,sz[1:M])
#   end
#   X
# end

# TODO:
# - consider automatically backing off to CSC for all SparseTensor{T,2}
function SparseArrays.sparse(t::SparseTensor{T,2}) where T
    i = _indarray(t)
    sparse(i[:,1],i[:,2],t.dict|>values|>collect,size(t)...)
end

Base.:*(t::SparseTensor{T,2},v::AbstractArray{T2,1}) where {T,T2} = begin
    SparseTensor(SparseArrays.sparse(t)*v)
end

Base.promote_rule(::Type{SparseTensor{T1,D}},::Type{Array{T2,D}}) where {T1,T2,D} = SparseTensor{promote_type(T1,T2),D}

Base.permutedims(t::SparseTensor,dims) = SparseTensor(Dict(k[vec(dims)] => v for (k,v) in t.dict), t.dims)

# heart of ttv for s(2,2,2) and v(2) is thus:
# julia> reshape(reshape(permutedims(s, [2 3 1]),2*2,2)*v,2,2)
#
# which uses dense arrays. Reshape has a special type itself.
#
# TS = promote_op(matprod, T, S) # where matprod = x*y + x*y
# julia> Base.promote_op(+,SparseTensor{Int,2},Array{Int,2})
# Array{Int64,2}
#
# Even though we have the promotion rule above. It uses julia> Core.Compiler.return_type(+,Tuple{T1,T2}) internally.
# Ideally we should be able to avoid it altoghter.

mutable struct ReshapedSparseTensor{T,D}  # <: AbstractSparseTensor
    # Will need custom get/set index
    underlying::SparseTensor{T,D}
    dims::Dims{D}
end

function Base.reshape(t::SparseTensor,dims::Dims)
end

# TODO: Support slices / Colon()

# TODO: stop implementing these myself. I should just write a TensorToolbox.jl compatible type and then optimise when I can be bothered. _headdesk_
# (bonus: it would let me check my maths)


text:
# This file investigates integration error with respect to free variables in hidden 1D integration

using FractalTools
using Makie

# Construct interpolation data
f(x) = sin(2π * x)  + 1
g(x) = cos(2π * x)  + 1
dx   = 0.1
xi   = 0.
xf   = 1.
x    = collect(xi : dx : xf)
y    = f.(x)
z    = g.(x)
pts  = collect.(zip(x, y, z))

# Compute errors
fvals = (cos(2π * xi) - cos(2π * xf)) / (2π) + (xf - xi)
freevars = 0.001 : 0.01 : 0.999
mse = map(freevars) do freevar
    ival = integrate(pts, HInterp1D(fill(freevar, 2, 2)))[1]
    abs(fval - ival)
end

# Plot mse
fig = Figure()
ax = fig[1, 1] = Axis(fig, xlabel="Free Variable", ylabel="MSE", title="1D Hidden Interpolation MSE")
stem!(ax, freevars, mse, color=:black)
save(joinpath(@__DIR__, "hinteg1d_error.png"), fig)
display(fig)


text:
function extract_cload!(
    cloads::Dict{TI,Vector{TF}}, file, ::Type{Val{dim}}
) where {TI,TF,dim}
    pattern = r"(\d+)\s*,\s*(\d)\s*,\s*(\-?\d+\.\d*E[\+\-]\d{2})"
    line = readline(file)
    m = match(stopping_pattern, line)
    while m isa Nothing
        m = match(pattern, line)
        if m != nothing
            nodeidx = parse(TI, m[1])
            dof = parse(Int, m[2])
            load = parse(TF, m[3])
            if haskey(cloads, nodeidx)
                cloads[nodeidx][dof] += load
            else
                cloads[nodeidx] = zeros(TF, dim)
                cloads[nodeidx][dof] = load
            end
        end
        line = readline(file)
        m = match(stopping_pattern, line)
    end
    return line
end


text:
module TurbulenceStatistics

using Oceananigans
using Oceananigans.AbstractOperations
using Oceananigans.Architectures
using Oceananigans.Diagnostics
using Oceananigans.Grids
using Oceananigans.Fields
using Oceananigans.Operators
using Oceanostics

ViscousDissipation(model; kw...) = Oceanostics.IsotropicViscousDissipationRate(model)

include("first_through_third_order.jl")
include("subfilter_fluxes.jl")
include("turbulent_kinetic_energy_budget.jl")

end # module

max_stars_repo_name: ali-ramadhan/LESbrary

max_issues_repo_name: ali-ramadhan/LESbrary

max_forks_repo_name: ali-ramadhan/LESbrary