Fixing Problems Integrating Function in MATLAB

[Sarah1287]

Hey I'm having problems trying to integrate my function func1 in MATLAB. Really don't know where I am going wrong so would really appreciate if anyone could help. The code I'm using is below, thanks:

Function Mfile:

function [ f1 ] = func1( E, Delta, kB, Temp )
%func1: calculates integrand for eqn 1
% Detailed explanation goes here
f1=(tanh(sqrt(((E*E)+(Delta*Delta))))/(2*kB*Temp))/(2*(sqrt(((E*E)+(Delta*Delta)))));
end

Command line code:

hbar=1.054E-34; Tdebye=420; h=6.626E-34; kB1= 1.38065E-23;
fdebye=(Tdebye*kB1)/h;
Temp1=0.5;
Delta1=0;
a=hbar*fdebye;
Q=quad(@(E1)func1(E1,0,1.38065E-23,0.5),-a,a)
? Undefined function or variable 'kB'.

? Error using ==> mtimes
Inner matrix dimensions must agree.

Error in ==> func1 at 4
f1=(tanh(sqrt(((E*E)+(Delta*Delta))))/(2*kB*Temp))/(2*(sqrt(((E*E)+(Delta*Delta)))));

Error in ==> @(E1)func1(E1,0,1.38065E-23,0.5)


Error in ==> quad at 77
y = f(x, varargin{:});


Thanks

 

[matonski]

quad needs the function it's integrating to accept a vector argument. In other words, in your anonymous function, E1 will be a vector. E1*E1 will cause an error because you can't multiply to row matrices together. If you want to do it element by element, you need to do E1.*E1 (note the period). Similarly with dividing. If you want to divide element by element, you need to use ./ instead of just /

 

[Sarah1287]

Thankyou. I have actually changed that now but I am getting a new error:

? Index exceeds matrix dimensions.

Error in ==> quad at 85
if ~isfinite(y(7))


new function code:

function [ f1 ] = func1( E, Delta, kB, Temp )
%func1: calculates integrand for eqn 1
% Detailed explanation goes here
f1=(tanh(sqrt(((E.*E)+(Delta.*Delta))))/(2.*kB.*Temp))/(2.*(sqrt(((E.*E)+(Delta.*Delta)))));
end


Sorry but you couldn't tell me where I've gone wrong this time could you?

Thankyou very much

 

[matonski]

You are still missing one period in dividing:

f1 = tanh(sqrt(E.*E+Delta*Delta))/(2.*kB.*Temp)./(2.*sqrt(E.*E + Delta*Delta));

I took out all those extra parenthesis you put in. Matlab follows standard order of operations. However, now you will get NaN's. To be sure, this is the function that Matlab thinks you are integrating:
$$\frac{\tanh \sqrt{E^2 + \delta^2}}{2k_BT} \frac{1}{2\sqrt{E^2+\delta^2}}$$
If that's not right, you need to fix the parenthesis and stuff. The expression you had and the one I put without the extra parenthesis are the same.

If that is right, then perhaps it's because you are dealing with numbers so small (a is of the order 10^-22) and Matlab doesn't really like dealing numbers smaller than 10^-16. It will be better to change units and then convert at the end to avoid round off errors.

 

[Sarah1287]

Thanks yes that does fix it. Also the NAN answers are only returned when delta is set to zero since the computer then thinks that when E is zero too it is dividing by zero. If i put an extra thing into avoid that it should work fine. Thanks a lot for your help.