[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.1 and 1.4

version 1.1, 1999/10/08 02:12:14 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/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 25  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 41  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 55  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 lenstr);  int phc_scan_for_string(FILE *fp, char str[]);
 struct phc_object phc_scan_solutions(FILE *fp, int npaths, int dim );  struct phc_object phc_scan_solutions(FILE *fp, int npaths, int dim );
 struct phc_object phc_scan_output_of_phc(char *fname);  struct phc_object phc_scan_output_of_phc(char *fname);
 struct phc_object phc_call_phc(char *sys);  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;
   for (i=1; i<argc; i++) {    for (i=1; i<argc; i++) {
     if (strcmp(argv[i],"-v") == 0) {      if (strcmp(argv[i],"-v") == 0) {
       phc_verbose = 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 98  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 (scanf("%d",&dim)<0) break;      if (fgets(input,INPUTSIZE,stdin) <= 0) break;
       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 128  main(int argc, char *argv[]) {
Line 144  main(int argc, char *argv[]) {
   }    }
 }  }
   
 int phc_scan_for_string(FILE *fp, char str[], int lenstr)  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[lenstr+1];    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
     char buf[BUF_SIZE];
   char ch;    char ch;
   int index,i,compare,npaths,dim,found;    int index,i,compare,npaths,dim,found;
     int lenstr;
     lenstr = strlen(str);
     if (lenstr >= BUF_SIZE-1) {
       fprintf(stderr,"Too long string in phc_scan_for_string\n");
       exit(-1);
     }
   index = -1;    index = -1;
   found = 0;    found = 0;
   while ((fscanf(fp,"%c",&ch)!=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; i<lenstr+1; 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;
Line 183  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;
   float 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;
   while (fscanf(fp,"%c",&ch)!=EOF)    realparts = (double *)sGC_malloc(sizeof(double)*npaths*dim+1);
     imagparts = (double *)sGC_malloc(sizeof(double)*npaths*dim+1);
     while ((ch = fgetc(fp)) != EOF)
     {      {
       fnd = phc_scan_for_string(fp,"start residual :",15);        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) nsols = nsols+1;            if (res < 1.0E-12) {
           fnd = phc_scan_for_string(fp,"the solution for t :",19);              nsols = nsols+1;
           for (i=0;i<dim;i++)              if (nsols > npaths) {
             {                fprintf(stderr,"Something is wrong in phc_scan_solutions\n");
               fnd = phc_scan_for_string(fp,":",0);                fprintf(stderr,"npaths=%d, nsols=%d \n",npaths,nsols);
               fscanf(fp,"%LE",&realpart);                exit(-1);
               fscanf(fp,"%LE",&imagpart);              }
               if (res < 1.0E-12)            }
                 {            fnd = phc_scan_for_string(fp,"the solution for t :");
                   realparts[nsols-1][i] = realpart;            for (i=0;i<dim;i++)
                   imagparts[nsols-1][i] = imagpart;              {
                 }                fnd = phc_scan_for_string(fp,":");
             }                            fgets(buf,BUFSIZE-1,fp);
         }                            sscanf(buf,"%le %le",&realpart,&imagpart);
                             if (phc_verbose) fprintf(stderr,"sols in string: %s\n",buf);
                 /*fscanf(fp,"%le",&realpart);
                 fscanf(fp,"%le",&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 222  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 244  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 :",14);    fnd = phc_scan_for_string(otp,"THE SOLUTIONS :");
   fscanf(otp,"%i",&npaths);    fscanf(otp,"%i",&npaths);
   if (phc_verbose) fprintf(stderr,"  number of paths traced = %i\n",npaths);    if (phc_verbose) fprintf(stderr,"  number of paths traced = %i\n",npaths);
   fscanf(otp,"%i",&dim);    fscanf(otp,"%i",&dim);
Line 290  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 323  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 355  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 363  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 377  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.1  
changed lines
  Added in v.1.4
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>