text:
# By default, Julia/LLVM does not use fused multiply-add operations (FMAs).
# Since these FMAs can increase the performance of many numerical algorithms,
# we need to opt-in explicitly.
# See https://ranocha.de/blog/Optimizing_EC_Trixi for further details.
@muladd begin

"""
    TimeSeriesCallback(semi, point_coordinates;
                       interval=1, solution_variables=cons2cons,
                       output_directory="out", filename="time_series.h5",
                       RealT=real(solver), uEltype=eltype(cache.elements))

Create a callback that records point-wise data at points given in `point_coordinates` every
`interval` time steps. The point coordinates are to be specified either as a vector of coordinate
tuples or as a two-dimensional array where the first dimension is the point number and the second
dimension is the coordinate dimension. By default, the conservative variables are recorded, but this
can be controlled by passing a different conversion function to `solution_variables`.

After the last time step, the results are stored in an HDF5 file `filename` in directory
`output_directory`.

The real data type `RealT` and data type for solution variables `uEltype` default to the respective
types used in the solver and the cache.

!!! warning "Experimental implementation"
    This is an experimental feature and may change in future releases.
"""
mutable struct TimeSeriesCallback{RealT<:Real, uEltype<:Real, SolutionVariables, VariableNames, Cache}
  interval::Int
  solution_variables::SolutionVariables
  variable_names::VariableNames
  output_directory::String
  filename::String
  point_coordinates::Array{RealT, 2}
  # Point data is stored as a vector of vectors of the solution data type:
  # * the "outer" `Vector` contains one vector for each point at which a time_series is recorded
  # * the "inner" `Vector` contains the actual time series for a single point,
  #   with each record  adding "n_vars" entries
  # The reason for using this data structure is that the length of the inner vectors needs to be
  # increased for each record, which can only be realized in Julia using ordinary `Vector`s.
  point_data::Vector{Vector{uEltype}}
  time::Vector{RealT}
  step::Vector{Int}
  time_series_cache::Cache
end

function Base.show(io::IO, cb::DiscreteCallback{<:Any, <:TimeSeriesCallback})
  @nospecialize cb # reduce precompilation time

  time_series_callback = cb.affect!
  @unpack interval, solution_variables, output_directory, filename = time_series_callback
  print(io, "TimeSeriesCallback(",
            "interval=", interval, ", ",
            "solution_variables=", interval, ", ",
            "output_directory=", "\"output_directory\"", ", ",
            "filename=", "\"filename\"",
            ")")
end

function Base.show(io::IO, ::MIME"text/plain", cb::DiscreteCallback{<:Any, <:TimeSeriesCallback})
  @nospecialize cb # reduce precompilation time

  if get(io, :compact, false)
    show(io, cb)
  else
    time_series_callback = cb.affect!

    setup = [
             "#points" => size(time_series_callback.point_coordinates, 2),
             "interval" => time_series_callback.interval,
             "solution_variables" => time_series_callback.solution_variables,
             "output_directory" => time_series_callback.output_directory,
             "filename" => time_series_callback.filename,
            ]
    summary_box(io, "TimeSeriesCallback", setup)
  end
end

# Main constructor
function TimeSeriesCallback(mesh, equations, solver, cache, point_coordinates;
                            interval::Integer=1,
                            solution_variables=cons2cons,
                            output_directory="out",
                            filename="time_series.h5",
                            RealT=real(solver),
                            uEltype=eltype(cache.elements))
  # check arguments
  if !(interval isa Integer && interval >= 0)
    throw(ArgumentError("`interval` must be a non-negative integer (provided `interval = $interval`)"))
  end

  if ndims(point_coordinates) != 2 || size(point_coordinates, 2) != ndims(mesh)
    throw(ArgumentError("`point_coordinates` must be a matrix of size n_points × ndims"))
  end

  # Transpose point_coordinates to our usual format [ndims, n_points]
  # Note: They are accepted in a different format to allow direct input from `readdlm`
  point_coordinates_ = permutedims(point_coordinates)

  # Invoke callback every `interval` time steps or after final step (for storing the data on disk)
  if interval > 0
    condition = (u, t, integrator) -> (integrator.iter % interval == 0 || isfinished(integrator))
  else # disable the callback for interval == 0
    condition = (u, t, integrator) -> false
  end

  # Create data structures that are to be filled by the callback
  variable_names = varnames(solution_variables, equations)
  n_points = size(point_coordinates_, 2)
  point_data = Vector{uEltype}[Vector{uEltype}() for _ in 1:n_points]
  time = Vector{RealT}()
  step = Vector{Int}()
  time_series_cache = create_cache_time_series(point_coordinates_, mesh, solver, cache)

  time_series_callback = TimeSeriesCallback(interval,
                                           solution_variables,
                                           variable_names,
                                           output_directory,
                                           filename,
                                           point_coordinates_,
                                           point_data,
                                           time,
                                           step,
                                           time_series_cache)

  return DiscreteCallback(condition, time_series_callback, save_positions=(false,false))
end

# Convenience constructor that unpacks the semidiscretization into mesh, equations, solver, cache
function TimeSeriesCallback(semi, point_coordinates; kwargs...)
  mesh, equations, solver, cache = mesh_equations_solver_cache(semi)

  return TimeSeriesCallback(mesh, equations, solver, cache, point_coordinates; kwargs...)
end

# Convenience constructor that converts a vector of points into a Trixi-style coordinate array
function TimeSeriesCallback(mesh, equations, solver, cache, point_coordinates::AbstractVector;
                            kwargs...)
  # Coordinates are usually stored in [ndims, n_points], but here as [n_points, ndims]
  n_points = length(point_coordinates)
  point_coordinates_ = Matrix{eltype(eltype(point_coordinates))}(undef, n_points, ndims(mesh))

  for p in 1:n_points
    for d in 1:ndims(mesh)
      point_coordinates_[p, d] = point_coordinates[p][d]
    end
  end

  return TimeSeriesCallback(mesh, equations, solver, cache, point_coordinates_; kwargs...)
end

# This method is called as callback during the time integration.
function (time_series_callback::TimeSeriesCallback)(integrator)
  # Ensure this is not accidentally used with AMR enabled
  if uses_amr(integrator.opts.callback)
    error("the TimeSeriesCallback does not work with AMR enabled")
  end

  @unpack iter = integrator
  @unpack interval = time_series_callback

  # Create record if in correct interval (needs to be checked since the callback is also called
  # after the final step for storing the data on disk, indepdendent of the current interval)
  if integrator.iter % interval == 0
    @trixi_timeit timer() "time series" begin
      # Store time and step
      push!(time_series_callback.time, integrator.t)
      push!(time_series_callback.step, iter)

      # Unpack data
      u_ode = integrator.u
      semi = integrator.p
      mesh, equations, solver, cache = mesh_equations_solver_cache(semi)
      u = wrap_array(u_ode, mesh, equations, solver, cache)

      @unpack (point_data, solution_variables,
              variable_names, time_series_cache) = time_series_callback

      # Record state at points (solver/mesh-dependent implementation)
      record_state_at_points!(point_data, u, solution_variables, length(variable_names), mesh,
                              equations, solver, time_series_cache)
    end
  end

  # Store time_series if this is the last time step
  if isfinished(integrator)
    semi = integrator.p
    mesh, equations, solver, _ = mesh_equations_solver_cache(semi)
    save_time_series_file(time_series_callback, mesh, equations, solver)
  end

  # avoid re-evaluating possible FSAL stages
  u_modified!(integrator, false)

  return nothing
end

include("time_series_dg.jl")
include("time_series_dg2d.jl")

end # @muladd


text:
function vec2LDTbyRow(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Integer,(sqrt(8n+1)-1)/2)
    s*(s+1)/2 == n || error("vec2LTbyRow: length of vector is not triangular")
    k=0
    [i<=j ? (k+=1; v[k]) : z for i=1:s, j=1:s]'
end

function vec2UDTbyRow(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Integer,(sqrt(8n+1)-1)/2)
    s*(s+1)/2 == n || error("vec2UTbyRow: length of vector is not triangular")
    k=0
    [i>=j ? (k+=1; v[k]) : z for i=1:s, j=1:s]'
end

function vec2LDTbyCol(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Integer,(sqrt(8n+1)-1)/2)
    s*(s+1)/2 == n || error("vec2LTbyCol: length of vector is not triangular")
    k=0
    [i>=j ? (k+=1; v[k]) : z for i=1:s, j=1:s]
end

function vec2UDTbyCol(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Integer,(sqrt(8n+1)-1)/2)
    s*(s+1)/2 == n || error("vec2UTbyCol: length of vector is not triangular")
    k=0
    [i<=j ? (k+=1; v[k]) : z for i=1:s, j=1:s]
end

function vec2UTbyRow(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Integer,((sqrt(8n+1)-1)/2)+1)
    (s*(s+1)/2)-s == n || error("vec2UTbyRow: length of vector is not triangular")
    k=0
    [i>j ? (k+=1; v[k]) : z for i=1:s, j=1:s]'
end

function vec2LTbyRow(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Integer,((sqrt(8n+1)-1)/2)+1)
    (s*(s+1)/2)-s == n || error("vec2LTbyRow: length of vector is not triangular")
    k=0
    [i<j ? (k+=1; v[k]) : z for i=1:s, j=1:s]'
end

function vec2FMbyRow(v::AbstractVector{T}, z::T=zero(T)) where T
    n = length(v)
    s = round(Int,sqrt(n))
    s^2 == n || error("vec2FMbyRow: length of vector is not square")
    k=0
    [(k+=1; v[k]) for i=1:s, j=1:s]
end

function findTSP(path::AbstractString)
  if isdir(path)
    syms = [Symbol(split(file,".")[1]) for file in readdir(path) if (split(file,".")[end] == "tsp")]
  else
    error("Not a valid directory")
  end
  return syms
end

ext: jl


text:
# run with julia --track-allocation=user

using Profile
using GeometricIntegrators

const Δt    = 0.1
const ntime = 10000

function f(t, x, v)
    v[1] = x[2]
    v[2] = sin(x[1])
    nothing
end

function profile_pendulum(tableau, prof=false)
    x0 = [acos(0.4), 0.0]
    ode = ODE(f, x0)
    int = Integrator(ode, tableau, Δt)
    sol = Solution(ode, Δt, ntime)

    integrate!(int, sol)
    set_initial_conditions!(sol, ode)

    @time integrate!(int, sol)

    if prof
        set_initial_conditions!(sol, ode)

        Profile.clear()
        Profile.clear_malloc_data()

        @profile integrate!(int, sol)
    end
end

profile_pendulum(getTableauExplicitEuler())
profile_pendulum(getTableauExplicitMidpoint())
profile_pendulum(getTableauKutta())
profile_pendulum(getTableauERK4())

profile_pendulum(getTableauCrouzeix())

profile_pendulum(getTableauGLRK(1))
profile_pendulum(getTableauGLRK(2))
profile_pendulum(getTableauGLRK(3))
profile_pendulum(getTableauGLRK(4))

profile_pendulum(getTableauGLRK(1), true)


text:
# Top level test script for Stan.jl
using CmdStan, Test, Statistics

println("\nRunning tests for CmdStan-j1-v5:\n")

code_tests = [
  "test_env.jl",
  "test_utilities.jl",
  "test_cmdtype.jl"
]

# Run execution_tests only if cmdstan is installed and CMDSTAN_HOME is set correctly.
execution_tests = [
  "test_bernoulli.jl",
  "test_bernoulli_diagnostics.jl",
  "test_bernoulli_optimize.jl",
  "test_bernoulli_diagnose.jl",
  "test_bernoulli_variational.jl",
  "test_bernoulli_namedarray.jl",
  "test_bernoulliinittheta.jl",
  "test_bernoulliinittheta2.jl",
  "test_bernoulliinittheta3.jl",
  "test_bernoulliscalar.jl",
  "test_binomial.jl",
  "test_binormal.jl",
  "test_schools8.jl",
  "test_dyes.jl",
  "test_kidscore.jl",
  "test_fixed_param.jl",
  "test_zerolengtharray.jl",
  "test_parse_interpolate.jl"
]

if CMDSTAN_HOME != ""
  println("CMDSTAN_HOME set. Try to run tests.")
  @testset "CmdStan.jl" begin
    for my_test in code_tests
        println("\n\n\n  * $(my_test) *")
        include(my_test)
    end

    for my_test in execution_tests
        println("\n\n  * $(my_test) *\n")
        include(my_test)
    end

    println("\n")
  end
else
  println("\n\nCMDSTAN_HOME not set or found.")
  println("Skipping all tests that depend on CmdStan!\n")
end

println("\n")


text:
using DashBootstrapComponents

form = dbc_form(
    dbc_row(
        [
            dbc_label("Email", width = "auto"),
            dbc_col(
                dbc_input(type = "email", placeholder = "Enter email"),
                className = "me-3",
            ),
            dbc_label("Password", width = "auto"),
            dbc_col(
                dbc_input(type = "password", placeholder = "Enter password"),
                className = "me-3",
            ),
            dbc_col(dbc_button("Submit", color = "primary"), width = "auto"),
        ],
        className = "g-2",
    ),
)

max_stars_repo_name: glsdown/dash-bootstrap-components

max_issues_repo_name: glsdown/dash-bootstrap-components

max_forks_repo_name: glsdown/dash-bootstrap-components


text:
quatErr(q_err::UnitQuaternion{TD}) where {TD} = (0 <= q_err.w) ? SVector( q_err.x,  q_err.y,  q_err.z) : SVector(-q_err.x, -q_err.y, -q_err.z)

function quatErr(q1::UnitQuaternion{TD}, qRef::UnitQuaternion{TD}) where {TD}
  q_err = q1 * inv(qRef)
  return quatErr(q_err)
end

cheapRV(q::UnitQuaternion) = 2 * quatErr(q)
cheapRV(spq::MRP) = cheapRV(UnitQuaternion(spq))

components(q::UnitQuaternion{T}) where {T} = SVector{4,T}(q.w, q.x, q.y, q.z)
components(spq::MRP{T}) where {T} = SVector{3,T}(spq.x, spq.y, spq.z)