function B = b (i,j,m)
min=0.0;max=1.0;
x=i*(max-min)/(m+1)+min;
y=j*(max-min)/(m+1)+min;

%with cos()
% -- eq=-(diff(cos(x)*diff(10*x*y*(1-x)*(1-y)*exp(x^4.5),x),x) + diff(cos(x)*diff(10*x*y*(1-x)*(1-y)*exp(x^4.5),y),y))
B=sin(x)*(10*y*(1-x)*(1-y)*exp(x^(9/2))-10*x*y*(1-y)*exp(x^(9/2))+45*x^(9/2)*y*(1-x)*(1-y)*exp(x^(9/2)))-cos(x)*(-20*y*(1-y)*exp(x^(9/2))+495/2*y*(1-x)*(1-y)*x^(7/2)*exp(x^(9/2))-90*x^(9/2)*y*(1-y)*exp(x^(9/2))+405/2*x^8*y*(1-x)*(1-y)*exp(x^(9/2)))+20*cos(x)*x*(1-x)*exp(x^(9/2));

%with sqrt(0.1+x)
% -- eq=-(diff(sqrt(0.1+x)*diff(10*x*y*(1-x)*(1-y)*exp(x^4.5),x),x) + diff(sqrt(0.1+x)*diff(10*x*y*(1-x)*(1-y)*exp(x^4.5),y),y))
%B=-5/(10+100*x)^(1/2)*(10*y*(1-x)*(1-y)*exp(x^(9/2))-10*x*y*(1-y)*exp(x^(9/2))+45*x^(9/2)*y*(1-x)*(1-y)*exp(x^(9/2)))-1/10*(10+100*x)^(1/2)*(-20*y*(1-y)*exp(x^(9/2))+495/2*y*(1-x)*(1-y)*x^(7/2)*exp(x^(9/2))-90*x^(9/2)*y*(1-y)*exp(x^(9/2))+405/2*x^8*y*(1-x)*(1-y)*exp(x^(9/2)))+2*(10+100*x)^(1/2)*x*(1-x)*exp(x^(9/2));

%xx and yy only
% -- eq=-(diff(diff(10*x*y*(1-x)*(1-y)*exp(x^4.5),x),x) + diff(diff(10*x*y*(1-x)*(1-y)*exp(x^4.5),y),y))
%B=20*y*(1-y)*exp(x^(9/2))-495/2*y*(1-x)*(1-y)*x^(7/2)*exp(x^(9/2))+90*x^(9/2)*y*(1-y)*exp(x^(9/2))-405/2*x^8*y*(1-x)*(1-y)*exp(x^(9/2))+20*x*(1-x)*exp(x^(9/2));

%eqn
%B=10*x*y*(1-x)*(1-y)*exp(x^(4.5));