/* 単精度での情報落ちの検査 */
/* C コンパイラ gcc が使える環境において, 次を実行してみよ. 
/* gcc rep1-2.c -o rep1-2
   rep1-2
*/
#include <stdio.h>
#include <stdlib.h>
void memoryimage(float x,char[12]);

int main(void) {
    unsigned long k,ke; /* 符号なし unsigned*/
    float r,t,s;
    double sd;
    char moji[12];
    k = 1;
    t = (float)0;
    do {
      r=s; s = t;
      t = t + (float)1/k;
      /*       memoryimage(t,moji);
	       printf("%lu:%f  %s\n",k,t,moji); */
      k++;
    } while (s != t);
    k-=2;
    ke=k;
    printf("Stabilized to the value %f at N = %lu\n",s,k); 
    memoryimage(r,moji);
    printf("%lu:%f  %s\n",k-1,r,moji);
    memoryimage(s,moji);
    printf("%lu:%f  %s\n",k,s,moji);
    memoryimage(t,moji);
    printf("%lu:%f  %s\n",k+1,t,moji);
    s = (float)0;
    while (k>0) {
        s = s + (float)1/k;
        k--;
    }
    printf("Result of backward summation = %f\n",s); 
    sd = (double)0;
    for (k=1;k<=ke;k++){
        sd=sd+(double)1/k;
    }
    printf("Cf. Forward summation of the same range via double = %18.15f\n",sd);
}

char hex(unsigned short x,char m[2]){
  unsigned int k;
  k=x/16;
  if (k<10) m[0]='0'+k;
  else m[0]='A'+k-10;
  k=x%16;
  if (k<10) m[1]='0'+k;
  else m[1]='A'+k-10;
}
void memoryimage(float x,char m[12]){
  unsigned char *p;
  p=(unsigned char *)&x;
  hex(*p,&m[0]);
  m[2]=' ';
  hex(*(p+1),&m[3]);
  m[5]=' ';
  hex(*(p+2),&m[6]);
  m[8]=' ';
  hex(*(p+3),&m[9]);
  m[11]=0;
}