MATLAB Answers

0

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

Asked by mk_ballav on 4 Jul 2018
Latest activity Commented on by Anton Semechko on 4 Jul 2018
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

  0 Comments

Sign in to comment.

1 Answer

Answer by Anton Semechko on 4 Jul 2018

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 Comments

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.
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.

Sign in to comment.