#include <iostream>
#include <complex>

extern "C"{
  void wofz_(double* xi, double* yi, double* ui, double* vi, int* flag);
}


std::complex<double> erfc(const std::complex<double>& z)
// C++ wrapper function for call to fortran routine to calculate erfc of a complex number
{
  int info;
  double rz = z.real(), iz = z.imag();
  double ru, iu;
  wofz_(&rz, &iz, &ru, &iu, &info);
  if (info!=0) std::cerr <<"Can't compute error function at "<<z<<" got "<<ru<<","<<iu<<" info="<<info<<std::endl;
  return std::complex<double>(ru,iu); // return value optimization!
}

int main()
{
  const int M=200;
  const double L=10;
  
  using namespace std;
 
  for (int i=0; i<M; i++){
    double x = -L/2 + (L*i)/(M-1);
    complex<double> z = erfc(complex<double>(x,0.01));
    cout<<x<<" "<<z.real()<<" "<<z.imag()<<endl;
  }
  return 0;
}
