Today I discovered a nice package to fit equivalent circuits: EquivalentCircuits.jl
Everything works fine, I am able to fit the electrical values of the components.
But I have not yet figured out, how to transfer the result back into impedance values.
Has someone experience how to achieve that?
By investigation of the source code, I figured out how to simulate the impedance of a given equivalent circuit (including the parameters). Here is my sample code, I hope it helps others to get started:
using EquivalentCircuits, PlotlyJS, Printf
# --- sample impedance data:
impedance_values_ = ComplexF64[5919.90084073586 - 15.794826681804063im, 5919.575521325405 - 32.677443741883025im, 5918.183674897797 - 67.57666460870544im,
5912.242152808868 - 139.49441081175667im, 5887.119965779756 - 285.73699600024963im, 5785.038233646888 - 566.878749499386im, 5428.935296370544 - 997.1881947423717im,
4640.2144606930815 - 1257.8277219098052im, 3871.8361085209845 - 978.9656717819273im, 3537.682636142598 - 564.9627167404748im, 3442.9419240480606 - 315.3996363823805im,
3418.140460121871 - 219.68986957046025im, 3405.513645508888 - 242.57272592660013im, 3373.904450003017 - 396.0671717029891im, 3249.673719526995 - 742.0301829777005im,
2808.423185495856 - 1305.924162464249im, 1779.4087896271944 - 1698.9660879948128im, 701.9588433822435 - 1361.4674351816855im, 208.28978681589047 - 777.6453690080142im,
65.93273498232111 - 392.50667235657im]
# --- corresponding frequencies:
frequency_values = [0.1, 0.20691380811147891, 0.42813323987193935, 0.8858667904100828, 1.8329807108324356, 3.792690190732248, 7.847599703514613, 16.237767391887218,
33.59818286283781, 69.51927961775601, 143.84498882876616, 297.63514416313194, 615.8482110660267, 1274.2749857031336, 2636.650898730358, 5455.5947811685255, 11288.378916846883,
23357.214690901215, 48329.30238571752, 100000.0]
# --- generate frequency vector with n_elements with the same range as given in the measurement:
n_elements = 100
frequ_vec = exp10.(LinRange(log10(frequency_values[1]), log10(frequency_values[end]), n_elements))
# --- examples of suitable equivalent circuits
nr_best_circuits = 1
if nr_best_circuits == 1
circuit_model_preset = "[C1,[C2-[R3,C4],R5]]"
circuit_params_preset = (C1 = 2.322248710116646e-9, C2 = 7.146778669252158e-7, R3 = 8015.389370331851, C4 = 1.6325663887245989e-9, R5 = 5918.9481528813185)
elseif nr_best_circuits == 2
circuit_model_preset = "[C1,R2-[R3,C4]]"
circuit_params_preset = (C1 = 0.025036871360843482, R2 = 396.73873944116787, R3 = 2178.061127814435, C4 = 1.1589755057609664e-5)
elseif nr_best_circuits == 3
circuit_model_preset = "R1-[C2,R3-[C4,R5]]"
circuit_params_preset = (R1 = 20.012355936798915, C2 = 4.000335253194046e-9, R3 = 3400.1181419153604, C4 = 4.0010158529883644e-6, R5 = 2499.9496058123705)
# circuit_par_preset = (R1 = 18.133936476549355, C2 = 3.967856543272228e-9, R3 = 3401.2083814207344, C4 = 3.9972681328104986e-6, R5 = 2500.4785300668846)
else
error(string("Choise nr_best_circuits = ", nr_best_circuits,"does not exist!"))
end
# --- simulate impedance based on suitable preset of a circuit model and its corresponding parameter-set:
circfunc_preset = EquivalentCircuits.circuitfunction(circuit_model_preset)
impedance_preset = EquivalentCircuits.simulateimpedance_noiseless(circfunc_preset, circuit_params_preset, frequ_vec)
# **************************************************************************************************************************
# --- function "circuitevolution()" to find a suitable equivalent circuit model and its parameters:
# **************************************************************************************************************************
# # Arguments
# - `generations::Integer=10`: the maximum number of iterations of the evolutionary algorithm.
# - `population_size::Integer=30`: the number of individuals in the population during each iteration.
# - `terminals::String="RCLP"`: the circuit components that are to be included in the circuit identification.
# - `head::Integer=8`: a hyperparameter than controls the maximum considered complexity of the circuits.
# - `cutoff::Float64=0.8`: a hyperparameter that controls the circuit complexity by removing redundant components.
# - `initial_population::Array{Circuit,1}=nothing`:the option to provide an initial population of circuits
# (obtained by using the loadpopulation function) with which the algorithm starts.
# -------------------------------------------------------------------------------------------------------------------------
equiv_circ_evo = circuitevolution(impedance_values_, frequency_values, terminals="RC", generations=200, population_size=50, head=6)
circfunc_evo = EquivalentCircuits.circuitfunction(equiv_circ_evo.Circuit)
impedance_evo = EquivalentCircuits.simulateimpedance_noiseless(circfunc_evo, equiv_circ_evo.Parameters, frequ_vec)
# --- plot measured against simulated impedance:
function plot_nyquist_comp_preset_evo()
s_pts_info = Vector{String}(undef, 0)
for i_ndx in eachindex(frequency_values)
push!(s_pts_info, @sprintf("#:%i, f=%.2fHz", i_ndx, frequency_values[i_ndx]))
end
# ---
trace_impedance = PlotlyJS.scatter(; x= real(impedance_values_), y= imag(impedance_values_), name = "measured impedance", text = s_pts_info, mode = "markers")
trace_simulated_preset = PlotlyJS.scatter(; x= real(impedance_preset), y= imag(impedance_preset), name = string("preset: ", circuit_model_preset))
trace_simulated_auto = PlotlyJS.scatter(; x= real(impedance_evo), y= imag(impedance_evo), name = string("evolution: ", equiv_circ_evo.Circuit))
# ---
plt_layout = PlotlyJS.Layout(;
title_text = "Sample <-> Preset <-> Evolution",
xaxis_title_text = "z<sub>Real</sub>",
xaxis_dtick = 1000,
yaxis_title_text = "z<sub>Imag</sub>",
yaxis_dtick = 1000,
# --- Fixed Ratio Axes Configuration
yaxis_scaleanchor = "x",
yaxis_scaleratio = 1,
)
return PlotlyJS.Plot([trace_impedance, trace_simulated_preset, trace_simulated_auto], plt_layout)
end
plot_nyquist_comp_preset_evo()