Info

此问题已关闭。 请重新打开它进行编辑或回答。

Array issue in numerical integration

1 次查看(过去 30 天)
heavygravitation
heavygravitation 2017-11-3
关闭: MATLAB Answer Bot 2021-8-20
Disclaimer: Brand New to Matlab; Apologies in advance for syntax or structure errors in my question.
I'm trying to solve the equation in the image. The integral needs to be solved numerically. Ultimately I'm trying to create a plot of Amin vs T for T between 1 and 1000. I keep getting dimensional errors, among other things, though. I'm not sure how to set it up so I can solve Amin for a range of T values to plot Amin vs T.
[For the image, note that I'm doing a 1D case so you can ignore the summation over i]
k=1.38e-23; %Boltzmann constant%
n=2.11e29; %Number density of atoms%
v=6134.8; %Speed of sound in YTaO4%
T=[1:1000]; %Temperature range in Kelvin%
hbar= 1.055e-34; %h-bar constant%
theta=v*(hbar/k)*(6*pi^2*n)^(1/3); %cutoff frequency for each polarization in Kelvin%
y = theta/T;
fun = @(x) (x^3*exp(x))/(e^x-1)^2;
I = integrate(fun,1,y);
Amin=(pi/6)^(1/3)*k*n^(2/3)*v*(T/theta)^2*I;
plot(T,Amin)
I know I need to get the sizes of T and Amin equal but don't know how to do this in Matlab. Thanks in advance for any help anyone can offer.

回答(1 个)

Star Strider
Star Strider 2017-11-3
This runs:
k=1.38e-23; %Boltzmann constant%
n=2.11e29; %Number density of atoms%
v=6134.8; %Speed of sound in YTaO4%
T=[1:1000]; %Temperature range in Kelvin%
hbar= 1.055e-34; %h-bar constant%
theta=v*(hbar/k)*(6*pi^2*n)^(1/3); %cutoff frequency for each polarization in Kelvin%
y = theta./T;
fun = @(x) (x.^3.*exp(x))./(exp(x)-1).^2;
I = arrayfun(@(z)integral(fun,1E-4,z),y);
Amin=(pi/6)^(1/3)*k*n^(2/3)*v*(T./theta).^2.*I;
plot(T,Amin)
You may need to tweak it to get the result you want.
I substituted ‘1E-4’ for 0, but it still gives:
Warning: Infinite or Not-a-Number value encountered.
I leave that to you to sort.

此问题已关闭。

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by