I want to plot a numerical integral in a mesh which is dependent of x and y. How can I do this?

조회 수: 5 (최근 30일)
I have a two dimensional function and I want to plot the integral function in 2D mesh. I want to plot the C2 in color mesh and also want to plot diagonal element of C2. I am not getting the correct result. Please see the code I am working.
if true
%
L = 20;
z1 = 0:0.01:L;
z2 = 0:0.01:L;
dz = z1(2)-z1(1);
[Z1,Z2] = meshgrid(z1,z2);
W1 = 1+Z1/L;
W2 = 2+Z2/L;
for j = 1:size(z1,2)
for k= 1:size(z2,2)
A(j,k) = sum(sum(1./W1(1:j,1:k))*dz)*dz;
B(j,k) = sum(sum(1./W2(1:j,1:k))*dz)*dz;
end
end
C2_h = (1-A).^2+(1-B);
figure(32);pcolor(Z1,Z2,C2_h);axis equal tight; shading flat;colorbar;colormap jet;
figure(5);
plot(z1/L, diag(C2_h) ,'r','Linewidth',2 );hold on;
end

답변 (1개)

Anton Semechko
Anton Semechko 2018년 7월 4일
The integrals in your equation can be evaluated analytically. Here how you can visualize C(z1,z2) and C(z1,z1):
L = 20;
dz = 0.01;
z = 0:dz:L;
[Z1,Z2] = meshgrid(z);
I1=@(z) log(1+z);
I2=@(z) log(2+z);
C=(1-I1(Z1)).^2 + (1+log(2)-I2(Z2));
figure('color','w')
ha=subplot(1,2,1);
z_lim=[z(1) z(end)];
imagesc(z_lim,z_lim,C)
axis equal
hold on
plot(z_lim(:),z_lim(:),'--k','LineWidth',2)
set(ha,'XLim',z_lim+[-1 1],'YLim',z_lim+[-1 1],'YDir','normal')
xlabel('z_1','Interpreter','tex','FontSize',20)
ylabel('z_2','Interpreter','tex','FontSize',20)
colorbar('Location','NorthOutside')
subplot(1,2,2)
Cd=diag(C);
plot(z(:),Cd(:),'--k','LineWidth',2)
xlabel('z_1','Interpreter','tex','FontSize',20)
ylabel('C(z_1,z_1)','Interpreter','tex','Interpreter','tex','FontSize',20)
  댓글 수: 2
mk_ballav
mk_ballav 2018년 7월 4일
편집: mk_ballav 2018년 7월 4일
This is one of the special case where W1 and W2 are linear functions but I have other cases where W1 and W2 are not linear and I can't directly evaluate integral anlytically,so I have to do numerical integration. Here C(z1,z2) is the whole matrix elements and C(z1,z1) is just the diagonal elements of C.
Anton Semechko
Anton Semechko 2018년 7월 4일
OK. Assuming W1 and W2 are univariate functions of z1 and z2, respectively, you can integrate 1/W1 and 1/W2 numerically using builtin 'integral' function.. Here is an example:
% Visualization grid
L = 20;
dz = 0.01;
z = 0:dz:L;
% Evaluate definite integrals of 1/W1 and 1/W2 beween z(i-1) and z(i)
W1=@(z) 1./(z+1);
W2=@(z) 1./(z+2);
[I1,I2]=deal(zeros(size(z)));
for i=2:numel(z)
I1(i)=integral(W1,z(i-1),z(i));
I2(i)=integral(W2,z(i-1),z(i));
end
% Definite intergral between z(1) and z(i) = sum(I(1:i))
I1=cumsum(I1);
I2=cumsum(I2);
% Function
C=bsxfun(@plus,(1-I1).^2,(1-I2(:)));
% Visualization
figure('color','w')
ha=subplot(1,2,1);
z_lim=[z(1) z(end)];
imagesc(z_lim,z_lim,C)
axis equal
hold on
plot(z_lim(:),z_lim(:),'--k','LineWidth',2)
set(ha,'XLim',z_lim+[-1 1],'YLim',z_lim+[-1 1],'YDir','normal')
xlabel('z_1','Interpreter','tex','FontSize',20)
ylabel('z_2','Interpreter','tex','FontSize',20)
colorbar('Location','NorthOutside')
subplot(1,2,2)
Cd=diag(C);
plot(z(:),Cd(:),'--k','LineWidth',2)
xlabel('z_1','Interpreter','tex','FontSize',20)
ylabel('C(z_1,z_1)','Interpreter','tex','Interpreter','tex','FontSize',20)
You can confirm that the final result is visually indistinguishable from the one obtained analytically.

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Surface and Mesh Plots에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by