#include <stdio.h>
#include <math.h>
#include <sys/time.h>
/*.proc files #include at end*/

#define TRUE 1 
#define FALSE 0
#define EMPTY (-1)       
     
#define LEVELS 20        
#define DLEVELS 10       
#define DRND 1.0/15.0
#define DALPHA 0.9524    
                   
#define NODESPERBLOCK 500                
#define MAXFILES 100                     
#define NODESINMEM 8000                 
#define BLOCKSINMEM 100                  
struct   pi {      
  float  Jmin,     
         denom[3][2][8],                 
         prob,     
         value[8], 
         cube;     
                   
  short  Opti,     
         Optj,     
         optflag,  
         child[3][2][8],
         left,
         right,
         parent;
};                 
typedef struct pi *list;                 
                   
short levels=LEVELS,                     
      jval,
      lrutimer=50; 
                   
int   next=1;      
      nnod[LEVELS+1],
      stop[LEVELS+1],                    
      virtmem[MAXFILES],                 
      physmem[BLOCKSINMEM],              
      used[BLOCKSINMEM],                 
      count=0, external=0,               
      space=0, maxdeep=0, made,                
      nodesinmem=NODESINMEM,             
      nodesperblock=NODESPERBLOCK,       
      blocksinmem=(NODESINMEM/NODESPERBLOCK);                  
                   
char datafile[20]; 
                   
float alpha=DALPHA,
      rnd=DRND,    
      sumprobj[2], sumprobi[3];
                   
                   
FILE  *fsens,      
      *fdata,      
      *ffact,      
      *fin,        
      *fout;       
                   
list  node[NODESINMEM];                  
                   
main(argc,argv)    
int argc;          
char *argv[]; {    
                   
  void  make_sens(),                    get_constants(),       
	sort_Jmin(),                    determine_optimum(),
	init_to_zero(),                 pick_mtrx(),
	set_p_mtrx(),                   get_make48(),
	open_files(),                   set_cu_mtrx(),
	print_tree(),                   set_ca_mtrx(),
	set_cx_mtrx(),                  unalloc(),
	set_q_mtrx(),                   write_pi(),                   
        graph_data(),
        get_pi_val(),                   close_files(),
        expand_tree(),                  data_head(),
        clock(),  cp_node_make(),
        init_tree(),                    make_fortyeight(),
        read_struct(),                  write_struct(),
        read_block(),                   write_block(),
        set_mem_loc(),                  set_file_ptr(),
        cp_node_find(),                 cp_node_sort(),
        init_node(),                    init_block(),
        dots();

  float assign_value(),
        round(), calc_denom(),
        p[3][8][8],
        q[2][8][8],
        ca[3][8], 
        cu[8],     
        cx[2][7];  
                   
  int   rndmax(), check_mem_size(),      
        get_insert_node(), menu(),       
        get(), lru(), insert_node(),     
        adjust_mem(), alloc(),           
        pchoice, qchoice,                
        start[LEVELS+1];                 
                   
  short time[3], z;
                   
  jval=menu();     
                   
  printf("\nNodesinmem=%d, Nodesperblock=%d, Blocksinmem=%d\n",nodesinmem,
          nodesperblock, blocksinmem);   
  init_to_zero(p, 3);                    
  init_to_zero(q, 2);                    
                   
  pick_mtrx(&pchoice, &qchoice);         
  get_constants(&rnd, &alpha, &levels);  
                   
  set_p_mtrx(p, pchoice);                
  set_q_mtrx(q, qchoice);                
  set_ca_mtrx(ca); 
  set_cu_mtrx(cu); 
  set_cx_mtrx(cx); 
                   
  init_tree(start);
  clock(time);     
                   
  expand_tree(start, p, q);              
                   
  determine_optimum(start, ca, cu, cx, p, qchoice);            
  sort_Jmin(p,q);  

  open_files(argc,argv);                 
  data_head(alpha, rnd, pchoice, qchoice, time, start);        
  for(z=0;z<levels;z++) {  
    printf(".");              
    make_sens (start[0], z, pchoice, qchoice);                 
    print_tree (start,z);                
  }
  close_files();

  printf("\n\nRun Time: %d:%02d:%02d",time[0],time[1],time[2]);
  printf("                   Jmin:%f.\n\n",                    
          node[get(start[0])]->Jmin);
  printf("Tree depth: %d\n",maxdeep);    
  printf("Total nodes: %d\n",next);      
  printf("Total calls to get: %d\n",count);                    
  printf("Calls to external locations: %d",external);          
  printf("   percentage: %.3f.\n",(100*((float) external)/((float) count)));
}                  
/********************************************************************************/
void expand_tree(start, p, q)            
int start[];       
float p[][8][8], q[][8][8]; {            
  int temp,z,overflow=0;                 
  int futnode=0, levnode=0;
  
  made=0;
  for(z=0;z<levels;z++) {
    nnod[z]=0;     
    temp=start[z]; 
    levnode=stop[z]-start[z];            
    if(overflow==1) {                    
      printf("Starting use of external files...\n");           
      overflow++;  
      space=alloc(space,nodesinmem);     
      if(space==(-1))                    
        break;     
    }              
    printf("Level %2d ", z+1);           
    if (overflow==2) {                   
      while(temp<stop[z]) {              
        get_make48(temp,start, p, q, z); 
        temp++;    
        dots(levnode,z);                 
      }            
      dots(-1,z);  
    }              
    else {         
      while(temp<stop[z]) {              
        make_fortyeight(temp,start, p, q, z);                  
        temp++;    
        dots(levnode,z);                 
      }            
      dots(-1,z);  
      if(z!=(levels-1)) { 
        if(!z)            
          futnode=next*next+100;   
        else
          futnode=next+(nnod[z]*nnod[z]/nnod[z-1])+100;
        if(futnode>nodesinmem) 
          overflow=1;
        else 
          if(futnode>space) { 
            if(futnode>(space+2000))
              space=alloc(space,futnode);
            else
              space=alloc(space,space+2000);
            if(space==(-1))
              break; 
          }
      }            
    }              
    stop[z+1]=next;
    printf("made: %06d  unique:%04d.\n",made,nnod[z]);                 
  }                
}                  
/*******************************************************************************/
void get_make48(dad,start, p, q, z)      
int start[],dad;   
float p[][8][8], q[][8][8];              
short z; {         
  short i, j, l, combined;               
  int kid, gdad, gkid;                   
  list tempnode;   
                   
  gdad=get(dad);   
  tempnode=(list) malloc(sizeof(struct pi));                   
  cp_node_make(node[gdad],tempnode);     
                   
  for(i=0;i<3;i++)                     
    for(j=jval;j<2;j++)                  
      for(l=0;l<8;l++){                  
        kid=next;  
        gkid=get(kid);                   
        node[gkid]->cube=assign_value(tempnode, node[gkid], i, j, l, p, q);
        if(node[gkid]->cube==0.0)        
          tempnode->child[i][j][l]=EMPTY;
        else {     
          tempnode->child[i][j][l]=kid;  
          combined=FALSE;                
          if(start[z+1]==EMPTY)          
            start[z+1]=kid;              
          else {   
            combined=get_insert_node(tempnode,start[z+1],i,j,l);
            gkid=get(kid);               
	  }        
          if(!combined) {                
            init_node(node[gkid]);       
            node[gkid]->parent=dad;      
            nnod[z]++;                   
            next++;
          }        
        }          
      }           
  cp_node_make(tempnode,node[get(dad)]);
  free((char *)tempnode); 
}                  
/*******************************************************************************/
int get_insert_node(ndad,slot,i,j,l)     
int slot;          
list ndad;         
short i,j,l; {     
  list temp;       
  float kcube;     
  int depth=1;     
                   
  kcube=node[get(ndad->child[i][j][l])]->cube;                 
                   
  while(TRUE) {    
    temp=node[get(slot)];                
    if(temp->cube==kcube) {              
        ndad->child[i][j][l]=slot;       
        return(TRUE);                    
    }              
    else {         
      if(kcube<temp->cube) {             
        if(temp->left==EMPTY) {          
          temp->left=ndad->child[i][j][l];                     
          if (depth>maxdeep)             
            maxdeep=depth;               
          return(FALSE);                 
        }          
        slot=temp->left;                 
      }            
      else {       
        if(temp->right==EMPTY) {
           temp->right=ndad->child[i][j][l];
           if (depth>maxdeep) 
             maxdeep=depth;
           return(FALSE);
        }
        slot=temp->right;
      }
    }
    depth++;
  }
}
/*******************************************************************************/
void make_fortyeight(dad,start, p, q, z)                   
int start[],dad;   
float p[][8][8], q[][8][8];              
short z; {         
  short i, j, l, combined;
  int kid;

  for(i=0;i<3;i++) 
    for(j=jval;j<2;j++)  
      for(l=0;l<8;l++){
        kid=next;
          node[kid]->cube=assign_value(node[dad],node[kid],i,j,l,p,q);
        if(node[kid]->cube==0.0)
          node[dad]->child[i][j][l]=EMPTY;
        else {
          node[dad]->child[i][j][l]=kid;
          combined=FALSE;
          if(start[z+1]==EMPTY)
            start[z+1]=kid;
          else {
            combined=insert_node(dad,start[z+1],i,j,l);
            made++;
          }
          if(!combined) {
            init_node(node[kid]);
            node[kid]->parent=dad;
            nnod[z]++;
            next++;
          }
        }
      }
}
/*******************************************************************************/
int insert_node(dad,slot,i,j,l)          
int slot,dad;      
short i,j,l; {     
  list temp,kid;   
  int depth=1;
  kid=node[node[dad]->child[i][j][l]];

  while(TRUE) {
    temp=node[slot];
    if(temp->cube==kid->cube) {
        node[dad]->child[i][j][l]=slot;
        return(TRUE);
    }
    else {
      if(kid->cube<temp->cube) {
        if(temp->left==EMPTY) {
          temp->left=node[dad]->child[i][j][l];
          if (depth>maxdeep) 
            maxdeep=depth;
          return(FALSE);
        }
        slot=temp->left;
      }
      else {
        if(temp->right==EMPTY) {
           temp->right=node[dad]->child[i][j][l];
           if (depth>maxdeep) 
             maxdeep=depth;
           return(FALSE);
        }
        slot=temp->right;
      }
    }
  depth++;
  }
}
/*******************************************************************************/
void init_node(side)                     
list side; {       
  int i,j,l;
  
  for(i=0;i<3;i++)
    for(j=jval;j<2;j++)
      for(l=0;l<8;l++) {
        side->child[i][j][l]=EMPTY;
        side->denom[i][j][l]=0.0;
      }
   side->prob=0.0;
   side->Jmin=0.0;
   side->left=EMPTY;
   side->right=EMPTY;
   side->parent=EMPTY;
   side->optflag=FALSE;
   side->Opti=side->Optj=EMPTY;
}
/*******************************************************************************/
void determine_optimum(start, ca, cu, cx, p, qchoice)          
int start[];       
int qchoice;       
float ca[][8], cu[], cx[2][7], p[][8][8]; {                    
  int temp;        
  float cofm, cofu, cofa, future, cost;  
  short z;         
  int i, j, l, m, k;                     
  list tempnode;   
  printf("Finding optimal tree");
  tempnode=(list) malloc(sizeof(struct pi));
  for(z=levels-1;z>(-1);z--) {
    printf(".");
    temp=start[z]; 
    while(temp<stop[z]) {                
      init_node(tempnode);               
      cp_node_find(node[get(temp)],tempnode);                  
      tempnode->Jmin=10000000.0;         
      for(i=0;i<3;i++) {                 
	cofa=cofu=0.0;                   
	for(k=0;k<8;k++) {               
	  cofa+=(tempnode->value[k]*ca[i][k]);                 
	  for(m=0;m<8;m++)               
	    cofu+=(tempnode->value[k]*p[i][k][m]*cu[m]);       
	}          
	cofu*=alpha;                     
	for(j=jval;j<2;j++) {            
	  cofm=alpha*cx[j][qchoice];     
	  future=0.0;                    
          cost=cofa+cofu+cofm;           
	  for(l=0;l<8;l++) {             
            if(tempnode->child[i][j][l]!=EMPTY)                
	      future+=(tempnode->denom[i][j][l]*               
		       node[get(tempnode->child[i][j][l])]->Jmin);
          }        
	  future=(alpha*future)+cost;    
	  if(future<tempnode->Jmin) {    
	    tempnode->Jmin=future;       
	    tempnode->Opti=i;            
	    tempnode->Optj=j;            
          }        
        }          
      }            
      cp_node_find(tempnode,node[get(temp)]);                  
      temp++;      
    }              
  }                
}                  
/*******************************************************************************/
void sort_Jmin(p, q) 
float p[][8][8], q[][8][8]; {
  short l, iopt, jopt, temp=0, gkid;
  list tempnode;
  float tdenom;

  printf("Marking tree...\n");
  tempnode=(list) malloc(sizeof(struct pi));

  while(temp<stop[levels]) {             
    if(node[get(temp)]->optflag==TRUE) { 
      cp_node_sort(node[get(temp)],tempnode);                  
      iopt=tempnode->Opti;               
      jopt=tempnode->Optj;               
      if((iopt!=EMPTY)&&(jopt!=EMPTY)) { 
        for(l=0;l<8;l++)                 
          if(tempnode->child[iopt][jopt][l]!=EMPTY) {          
            gkid=get(tempnode->child[iopt][jopt][l]);          
            node[gkid]->optflag=TRUE;    
            tdenom=calc_denom(tempnode->value,iopt,jopt,l,p,q);
            node[gkid]->prob+=tempnode->prob*tdenom;           
          }        
      }            
    }
    temp++;
  }
}
/*******************************************************************************/
float calc_denom(pi,i,j,l,p,q)           
float pi[8],p[][8][8], q[][8][8];        
int i,j,l; {       
  float dummy, dummy2;
  int m,n;

  dummy2=0.0;
  for(m=0;m<8;m++) {
    dummy=0.0;
    for(n=0;n<8;n++)
      dummy+=(p[i][n][m]*pi[n]);
    dummy2+=(q[j][m][l]*dummy);
  }
  return(dummy2);
}
/******************************************************************************/
void init_to_zero(a,i)                   
short i;
float a[][8][8];{
  short loop1, loop2, loop3;

  for(loop1=0;loop1<i;loop1++)  
    for(loop2=0;loop2<8;loop2++)
	for(loop3=0;loop3<8;loop3++)
	  a[loop1][loop2][loop3]=0.0;
}
/******************************************************************************/
void pick_mtrx(pchoice, qchoice)         
int *pchoice,*qchoice; {
  printf("\nEnter choices of matrixes...\n\n");
  printf("  Choice of p matrix:");
  scanf("%d", pchoice);
  printf("  Choice of q matrix:");
  scanf("%d", qchoice);
  printf("\n");
}
/******************************************************************************/
void get_pi_val() {
  int k;

  printf("\nEnter intial pi values...\n\n");
  for(k=0; k<8; k++){
    printf("  Pi value[%d]=",k);
    scanf("%f,",&node[0]->value[k]);
  }
  printf("\n");
}
/******************************************************************************/
void init_tree(start)
int start[]; {
  int i;

  space=alloc(0,2000);
                   
  for(i=1;i<LEVELS+1;i++)                
    start[i]=stop[i]=(-1);               
  start[0]=0;      
  stop[0]=1;       
                   
  init_node(node[0]);                    
                   
  node[0]->prob=1.0;                     
  node[0]->optflag=TRUE;                 
  get_pi_val();    
                   
  node[0]->cube=round(node[0]->value,rnd);                     
                   
  setbuf(stdout,NULL);                   
  set_mem_loc();   
}                  
/******************************************************************************/
void get_constants(rnd,alpha,levels)     
short *levels;     
float *rnd, *alpha; {                    
  printf("\nTo assign the default value to the constants, enter a zero.\n\n");
  printf("  Enter value for alpha (default is %.4f):",DALPHA);
  scanf("%f",alpha);
  if (*alpha==0.0)
    *alpha=DALPHA;
  printf("  Enter value for levels (default is %d, maximum is %d):",
         DLEVELS,LEVELS);
  scanf("%d",levels);
  if (*levels==0)
    *levels=DLEVELS;
  printf("  Enter value for the interval (default is %.3f, minimum is 0.005):",
         DRND);
  scanf("%f",rnd);
  if (*rnd==0.0)
    *rnd=DRND;
}
/******************************************************************************/
int alloc(bottom,top)                    
int bottom,top; {  
  int i,pisize;    

  if(top>nodesinmem)
    top=nodesinmem;
  pisize=sizeof(struct pi);
  for(i=bottom;i<top;i++) {
    node[i]=(list)malloc(pisize);
    if (node[i]==NULL) {
      printf("\nError in malloc for node %d. ",i);
      printf("Stopping malloc procedure...\n");
      i=top;
    }
  }
  return(top);
}
/******************************************************************************/
int menu() {       
  int choice;      
                   
  printf("\n\n\nSelect:\n\n");           
  printf("  Standard run:         0\n"); 
  printf("  Constrain J to 1:     1\n"); 
  printf("  Adjust Memory Bounds: 2\n");
  printf("  Check memory size:    3\n");
  printf("  Choice:");
  scanf("%d",&choice);
  
  if(choice==3)
    return(check_mem_size());
  if(choice==2)
    return(adjust_mem());
  if(choice==1)
    return(1);
  else {
    if(choice!=0)
      printf("Incorrect selection. Executing standard run.\n\n");
    return(0);
  }
}
/******************************************************************************/
int check_mem_size() {                   
  int i,j,pisize,blocksize=512;          
                   
  pisize=sizeof(struct pi);

  printf("\n\nElements in Node Array: %d.\n",nodesinmem);
  printf("Size of Node Structure: %d bytes.\n",pisize);
  printf("Enter Blocksize: ");
    scanf("%d",&blocksize);

  setbuf(stdout,NULL);
  printf("\nchecking...");
  for(i=0;i<50000;i++) 
    if ((node[i]=(list)malloc(blocksize))==0) 
      break;
  printf("freeing...");
  for(j=0;j<(i+1);j++)
    free((char *) node[j]);
  printf("done.\n\n");

  printf("Blocks Allocated: %d.  Block Size: %d.\n",i,blocksize);
  printf("Amount of Memory Used: %d bytes.\n",blocksize*i);
  printf("Maximum Possible Elements in Node Array: %d.\n",
         (blocksize*(i+1))/pisize);
  printf("Recommended Maximum:%d.\n",((blocksize*(i+1))/512)-1000);
  return(menu());
}
/******************************************************************************/
int adjust_mem() { 
                   
  printf("\n\n  Enter number of internal nodes (less than %d):",NODESINMEM);
    scanf("%d",&nodesinmem);
  printf("  Enter memory blocksize (factor of above number):");
    scanf("%d",&nodesperblock);
  if((nodesinmem/nodesperblock)>BLOCKSINMEM)
    printf("    Err: Blocksize too small\n\n");
  if(nodesperblock>nodesinmem)
    printf("    Err: Blocksize to large.");
  printf("  Enter lru time interval (between 50 & 50000):");
    scanf("%d",&lrutimer);
  blocksinmem=(nodesinmem/nodesperblock);
  return(menu());
}
/******************************************************************************/
void dots(levnode,z)                     
int levnode, z; {  
  static float dots;                     
  static int line;
  float interval;
  float length=40.0;

  if(levnode==(-1)) {
    while(line<((int) length)) {
      printf(".");
      line++;
    }
    line=0;
    dots=0.0;
  }
  else if(z>0){
    interval=(float)levnode/length;
    dots+=1.0;
    while ((dots>interval)&&(line<((int)length))) {
      printf(".");
      dots-=interval;
      line++;
    }
  }
}
/******************************************************************************/
#include "proc/assign_value.proc"
#include "proc/set_mtrx.proc~~"
#include "proc/cp_node.proc"
#include "proc/outfiles.fix"
#include "proc/clock.proc"
#include "proc/get.proc"
