function z_1 = integrand_1(x, y)
  global  h x_min x_max y_min y_max
  
%%    For Q = A^(-1) case
  %z_1 = (y-y_min).*(x-x_min).*Green(x,y);  
  
%%    For Gaussian case  
  Sigma = 1/2^10;
  z_1 = (y-y_min).*(x-x_min).*exp(-(x-y).*(x-y)./(2*Sigma^2));