Conversion to double from function_handle is not possible message

Hello everyone. I am trying to calculate a double integral in matlab. However my program gave me "Conversion to double from function_handle is not possible" message. Please help me how to fix it. I am a beginner and I do not know how to do. Thank you.

clear all;

clc;

format long

syms x

syms y

syms ome1

syms ome2

syms ome12

syms z1

syms z2

syms z12

A = 1.85; %alpha

a = 2.656; % r_0

D = 0.269; %D_0

M_Fe = 56;

M_B = 27;

e = 2.71828;

h = 6.625*10^(-34);

h_bar = h/(2*pi);

k_B = 1.3824*10^(-23);

muy_Fe = M_Fe/(M_Fe+M_B);

muy_B = M_B/(M_Fe+M_B);

C_Fe = 0.8; %Ti le phan tram cua Fe

C_B = 0.2;

k_eff = (1+(5/3)*(muy_Fe^2+muy_B^2))*D*A^2;

k3 = -(1-muy_Fe^3+muy_B^3)*D*A^3;

k4 = (7/12+(203/324)*(muy_Fe^4+muy_B^4))*D*A^4;

M_FeB = (M_Fe+M_B)/2;

M = 2*C_B*M_FeB+(C_Fe-C_B)*M_Fe;

omega_D = 2*sqrt(k_eff/M);

theta_D = h_bar*omega_D/k_B;

sigma_01 = (3*h_bar*a/(2*pi))*(2*C_B*k3+(C_Fe-C_B)*k3)/(2*C_B*k_eff+(C_Fe-C_B)*k_eff);

sigma_02 = (h_bar*a/(2*pi))*(1/(2*C_B*k_eff+(C_Fe-C_B)*k_eff));

sigma_03 = (3*h_bar^2*a^2/(4*pi))*(2*C_B*k3+(C_Fe-C_B)*k3)/(2*C_Fe*k_eff+(C_Fe-C_B)*k_eff);

k0_eff = 0.25*M*omega_D^2;

p = 10;

Tmax = 1000;

dT = Tmax/p;

T = 0.000001;

C3 = zeros(1,p);

Tem = (0:dT:Tmax-dT)

for n = 1:p

beta = 1/(k_B*T);

ome1 = @(x) 2.*sqrt(k_eff/M).*abs(sin(a.*x./2));

ome2 = @(y) 2.*sqrt(k_eff/M).*abs(sin(a.*y./2));

ome12 = @(x,y) 2.*sqrt(k_eff/M).*abs(sin(a.*(x+y)./2));

z1 = exp(beta.*h_bar.*ome1(x));

z2 = exp(beta.*h_bar.*ome2(y));

z12 = exp(beta.*h_bar.*(ome1(x)+ome2(y)));

fun = @(x,y) (ome1(x).*ome2(y).*ome12(x,y)/(ome1(x)+ome2(y)+ome12(x,y)))*(1+(6.*(ome1(x)+ome2(y))./(ome1(x)+ome2(y)-ome12(x,y)))*((z1(ome1).*z2(ome2)-2.*z12(ome1,ome2))./((z1(ome1)-1).*(z2(ome2)-1).*(z12(ome1,ome2)-1))));

ymax = @(x) pi/a - x;

C3 = integral2(fun,0,1,0,ymax)

T = T+dT;

waitbar(n/p)

end

Answer 1

[Stephen23]

You define ome1, ome2, and ome12 as functions of 1 or 2 input arguments, but in some cases you do not call their function handles with any inputs:

fun = @(x,y) (ome1(x).*ome2(y).*ome12(x,y)/(ome1(x)+ome2(y)+ome12(x,y)))*...
    (1+(6.*(ome1(x)+ome2(y))./(ome1(x)+ome2(y)-ome12(x,y)))*...
    ((z1(ome1).*z2(ome2)-2.*z12(ome1,ome2))./((z1(ome1)-1).*... more here
%           ^         ^            ^    ^            ^      ... more here

You need to call those functions with appropriate input arguments.

[Vinh Le Duc]

Thank you so much