clear all;
clc;

syms x y z a1 a2 a3 b1 b2 b3 N
A = [a1;a2;a3];
B = [b1;b2;b3];
R = [x;y;z];
E = R/sqrt(x^2+y^2+z^2);

V = cross(A,cross(B,E))/dot(A',E);

div = divergence(b1*diff(V,x)+b2*diff(V,y)+b3*diff(V,z), R);
div_szamolt = 2*(dot(A',B)/dot(A',R))^2;

simplify(div-div_szamolt) % örülünk, ez nulla