Issue in sloving numerical integration

Running this code I get "ImV = function_handle with value: @(r)-ft1-ft2" instead of number. I want to plot ImV Vs rho but I am not getting the numerical values. I don't know where I am making mistake, please help out. Any help shall be highly appreciated. Thanks in advance

clear all; close all; clc
    
    fmGeV_scale = 5.0574;
    sigma = 0.15; 			
    Nc = 3;
    Nf = 3;
    x = 0.15:0.1:0.5;
    s = size(x,2);
    
    rho = 0.3:0.1:2.0;
    
    for i = 1:s
        T = x(i);
        Lambda_MS = 0.176; %GeV
        %alpha = 0.48;
        
        al1 = 33 - 2*Nf;
        al2 = log(2*pi*T/Lambda_MS);
        al3 = al1*al2;
        alpha = 6*pi/al3;
        
        m_pi2 = 0.0196; %in GeV^2 and 19600 in MeV^2
        
        b1 = 5*m_pi2;% When B = 5mpi2
        b2 = 15*m_pi2;% When B = 15mpi2
        b3 = 25*m_pi2;% When B = 25mpi2
        
        n1 = Nc/3;
        n2 = Nf/6;
        n3 = 4*pi*alpha;
        n4 = n3*(n1 + n2);
        
        d1 = 3*b1*alpha/(2*pi*T.^2);% When B = 5mpi2
        d2 = 3*b2*alpha/(2*pi*T.^2);% When B = 15mpi2
        d3 = 3*b3*alpha/(2*pi*T.^2);% When B = 25mpi2
        
        
        md10(:,i) = T.*sqrt(n4); % Our old form of Debye Mass, B = 0
        md11(:,i) = T.*sqrt(n4 + d1);% When B = 5mpi2
        md12(:,i) = T.*sqrt(n4 + d2);% When B = 15mpi2
        md13(:,i) = T.*sqrt(n4 + d3);% When B = 25mpi2
        
        
        n5 = n3.*n1;
        md20(:,i) = T.*sqrt(n5); % New form of Debye mass, B = 0
        md21(:,i) = T.*sqrt(n5 + d1);% When B = 5mpi2
        md22(:,i) = T.*sqrt(n5 + d2);% When B = 15mpi2
        md23(:,i) = T.*sqrt(n5 + d3);% When B = 25mpi2
        
        
       if T == x(2)
           
           TT = T
           
        r = rho*fmGeV_scale;
             rr = rho;
             
             %first term of the potential
             v20 = alpha*exp(-md20(:,i).*r)./r;% When B = 0
             v21 = alpha*exp(-md21(:,i).*r)./r;% When B = 5mpi2
             v22 = alpha*exp(-md22(:,i).*r)./r;% When B = 15mpi2
             v23 = alpha*exp(-md23(:,i).*r)./r;% When B = 25mpi2
             
    
      
             %second term of the potential
             dc20 = alpha*md20(:,i);% When B = 0
             dc21 = alpha*md21(:,i);% When B = 5mpi2
             dc22 = alpha*md22(:,i);% When B = 15mpi2
             dc23 = alpha*md23(:,i);% When B = 25mpi2
             
    
             
             %third term of the potential
             g2 = 4*pi*alpha;
             mdg2 = 0.33*g2*Nc*(T^2); 
             mu0 = mdg2*sigma/alpha;
             mu = mu0^(1/4);
             var = sqrt(2).*mu.*r;
             bes = besselk(-1/2,var); % Bessel Function
             gm1 = gamma(1/4); %Gamma(1\4)
             nom = gm1*sigma*bes;%nominator
             dnom = sqrt(pi)*mu*(2^(3/4));%denominator
             v30 = nom/dnom;
             
    
             
    
             gm2 = gamma(3/4); %Gamma(3\4)
             v40 = gm1*sigma/(2*gm2*mu);
             
     
     
     syms x5 %calculation for g(mDr)
    
      h1 = @(x5) (x5./(x5.*x5 + 1));
      h2 = @(x5,r) (1 - sin(md21.*r.*x5)./(md21.*r.*x5));
      h   = @(x5,r) h1(x5) .* h2(x5,r);
    
      g_mDr = @(r) integral(@(x5) h(x5), 0, inf);
      
      
     
     
     syms x1 r1 %calculation for g(mDx)
    
      h1 = @(x1) (x1./(x1.*x1 + 1));
      h2 = @(x1,r1) (1 - sin(md21.*r1.*x1)./(md21.*r1.*x1));
      h   = @(x1,r1) h1(x1) .* h2(x1,r1);
    
      g_mD = @(r1) integral(@(x1) h(x1,r1), 0, inf);
      
      mr = sqrt(2)*mu.*r;
      Dj1 = besselk(-1/2,mr);
      im = 1i;
      
      Dj2 = real(besselk(-1/2,im*mr));% real part of imaginary bessel function
      
      Dj3 = besselk(-1/2,0);
      
      
     syms x2
     
       mmu = sqrt(2).*mu
       Dx1 =@(x2) besselk(-1/2,mmu.*x2);
       Dx2 =@(x2) real(besselk(-1/2,im*mmu.*x2));
    
    % first integrand and integral
    
    Im1 = @(x2) Dx2(x2).*x2.*x2.*g_mD(x2);
    ph1 =@(r) Dj1.*integral(Im1(x2),0,r);
    
    Im2 = @(x2) Dx1(x2).*x2.*x2.*g_mD(x2);
    ph2 =@(r) Dj2.*integral(Im1(x2),r,inf);
     
    
    Im3 = @(x2) Dx1(x2).*x2.*x2.*g_mD(x2);
    ph3 =@(r) Dj3.*integral(Im1(x2),0,inf);
    
    % r = 0.1:0.1:2.0;
    
    phi =@(r) ph1(r) + ph2(r) - ph3(r);
    
    
    ft1 =@(r) 2*alpha*T*g_mDr;
    ft2 =@(r) mdg2*sigma*T*phi/mu;
    
    ImV =@(r) -ft1-ft2
             
       end
        
    end   

TT = 0.2500

mmu = 0.8265

ImV = function_handle with value:
    @(r)-ft1-ft2

        
    
    
     figure;
    
     imp1 = plot(rho,phi,'Color','blue','LineStyle','-');

Error using plot
Invalid data argument.

[Jan]

When I try:

imp1 = plot(rho,phi(rho),'Color','blue','LineStyle','-');
%                  ^^^^^

I get an error message about not matching array sizes.

[Captain Rituraj Singh]

my main concern was hw to convert "ImV = function_handle with value: @(r)-ft1-ft2" into number. Also I gave full code to see why am I getting any error. Let me know if you have any solution

[Torsten]

ImV =@(r) -ft1(r)-ft2(r)

instead of

ImV =@(r) -ft1-ft2

I wonder why you define ft1 and ft2 as function handles dependent on r if they don't depend on r:

ft1 =@(r) 2*alpha*T*g_mDr;
ft2 =@(r) mdg2*sigma*T*phi/mu;

A, I see: Here you have made the same mistake:

ft1 =@(r) 2*alpha*T*g_mDr(r);
ft2 =@(r) mdg2*sigma*T*phi(r)/mu;

[Captain Rituraj Singh]

Even after making the suggested changes I am still getting same output:

ImV =  function_handle with value: @(r)-ft1(r)-ft2(r)

while I need, proper solution in numbers not in algebric format, please help me in that.

Thanks

[Torsten]

A function handle is a function. It is necessary to give numerial input to a function handle to get numerical output.

f = @(x) x.^2;
x = 0:0.1:1;
f(x)

ans = 1×11

0 0.0100 0.0400 0.0900 0.1600 0.2500 0.3600 0.4900 0.6400 0.8100 1.0000

It's impossible for us to look through your code and decide which function handle you want to evaluate with which numerical input.