MATLAB: Problems using diff() to compute function handle derivative

diff()function derivativefunction handleMATLABmatlab function

Hi everybody,
First of all, sorry for the large code-insert. But as explained below, I wasn't able to reproduce the error within a simple example.
I try to compute the derivative of the function handles "B_Teeth" and "B_Yoke" with respect to the only variable t (see last few rows of the code). As I need to transform the output of this operation further, I need to have it as a function of t itself.
As I read here, this is in principle feasable but although I've been trying for a while it does not work. What I'm wondering about are two things:
1) The error code states that the input argument of diff() is sym. But when I doublecheck this via class(), it says that "B_Teeth" is a function handle (as it should be).
2) I tried to simplify the problem and reproduce the error code (see second code) but in this case, everything works well.
I highly appreciate your help!
%% Parameter
p=4;N=32;a_p=0.742;b=0.886;B_r=1.25;mu_r=1.024;Rs=0.091;Rr=0.090;t_b=0.001;w_m=0.020;w_r=0.008;mu_0=4*pi*10^(-7);   
l_m=0.006;l_bar1=0.005;l_bar2=0.0059;w_bar=0.005;w_b=w_bar+2*t_b; Q=48;L=135/149*Rr;Wtau_p=5/6;m=3;q=Q/(2*p*m);J=10;Omega=100;a_PM=140/360*2*pi;
tau_p=pi*Rs/p;tau_s=pi*2*Rs/Q;b_StZa=0.004;a_s=pi/(m*q);mu_b=2.3696e-5;tau_2=pi*(Rs+0.004)/p; B_b=1.6178;
phi=pi; I=100; k_v=[-1;1;-1;1;-1;1]; v=[1;5;7;11;13;17];
%% Reaktances
P1=mu_0./(p.*(Rs-Rr)).*((Rs+Rr)./2).*L; P2=P1.*a_p.*pi; P3=2.*mu_0.*mu_r.*w_m.*L./l_m; P4=2*mu_0*(l_bar1+l_bar2)*L/w_bar; P5=2*mu_b*t_b*L/w_b;
P6=mu_0/(p*(Rs-Rr))*(Rs+Rr)/2*L*pi*(b-a_p)/2;
syms nu theta j t k_vv
k_wv=sin(Wtau_p.*nu.*pi./2).*(sin(nu.*pi./(2.*m))./(q.*sin(nu.*pi./(2.*m.*q))));
F_rdModHiVarUd1=4./pi.*symsum(sin((2*j+1).*a_p.*pi./2)./(2*j+1).*cos((2*j+1).*(theta-Omega.*t)),j,0,J);
F_rdModHiVarUd2=4./pi.*symsum((sin((2*j+1).*b.*pi./2)-sin((2*j+1).*a_p.*pi./2))./(2*j+1).*cos((2.*j+1).*(theta-Omega.*t)),j,0,J); 
Fr_qModHiVar=4./pi.*symsum((cos((2*j+1).*a_p.*pi./2)-cos((2*j+1).*b.*pi./2))./(2*j+1).*sin((2.*j+1).*(theta-Omega.*t)),j,0,J);
B_airgapPM=B_r.*4./pi.*symsum(sin((2*j+1).*a_p.*pi./2)./(2.*j+1).*cos((2.*j+1).*(theta-Omega.*t)),j,0,J);
F=3.*N.*I.*k_wv./(nu.*p.*pi); 
U_d1Mod_HiVar=2.*P1.*F.*sin(nu.*a_p.*pi./2)./(nu.*(P2+P3+P4+P5)).*cos(phi).*cos((k_vv+nu).*Omega.*t);
U_d2Mod_HiVar=2.*P1.*F.*sin(nu.*pi./4.*(b-a_p)).*cos(phi)./(nu.*(0.5.*P5+P6)).*cos((nu+k_vv).*Omega.*t+nu.*(b+a_p).*pi./4);
U_qMod_HiVar=-2.*k_vv.*P1.*F.*sin(phi).*sin(nu.*pi./4.*(b-a_p)).*sin(Omega.*t.*(k_vv+nu).*nu+pi./4.*(a_p+b))./(nu.*(2.*P5+P6));
for i=1:length(v)
    U_d1Mod(i)=subs(U_d1Mod_HiVar,[nu,k_vv],[v(i),k_v(i)]);
    U_d2Mod(i)=subs(U_d2Mod_HiVar,[nu,k_vv],[v(i),k_v(i)]);
    U_qMod(i)=subs(U_qMod_HiVar,[nu,k_vv],[v(i),k_v(i)]);
end
U_d1Mod=sum(U_d1Mod,'all');
U_d2Mod=sum(U_d2Mod,'all');
U_qMod=sum(U_qMod,'all');
Fs_dMod=F.*cos(nu.*theta+k_vv.*Omega.*t).*cos(phi); 
Fs_qMod=F.*k_vv.*sin(-k_vv*Omega.*t-nu.*theta).*sin(phi);
Fr_dMod=U_d1Mod.*F_rdModHiVarUd1+U_d2Mod.*F_rdModHiVarUd2;
Fr_qMod=U_qMod.*Fr_qModHiVar;
for i=1:length(v)
Fs_DMod(i)=subs(Fs_dMod,[nu,k_vv],[v(i),k_v(i)]);
Fs_QMod(i)=subs(Fs_qMod,[nu,k_vv],[v(i),k_v(i)]);
end
Fs_DMod=(sum(Fs_DMod,'all'));
Fs_QMod=(sum(Fs_QMod,'all'));
B_dMod=(Fs_DMod-Fr_dMod).*mu_0./(Rs-Rr);
B_qMod=(Fs_QMod-Fr_qMod).*mu_0./(Rs-Rr);
%% Airgap Flux density in yoke and teeth 
B_airgapMod=matlabFunction(B_dMod+B_qMod+B_airgapPM);
syms t 
B_Teeth=@(t)tau_s./(b_StZa.*a_s).*integral(@(theta)B_airgapMod(t,theta),-a_s/2,a_s/2);
B_Yoke=@(t)tau_p./(3.*tau_2.*pi).*sqrt(2./(3.*q)).*integral(@(theta)B_airgapMod(t,theta),-pi/2,pi/2);
DeltaB_Yoke=matlabFunction(diff(B_Yoke(t)));
DeltaB_Teeth=matlabFunction(diff(B_Teeth(t)));
The good-working simple trial is:
syms x
Fun=@(x)cos(x);
DiffFun=diff(Fun(x));
FunDiffFun=matlabFunction(DiffFun);
FunDiffFun(2)

Best Answer

You do not have any choice. You need a symbolic expression of B_teeth and you cannot get one by doing numeric integration. You will have to do symbolic integration.
vpaintegral() will only evaluate to a placeholder if there are any free variables not being integrated over. However, diff() is able to reason about those placeholders, and so there is a possibility that you might be able to get somewhere. I would not count on it though.
Note: if you do use int() on a summation, in some cases it is faster to take children() of the summation to get the individual terms, int() those, and add the results together. This will keep int() from trying to pattern match each combination of terms to find something in its library that it has a rule for, which can be a waste of time if you have reasons to believe that there is no such useful combination.
Related Question