#This file finds the complex variable form of a real harmonic function, in the form f=g+bar{g} or f=2*Re{g} (see equation (18.24)).

#Replace the next line by the required real harmonic function in Cartesian coordinates, or comment it, specify the function f(x,y) in the exterior part of the program and then call this file as a subroutine.

#f:=ln(x^2+y^2);

f:=expand(subs({x=(z1+z2)/2,y=I*(z2-z1)/2},f));

#Evaluate g^prime from equation (18.26).

gp:=diff(f,z1);

#Integrate to recover g.

g:=int(gp,z1)+A;

#Complex form of f from (18.24).

ff:=simplify(g+subs({z1=z2,I=-I},g));

err:=simplify(ff-f);

solution:=solve({err=0},{A});

ff:=simplify(subs(solution,ff));