Annotation of OpenXM_contrib/PHC/C/call_phc.c, Revision 1.1.1.1
1.1 maekawa 1: /* This is a simple C interface to the black-box solver of phc,
2: ** written by Nobuki Takayama.
3: ** Requirements:
4: ** 1) executable version of phc must be on /tmp;
5: ** 2) user of this program has write permissions to create
6: ** the files "input" and "output" in the directory where
7: ** this program is executed.
8: */
9:
10: #include <stdio.h>
11: #include <sys/stat.h>
12: #include <unistd.h>
13:
14: /* Definition of class identifiers. */
15: #define SstringObject 5
16: #define Sarray 6
17:
18: /* Definition of Object */
19: union cell {
20: int ival;
21: char *str;
22: struct object *op;
23: };
24: struct object{
25: int tag; /* class identifier */
26: union cell lc; /* left cell */
27: union cell rc; /* right cell */
28: };
29:
30: /* Memory allocation function.
31: I recommend not to use malloc and to use gc4.14 for large applications. */
32: #define sGC_malloc(n) malloc(n)
33:
34: /********** macros to use Sarray **************/
35: /* put to Object Array */
36: #define phc_putoa(ob,i,cc) {\
37: if ((ob).tag != Sarray) {fprintf(stderr,"Warning: PUTOA is for an array of objects\n");} else \
38: {if ((0 <= (i)) && ((i) < (ob).lc.ival)) {\
39: (ob.rc.op)[i] = cc;\
40: }else{\
41: fprintf(stderr,"Warning: PUTOA, the size is %d.\n",(ob).lc.ival);\
42: }}}
43: #define phc_getoa(ob,i) ((ob.rc.op)[i])
44: #define phc_getoaSize(ob) ((ob).lc.ival)
45:
46:
47: /* prototypes */
48: struct object phc_newObjectArray(int size);
49: void phc_printObject(FILE *fp,struct object ob);
50: char *phc_generateUniqueFileName(char *s);
51: struct object phc_complexTo(long double r, long double i);
52: struct object phc_longdoubleToStringObject(long double r);
53:
54: int phc_scan_for_string(FILE *fp, char str[], int lenstr);
55: struct object phc_scan_solutions(FILE *fp, int npaths, int dim );
56: struct object phc_scan_output_of_phc(char *fname);
57: struct object phc_call_phc(char *sys);
58:
59: main() {
60: struct object ob;
61: ob = phc_call_phc("2\n x**2 + y**2 - 1;\n x**2 + y**2 - 8*x - 3;\n");
62: printf("-----------------------------------------------------------\n");
63: phc_printObject(stdout,ob);
64: printf("\n");
65: }
66:
67: int phc_scan_for_string(FILE *fp, char str[], int lenstr)
68: /*
69: ** Scans the file fp for a certain string str of length lenstr+1.
70: ** Reading stops when the string has been found, then the variable
71: ** on return equals 1, otherwise 0 is returned.
72: */
73: {
74: char buf[lenstr+1];
75: char ch;
76: int index,i,compare,npaths,dim,found;
77: index = -1;
78: found = 0;
79: while ((fscanf(fp,"%c",&ch)!=EOF) && found == 0)
80: {
81: if (index == -1 && ch == str[0])
82: {
83: index = 0;
84: buf[index] = ch;
85: }
86: else
87: {
88: if (index == lenstr)
89: {
90: compare = 0;
91: for (i=0; i<lenstr+1; i++)
92: {
93: if (buf[i]!=str[i])
94: {
95: compare = compare+1;
96: }
97: }
98: if (compare == 0)
99: {
100: found = 1;
101: }
102: index = -1;
103: }
104: else
105: if (index > -1 && index < lenstr)
106: {
107: index = index+1;
108: buf[index] = ch;
109: }
110: }
111: if (found == 1) break;
112: }
113: return found;
114: }
115: struct object phc_scan_solutions(FILE *fp, int npaths, int dim )
116: /*
117: ** Scans the file for the solutions, from a list of length npaths,
118: ** of complex vectors with dim entries.
119: ** The tolerance for the residual to a solution is set to 1.0E-12.
120: ** Returns the number of solutions.
121: */
122: {
123: struct object rob,sob;
124: char ch;
125: int fnd,i,j,nsols;
126: float res;
127: long double realpart;
128: long double imagpart;
129: long double realparts[npaths][dim];
130: long double imagparts[npaths][dim];
131: nsols = 0;
132: while (fscanf(fp,"%c",&ch)!=EOF)
133: {
134: fnd = phc_scan_for_string(fp,"start residual :",15);
135: if (fnd==1)
136: {
137: fscanf(fp,"%E",&res);
138: /* printf(" residual = "); printf("%E\n",res); */
139: if (res < 1.0E-12) nsols = nsols+1;
140: fnd = phc_scan_for_string(fp,"the solution for t :",19);
141: for (i=0;i<dim;i++)
142: {
143: fnd = phc_scan_for_string(fp,":",0);
144: fscanf(fp,"%LE",&realpart);
145: fscanf(fp,"%LE",&imagpart);
146: if (res < 1.0E-12)
147: {
148: realparts[nsols-1][i] = realpart;
149: imagparts[nsols-1][i] = imagpart;
150: }
151: }
152: }
153: }
154: /* fprintf(stderr," number of solutions = %i\n",nsols); */
155: rob = phc_newObjectArray(nsols);
156: for (i=0;i<nsols;i++)
157: {
158: /* fprintf(stderr,"Solution %i :\n",i+1); */
159: sob = phc_newObjectArray(dim);
160: for (j=0;j<dim;j++)
161: {
162: /*
163: printf("%20.14LE",realparts[i][j]); printf(" ");
164: printf("%20.14LE",imagparts[i][j]); printf("\n");
165: */
166: phc_putoa(sob,j,phc_complexTo(realparts[i][j],imagparts[i][j]));
167: }
168: phc_putoa(rob,i,sob);
169: }
170: return(rob);
171: }
172: struct object phc_scan_output_of_phc(char *fname)
173: /*
174: ** Scans the file "output" of phc in two stages to get
175: ** 1) the number of paths and the dimension;
176: ** 2) the solutions, vectors with residuals < 1.0E-12.
177: */
178: {
179: FILE *otp;
180: char ch;
181: int fnd,npaths,dim,i,nsols;
182: otp = fopen(fname,"r");
183: fprintf(stderr,"Scanning the %s of phc.\n",fname);
184: fnd = phc_scan_for_string(otp,"THE SOLUTIONS :",14);
185: fscanf(otp,"%i",&npaths);
186: fprintf(stderr," number of paths traced = %i\n",npaths);
187: fscanf(otp,"%i",&dim);
188: fprintf(stderr," dimension of solutions = %i\n",dim);
189: return(phc_scan_solutions(otp,npaths,dim));
190: }
191: struct object phc_call_phc(char *sys) /* call phc, system as string */
192: {
193: FILE *inp;
194: char *f,*outf;
195: char cmd[1024];
196: f = phc_generateUniqueFileName("tmp.input");
197: outf = phc_generateUniqueFileName("tmp.output");
198: fprintf(stderr,"Creating file with name %s.\n",f);
199: inp = fopen(f,"w");
200: fprintf(inp,sys);
201: fclose(inp);
202: sprintf(cmd,"/tmp/phc -b %s %s",f,outf);
203: fprintf(stderr,"Calling %s, black-box solver of phc.\n",cmd);
204: system(cmd);
205: fprintf(stderr,"See the file %s for results.\n",outf);
206: return(phc_scan_output_of_phc(outf));
207: }
208:
209:
210: struct object phc_newObjectArray(size)
211: int size;
212: {
213: struct object rob;
214: struct object *op;
215: if (size > 0) {
216: op = (struct object *)sGC_malloc(size*sizeof(struct object));
217: if (op == (struct object *)NULL) {fprintf(stderr,"No memory\n");exit(1);}
218: }else{
219: op = (struct object *)NULL;
220: }
221: rob.tag = Sarray;
222: rob.lc.ival = size;
223: rob.rc.op = op;
224: return(rob);
225: }
226:
227: void phc_printObject(FILE *fp,struct object ob)
228: {
229: int n,i;
230: if (ob.tag == SstringObject) {
231: fprintf(fp,"%s",ob.lc.str);
232: }else if (ob.tag == Sarray) {
233: n = phc_getoaSize(ob);
234: fprintf(fp,"[");
235: for (i=0; i<n; i++) {
236: phc_printObject(fp,phc_getoa(ob,i));
237: if (i <n-1) fprintf(fp,", ");
238: }
239: fprintf(fp,"]");
240: }else{
241: fprintf(stderr,"Unknown object tag %d",ob.tag);
242: }
243: }
244:
245:
246: char *phc_generateUniqueFileName(char *s)
247: {
248: char *t;
249: int i;
250: struct stat statbuf;
251: t = (char *)sGC_malloc(sizeof(char)*strlen(s)+4+2);
252: for (i=0; i<1000; i++) {
253: /* Give up if we failed for 1000 names. */
254: sprintf(t,"%s.%d",s,i);
255: if (stat(t,&statbuf) < 0) {
256: return(t);
257: }
258: }
259: fprintf(stderr,"Could not generate a unique file name.");
260: return(NULL);
261: }
262:
263: struct object phc_longdoubleToStringObject(long double r)
264: {
265: struct object rob;
266: rob.tag = SstringObject;
267: rob.lc.str = (char *) sGC_malloc(sizeof(char)*50); /*??*/
268: sprintf(rob.lc.str,"%20.14LE",r);
269: return(rob);
270: }
271:
272: struct object phc_complexTo(long double r, long double i)
273: {
274: struct object rob;
275: rob = phc_newObjectArray(2);
276: phc_putoa(rob,0,phc_longdoubleToStringObject(r));
277: phc_putoa(rob,1,phc_longdoubleToStringObject(i));
278: return(rob);
279: }
280:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>