How to solve: Error using symengine 'Unable to convert expression containing symbolic variables into double array.'

Hello together,

I'm still quite new to Matlab, but as i try to solve a thermal model for my studies I tried quite a bit to solve on my own without success.

For an explanation:

First I'm having two matrices for local distance to a given radius and the angle respective of a vertical through the center of the circle.

r=2.96; %example radius
lc=14; %example arc of contact
ap=0.12;%example depth of cut
a=3.29e-6;%m^2/s
lambda=12.47;%W/mK

K=512;
L=512;
T_Mi=ones(K,L);
distance=ones(K,L);
rotangle=ones(K,L);

dpi=0.1;
ri=r/dpi;
lci=lc/dpi;
mp_x=-ri+fz/dpi;
mp_y=K/2;
M=[mp_y mp_x];
Q=[1 mp_x];
for i=1:size(distance,1)
    for j = 1:size(distance,2)
        distance(i,j) = norm([i j]-[mp_y,mp_x],2);
        if distance(i,j)<=ri
            distance(i,j)=0;
            rotangle(i,j)=0;
        elseif distance(i,j)>ri
            dR=sqrt((mp_x-j)^2+(mp_y-i)^2);
            u=M-Q;
            v=M-[i j];
            rotangle(i,j)=acos(dot(v,u)/(norm(u)*norm(v)));
            distance(i,j)=dR-ri;
        else
            distance(i,j)=0;
        end
       
    end
end

And in the following I give symbolic equations, which should by my understanding fill the matrix T_Mi with values. But the error:

'Unable to perform assignment because value of type 'sym' is not convertible to 'double'. Caused by: Error using symengine

Unable to convert expression containing symbolic variables into double array. Apply 'subs' function first to substitute values for variables.'

keeps coming, though I tried any variation of assumptions and substitutions.

syms phi R;
vf=fz/dpi*Z_total*N;
vc=2*ri*pi*N/1000;

p00=0.01949;
p10=0.5988;
p01=-0.0001006;
p20=-4.164;
p11=4.871e-05;
p02=-2.923e-08;

h=fz/dpi*sin(phi);
qwp=p00+p10*h+p01*vc+p20*h^2+p11*h*vc+p02*vc^2;
qavg=ap*ri*vpaintegral(qwp,phi,phi_st,phi_ex)*n/lci/ap/dpi*2;
qmax=qavg*((phi_c*sin(phi_c))/(1-cos(phi_c)));
F=(vf/(2*a)*sqrt(((R*cos(phi))+ri*(1-cos(phi))).^2+((R*sin(phi))-ri*sin(phi)).^2));
K_0=bessely(0,F);
G=(qmax*sin(phi))/sin(phi_c)*exp(-(vf*((R*sin(phi))-ri*sin(phi)))/(2*a)).*K_0*ri*cos(phi);
H=vpaintegral(G,phi,phi_st,phi_ex);
T_M=1/pi/lambda*H;

for i=1:size(T_Mi,1)
    for j = 1:size(T_Mi,2)
        subs(T_M,R,distance(i,j));
        T_Mi(i,j)=(T_M); %%Error created in this line
    end
end

I hope it's not too much code, and more so that someone can help me solve this.

[Torsten]

for i=1:size(T_Mi,1)
    for j = 1:size(T_Mi,2)
        T_Mi(i,j) = double(subs(T_M,R,distance(i,j)));     
    end
end

might work if T_M does not contain the symbolic variable "phi" or other symbolic variables apart from R.

Note that several variables are undefined in the code you supplied

fz=1.0;   %added
Z_total = 1.0;   %added
N=1.0;  %added
phi_st = 1.0;   % added
phi_ex = 2.0;   % added
n=N;    % added
phi_c = 1.0;  % added

[Niklas M.]

Unfortunately it does not work, looks like the problem is created by:

H=vpaintegral(G,phi,phi_st,phi_ex);

which returns the foregoing equations with 'phi' still in it, instead of solving the equation numerically.

phi_st and phi_ex are defined as 0 and pi, respectively, so I was hoping it would calculate a real value instead of staying symbolic.

[Torsten]

phi_st and phi_ex are defined as 0 and pi, respectively, so I was hoping it would calculate a real value instead of staying symbolic.

G does not depend on R for that

H=vpaintegral(G,phi,phi_st,phi_ex); 

can make sense ?

If G still depends on other symbolic variables apart from phi, you cannot expect a numerical answer. Then you would have to use "int" instead of "vpaintegral". But "int" won't most probably succeed because your integrand is too complex for an analytical antiderivative.

I suggest you do all your computations numerically using "integral" and without symbolic variables.

[Niklas M.]

G depends on R and phi, as shown above the function T_M is expected to be symbolic until it is substituted in the for loop at the end, where R is replaced by the numerical distances from the 'distance'-matrix.

But with that I expect it to be converted to a numerical value.

I also tried the whole thing with "integral" and function handles, but it creates other errors I could not resolve.

[Torsten]

Substitute distance(i,j) for R in G, not in T_M. Then G will only depend on phi, and you'll get numerical value for H(i,j) (hopefully).

[Niklas M.]

I will definitely try that in an instance.

But here now my numerical proceeding

h=@(phi) fz/dpi*sin(phi);
qwp=@(h) p00+p10*h+p01*vc+p20*h^2+p11*h*vc+p02*vc^2;
qavg=@(phi) ap*ri*vpaintegral(qwp,phi,phi_st,phi_ex)*n/lci/ap/dpi*2;
qmax=@(qavg) qavg*((phi_c*sin(phi_c))/(1-cos(phi_c)));
F=@(R,phi) (vf/(2*a)*sqrt(((R*cos(phi))+ri*(1-cos(phi))).^2+((R*sin(phi))-ri*sin(phi)).^2));
K_0=@(F) bessely(0,F);
G=@(R,phi) (qmax*sin(phi))/sin(phi_c)*exp(-(vf*((R*sin(phi))-ri*sin(phi)))/(2*a)).*K_0*ri*cos(phi);
H=@(R) integral(G,phi,phi_st,phi_ex);
T_M=@(R) 1/pi/lambda*H;
%functionhandle
for i=1:size(T_Mi,1)
    for j = 1:size(T_Mi,2)
        T_Mi(i,j) = T_M (distance(i,j));     %% error output created here
    end
end

which also leads to an error : 'Operator '*' is not supported for operands of type 'function_handle'.

Error in live_tmod_handle (line 116)

T_M=@(R) 1/pi/lambda*H;'

Does the syntax I used make sense?

Answer 1

[Torsten]

vf=fz/dpi*Z_total*N;
vc=2*ri*pi*N/1000;
p00=0.01949;
p10=0.5988;
p01=-0.0001006;
p20=-4.164;
p11=4.871e-05;
p02=-2.923e-08;
h=@(phi)fz/dpi*sin(phi);
qwp=@(phi)p00+p10*h(phi)+p01*vc+p20*h(phi).^2+p11*h(phi)*vc+p02*vc^2;
qavg=ap*ri*integral(qwp,phi_st,phi_ex)*n/lci/ap/dpi*2
qmax=qavg*((phi_c*sin(phi_c))/(1-cos(phi_c)));
F=@(R,phi)(vf/(2*a)*sqrt(((R.*cos(phi))+ri*(1-cos(phi))).^2+((R.*sin(phi))-ri*sin(phi)).^2));
K_0=@(R,phi)bessely(0,F(R,phi));
G=@(R,phi)(qmax*sin(phi))/sin(phi_c).*exp(-(vf*(R*sin(phi)-ri*sin(phi)))/(2*a)).*K_0(R,phi)*ri.*cos(phi);
for i=1:size(T_Mi,1)
    for j=1:size(T_Mi,2)
        R = distance(i,j);
        fun = @(phi)G(R,phi);
        H = integral(fun,phi_st,phi_ex);
        T_Mi(i,j) = 1/pi/lambda*H;
    end
end

And have in mind that MATLAB is case-sensitive: "n" is not the same as "N".

[Niklas M.]

Thank you so much, really appreciate your help !!