/* 
complex.c : Komplexe Zahlen, zu uebersetzen mit :
    gcc complex.c -lm 
*/

#define PI 3.14159265358979323846L

#include <stdio.h>
#include <math.h>

struct complex { 
  double r; double i; 
};

typedef struct complex comp;   /* 'comp' wird Synomym fuer 'struct complex' */

comp   makecomp(double, double);
comp   sum(comp, comp);
comp   product(comp*, comp*);  /* Zeigervariable zur Demonstration  */
comp   power(comp, int);
void   compprint(comp);
double ab(comp);              /* Absolutbetrag */
double arg(comp);              /* Argument im Bogenmass */
double warg(comp);             /* Argument im Gradmass (als "Winkel") */


#define N 16

int main() {                       /* Einige Rechnungen mit komplexen Zahlen */
  int k;
  comp arr[N];                 /* Array komplexer Zahlen */  
  comp basis, result;
  double w = 1/sqrt(2);
  /* basis    = makecomp(w,w);          8-te Einheitswurzel  */
  basis    = makecomp(1,1);       
  for (k = 0; k < N; k++) {         
    arr[k] = power( basis, k+1 ); 
    compprint(arr[k]);
    printf(" .. Betrag: %6.2f : Argument: %6.2f : Winkel: %6.2f\n", 
	   ab(arr[k]), arg(arr[k]), warg(arr[k]));
  }
  result = makecomp(0,0);
  for (k = 0; k < N; k++) 
    result = sum(result, arr[k]);
  compprint(result);
  printf("\n");
}

comp makecomp(double r, double i) { /*  baut komplexe Zahl mit      */
  comp u = {r,i};                   /*  Realteil r, Imaginaerteil i */
  return u;
}

comp sum(comp a, comp b) {           /* Summe zweier komplexer Zahlen */
  comp u;
  u.r = a.r + b.r;
  u.i = a.i + b.i;
  return u;
}

double ab(comp x) {
  return sqrt(x.r*x.r + x.i*x.i);
}

double arg(comp x) {
  if (x.r == 0) return (x.i == 0) ? 0 : ( x.i > 0 ? 0.5*PI : 1.5*PI);
  if (x.r <  0) return atan(x.i/x.r) + PI;
  return  (x.i >= 0) ? atan(x.i/x.r) : atan(x.i/x.r) + 2*PI ;
}

double warg(comp x) {
  return arg(x)/PI*180; 
}

/* Hier werden zur Demonstration Zeigervariablen fuer die Argumente verwendet */

comp product(comp *x, comp *y) {     /* Produkt zweier komplexer Zahlen */
  comp u;
  u.r = x->r * y->r - x->i * y->i;   /* -> als Zugriff auf Komponente */
  u.i = x->r * y->i + x->i * y->r;
  return u;
}

comp power( comp b, int n) {  /* b hoch n, b komplex */
  comp u = {1, 0};                   /* 1 komplex */
  while (n > 0) {
    if (n % 2) {
      n--; u = product(&b, &u); /* Adressen werden uebergeben */
    }
    else {
       n /= 2; b = product(&b, &b);
    }
  }                    /* beachte: b^{2 k + 1} = (b^2)^k * b */
  return u;
}


void compprint(comp z) {            /* Druckt eine komplexe Zahl  */ 
  if (z.i >= 0)
    printf("%6.2f + %6.2f i", z.r,  z.i);
  else 
    printf("%6.2f - %6.2f i", z.r, -z.i);
}