clear all
close all

Ix  = [-1 1];
Iy  = [-1 1];
DIM = [20 20];
hx  = (Ix(2)-Ix(1))/(DIM(1)-1);
hy  = (Iy(2)-Iy(1))/(DIM(2)-1);

[x1, x2] = meshgrid(Ix(1):hx:Ix(2),Iy(1):hy:Iy(2));

u1 = 2*x1.^2.*x2;
u2 = 3*x1-x2;

Du11 = @(x1,x2) 4*x1.*x2;
Du12 = @(x1,x2) 2*x1.^2;
Du21 = @(x1,x2) 3*ones(length(x1));
Du22 = @(x1,x2) -1*ones(length(x1));

div = @(x1,x2) Du11(x1,x2)+Du22(x1,x2);
rot = @(x1,x2) Du21(x1,x2)-Du12(x1,x2);
bet = @(x1,x2) sqrt(u1.^2+u2.^2);

subplot(2,2,1)
% plot vector field
quiver(x1, x2, u1, u2)
title('Vektorfeld')
grid
xlabel('x1');
ylabel('x2');

x3=div(x1,x2);
x4=rot(x1,x2);
x5=bet(x1,x2);

subplot(2,2,2)
surfc(x1,x2,x3);
title('Divergenz')
xlabel('x1');
ylabel('x2');
subplot(2,2,3)
surfc(x1,x2,x4);
title('Rotation')
xlabel('x1');
ylabel('x2');
subplot(2,2,4)
surfc(x1,x2,x5);
title('||u||')
xlabel('x1');
ylabel('x2');

