Implement a multiplier to rate constants

[jiweiqi]

Pass a vector that contains factors to create the solution, and multiply it in wdot_func

[jiweiqi]

This can be realized by reload the relevant functions and pass the parameters there, for example

@inbounds function dudt!(du, u, p, t)
    Y = @view(u[1:ns])
    T = u[end]
    mean_MW = 1.0 / dot(Y, 1 ./ gas.MW)
    ρ_mass = P / R / T * mean_MW
    X = Y2X(gas, Y, mean_MW)
    C = Y2C(gas, Y, ρ_mass)
    cp_mole, cp_mass = get_cp(gas, T, X, mean_MW)
    h_mole = get_H(gas, T, Y, X)
    S0 = get_S(gas, T, P, X)
    # _p = reshape(p, nr, 3)
    # kp = @. @views(exp(_p[:, 1] + _p[:, 2] * log(T) - _p[:, 3] * 4184.0 / R / T))
    kp = exp.(p)
    qdot = wdot_func(gas.reaction, T, C, S0, h_mole; get_qdot = true) .* kp
    qdot[[2,3,4,6,7]] .*= 0.0
    wdot = gas.reaction.vk * qdot
    Ydot = wdot / ρ_mass .* gas.MW
    Tdot = -dot(h_mole, wdot) / ρ_mass / cp_mass
    du .= vcat(Ydot, Tdot)
    return du
end

and

@inbounds function dudt!(du, u, p, t)
    w_in, w_b, w_out = p2vec(p)
    Y = clamp.(@view(u[1:ns]), -1.e-3, 1.0)
    T = clamp(u[end], 600.0, 3500.0) * ones(eltype(p), 1)[1]
    mean_MW = 1.0 / dot(Y, 1 ./ gas.MW)
    ρ_mass = P / R / T * mean_MW
    X = Y2X(gas, Y, mean_MW)
    C = Y2C(gas, Y, ρ_mass)
    cp_mole, cp_mass = get_cp(gas, T, X, mean_MW)
    h_mole = get_H(gas, T, Y, X)
    S0 = get_S(gas, T, P, X)
    qdot = wdot_func(gas.reaction, T, C, S0, h_mole; get_qdot=true)
    qdot[2:7] .*= 0.0
    wdot = gas.reaction.vk * qdot
    crnn_in = vcat(log.(clamp.(@view(C[ind_crnn]), lb, 1.0)), -1.0 / T, log(T))
    wdot[ind_crnn] .+= w_out * exp.(w_in' * crnn_in + w_b)
    Ydot = wdot / ρ_mass .* gas.MW
    Tdot = -dot(h_mole, wdot) / ρ_mass / cp_mass
    du .= vcat(Ydot, Tdot)
    return du
end

@inbounds function p2vec(p)
    slope = p[end] .* 10.0
    w_b = @view(p[1:crnn_nr]) .* slope
    w_in_b = @view(p[crnn_nr + 1:crnn_nr * 2])
    w_in_Ea = @view(p[crnn_nr * 2 + 1:crnn_nr * 3]) .* slope
    w_out = E_null' *
            reshape(@view(p[crnn_nr * 3 + 1:crnn_nr * (nse + 3)]),
                    nse, crnn_nr)
    w_in = vcat(clamp.(-w_out, 0.0, 2.5), w_in_Ea', w_in_b')
    return w_in, w_b, w_out
end