How to convert a symbolic output to its numerical form ( for vector integration) ?
조회 수: 5 (최근 30일)
이전 댓글 표시
I am stuck in a resilient difficulty due to my inexperience with vector computation.
The code is all about calculating outward vector integrals on a parameterized surface.
The parametric form I have used here is a superellipsoid, which has a general form
(x^eps1)/a1+(y^eps2)/a2+(z^eps3)/a3=1.
With suitable parameterization, I have formed a matrix (x, y, z) which returns the geometry for superellipsoidal shapes, as,
function [x,y,z]=superellipse(epsilon1,epsilon2,epsilon3,a1,a2,a3)
n=100;
etamax=pi/2;
etamin=-pi/2;
wmax=pi;
wmin=-pi;
deta=(etamax-etamin)/n;
dw=(wmax-wmin)/n;
[i,j] = meshgrid(1:n+1,1:n+1);
eta = etamin + (i-1) * deta;
w = wmin + (j-1) * dw;
x = a1.*sign(cos(eta)).*abs(cos(eta)).^epsilon1.*sign(cos(w)).*abs(cos(w)).^epsilon1;
y = a2.*sign(cos(eta)).*abs(cos(eta)).^epsilon2.* sign(sin(w)).*abs(sin(w)).^epsilon2;
z = a3.*sign(sin(eta)).*abs(sin(eta)).^epsilon3;
surf(x,y,z);
surf2stl('cube3D1.stl',x,y,z,'binary');
end
now for the second part, we need to compute the pressure tensor on the surface points. I follow the procedure of integral mentioned here. Please follow the link.
syms eta w real
ellip=[(a1.*sign(cos(eta)).*abs(cos(eta)).^epsilon1.*sign(cos(w)).*abs(cos(w)).^epsilon1),(a2.*sign(cos(eta)).*abs(cos(eta)).^epsilon2.*sign(sin(w)).*abs(sin(w)).^epsilon2),(a3.*sign(sin(eta)).*abs(sin(eta)).^epsilon3)];
F=[(a1.*sin(eta).*cos(w)),(a2.*sin(eta).*sin(w)),(a3.*cos(eta))];
nds=simplify(cross(diff(ellip,eta),diff(ellip,w)));hold on
Fpar=subs(F,[(a1.*sin(eta).*cos(w)),(a2.*sin(eta).*sin(w)),(a3.*cos(eta))],ellip);
Fj=@(Fpar,nds)Fpar.*transpose(nds);
flux=(symint2(Fj,eta,-pi,pi,w,-(pi/2),(pi/2)));
I have solved the integral analytically with Matlab syms function. But, after computing the integral,
the values have been returned in 'char' form. I need to have the values in 'double'.
I have used a function symint2, which is attached herewith.
Please take a look at the code.
Hirak
댓글 수: 0
답변 (3개)
madhan ravi
2018년 12월 31일
편집: madhan ravi
2018년 12월 31일
Simply use double() or vpa() to convert symbolic answer to numerical also have a look at matlabFunction() to convert symbolic operations to numeric ones thereby you can use numerical integration (integral()) too.
댓글 수: 7
madhan ravi
2019년 1월 5일
You can unaccept my answer so that someone else can help , I am not able to figure it out.
Walter Roberson
2019년 1월 6일
flux=vpaintegral(Fj(eta,w),eta,-pi,pi,w,-(pi/2),(pi/2));
vpaintegral() does not have a syntax for 2D integrals. You need
flux=vpaintegral(vpaintegral(Fj(eta,w),eta,-pi,pi),w,-(pi/2),(pi/2));
This will give you an answer that is 0 to within round-off.
This is because...
int(Fj(eta,w),eta,-pi,pi)
is 0 becaus Fj(eta,w) is simply eta*w
You need to look more closely at your lines
nds=cross(jacobian(ellip_par,eta),jacobian(ellip_par,w));
Fj=@(Fpar,nds)Fpar.*transpose(nds);
%Fparam=matlabFunction(Fj);
flux=vpaintegral(Fj(eta,w),eta,-pi,pi,w,-(pi/2),(pi/2));
The first of those lines defines a variable named nds . The second line defines an anonymous function with dummy parameter name nds that does a multiplication involving the dummy variable nds . With nds being a dummy variable name, the nds after the @(Fpar,nds) has nothing to do with the variable nds that was define in the previous line.
Now, the third line invokes Fj(eta,w) with the symbolic variables eta and w. Those get passed into Fj(), where they become locally known as Fpar and nds . transpose(nds) is then transpose(w) which is just w. Fpar is the passed in eta. So Fj(eta,w) is going to be just eta .* w .
I have modified your code to get closer, but I am not convinced that this is what you are looking for.
댓글 수: 0
Hirak
2019년 1월 7일
댓글 수: 2
Walter Roberson
2019년 1월 7일
We do not have a value for thetamin or dtheta or phimin so we cannot execute
theta = thetamin + (i-1) * dtheta;
phi= phimin + (j-1) * dphi;
On the other hand your loops after that go
for theta=i
i=i+1;
for phi=j
j=j+1;
superquad2(a1,a2,a3,epsilon1,epsilon2,epsilon3,theta,phi);
end
end
which overwrite theta and phi with entries from the results of the meshgrid() . After the meshgrid, i and j will be 2D arrays of results. When you use for theta=i and i is 2D then the result is defined as assigning theta to be successive columns of i, so the first iteration through theta would be the entire column vector i(:,1), the second iteration it would be the entire column vector i(:,2) and so on. It is unlilkely that is the result you want. Also, you are not assigning the results of the superquad2 call to anything.
참고 항목
카테고리
Help Center 및 File Exchange에서 Calculus에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!