How to integrate an array properly in matlab

2013 6 6
1 个回答
回答已采纳
16 次查看（30 天）

0 个投票

Hello, I am trying to finish a m-file to find the inductance in 2 coils. I have finished the program to the point of integration. I believe the problem is the integration of an array. I have tried various methods such as int, trapz, and quad but all of these seem to be returning an error. I am not sure as to whether I am implementing the commands wrong or rather I have a bad equation. Here is my code
% This program finds mutual inductance for TET coil system %
% Using Neumann's definition %
% M = sqrt(ap*as)*(1/(2*pi))*(int(((cos(theta)-(d/as)*((cos(psi)*cos(phi))-(sin(psi)*sin(phi)*cos(theta))))/(R^(3/2)))*f,phi,0,2*pi))
%
% Increment angles
%psi = [0:1:90];
%theta = 0:1:90;
%phi = 0:1:90;
psi = input('Enter psi value in degrees \n')
theta = input('Enter theta value in degrees \n')
phi = input('Enter phi value in degrees \n')
% Arguments
%dim = 0:2*pi
h = input('Enter h value in mm \n')
ap = input('Enter ap value in mm \n')
as = sqrt((ap^2)-(h^2))
delta = h/ap
alpha = as/ap
d = h.*tan(phi)
Ra = (1-(cos(phi).*cos(phi).*sin(theta).*sin(theta)))
Rb = ((2)*(d/as)).*((sin(psi).*sin(phi))-(cos(psi).*cos(phi).*cos(theta)))
Rc = ((d.^2)/(as.^2))
R = sqrt(Ra+Rb+Rc)
z = delta-(alpha*sin(theta)*cos(phi))
%kprime_2 = (((1-(alpha.*R)).^2)+z.^2)/(((1+(alpha.*R)).^2+z.^2))
kprime_2a = ((1-(alpha.*R)).^2)+z.^2
kprime_2b = ((1+(alpha.*R)).^2)+z.^2
kprime_2 = kprime_2a./kprime_2b
f = -0.011*(log(kprime_2))-0.0021
integrand = ((cos(theta)-(d/as).*((cos(psi).*cos(phi))-
(sin(psi).*sin(phi).*cos(theta))))/(R.^(3/2))).*f
%integrand1 = double(int(((cos(theta)-(d/as).*((cos(psi).*cos(phi))-(sin(psi).*sin(phi).*cos(theta))))./(R.^(3/2))).*f,phi,0,2.*pi))
%integrand = double(integrand);
stuff = trapz(phi,integrand)
M = sqrt(ap.*as).*(1/(2*pi)).*stuff
I am setting psi and phi to 0 and setting theta to 0:10:90. H is usually 3 and ap is usuall 6. This gives me different error messages for each method of integration I use.
Any help would be appreciated. Thanks.

5 个评论
显示 3更早的评论 隐藏 3更早的评论

Andrew Newell
Andrew Newell 2013-6-6
Are you trying to integrate over theta or phi?
Alexander Lampe
Alexander Lampe 2013-6-6
Hey, I tested your script with your values psi, phi... But then "phi" and "integrand" are both scalar. Then MATLAB expects, that the second Input off trapz() is the dimension, so it needs to be an integer. Write trapz(phi,integrand,1) or rather trapz(phi,integrand,dim) and it will work. But why you use trapz() for scalar values?
Alex
Andrew Matthews
Andrew Matthews 2013-6-7
Hello, thanks for the help. When I tried that earlier it did not work. Thank you for your time.
Andrew Matthews
Andrew Matthews 2013-6-7
编辑：Andrew Matthews 2013-6-7
Whenever I use 1 my integral becomes 0. It was not supposed to be 0. Would you know the reasoning? Thank you.
Andrew Matthews
Andrew Matthews 2013-6-7
I just realized, it is because my program is filling in the angles before integration. Because of this I am getting 0.

请先登录，再进行评论。

请先登录，再回答此问题。

 采纳的回答

Andrew Newell
Andrew Newell 2013-6-6
编辑：Andrew Newell 2013-6-6

0 个投票

Assuming you are trying to integrate over theta, a good approach is to create a function for your integrand like this:
function y = inductanceIntegrand(theta,psi,phi,h,ap)
as = sqrt((ap^2)-(h^2));
delta = h/ap;
alpha = as/ap;
d = h.*tan(phi);
Ra = (1-(cos(phi).*cos(phi).*sin(theta).*sin(theta)));
Rb = ((2)*(d/as)).*((sin(psi).*sin(phi))-(cos(psi).*cos(phi).*cos(theta)));
Rc = ((d.^2)/(as.^2));
R = sqrt(Ra+Rb+Rc);
z = delta-(alpha*sin(theta).*cos(phi));
kprime_2a = ((1-(alpha.*R)).^2)+z.^2;
kprime_2b = ((1+(alpha.*R)).^2)+z.^2;
kprime_2 = kprime_2a./kprime_2b;
f = -0.011*(log(kprime_2))-0.0021;
y = sqrt(ap.*as).*(1/(2*pi)).*((cos(theta)-(d/as).*((cos(psi).*cos(phi))-(sin(psi).*sin(phi).*cos(theta))))/(R.^(3/2))).*f;
Then you can create an anonymous function that is a parameter of theta only, and integrate:
psi = 0;
phi = 0;
h = 3;
ap = 6;
f = @(theta) inductanceIntegrand(theta,psi,phi,h,ap);
M = quadl(f,0,2*pi);
disp(M)
-4.5074e-04
If you're integrating over phi, just change the argument of f to phi and provide a value for theta.

3 个评论
显示 1更早的评论 隐藏 1更早的评论

Andrew Matthews
Andrew Matthews 2013-6-7
编辑：Andrew Matthews 2013-6-7
Hello again. This makes sense, and I realized that is why my code was not working, it was integrating an answered problem. However, I tried the way you have it specified but now it is saying that I have stuff undefined, such as the variable as. I see you actually had an answer for M. How was that obtained? I am just getting errors. Thank you for your time. Here is what I used:
function y = inductanceIntegrand(phi,psi,theta,h,ap)
% Arguments
as = sqrt((ap^2)-(h^2));
delta = h/ap;
alpha = as/ap;
d = h.*tan(phi);
Ra = (1-(cos(phi).*cos(phi).*sin(theta).*sin(theta)));
Rb = ((2)*(d/as)).*((sin(psi).*sin(phi))(cos(psi).*cos(phi).*cos(theta))); Rc = ((d.^2)/(as.^2));
R = sqrt(Ra+Rb+Rc);
z = delta-(alpha*sin(theta).*cos(phi));
kprime_2a = ((1-(alpha.*R)).^2)+z.^2;
kprime_2b = ((1+(alpha.*R)).^2)+z.^2;
kprime_2 = kprime_2a./kprime_2b;
f = -0.011*(log(kprime_2))-0.0021;
y = sqrt(ap.*as).*(1/(2*pi)).*((cos(theta)-(d/as).*((cos(psi).*cos(phi))-(sin(psi).*sin(phi).*cos(theta))))/(R.^(3/2))).*f;
and then
psi = 0;
theta = 0:10:90;
h = 3;
ap = 6;
f = @(phi) inductanceIntegrand(phi,psi,theta,h,ap);
M = quadl(f,0,2*pi);
disp(M)
Thanks, any help would be appreciated.
Andrew Matthews
Andrew Matthews 2013-6-7
编辑：Andrew Matthews 2013-6-7
When I run the code the way you have it written it says 'requires more input arguments to run'. And when I have variables for input when I run the second m file it says the output and input vectors must be the same length. What do I do from here? Thank you.
Andrew Matthews
Andrew Matthews 2013-6-7
I just realize my error. I was trying to run to man variables at a time. Thank you all, you have been a great help.

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Numerical Integration and Differentiation 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by