/* Farlie-Gumbel-Morgenstern Copula ( C(u,v)=u.*v.*(1+theta*(1-u).*(1-v)) ) */
/* Solve q=v.*(1+theta-2*theta*u-theta*v+2*theta*u.*v) (=Cu) in terms of v. */
new; cls;
theta=-0.99;
n=1000;
U=Cfgm(theta,n);

library pgraph;
graphset;
_plctrl=-1;
xy(U[.,1],U[.,2]);



proc Cfgm(theta,n);
   local u,q,a,b,c,v;
   if theta<=-1 or theta>=1;
      errorlog "ERROR: Parameter theta must be -1<theta<1.";
      retp(".");
   endif;
   u=rndu(n,1);
   q=rndu(n,1);
   if theta==0;
      v=q;
   else;
      a=2*theta*u-theta; b=1+theta-2*theta*u; c=-q;
      v=1./(2*a).*(-b+sqrt(b^2-4*a.*c));
   endif;
   retp(u~v);
endp;