EIS with SPME+SR PyBaMM Model

[Akila1993]

Hello,

I am trying to obtaing the EIS from simscape battery model Dr Widanage based on the SPMe+SR model proposed by @brosaplanella.

The observation is I am getting an EIS quite different from the EIS plot we obtained previously in this discussion. The parameter set is Chen2020 for both cases. I am attaching my Matlab code and the EIS plot obtained below. If anyone can guide or propose some suggestions to this will be highly appreciated.

Thank you.

MATLAB Code

simuLinkMdlStr = 'TSPMeA_EIS_New1'; % Specify the simulink file name (without the file extension)
load_system(simuLinkMdlStr+".slx")
set_param(simuLinkMdlStr,'SimscapeUseOperatingPoints','off');
Vmax = 4.2; % Set upper cut-off voltage [V]
Vmin = 2.5; % Set lower cut-off voltage for terminating discharge [V]

modelProp = Simulink.SimulationInput(simuLinkMdlStr);

frequencies = logspace(-4,2,30);
Z_real_array = [];
Z_imag_array = [];

for FF = frequencies
Freq = 2piFF; % rad/s
Ts = 2*pi/Freq;
Amp = 0.01; % A

%Changing any intial conditions if needed:
modelProp = modelProp.setBlockParameter([simuLinkMdlStr,'/TSPMeA'],'z0',num2str(0.9)); % Initial cell SoC, let's assume 70%
modelProp = modelProp.setBlockParameter([simuLinkMdlStr,'/Temperature Source'],'temperature',num2str(298)); % Ambient temperature [K], let's assume room temperature
modelProp = modelProp.setBlockParameter([simuLinkMdlStr,'/TSPMeA'],'T0',num2str(1500)); % Initial cell temperature [K], let's assume 283K (10degC)

modelProp = modelProp.setModelParameter('StopTime','Ts*20','SaveState','on','SaveOutput','on','SaveFinalState','on','FinalStateName','OperPoint'); % Set simuation options

simOut_P = sim(modelProp); % Simulate the "modelProp" object with the above specifications

%period for the desired frequency 10Hz
period = 1 / FF;

% Calculate the number of samples per period
samples_per_period = round(simOut_P.tout(end) * FF);

Current_new = simOut_P.logsout{1}.Values.Data; %(end - 2 * samples_per_period:end);
Voltage_new = simOut_P.logsout{2}.Values.Data; %(end - 2 * samples_per_period:end);

save_system(simuLinkMdlStr+".slx")

% Extract time vector from simulation results
time_vector = simOut_P.tout;

% Calculate the sampling interval
sampling_interval = mean(diff(time_vector));

% Calculate impedance using Fourier transform
Fs = 1 / sampling_interval; % Use the determined sampling interval

N = length(Current_new); % Number of samples

current_fft = fft(Current_new)/Fs;
voltage_fft = fft(Voltage_new)/Fs;

%voltage_fft = fft(Voltage_new-mean(Voltage_new))/Fs;

% Find the index of the dominant frequency in the current signal
[~, idx] = max(abs(current_fft));

impedance_at_dominant_freq = voltage_fft(idx) / current_fft(idx);

% Store impedance results
Z_real_array = [Z_real_array, real(impedance_at_dominant_freq)];
Z_imag_array = [Z_imag_array, imag(impedance_at_dominant_freq)];

end

[DrSOKane]

That doesn't look physically plausible, but I don't know much about EIS so not sure I can help!