Speed up MatlabFunction or use alternatives

Hello,

I am currently processing a big equation (in fact a 36 - Vector with ~15.000 characters each) and am trying to get a function using MatlabFunction. The goal is to receive a function which can be evaluated as quickly as possible, as it's part of a system solved by an ODE solver. The equation itself only contains polynomials to the degree of 3 and square-roots (^3/2 e.g.) but no sin, cos, etc. . Calling MatlabFunction took longer than 18 Hours, and now I had to change the equation and I really don't want to wait for the full 18 Hours. Is setting optimize to false advisable, or will this make my evaluation time worse? Are there other options beside MatlabFunction? I am not very familiar with this kind of stuff as I never had to deal with calculations taking so long.

Any help is greatly appreciated!

[darova]

Can you replace your equations? Did you try fitting (fit function)?

[Finn Busch]

Unfortunately not :/

[Stephen23]

Is it required to use symbolic expressions?

If you could stick to numeric operations then you would likely get a significant speed up.

[Finn Busch]

It is needed to use symbolic expressions as it is a calculation which is only to be done once. Precalculating it and using MatlabFunction will boost the speed of the actual solving process (it's part of a flexible multibody system simulation). I can attach the script, if that helps. The goal is to get a function that gives me the solution for actual values as fast as possible, and I found that precalculating it symbolically helps a lot..

[Finn Busch]

and the calculation contains derivatives, which would not work numerically.

[darova]

Can show your formulas or code?

[John D'Errico]

Far too often, these huge things with endless equations result in something mathematically and numerically useless - not to be trusted. It looks impressive, but so what? There may be numerical problems just evaluating that mess. Massive internal subtractive cancellation is not uncommon, so that even while a given computation is mathematically correct, it is numerical garbage.

Can you usea GPU or parallel computations to do something? Sadly not. MatlabFunction is a symbolic tool, it does not work well in parallel. GPUs & parallel tools are GREAT when you need to do the same thing over and over again. Not here.

As I said, it looks impressive. It sounds impressive - taking 18 hours just to process. Computers are great, right?

Good luck. Get some coffee while you wait, and read a good, LONG book.

[Finn Busch]

John I think you are certainly misunderstanding something here. I'm not bragging about the process time it takes, in fact I'd love to reduce it. Other than that it's kind of bold of you to just assume that what I am doing is useless or garbage based on the fact that it takes long to solve.

Thanks for pointing out what's possible and what not though, even though I'd highly appreciate if you could keep your opinion based on random guesses by yourself.

[darova]

I don't want it if it's too big (more 1mb)

[Finn Busch]

dont worry. 2 files, 2kb each

[darova]

Please attach them here!

[Finn Busch]

Run the elastis forces gen.m

[darova]

How long does it take to evaluate? 18h?

[darova]

I could'n handle that. I'm sorry

And what is wrong with numerical calculations? Can i help you with them?

[Finn Busch]

Nono, I don't ask you to solve that for me, sorry! But maybe you can spot an error or something that slows the process down significantly without changing the outcome?

[darova]

I think there are no significant errors or mistakes. I think those calculations should be numerical. It's madness for symbolic evaluations. I can't even open dude variable!

[Finn Busch]

I agree that usually those calculations should be numerical. However these calculations need to be done each time step while solving an ODE. I want to minimize the time it takes to get these calculations done, so I try to precalculate them symbolically, get a MatlabFunction out of it and just call that on each time step (which is quite fast to my experience). Do you think doing these calculations numerically for each time step will be faster after all? Just to be clear, the PRECALCULATION (the script I sent you) only has to be done once.

EDIT: In addition: it's not the calculations that take so long symbolically, it's the MatlabFunction call..

[John D'Errico]

Finn - you misunderstand. Nobody ever said you were bragging. My point is it is too easy to get hung up on the idea that a computer produced these equations, therefore they are correct, and anything that comes from them is meaningful. Remember that any model is merely an approximation of some process. A simpler model may be far more useful, because it is possible to make vaild predictions using it. The art then comes down to which terms are safely excluded.

Anyway, when your code eventually runs, it will produce a result, SOME prediction. You need to heavily validate what you get. Its that result meaningful in any way? Is it consistent with reality? All code like this is a model of some process. But did your model predict something that seems reasonable? The first filter is your eye, as you should know what is realistic. That is not a proof of validity, but just an aid to intuition that you may have done something reasonable.

If you just accept the result as "truth" becuase a computer produced it, then you are asking for trouble. Computers are not infallible. They just do what you ask them to do. As far as you can, then make more conclusive tests to verify the result. It is often true that you can find ways to drive even a complex system to produce an expected simple result. Does your model verify that result?

As I said, don't trust the results from such a complex computation, not until you can verify the computations.

[darova]

What about subs? How much does matlabFunction take for calculations?

vars1 = [e(:)', [a, b, h], [x, y, z],lam,mu];
tic
A = subs(dude, vars1, 2*rand(size(vars1)));
double(A)
toc
% Elapsed time is 30.214898 seconds.

[Finn Busch]

MatlabFunction takes <<1ms for one calculation.

The model I am using is based on continuum mechanics and mostly validated. I cannot change the model itself as it's part of a thesis in which I compare different models. It is expected to be consistent with reality and I have validated it for other models.

[darova]

Im out of ideas. I give up

[Finn Busch]

thanks for trying :)

[Franz Enkelmann]

Have you just found a solution?

Answer 1

[Steven Lord]

Have you tried telling matlabFunction to skip trying to Optimize the code?

cd(tempdir)
syms a b c d positive
syms x
s = solve(a*x^3+b*x^2+c*x+d == 0, x, MaxDegree=3);

s is a pretty complicated expression, but the final answer can be simplified through the use of temporary sub-expressions. But that simplification can result in longer files (albeit with simpler expressions and shorter lines) and longer time creating the files.

tic
matlabFunction(s, File='cubicSolver1.m', Optimize=true);
toc

Elapsed time is 0.477392 seconds.

dbtype cubicSolver1.m

1     function s = cubicSolver1(a,b,c,d)
2     %cubicSolver1
3     %    S = cubicSolver1(A,B,C,D)
4     
5     %    This function was generated by the Symbolic Math Toolbox version 9.2.
6     %    28-Nov-2022 14:57:43
7     
8     t2 = b.^2;
9     t3 = b.^3;
10    t4 = 1.0./a;
11    t7 = sqrt(3.0);
12    t5 = t4.^2;
13    t6 = t4.^3;
14    t8 = (b.*t4)./3.0;
15    t9 = (c.*t4)./3.0;
16    t10 = (d.*t4)./2.0;
17    t11 = -t8;
18    t12 = -t10;
19    t13 = (b.*c.*t5)./6.0;
20    t15 = (t2.*t5)./9.0;
21    t16 = (t3.*t6)./2.7e+1;
22    t14 = -t13;
23    t17 = -t15;
24    t18 = -t16;
25    t19 = t9+t17;
26    t21 = t10+t14+t16;
27    t20 = t19.^3;
28    t22 = t21.^2;
29    t23 = t20+t22;
30    t24 = sqrt(t23);
31    t25 = t12+t13+t18+t24;
32    t26 = t25.^(1.0./3.0);
33    t27 = 1.0./t26;
34    t28 = t26./2.0;
35    t29 = -t28;
36    t30 = t19.*t27;
37    t31 = t30./2.0;
38    t32 = t26+t30;
39    t33 = t7.*t32.*5.0e-1i;
40    s = [t11+t26-t30;t11+t29+t31-t33;t11+t29+t31+t33];

tic
matlabFunction(s, File='cubicSolver2.m', Optimize=false);
toc

Elapsed time is 0.121004 seconds.

dbtype cubicSolver2.m

1     function s = cubicSolver2(a,b,c,d)
2     %cubicSolver2
3     %    S = cubicSolver2(A,B,C,D)
4     
5     %    This function was generated by the Symbolic Math Toolbox version 9.2.
6     %    28-Nov-2022 14:57:43
7     
8     et1 = (c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).*1.0./((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0);
9     et2 = ((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0);
10    et3 = ((c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).*1.0./((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0))./2.0+sqrt(3.0).*(et1+et2).*5.0e-1i-b./(a.*3.0);
11    et4 = ((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0).*(-1.0./2.0);
12    et5 = (c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).*1.0./((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0);
13    et6 = ((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0);
14    et7 = ((c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).*1.0./((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0))./2.0-sqrt(3.0).*(et5+et6).*5.0e-1i-b./(a.*3.0);
15    et8 = ((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0).*(-1.0./2.0);
16    et9 = -(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).*1.0./((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0)-b./(a.*3.0);
17    et10 = ((d.*(-1.0./2.0))./a-(1.0./a.^3.*b.^3)./2.7e+1+sqrt((d./(a.*2.0)+(1.0./a.^3.*b.^3)./2.7e+1-(1.0./a.^2.*b.*c)./6.0).^2+(c./(a.*3.0)-(1.0./a.^2.*b.^2)./9.0).^3)+(1.0./a.^2.*b.*c)./6.0).^(1.0./3.0);
18    s = [et9+et10;et7+et8;et3+et4];