How do I make my code faster (currently sym.vpasolve takes 43,434 seconds for 30,000 calls)
13 views (last 30 days)
Show older comments
% Universal Constants
F = 96485; % farads
R = 8.31; % gas
NA = 6.02e23; % avagadros
T = 298; % in kelvin
% Cell constants
Length = 2.95 * 1e-6;
Width = 1.07 * 1e-6;
Radius = Width / 2;
% Definition of surface and volume (for a spherocylinder) = π•r²•(l-2r) [volume of cylindrical section whose length is l-2r] + 4/3•π•r³ [volume of speherical section]
S = pi * Width * Length;
Vol = pi * (Width / 2)^2 * (Length - Width / 3);
% Membrane specific capacitance
Cm = 6.5e-3;
constant = (F * Vol) / (S * Cm); % FV/SCm constant in memb potential eqn
eta = F / (R * T); % F/RT
pKa_w = 14; % pKa of water
% External environments
pHe = 7; % Neutral environment
% Concentration of protons [mM]
H_c_e = 10^(pHe) * 1e3;
% Concentration of hydroxides [mM]
OH_c_e = 10^(pHe - pKa_w) * 1e3;
% intracellular trapped ions conc in mM
z_y = [-10, 0, 10];
% In the case of pHe= 7= pHi, PMF= psi (membrane voltage)
pmf = linspace(-0.24, 0, 100); % setting 101 regularly spaced values for PMF
CA_e = linspace(0.1, 100, 100); %0.1 mM to 100mM in total contribution from [Na] and [Cl] ions
total_iterations = length(CA_e) * length(z_y) * length(pmf);
count = 0;
% Preallocate cell arrays to store data
%variables_values = cell(total_iterations, 1);
%solutions = cell(total_iterations, 1);
syms x;
for CA_e_i = CA_e
%
C_c_e = 0.5*CA_e_i; % in mM
% fractional ionic composition
if CA_e_i == 0
continue;
else
alpha_c = C_c_e / CA_e_i;
end
% Concentration of anions
A_c_e = 0.5*CA_e_i; % in mM
if CA_e_i == 0
continue;
else
alpha_a = A_c_e / CA_e_i;
end
for z_y_i = z_y
if CA_e_i == 0
continue;
else
alpha_y = z_y_i / CA_e_i;
end
for pmf_i = pmf
psi = pmf_i;
count = count + 1;
% Define the equation, Eq is an equation object, first term is the rhs and the second term is the LHS (normally what it is equal to)
if CA_e_i == 0
continue;
else
lhs = (psi) / constant * CA_e_i; % from eqn S42 in SI
end
% assuming only cation exporting antiporters work
equation = alpha_y + alpha_c * exp(eta * (x - psi) - alpha_a * exp(eta * psi)) == lhs;
solution = vpasolve(equation);
solution_double = double(solution);
% Store variable values and solution
variables_values{count} = struct('CA_e_i', CA_e_i, 'z_y_i', z_y_i, 'pmf_i', pmf_i);
%this has empty values for places where the solution is not a
%real number (since the equation represents a physical system)
solutions{count} = solution_double;
end
end
end
0 Comments
Answers (1)
Ronit
on 25 Apr 2024
Hi Keshav,
I see that you are looking for ways to optimize MATLAB code that currently takes an excessively long time to execute, specifically due to repeated calls to ‘vpasolve’ within nested loops.
I can suggest some strategies to make your code run faster:
1. Reduce Precision When Possible: ‘vpasolve’ allows you to specify the precision of the solution. If you do not need very high precision, reducing it can significantly speed up the computation. Use the 'Digits' option to lower the precision. You can refer to the documentation: https://www.mathworks.com/help/symbolic/sym.vpasolve.html
solution = vpasolve(equation, x, 'Random', true, 'Digits', 4); % Example with reduced precision
2. Preallocate memory for your variables outside the loop. It seems you have commented out the preallocation lines; ensure they are active and correctly sized to avoid dynamic resizing during the loop.
3. Profiling: Use MATLAB’s profiler (profile on; profile viewer) to identify the exact bottlenecks in your code.
4. Parallel Computing: You can use the Parallel Computing Toolbox to distribute the iterations across multiple cores. This is particularly effective for loop structures that are independent of each iteration. More about the Parallel Computing Toolbox can be gathered here: https://www.mathworks.com/products/parallel-computing.html
I hope these strategies will aid in making your code run faster.
0 Comments
See Also
Categories
Find more on Particle & Nuclear Physics in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!