Green's function, "double mesh", and weak singularities.

Hello, 

I am trying to compute a first order corrector for a perturbed resonance problem. To do this I need to calculate 
Where B and B' are the unit squares
x = [x_1, x_2]
x' = [x_1',x_2']
G is the Hankel function of the first kind of order 0 (it has a log singularity when x=x')
u and u* are found using an Arnoldi solver (similar to 1.7.1 in the tutorials) and are a resonance function and its conjugate. 

I can do this using base python but the estimated runtime is very large and I know NGSolve is much smarter about doing these types of problems than my basic python for loops.

My questions are as follows: 
First, to evaluate the integral I need u, and u* to be defined on duplicate mesh (same mesh twice once for u and one for u*) but with different variables and I'm unsure how to do such a thing. Do I need to define a new mesh (different from the one used in the Arnoldi) or can I just define functions such that u(x_1, x_2) and u*(x_1', x_2') and then use the same mesh twice somehow. 

Second, can the integral (or assemble or apply) handle the weak singularity since or is there some built-in feature that deals with this issue. (k is the resonance produced from the Arnoldi solver). If not how would I go about modifying integrate, or assemble, or apply such that if NaN occurs do this instead while still retaining the speed of NGSolve.

Any help or ideas would be greatly appreciated
Alex