[BACK]Return to phc6.c CVS log [TXT][DIR] Up to [local] / OpenXM / src / phc

Diff for /OpenXM/src/phc/phc6.c between version 1.3 and 1.4

version 1.3, 2000/10/31 12:22:05 version 1.4, 2002/05/02 02:51:37
Line 1 
Line 1 
 /*  phc6.c ,  yama:1999/sm1-prog/phc6.c */  /*  phc6.c ,  yama:1999/sm1-prog/phc6.c */
 /* $OpenXM$ */  /* $OpenXM: OpenXM/src/phc/phc6.c,v 1.3 2000/10/31 12:22:05 takayama Exp $ */
 /* This is a simple C interface to the black-box solver of phc.  /* This is a simple C interface to the black-box solver of phc.
 ** Requirements:  ** Requirements:
 **  1) executable version of phc will be searched in the following order:  **  1) executable version of phc will be searched in the following order:
Line 26  union cell {
Line 26  union cell {
   int ival;    int ival;
   char *str;    char *str;
   struct phc_object *op;    struct phc_object *op;
   long double longdouble;    double longdouble;
 };  };
 struct phc_object{  struct phc_object{
   int tag;                /* class identifier */    int tag;                /* class identifier */
Line 42  struct phc_object{
Line 42  struct phc_object{
 /********** macros to use Sarray **************/  /********** macros to use Sarray **************/
 /* put to Object Array */  /* put to Object Array */
 #define phc_putoa(ob,i,cc) {\  #define phc_putoa(ob,i,cc) {\
 if ((ob).tag != Sarray) {fprintf(stderr,"Warning: PUTOA is for an array of objects\n");} else \                     if ((ob).tag != Sarray) {fprintf(stderr,"Warning: PUTOA is for an array of objects\n");} else \
 {if ((0 <= (i)) && ((i) < (ob).lc.ival)) {\                       {if ((0 <= (i)) && ((i) < (ob).lc.ival)) {\
   (ob.rc.op)[i] = cc;\                        (ob.rc.op)[i] = cc;\
 }else{\                      }else{\
   fprintf(stderr,"Warning: PUTOA, the size is %d.\n",(ob).lc.ival);\                         fprintf(stderr,"Warning: PUTOA, the size is %d.\n",(ob).lc.ival);\
 }}}                    }}}
 #define phc_getoa(ob,i) ((ob.rc.op)[i])  #define phc_getoa(ob,i) ((ob.rc.op)[i])
 #define phc_getoaSize(ob) ((ob).lc.ival)  #define phc_getoaSize(ob) ((ob).lc.ival)
   
Line 56  struct phc_object phc_newObjectArray(int size);
Line 56  struct phc_object phc_newObjectArray(int size);
 void phc_printObject(FILE *fp,struct phc_object ob);  void phc_printObject(FILE *fp,struct phc_object ob);
 char *phc_generateUniqueFileName(char *s);  char *phc_generateUniqueFileName(char *s);
 char *phc_which(char *s);   /* search a path for the file s */  char *phc_which(char *s);   /* search a path for the file s */
 struct phc_object phc_complexTo(long double r, long double i);  struct phc_object phc_complexTo(double r, double i);
   
   
 int phc_scan_for_string(FILE *fp, char str[]);  int phc_scan_for_string(FILE *fp, char str[]);
Line 66  struct phc_object phc_call_phc(char *sys); 
Line 66  struct phc_object phc_call_phc(char *sys); 
   
 int phc_verbose = 0;  int phc_verbose = 0;
 int phc_overwrite = 1;  /* Always use tmp.input.0 and tmp.output.0  int phc_overwrite = 1;  /* Always use tmp.input.0 and tmp.output.0
                            for work files. */                             for work files. */
   
 main(int argc, char *argv[]) {  main(int argc, char *argv[]) {
   struct phc_object ob;    struct phc_object ob;
   int n,i,dim;    int n,i,dim;
 #define INPUTSIZE 4096  #define INPUTSIZE 4096
   char input[INPUTSIZE];    char input[INPUTSIZE];
     char fname[INPUTSIZE];
 #define A_SIZE 1024  #define A_SIZE 1024
   char a[A_SIZE];    char a[A_SIZE];
   int message = 0;    int message = 0;
Line 81  main(int argc, char *argv[]) {
Line 82  main(int argc, char *argv[]) {
       phc_verbose = 1; message=1;        phc_verbose = 1; message=1;
     }else if (strcmp(argv[i],"-g") == 0) {      }else if (strcmp(argv[i],"-g") == 0) {
       phc_overwrite = 0;        phc_overwrite = 0;
       }else if (strcmp(argv[i],"-file") == 0) {
         /* For debugging. */
         i++;
         strncpy(fname,argv[i],INPUTSIZE-1);
         ob = phc_scan_output_of_phc(fname);
         n = phc_getoaSize(ob);
         printf("[\n");
         for (i=0; i<n; i++) {
           phc_printObject(stdout,phc_getoa(ob,i));
           if (i != n-1) printf(" ,\n"); else printf(" \n");
         }
         printf("]\n");
         exit(0);
     }else if (strcmp(argv[i],"-i") == 0) {      }else if (strcmp(argv[i],"-i") == 0) {
       ob = phc_call_phc(argv[i+1]);        ob = phc_call_phc(argv[i+1]);
       n = phc_getoaSize(ob);        n = phc_getoaSize(ob);
       printf("[\n");        printf("[\n");
       for (i=0; i<n; i++) {        for (i=0; i<n; i++) {
         phc_printObject(stdout,phc_getoa(ob,i));          phc_printObject(stdout,phc_getoa(ob,i));
         if (i != n-1) printf(" ,\n"); else printf(" \n");          if (i != n-1) printf(" ,\n"); else printf(" \n");
       }        }
       printf("]\n");        printf("]\n");
       exit(0);        exit(0);
Line 99  main(int argc, char *argv[]) {
Line 113  main(int argc, char *argv[]) {
   }    }
   while (1) {    while (1) {
     if (message) printf("dim= ");      if (message) printf("dim= ");
         if (fgets(input,INPUTSIZE,stdin) <= 0) break;      if (fgets(input,INPUTSIZE,stdin) <= 0) break;
     sscanf(input,"%d",&dim);      sscanf(input,"%d",&dim);
     sprintf(input,"%d\n",dim);      sprintf(input,"%d\n",dim);
     if (message) printf("Input %d equations please.\n",dim);      if (message) printf("Input %d equations please.\n",dim);
     for (i=0; i<dim; i++) {      for (i=0; i<dim; i++) {
      if (message) {printf("eq[%d] = ",i);  fflush(stdout);}        if (message) {printf("eq[%d] = ",i);  fflush(stdout);}
      do {        do {
        fgets(a,A_SIZE-1, stdin);          fgets(a,A_SIZE-1, stdin);
      } while (strlen(a) == 0);        } while (strlen(a) == 0);
      if (strlen(a) >= A_SIZE-3) {        if (strlen(a) >= A_SIZE-3) {
        fprintf(stderr,"Too big input for the input buffer a.\n"); exit(10);          fprintf(stderr,"Too big input for the input buffer a.\n"); exit(10);
      }        }
      if (strlen(input)+strlen(a) >= INPUTSIZE) {        if (strlen(input)+strlen(a) >= INPUTSIZE) {
        fprintf(stderr,"Too big input for the input buffer input.\n"); exit(10);          fprintf(stderr,"Too big input for the input buffer input.\n"); exit(10);
      }        }
      sprintf(input+strlen(input),"%s\n",a);        sprintf(input+strlen(input),"%s\n",a);
     }      }
     ob = phc_call_phc(input);      ob = phc_call_phc(input);
     if (message) {      if (message) {
Line 130  main(int argc, char *argv[]) {
Line 144  main(int argc, char *argv[]) {
   }    }
 }  }
   
 int phc_scan_for_string(FILE *fp, char str[])  int phc_scan_for_string(FILE *fp, char s[])
      /*       /*
   **  Scans the file fp for a certain string str of length lenstr+1.    **  Scans the file fp for a certain string str of length lenstr+1.
   **  Reading stops when the string has been found, then the variable    **  Reading stops when the string has been found, then the variable
   **  on return equals 1, otherwise 0 is returned.    **  on return equals 1, otherwise 0 is returned.
   */    */
   #define BUF_SIZE 1024
 {  {
     char buf[BUF_SIZE];
     int pt,n,c,i;
     pt = 0;
     n=strlen(s);
     if (n > BUF_SIZE-2) {
           fprintf(stderr,"Too long string for scan_for_string.\n");
           exit(1);
     }
     for (i=0; i<n; i++) {
           buf[i] = fgetc(fp); buf[i+1]='\0';
           if (buf[i] == EOF) return 0;
     }
     if (strcmp(s,buf) == 0) return 0;
     while ((c=fgetc(fp)) != EOF) {
           for (i=0; i<n; i++) {
             buf[i] = buf[i+1];
           }
           buf[n-1]=c; buf[n] = '\0';
           if (strcmp(s,buf) == 0) return 1;
     }
     return 0;
   }
   
   int phc_scan_for_string_old(FILE *fp, char str[])
        /*
     **  Scans the file fp for a certain string str of length lenstr+1.
     **  Reading stops when the string has been found, then the variable
     **  on return equals 1, otherwise 0 is returned.
     */
   {
 #define BUF_SIZE 1024  #define BUF_SIZE 1024
   char buf[BUF_SIZE];    char buf[BUF_SIZE];
   char ch;    char ch;
Line 144  int phc_scan_for_string(FILE *fp, char str[])
Line 189  int phc_scan_for_string(FILE *fp, char str[])
   int lenstr;    int lenstr;
   lenstr = strlen(str);    lenstr = strlen(str);
   if (lenstr >= BUF_SIZE-1) {    if (lenstr >= BUF_SIZE-1) {
         fprintf(stderr,"Too long string in phc_scan_for_string\n");      fprintf(stderr,"Too long string in phc_scan_for_string\n");
         exit(-1);      exit(-1);
   }    }
   index = -1;    index = -1;
   found = 0;    found = 0;
   while (((ch = fgetc(fp))!=EOF) && found == 0)    while (((ch = fgetc(fp))!=EOF) && found == 0)
         {      {
           if (index == -1 && ch == str[0])        if (index == -1 && ch == str[0])
                 {          {
                   index = 0;            index = 0;
                   buf[index] = ch;            buf[index] = ch;
                 }          }
           else        else
                 {          {
                   if (index == lenstr)            if (index == lenstr)
                         {              {
                           compare = 0;                compare = 0;
                           for (i=0; str[i] != '\0'; i++)                for (i=0; str[i] != '\0'; i++)
                                 {                  {
                                   if (buf[i]!=str[i])                    if (buf[i]!=str[i])
                                         {                      {
                                           compare = compare+1;                        compare = compare+1;
                                         }                      }
                                 }                  }
                           if (compare == 0)                if (compare == 0)
                                 {                  {
                                   found = 1;                    found = 1;
                                 }                  }
                           index = -1;                index = -1;
                         }              }
                   else            else
                         if (index > -1 && index < lenstr)              if (index > -1 && index < lenstr)
                           {                {
                                 index = index+1;                  index = index+1;
                                 buf[index] = ch;                  buf[index] = ch;
                           }                }
                 }          }
           if (found == 1) break;        if (found == 1) break;
         }      }
   return found;    return found;
 }  }
 struct phc_object phc_scan_solutions(FILE *fp, int npaths, int dim )  struct phc_object phc_scan_solutions(FILE *fp, int npaths, int dim )
Line 192  struct phc_object phc_scan_solutions(FILE *fp, int npa
Line 237  struct phc_object phc_scan_solutions(FILE *fp, int npa
   **  The tolerance for the residual to a solution is set to 1.0E-12.    **  The tolerance for the residual to a solution is set to 1.0E-12.
   **  Returns solutions.    **  Returns solutions.
   */    */
   #define BUFSIZE 1024
 {  {
   struct phc_object rob,sob;    struct phc_object rob,sob;
   char ch;    char ch;
   int fnd,i,j,nsols;    int fnd,i,j,nsols;
   double res;    double res;
   long double realpart;    double realpart;
   long double imagpart;    double imagpart;
   long double realparts[npaths][dim];    /*  double realparts[npaths][dim];
   long double imagparts[npaths][dim];    double imagparts[npaths][dim]; */
     double *realparts;
     double *imagparts;
     char buf[BUFSIZE];
   nsols = 0;    nsols = 0;
     realparts = (double *)sGC_malloc(sizeof(double)*npaths*dim+1);
     imagparts = (double *)sGC_malloc(sizeof(double)*npaths*dim+1);
   while ((ch = fgetc(fp)) != EOF)    while ((ch = fgetc(fp)) != EOF)
     {      {
       fnd = phc_scan_for_string(fp,"start residual :");        fnd = phc_scan_for_string(fp,"start residual :");
       if (fnd==1)        if (fnd==1)
         {          {
           fscanf(fp,"%E",&res);            fscanf(fp,"%le",&res);
           /* printf(" residual = "); printf("%E\n",res);  */            /* printf(" residual = "); printf("%E\n",res);  */
           if (res < 1.0E-12) {            if (res < 1.0E-12) {
                 nsols = nsols+1;              nsols = nsols+1;
                 if (nsols > npaths) {              if (nsols > npaths) {
                   fprintf(stderr,"Something is wrong in phc_scan_solutions\n");                fprintf(stderr,"Something is wrong in phc_scan_solutions\n");
                   fprintf(stderr,"npaths=%d, nsols=%d \n",npaths,nsols);                fprintf(stderr,"npaths=%d, nsols=%d \n",npaths,nsols);
                   exit(-1);                exit(-1);
                 }              }
           }            }
           fnd = phc_scan_for_string(fp,"the solution for t :");            fnd = phc_scan_for_string(fp,"the solution for t :");
           for (i=0;i<dim;i++)            for (i=0;i<dim;i++)
             {              {
               fnd = phc_scan_for_string(fp,":");                fnd = phc_scan_for_string(fp,":");
               fscanf(fp,"%LE",&realpart);                            fgets(buf,BUFSIZE-1,fp);
               fscanf(fp,"%LE",&imagpart);                            sscanf(buf,"%le %le",&realpart,&imagpart);
               if (res < 1.0E-12)                            if (phc_verbose) fprintf(stderr,"sols in string: %s\n",buf);
                 {                /*fscanf(fp,"%le",&realpart);
                   realparts[nsols-1][i] = realpart;                fscanf(fp,"%le",&imagpart);*/
                   imagparts[nsols-1][i] = imagpart;                if (res < 1.0E-12)
                 }                  {
             }                    /*realparts[nsols-1][i] = realpart;
         }                    imagparts[nsols-1][i] = imagpart;*/
                     realparts[(nsols-1)*npaths+i] = realpart;
                     imagparts[(nsols-1)*npaths+i] = imagpart;
                   }
               }
           }
     }      }
   if(phc_verbose) fprintf(stderr,"  number of solutions = %i\n",nsols);    if(phc_verbose) fprintf(stderr,"  number of solutions = %i\n",nsols);
   rob = phc_newObjectArray(nsols);    rob = phc_newObjectArray(nsols);
Line 238  struct phc_object phc_scan_solutions(FILE *fp, int npa
Line 294  struct phc_object phc_scan_solutions(FILE *fp, int npa
       /* fprintf(stderr,"Solution %i :\n",i+1); */        /* fprintf(stderr,"Solution %i :\n",i+1); */
       sob = phc_newObjectArray(dim);        sob = phc_newObjectArray(dim);
       for (j=0;j<dim;j++)        for (j=0;j<dim;j++)
         {          {
           /*            /*
           printf("%20.14LE",realparts[i][j]); printf("  ");              printf("%20.14le",realparts[i][j]); printf("  ");
           printf("%20.14LE",imagparts[i][j]); printf("\n");              printf("%20.14le",imagparts[i][j]); printf("\n");
           */              */
           phc_putoa(sob,j,phc_complexTo(realparts[i][j],imagparts[i][j]));  
         }              printf("%20.14le",realparts[i*npaths+j]); printf("  ");
               printf("%20.14le",imagparts[i*npaths+j]); printf("\n");
   
             /*phc_putoa(sob,j,phc_complexTo(realparts[i][j],imagparts[i][j]));*/
             phc_putoa(sob,j,phc_complexTo(realparts[i*npaths+j],
                                           imagparts[i*npaths+j]));
           }
       phc_putoa(rob,i,sob);        phc_putoa(rob,i,sob);
     }      }
   return(rob);    return(rob);
Line 260  struct phc_object phc_scan_output_of_phc(char *fname)
Line 322  struct phc_object phc_scan_output_of_phc(char *fname)
   char ch;    char ch;
   int fnd,npaths,dim,i,nsols;    int fnd,npaths,dim,i,nsols;
   otp = fopen(fname,"r");    otp = fopen(fname,"r");
     if (otp == NULL) {
       fprintf(stderr,"The file %d is not found.\n",fname);
       exit(0);
     }
   if (phc_verbose) fprintf(stderr,"Scanning the %s of phc.\n",fname);    if (phc_verbose) fprintf(stderr,"Scanning the %s of phc.\n",fname);
   fnd = phc_scan_for_string(otp,"THE SOLUTIONS :");    fnd = phc_scan_for_string(otp,"THE SOLUTIONS :");
   fscanf(otp,"%i",&npaths);    fscanf(otp,"%i",&npaths);
Line 306  struct phc_object phc_call_phc(char *sys)  /* call phc
Line 372  struct phc_object phc_call_phc(char *sys)  /* call phc
   
   
 struct phc_object phc_newObjectArray(size)  struct phc_object phc_newObjectArray(size)
 int size;       int size;
 {  {
   struct phc_object rob;    struct phc_object rob;
   struct phc_object *op;    struct phc_object *op;
Line 339  void phc_printObject(FILE *fp,struct phc_object ob)
Line 405  void phc_printObject(FILE *fp,struct phc_object ob)
     fprintf(fp,"]");      fprintf(fp,"]");
   }else if (ob.tag == SlongdoubleComplex) {    }else if (ob.tag == SlongdoubleComplex) {
     /* Try your favorite way to print complex numbers. */      /* Try your favorite way to print complex numbers. */
     /*fprintf(fp,"(%20.14LE)+I*(%20.14LE)",ob.lc.longdouble,ob.rc.longdouble);*/      /*fprintf(fp,"(%20.14le)+I*(%20.14le)",ob.lc.longdouble,ob.rc.longdouble);*/
     fprintf(fp,"[%Lf , %Lf]",ob.lc.longdouble,ob.rc.longdouble);      fprintf(fp,"[%lf , %lf]",ob.lc.longdouble,ob.rc.longdouble);
   }else{    }else{
     fprintf(stderr,"Unknown phc_object tag %d",ob.tag);      fprintf(stderr,"Unknown phc_object tag %d",ob.tag);
   }    }
Line 371  char *phc_which(char *s) {
Line 437  char *phc_which(char *s) {
   a = getenv("OpenXM_HOME");    a = getenv("OpenXM_HOME");
   if (a != NULL) {    if (a != NULL) {
     cmd = (char *) sGC_malloc(sizeof(char)*(strlen(s)+strlen(a)      cmd = (char *) sGC_malloc(sizeof(char)*(strlen(s)+strlen(a)
                                             +strlen("/usr/local/bin/")+1));                                              +strlen("/usr/local/bin/")+1));
     if (cmd == NULL) {fprintf(stderr,"No more memory.\n"); exit(1);}      if (cmd == NULL) {fprintf(stderr,"No more memory.\n"); exit(1);}
     strcpy(cmd,a); strcat(cmd,"/bin/"); strcat(cmd,s);      strcpy(cmd,a); strcat(cmd,"/bin/"); strcat(cmd,s);
     if (stat(cmd,&statbuf) == 0) {      if (stat(cmd,&statbuf) == 0) {
Line 379  char *phc_which(char *s) {
Line 445  char *phc_which(char *s) {
     }      }
   }    }
   cmd = (char *) sGC_malloc(sizeof(char)*(strlen(s)    cmd = (char *) sGC_malloc(sizeof(char)*(strlen(s)
                                           +strlen("/usr/local/bin/")+1));                                            +strlen("/usr/local/bin/")+1));
   if (cmd == NULL) {fprintf(stderr,"No more memory.\n"); exit(1);}    if (cmd == NULL) {fprintf(stderr,"No more memory.\n"); exit(1);}
   strcpy(cmd,"/usr/local/bin/"); strcat(cmd,s);    strcpy(cmd,"/usr/local/bin/"); strcat(cmd,s);
   if (stat(cmd,&statbuf) == 0) {    if (stat(cmd,&statbuf) == 0) {
Line 393  char *phc_which(char *s) {
Line 459  char *phc_which(char *s) {
 }  }
   
   
 struct phc_object phc_complexTo(long double r, long double i)  struct phc_object phc_complexTo(double r, double i)
 {  {
   struct phc_object rob;    struct phc_object rob;
   rob.tag = SlongdoubleComplex;    rob.tag = SlongdoubleComplex;
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.4
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>