f<-function(x,y,p){
if (p==0){q=c(x,y)}
else{q=c((f(x,y,p-1)[1])^2-(f(x,y,p-1)[2])^2+x,2*(f(x,y,p-1)[1])*(f(x,y,p-1)[2])+y)}
return(q)
}

f(1,2,0)
f(1,2,1)

mf<-function(x,y,p){(f(x,y,p)[1])^2+(f(x,y,p)[2])^2}
mf(1,2,0)

p=4
m=100
my.fun <- function(x,y){mf(x,y,p)-4}
x<-seq(-2,1,length=m)
y<-seq(-1.3,1.3,length=m)
z=matrix(0,m,m)

for (i in 1:m){ for (j in 1:m){z[i,j]=my.fun(x[i],y[j])}}
contour(x,y,z,level=0)