- the amount of time spent parsing the unrolled code exceeds the amount of overhead of indexing and looping; or
- the Execution Engine is able to find code patterns that it can optimize native code for

# Can it be possible to use a for loop to precise the code

1 view (last 30 days)

Show older comments

Sushanth
on 19 Mar 2023

Commented: Walter Roberson
on 19 Mar 2023

M = 0.2; k = 1; L = 0.1; Pr = 2; G = 2; S1 = 0.1; Ec = 0.1; fw = 1; R = 0.3; a5 = 2; a6 = 3; a7 = 4;

a8 = 10; a9 = 1.72; a10 = 8; a11 = 2; a12 = 3;

xb = 5; %%% xb = Inf;

syms x f0(x) g0(x) h0(x) f1(x) g1(x) h1(x) f2(x) g2(x) h2(x) f3(x) g3(x) h3(x) f4(x) g4(x) h4(x)

eqn0 = [ diff(f0,3) == 0, diff(g0,2) == 0, diff(h0,2) == 0 ];

cond0 = [ f0(0) == fw, subs(diff(f0),0) == 1, g0(0) == 0, h0(0) == 1, subs(diff(f0),xb) == 0, g0(xb) == 0, h0(xb) == 0 ];

F0 = dsolve(eqn0,cond0); f0 = F0.f0; g0 = F0.g0; h0 = F0.h0;

%%% The following lines need to be in squeezed so that it will save time (using for loop)

eqn1 = [ diff(f1,3) + a5*( f0*diff(f0,2) - diff(f0)^2 + k*diff(g0) ) - a6*M^2*diff(f0) + a7*L*h0 == 0, G*diff(g1,2) - 2*g0 - diff(f0,2) == 0, (a9+R)*diff(h1,2) + Pr*( a10*f0*diff(h0) + S1*h0 + Ec*( a11*diff(f0,2)^2 + a12*M^2*diff(f0)^2 ) ) == 0 ];

cond1 = [ f1(0) == 0, subs(diff(f1),0) == 0, g1(0) == 0, h1(0) == 0, subs(diff(f1),xb) == 0, g1(xb) == 0, h1(xb) == 0 ];

F1 = dsolve(eqn1,cond1); f1 = F1.f1; g1 = F1.g1; h1 = F1.h1;

eqn2 = [ diff(f2,3) + a5*( f1*diff(f1,2) - diff(f1)^2 + k*diff(g1) ) - a6*M^2*diff(f1) + a7*L*h1 == 0,G*diff(g2,2) - 2*g1 - diff(f1,2) == 0, (a9+R)*diff(h2,2) + Pr*( a10*f1*diff(h1) + S1*h1 + Ec*( a11*diff(f1,2)^2 + a12*M^2*diff(f1)^2 ) ) == 0 ];

cond2 = [ f2(0) == 0, subs(diff(f2),0) == 0, g2(0) == 0, h2(0) == 0, subs(diff(f2),xb) == 0, g2(xb) == 0, h2(xb) == 0 ];

F2 = dsolve(eqn2,cond2); f2 = F2.f2; g2 = F2.g2; h2 = F2.h2;

eqn3 = [ diff(f3,3) + a5*( f2*diff(f2,2) - diff(f2)^2 + k*diff(g2) ) - a6*M^2*diff(f2) + a7*L*h2 == 0,G*diff(g3,2) - 2*g2 - diff(f2,2) == 0, (a9+R)*diff(h3,2) + Pr*( a10*f2*diff(h2) + S1*h2 + Ec*( a11*diff(f2,2)^2 + a12*M^2*diff(f2)^2 ) ) == 0 ];

cond3 = [ f3(0) == 0, subs(diff(f3),0) == 0, g3(0) == 0, h3(0) == 0, subs(diff(f3),xb) == 0, g3(xb) == 0, h3(xb) == 0 ];

F3 = dsolve(eqn3,cond3); f3 = F3.f3; g3 = F3.g3; h3 = F3.h3;

eqn4 = [ diff(f4,3) + a5*( f3*diff(f3,2) - diff(f3)^2 + k*diff(g3) ) - a6*M^2*diff(f3) + a7*L*h3 == 0, G*diff(g4,2) - 2*g3 - diff(f3,2) == 0, (a9+R)*diff(h4,2) + Pr*( a10*f3*diff(h3) + S1*h3 + Ec*( a11*diff(f3,2)^2 + a12*M^2*diff(f3)^2 ) ) == 0 ];

cond4 = [ f4(0) == 0, subs(diff(f4),0) == 0, g4(0) == 0, h4(0) == 0, subs(diff(f4),xb) == 0, g4(xb) == 0, h4(xb) == 0 ];

F4 = dsolve(eqn4,cond4); f4 = F4.f4; g4 = F4.g4; h4 = F4.h4;

f = f0 + f1 + f2 + f3 + f4; g = g0 + g1 + g2 + g3 + g4; h = h0 + h1 + h2 + h3 +h4;

vpa(f)

##### 0 Comments

### Accepted Answer

Walter Roberson
on 19 Mar 2023

Edited: Walter Roberson
on 19 Mar 2023

%%% The following lines need to be in squeezed so that it will save time (using for loop)

Using a for loop is only faster than "unrolling" to individual statements if

With the small number of lines, and the small number of iterations implied (5 at most), and taking into account that you cannot put symbolic functions into normal vectors (you need to put them into cell arrays instead), I suspect that you cannot save enough execution time by looping to make a notable difference. And the Execution Engine cannot do much to optimize symbolic operations beyond what you could already do by using, for example,

S0 = sym(0);

%...

cond2 = [ f2(S0) == S0, subs(diff(f2),S0) == S0, g2(S0) == S0, h2(S0) == S0, subs(diff(f2),xb) == S0, g2(xb) == S0, h2(xb) == S0 ];

to reduce the time needed to detect the numeric 0 and convert it to symbolic 0 repeatedly.

Symbolic calculations are done inside a separate process, inside the MuPAD-based Symbolic Engine, that does not use the Execution Engine (and cannot take advantage of multiple cores or of GPUs). It is those operations that are slow; looping or not in this code is timing noise compared to the time spent in the Symbolic Engine.

##### 5 Comments

Walter Roberson
on 19 Mar 2023

Edited: Walter Roberson
on 19 Mar 2023

F{k} = sol.("f" + k);

ummm, I will have to come up with a correction for the sum(). I will post that later.

Walter Roberson
on 19 Mar 2023

There might still need to be some tweaks to how the solutions are extracted and to exactly where in F is written to, and to forming the sum.

But in the meantime, the reason that dsolve() are failing complaining about extra indeterminates, is because you are trying to give two boundary conditions for the same function, g1(0) == 0 and g1(xb) == 0, and h1(0) == 0 and h1(xb) == 0. dsolve() can never handle more than one boundary condition for the same derivative of the same function.

M = 0.2; k = 1; L = 0.1; Pr = 2; G = 2; S1 = 0.1; Ec = 0.1; fw = 1; R = 0.3; a5 = 2; a6 = 3; a7 = 4;

a8 = 10; a9 = 1.72; a10 = 8; a11 = 2; a12 = 3;

xb = 5; %%% xb = Inf;

syms x f0(x) g0(x) h0(x) f1(x) g1(x) h1(x) f2(x) g2(x) h2(x) f3(x) g3(x) h3(x) f4(x) g4(x) h4(x)

eqn0 = [ diff(f0,3) == 0, diff(g0,2) == 0, diff(h0,2) == 0 ];

cond0 = [ f0(0) == fw, subs(diff(f0),0) == 1, g0(0) == 0, h0(0) == 1, subs(diff(f0),xb) == 0, g0(xb) == 0, h0(xb) == 0 ];

F0 = dsolve(eqn0,cond0); f0 = F0.f0; g0 = F0.g0; h0 = F0.h0;

%this might need some tweaks

S0 = sym(0); xb = sym(5);

f = {f0, f1, f2, f3, f4};

g = {g0, g1, g2, g3, g4};

h = {h0, h1, h2, h3, h4};

N = length(f);

F = cell(1,N);

for k = 1:N-1

df0 = diff(f{k}); d2f0 = diff(df0);

df1 = diff(f{k+1}); d3f1 = diff(df1,2);

dg0 = diff(g{k});

dh0 = diff(h{k});

d2g1 = diff(g{k+1},2);

d2h1 = diff(h{k+1},2);

eqn = [ d3f1 + a5*( f{k}*d2f0 - df0^2 + k*dg0 ) - a6*M^2*df0 + a7*L*h0 == S0,

G*d2g1 - 2*g{k} - d2f0 == S0,

(a9+R)*d2h1 + Pr*( a10*f{k}*dh0 + S1*h{k} + Ec*( a11*d2f0^2 + a12*M^2*df0^2 ) ) == 0 ]

cond = [ f{k+1}(S0) == S0,

df1(S0) == S0,

g{k+1}(S0) == S0,

h{k+1}(S0) == S0,

df1(xb) == S0,

g{k+1}(xb) == S0,

h{k+1}(xb) == S0 ]

sol = dsolve(eqn,cond);

if isempty(sol)

F{k} = sym(nan);

else

F{k} = sol.("f" + k);

end

end

celldisp(F)

%f = sum(F)

### More Answers (0)

### See Also

### Categories

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!