1 /****************************************************************
2 Copyright (C) 1997-2001 Lucent Technologies
3 All Rights Reserved
4
5 Permission to use, copy, modify, and distribute this software and
6 its documentation for any purpose and without fee is hereby
7 granted, provided that the above copyright notice appear in all
8 copies and that both that the copyright notice and this
9 permission notice and warranty disclaimer appear in supporting
10 documentation, and that the name of Lucent or any of its entities
11 not be used in advertising or publicity pertaining to
12 distribution of the software without specific, written prior
13 permission.
14
15 LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
16 INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.
17 IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY
18 SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
19 WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
20 IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
21 ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
22 THIS SOFTWARE.
23 ****************************************************************/
24
25 #include "jac2dim.h"
26 #include "opnos.hd"
27
28 #ifdef __cplusplus
29 extern "C" {
30 #endif
31
32 #define Egulp 400
33
34 static int com11, n_com1, ncom_togo, nvar0, nvinc, nvref;
35 static int *vrefnext, *vrefx;
36 static expr_if *if2list_end;
37 static expr_va *varg2list_end;
38
39 #define NFHASH 23
40
41 static derp *last_d;
42 static expr *last_e;
43
44 #define efunc efunc2
45
46 static relo *relolist, *relo2list;
47 static expr_if *iflist, *if2list;
48 static expr_va *varglist, *varg2list;
49 static real one = 1.;
50 static int nderp;
51 #undef nzc
52 static int amax1, imap_len, k_seen, lasta, lasta0, lasta00, lastj,
53 max_var, nocopy, nv0, nv01, nv0b, nv0c, nv1, nzc, nzclim;
54 static int co_first = 1;
55 static int *imap, *zc, *zci;
56
57 static expr *(*holread) ANSI((EdRead*));
58 extern real f_OPCPOW ANSI((expr* A_ASL));
59
60 static ASL_fgh *asl;
61
62 static void
63 #ifdef KR_headers
ed_reset(a)64 ed_reset(a) ASL *a;
65 #else
66 ed_reset(ASL *a)
67 #endif
68 {
69
70 a->i.memLast = a->i.memNext = 0;
71 k_seen = 0;
72 nvref = 0;
73 vrefnext = vrefx = 0;
74 relolist = relo2list = 0;
75 last_d = 0;
76 iflist = if2list = if2list_end = 0;
77 varglist = varg2list = varg2list_end = 0;
78 imap = 0;
79 amax1 = imap_len = lastj = nocopy = 0;
80 com11 = lasta = lasta0 = lasta00 = n_com1 = nderp = 0;
81 last_e = 0;
82 co_first = 1;
83 }
84
85 static void
86 #ifdef KR_headers
fscream(R,name,nargs,kind)87 fscream(R, name, nargs, kind) EdRead *R; char *name; int nargs; char *kind;
88 #else
89 fscream(EdRead *R, const char *name, int nargs, char *kind)
90 #endif
91 {
92 scream(R, ASL_readerr_argerr,
93 "line %ld: attempt to call %s with %d %sargs\n",
94 R->Line, name, nargs, kind);
95 }
96
97 static void
98 #ifdef KR_headers
sorry_CLP(R,what)99 sorry_CLP(R, what) EdRead *R; char *what;
100 #else
101 sorry_CLP(EdRead *R, char *what)
102 #endif
103 {
104 fprintf(Stderr,
105 "Sorry, %s cannot handle %s.\n",
106 progname ? progname : "", what);
107 exit_ASL(R,ASL_readerr_CLP);
108 }
109
110 #ifdef Double_Align
111 #define memadj(x) x
112 #else
113 #define memadj(x) (((x) + (sizeof(long)-1)) & ~(sizeof(long)-1))
114 #endif
115
116 extern efunc f_OPPLTERM, f_OPHOL, f_OPVARVAL, f_OPFUNCALL;
117 extern sfunc f_OPIFSYM;
118
119 static void
120 #ifdef KR_headers
new_derp(a,b,c)121 new_derp(a, b, c) int a; int b; real *c;
122 #else
123 new_derp(int a, int b, real *c)
124 #endif
125 {
126 derp *d;
127 if (a == nv1)
128 return;
129 nderp++;
130 d = (derp *)mem(sizeof(derp));
131 d->next = last_d;
132 last_d = d;
133 d->a.i = a;
134 d->b.i = b;
135 d->c.rp = c;
136 }
137
138 static derp *
139 #ifdef KR_headers
new_relo(e,Dnext,ap)140 new_relo(e, Dnext, ap) expr *e; derp *Dnext; int *ap;
141 #else
142 new_relo(expr *e, derp *Dnext, int *ap)
143 #endif
144 {
145 relo *r;
146 derp *d;
147
148 if (last_d != Dnext) {
149 *ap = e->a;
150 for(d = last_d; d->next != Dnext; d = d->next);
151 d->next = 0;
152 }
153 else {
154 last_d = 0;
155 new_derp(e->a, *ap = lasta++, &one);
156 }
157 if (!last_d)
158 return 0;
159 r = (relo *)mem(sizeof(relo));
160 r->next = relolist;
161 r->next2 = relo2list;
162 relo2list = relolist = r;
163 r->D = r->Dcond = last_d;
164 r->Dnext = Dnext;
165 return r->D;
166 }
167
168 static relo *
169 #ifdef KR_headers
new_relo1(Dnext)170 new_relo1(Dnext) derp *Dnext;
171 #else
172 new_relo1(derp *Dnext)
173 #endif
174 {
175 relo *r;
176
177 r = (relo *)mem(sizeof(relo));
178 r->next = relolist;
179 relolist = r;
180 r->D = 0;
181 r->Dnext = Dnext;
182 return r;
183 }
184
185 static expr *
186 #ifdef KR_headers
new_expr(opcode,L,R,deriv)187 new_expr(opcode, L, R, deriv) int opcode; expr *L; expr *R; int deriv;
188 #else
189 new_expr(int opcode, expr *L, expr *R, int deriv)
190 #endif
191 {
192 expr *rv;
193 extern efunc f_OP1POW, f_OP2POW, f_OPPOW;
194 efunc *o;
195 int L1, R1, i;
196
197 o = r_ops[opcode];
198 if (o == f_OPPOW)
199 if (R->op == f_OPNUM)
200 if (((expr_n *)R)->v == 2.) {
201 o = f_OP2POW;
202 R = 0;
203 }
204 else
205 o = f_OP1POW;
206 else if (L->op == f_OPNUM)
207 o = f_OPCPOW;
208 rv = (expr *)mem(sizeof(expr));
209 rv->op = o;
210 rv->L.e = L;
211 rv->R.e = R;
212 if (deriv) {
213 L1 = L && L->op != f_OPNUM;
214 R1 = R && R->op != f_OPNUM;
215 if (L1 | R1) {
216 rv->a = lasta++;
217 if (L1)
218 new_derp(L->a, rv->a, &rv->dL);
219 if (R1)
220 new_derp(R->a, rv->a, &rv->dR);
221 rv->bak = last_e;
222 last_e = rv;
223 if (R)
224 rv->dL2 = rv->dLR = rv->dR2 = 0;
225 else if (o == f_OP2POW)
226 rv->dL2 = 2;
227 else
228 rv->dL2 = 0;
229 if (L1)
230 if (R1)
231 switch(opcode) {
232 case PLUS: i = Hv_plusLR; break;
233 case MINUS: i = Hv_minusLR; break;
234 case MULT: i = Hv_timesLR; break;
235 default: i = Hv_binaryLR;
236 }
237 else switch(opcode) {
238 case PLUS:
239 case MINUS: i = Hv_plusL; break;
240 case MULT: i = Hv_timesL; break;
241 case UMINUS: i = Hv_negate; break;
242 default: i = Hv_unary;
243 }
244 else
245 switch(opcode) {
246 case PLUS: i = Hv_plusR; break;
247 case MINUS: i = Hv_minusR; break;
248 case MULT: i = Hv_timesR; break;
249 default: i = Hv_binaryR;
250 }
251 rv->dO.i = i;
252 }
253 }
254 return rv;
255 }
256
257 static char op_type[] = {
258 #include "op_type.hd"
259 };
260
261 static expr *
262 #ifdef KR_headers
eread(R,deriv)263 eread(R, deriv) EdRead *R; int deriv;
264 #else
265 eread(EdRead *R, int deriv)
266 #endif
267 {
268 char **sa;
269 int a0, a1, i, i1, j, j1, k, kd, kd2, ks, numargs, symargs;
270 int *at, *nn, *nn0;
271 real *b, **fh, *hes, r, *ra;
272 expr *L, *arg, **args, **args1, **argse, *rv;
273 expr_n *rvn;
274 expr_va *rva;
275 plterm *p;
276 de *d;
277 derp *dsave;
278 efunc *op;
279 expr_if *rvif;
280 expr_f *rvf;
281 func_info *fi;
282 arglist *al;
283 argpair *ap, *da, *sap;
284 char *dig;
285 short sh;
286 unsigned int Ls;
287 long L1;
288 static real dvalue[] = {
289 #include "dvalue.hd"
290 };
291
292 switch(edag_peek(R)) {
293 case 'f':
294 if (xscanf(R, "%d %d", &i, &j) != 2)
295 badline(R);
296 fi = funcs[i];
297 if (fi->nargs >= 0) {
298 if (fi->nargs != j) {
299 bad_nargs:
300 fscream(R, fi->name, j, "");
301 }
302 }
303 else if (-(1+fi->nargs) > j)
304 goto bad_nargs;
305 rvf = (expr_f *)mem(sizeof(expr_f)
306 + (j-1)*sizeof(expr *));
307 rvf->op = f_OPFUNCALL;
308 rvf->fi = fi;
309 rvf->dO.i = Hv_func;
310 args = rvf->args;
311 argse = args + j;
312 k = ks = symargs = numargs = 0;
313 while(args < argse) {
314 arg = *args++ = eread(R, deriv);
315 if ((op = arg->op) == f_OPHOL)
316 symargs++;
317 else if (op == (efunc*)f_OPIFSYM)
318 ks++;
319 else {
320 numargs++;
321 if (op != f_OPNUM)
322 k++;
323 }
324 }
325 symargs += ks;
326 if (symargs && !(fi->ftype & 1))
327 fscream(R, fi->name, symargs, "symbolic ");
328 if (deriv) {
329 kd2 = (kd = k) ? numargs*(numargs+1) >> 1 : 0;
330 rvf->a = lasta++;
331 }
332 else {
333 kd = kd2 = 0;
334 rvf->a = nv1;
335 }
336 Ls = sizeof(arglist)
337 + (k + kd + ks)*sizeof(argpair)
338 + kd*kd*sizeof(real *)
339 + (numargs+kd+kd2)*sizeof(real)
340 + symargs*sizeof(char *)
341 + j*sizeof(int);
342 if (kd)
343 Ls += numargs*sizeof(real);
344 ra = (real *)mem(Ls);
345 dig = 0;
346 if (kd < numargs && kd)
347 dig = (char*)mem(numargs);
348 b = kd ? ra + numargs : ra;
349 hes = b + numargs;
350 al = rvf->al = (arglist *)(hes + kd2);
351 al->n = numargs + symargs;
352 al->nr = numargs;
353 al->ra = ra;
354 if (kd) {
355 al->derivs = b;
356 al->hes = hes;
357 memset(b, 0, (numargs + kd2)*sizeof(real));
358 }
359 else
360 al->derivs = al->hes = 0;
361 al->dig = dig;
362 al->funcinfo = fi->funcinfo;
363 al->AE = asl->i.ae;
364 al->sa = (Const char**)(sa = (char **)(al + 1));
365 ap = rvf->ap = (argpair *)(sa + symargs);
366 rvf->da = da = ap + k;
367 sap = rvf->sap = da + kd;
368 rvf->fh = fh = (real **)(sap + ks);
369 at = al->at = (int *)(fh + kd*kd);
370 symargs = numargs = i = 0;
371 nn = nn0 = (int *)b;
372 for(args = rvf->args; args < argse; at++) {
373 arg = *args++;
374 if ((op = arg->op) == f_OPHOL) {
375 *at = --symargs;
376 *sa++ = ((expr_h *)arg)->sym;
377 }
378 else if (op == (efunc*)f_OPIFSYM) {
379 *at = --symargs;
380 sap->e = arg;
381 (sap++)->u.s = sa++;
382 }
383 else {
384 *at = numargs++;
385 j = 1;
386 if (op == f_OPNUM)
387 *ra = ((expr_n *)arg)->v;
388 else {
389 ap->e = arg;
390 (ap++)->u.v = ra;
391 if (kd) {
392 j = 0;
393 new_derp(arg->a,
394 rvf->a, b);
395 *b = 0;
396 da->e = arg;
397 (da++)->u.v = b;
398 *nn++ = i;
399 }
400 }
401 if (dig)
402 *dig++ = j;
403 b++;
404 ra++;
405 i++;
406 }
407 }
408 rvf->ape = ap;
409 rvf->sape = sap;
410 rvf->dae = da;
411 if (kd) {
412 rvf->bak = last_e;
413 last_e = (expr *)rvf;
414 rvf->dO.i = Hv_func;
415 for(i1 = 0; i1 < kd; i1++) {
416 i = nn0[i1];
417 for(j1 = 0; j1 < kd; j1++) {
418 j = nn0[j1];
419 *fh++ = &hes[i >= j ? (i*(i+1)>>1)+j
420 : (j*(j+1)>>1)+i];
421 }
422 }
423 }
424 return (expr *)rvf;
425
426 case 'h':
427 return holread(R);
428
429 case 's':
430 if (xscanf(R, "%hd", &sh) != 1)
431 badline(R);
432 r = sh;
433 goto have_r;
434
435 case 'l':
436 if (xscanf(R, "%ld", &L1) != 1)
437 badline(R);
438 r = L1;
439 goto have_r;
440
441 case 'n':
442 if (xscanf(R, "%lf", &r) != 1)
443 badline(R);
444 have_r:
445 rvn = (expr_n *)mem(sizeof(expr_n));
446 rvn->op = f_OPNUM_ASL;
447 rvn->v = r;
448 return (expr *)rvn;
449
450 case 'o':
451 break;
452
453 case 'v':
454 if (xscanf(R,"%d",&k) != 1 || k < 0)
455 badline(R);
456 if (k >= nvar0)
457 k += nvinc;
458 if (k > max_var)
459 badline(R);
460 if (k < nv01 && deriv && !zc[k]++)
461 zci[nzc++] = k;
462 return (expr *)(var_e + k);
463
464 default:
465 badline(R);
466 }
467
468 if (xscanf(R, "%d", &k) != 1 || k < 0 || k >= sizeof(op_type))
469 badline(R);
470 switch(op_type[k]) {
471
472 case 1: /* unary */
473 rv = new_expr(k, eread(R, deriv), 0, deriv);
474 rv->dL = dvalue[k]; /* for UMINUS, FLOOR, CEIL */
475 return rv;
476
477 case 2: /* binary */
478 if (dvalue[k] == 11)
479 deriv = 0;
480 L = eread(R, deriv);
481 rv = new_expr(k, L, eread(R, deriv), deriv);
482 rv->dL = 1.;
483 rv->dR = dvalue[k]; /* for PLUS, MINUS, REM */
484 return rv;
485
486 case 3: /* vararg (min, max) */
487 i = -1;
488 xscanf(R, "%d", &i);
489 if (i <= 0)
490 badline(R);
491 rva = (expr_va *)mem(sizeof(expr_va));
492 rva->op = r_ops[k];
493 rva->L.d = d =
494 (de *)mem(i*sizeof(de) + sizeof(expr *));
495 rva->next = varglist;
496 varglist = varg2list = rva;
497 if (!last_d) {
498 new_derp(lasta, lasta, &one);
499 lasta++;
500 }
501 rva->d0 = dsave = last_d;
502 rva->bak = last_e;
503 a0 = a1 = lasta;
504 for(j = 0; i > 0; i--, d++) {
505 last_d = dsave;
506 last_e = 0;
507 d->e = L = eread(R, deriv);
508 d->ee = last_e;
509 if (L->op == f_OPNUM || L->a == nv1) {
510 d->d = dsave;
511 d->dv.i = nv1;
512 }
513 else {
514 if (deriv)
515 d->d = new_relo(L, dsave,
516 &d->dv.i);
517 if (a1 < lasta)
518 a1 = lasta;
519 lasta = a0;
520 }
521 }
522 d->e = 0; /* sentinnel expr * */
523 rva->a = lasta = a1;
524 last_d = dsave;
525 if (deriv) {
526 new_derp(0, lasta++, &one);
527 /* f_MINLIST or f_MAXLIST will replace the 0 */
528 rva->R.D = last_d;
529 nocopy = 1;
530 last_e = (expr *)rva;
531 rva->dO.i = Hv_vararg;
532 }
533 else {
534 rva->R.D = 0;
535 last_e = rva->bak;
536 }
537 return (expr *)rva;
538
539 case 4: /* piece-wise linear */
540 i = -1;
541 xscanf(R, "%d", &i);
542 if (i <= 1)
543 badline(R);
544 j = 2*i - 1;
545 p = (plterm *)mem(sizeof(plterm) + (j-1)*sizeof(real));
546 p->n = i;
547 b = p->bs;
548 do {
549 switch(edag_peek(R)) {
550 case 's':
551 if (xscanf(R,"%hd",&sh) != 1)
552 badline(R);
553 r = sh;
554 break;
555 case 'l':
556 if (xscanf(R,"%ld",&L1) != 1)
557 badline(R);
558 r = L1;
559 break;
560 case 'n':
561 if (xscanf(R,"%lf",&r) == 1)
562 break;
563 default:
564 badline(R);
565 }
566 *b++ = r;
567 }
568 while(--j > 0);
569 if (edag_peek(R) != 'v'
570 || xscanf(R, "%d", &k) != 1
571 || k < 0 || k >= max_var)
572 badline(R);
573 if (k >= nvar0)
574 k += nvinc;
575 rv = (expr *)mem(sizeof(expr));
576 rv->op = f_OPPLTERM;
577 rv->L.p = p;
578 rv->R.e = (expr *)(var_e + k);
579 if (deriv) {
580 new_derp(k, rv->a = lasta++, &rv->dL);
581 rv->bak = last_e;
582 last_e = rv;
583 rv->dO.i = Hv_plterm;
584 }
585 return rv;
586
587 case 5: /* if */
588 rvif = (expr_if *)mem(sizeof(expr_if));
589 rvif->op = r_ops[k];
590 rvif->next = iflist;
591 iflist = if2list = rvif;
592 if (!last_d) {
593 new_derp(lasta, lasta, &one);
594 lasta++;
595 }
596 rvif->d0 = dsave = last_d;
597 rvif->bak = last_e;
598 rvif->e = eread(R, 0);
599 last_e = 0;
600 a0 = lasta;
601 rvif->T = L = eread(R, deriv);
602 rvif->Te = last_e;
603 j = 0;
604 if (L->op == f_OPNUM) {
605 rvif->dT = dsave;
606 rvif->Tv.i = nv1;
607 }
608 else if (j = deriv)
609 rvif->dT = new_relo(L, dsave, &rvif->Tv.i);
610 a1 = lasta;
611 lasta = a0;
612 last_d = dsave;
613 last_e = 0;
614 rvif->F = L = eread(R, deriv);
615 rvif->Fe = last_e;
616 if (L->op == f_OPNUM) {
617 rvif->dF = dsave;
618 rvif->Fv.i = nv1;
619 }
620 else if (j = deriv)
621 rvif->dF = new_relo(L, dsave, &rvif->Fv.i);
622 if (lasta < a1)
623 lasta = a1;
624 last_d = dsave;
625 if (j) {
626 new_derp(0, rvif->a = lasta++, &one);
627 rvif->D = last_d;
628 nocopy = 1;
629 last_e = (expr *)rvif;
630 rvif->dO.i = Hv_if;
631 }
632 else {
633 rvif->a = nv1;
634 rvif->D = 0;
635 last_e = rvif->bak;
636 }
637 return (expr *)rvif;
638
639 case 11: /* OPCOUNT */
640 deriv = 0;
641 /* no break */
642 case 6: /* sumlist */
643 i = 0;
644 xscanf(R, "%d", &i);
645 if (i <= 2 && (op_type[k] == 6 || i < 1))
646 badline(R);
647 rv = (expr *)mem(sizeof(expr) - sizeof(real)
648 + (deriv ? i+i+1 : i)*sizeof(expr *));
649 rv->op = r_ops[k];
650 rv->a = deriv ? lasta++ : nv1;
651 rv->L.ep = args = (expr **)&rv->dR;
652 if (deriv) {
653 j = 0;
654 do {
655 *args++ = L = eread(R, deriv);
656 if (L->op != f_OPNUM) {
657 new_derp(L->a, rv->a, &one);
658 j++;
659 }
660 }
661 while(--i > 0);
662 if (j) {
663 rv->bak = last_e;
664 last_e = rv;
665 rv->dO.i = Hv_sumlist;
666 argse = args;
667 args1 = rv->L.ep;
668 while(args1 < args) {
669 arg = *args1++;
670 if (arg->op != f_OPNUM)
671 *argse++ = arg;
672 }
673 *argse = 0;
674 }
675 }
676 else do
677 *args++ = eread(R, deriv);
678 while(--i > 0);
679 rv->R.ep = args;
680 return rv;
681 }
682 badline(R);
683 return 0;
684 }
685
686 static list *
687 #ifdef KR_headers
new_list(nxt)688 new_list(nxt) list *nxt;
689 #else
690 new_list(list *nxt)
691 #endif
692 {
693 list *rv = (list *)mem(sizeof(list));
694 rv->next = nxt;
695 return rv;
696 }
697
698 static list *
crefs(VOID)699 crefs(VOID)
700 {
701 int i;
702 list *rv = 0;
703
704 while(nzc > 0) {
705 if ((i = zci[--nzc]) >= nv0) {
706 rv = new_list(rv);
707 rv->item.i = i;
708 }
709 zc[i] = 0;
710 }
711 return rv;
712 }
713
714 static funnel *
715 #ifdef KR_headers
funnelfix(f)716 funnelfix(f) funnel *f;
717 #else
718 funnelfix(funnel *f)
719 #endif
720 {
721 cexp *ce;
722 funnel *fnext, *fprev;
723
724 for(fprev = 0; f; f = fnext) {
725 fnext = f->next;
726 f->next = fprev;
727 fprev = f;
728 ce = f->ce;
729 ce->z.i = ce->d->b.i;
730 }
731 return fprev;
732 }
733
734 static derp *
735 #ifdef KR_headers
derpadjust(d0,a,dnext)736 derpadjust(d0, a, dnext) derp *d0; int a; derp *dnext;
737 #else
738 derpadjust(derp *d0, int a, derp *dnext)
739 #endif
740 {
741 derp *d, *d1;
742 int *r, *re;
743 relo *rl;
744 expr_if *il, *ile;
745 expr_va *vl, *vle;
746 de *de1;
747
748 if (!(d = d0))
749 return dnext;
750 r = imap + lasta0;
751 re = imap + lasta;
752 while(r < re)
753 *r++ = a++;
754 if (amax < a)
755 amax = a;
756 r = imap;
757 for(;; d = d1) {
758 d->a.i = r[d->a.i];
759 d->b.i = r[d->b.i];
760 if (!(d1 = d->next))
761 break;
762 }
763 d->next = dnext;
764 if (rl = relo2list) {
765 relo2list = 0;
766 do {
767 d = rl->Dcond;
768 do {
769 d->a.i = r[d->a.i];
770 d->b.i = r[d->b.i];
771 }
772 while(d = d->next);
773 }
774 while(rl = rl->next2);
775 }
776 if (if2list != if2list_end) {
777 ile = if2list_end;
778 if2list_end = il = if2list;
779 do {
780 il->Tv.i = r[il->Tv.i];
781 il->Fv.i = r[il->Fv.i];
782 }
783 while((il = il->next) != ile);
784 }
785 if (varg2list != varg2list_end) {
786 vle = varg2list_end;
787 varg2list_end = vl = varg2list;
788 do {
789 for(de1 = vl->L.d; de1->e; de1++)
790 de1->dv.i = r[de1->dv.i];
791 }
792 while((vl = vl->next) != vle);
793 }
794 return d0;
795 }
796
797 static derp *
798 #ifdef KR_headers
derpcopy(ce,dnext)799 derpcopy(ce, dnext) cexp *ce; derp *dnext;
800 #else
801 derpcopy(cexp *ce, derp *dnext)
802 #endif
803 {
804 derp *d, *dprev;
805 int *map;
806 derp d00;
807
808 if (!(d = ce->d))
809 return dnext;
810 map = imap;
811 for(dprev = &d00; d; d = d->next) {
812 new_derp(map[d->a.i], map[d->b.i], d->c.rp);
813 dprev = dprev->next = last_d;
814 }
815 dprev->next = dnext;
816 return d00.next;
817 }
818
819 static void
imap_alloc(VOID)820 imap_alloc(VOID)
821 {
822 int i, *r, *re;
823
824 if (imap) {
825 imap_len += lasta;
826 imap = (int *)Realloc(imap, imap_len*Sizeof(int));
827 return;
828 }
829 imap_len = amax1 > lasta ? amax1 : lasta;
830 imap_len += 100;
831 r = imap = (int *)Malloc(imap_len*Sizeof(int));
832 for(i = 0, re = r + nv1+1; r < re;)
833 *r++ = i++;
834 }
835
836 static int
837 #ifdef KR_headers
compar(a,b,v)838 compar(a, b, v) char *a, *b, *v;
839 #else
840 compar(const void *a, const void *b, void *v)
841 #endif
842 {
843 Not_Used(v);
844 return *(int*)a - *(int*)b;
845 }
846
847 static void
848 #ifdef KR_headers
comsubs(alen,d)849 comsubs(alen, d) int alen; cde *d;
850 #else
851 comsubs(int alen, cde *d)
852 #endif
853 {
854 list *L;
855 int a, i, j, k;
856 int *r, *re;
857 cexp *ce;
858 derp *D, *dnext;
859 relo *R;
860
861 D = last_d;
862 a = lasta00;
863 dnext = 0;
864 R = 0;
865 for(i = j = 0; i < nzc; i++)
866 if ((k = zci[i]) >= nv0)
867 zci[j++] = k;
868 else
869 zc[k] = 0;
870 if (nzc = j) {
871 for(i = 0; i < nzc; i++)
872 for(L = cexps[zci[i]-nv0].cref; L; L = L->next)
873 if (!zc[L->item.i]++)
874 zci[nzc++] = L->item.i;
875 if (nzc > 1)
876 if (nzc < nzclim)
877 qsortv(zci, nzc, sizeof(int), compar, NULL);
878 else for(i = nv0, j = 0; i < max_var; i++)
879 if (zc[i])
880 zci[j++] = i;
881 if (nzc > 0) {
882 R = new_relo1(dnext);
883 i = 0;
884 do {
885 j = zci[i];
886 zc[j] = 0;
887 ce = &cexps[j - nv0];
888 if (ce->funneled)
889 imap[var_e[j].a] = a++;
890 else {
891 r = imap + ce->z.i;
892 re = r + ce->zlen;
893 while(r < re)
894 *r++ = a++;
895 }
896 dnext = R->D = derpcopy(ce, R->D);
897 }
898 while(++i < nzc);
899 nzc = 0;
900 }
901 }
902 if (D || R) {
903 if (!R)
904 R = new_relo1(dnext);
905 D = R->D = derpadjust(D, a, R->D);
906 if (d->e->op != f_OPVARVAL)
907 d->e->a = imap[d->e->a];
908 }
909 d->d = D;
910 a += alen;
911 d->zaplen = (a > lasta00 ? a - nv1 : 0)*sizeof(real);
912 if (amax < a)
913 amax = a;
914 }
915
916 static void
917 #ifdef KR_headers
co_read(R,d,wd)918 co_read(R, d, wd) EdRead *R; cde *d; int wd;
919 #else
920 co_read(EdRead *R, cde *d, int wd)
921 #endif
922 {
923 int alen;
924
925 d->com11 = com11;
926 d->n_com1 = n_com1;
927 com11 += n_com1;
928 n_com1 = 0;
929
930 if (amax1 < lasta)
931 amax1 = lasta;
932 if (co_first) {
933 co_first = 0;
934 if (imap_len < lasta)
935 imap_alloc();
936 f_b = funnelfix(f_b);
937 f_c = funnelfix(f_c);
938 f_o = funnelfix(f_o);
939 }
940 if (!lastj) {
941 lasta = lasta0;
942 last_d = 0;
943 }
944 lastj = 0;
945 last_e = 0;
946 d->e = eread(R, wd);
947 d->ee = last_e;
948 alen = lasta - lasta0;
949 if (imap_len < lasta)
950 imap_alloc();
951 comsubs(alen, d);
952 }
953
954 static linpart *
955 #ifdef KR_headers
linpt_read(R,nlin)956 linpt_read(R, nlin) EdRead *R; int nlin;
957 #else
958 linpt_read(EdRead *R, int nlin)
959 #endif
960 {
961 linpart *L, *rv;
962
963 if (nlin <= 0)
964 return 0;
965 L = rv = (linpart *)mem(nlin*sizeof(linpart));
966 do {
967 if (xscanf(R, "%d %lf", &L->v.i, &L->fac) != 2)
968 badline(R);
969 L++;
970 }
971 while(--nlin > 0);
972 return rv;
973 }
974
975 static int
976 #ifdef KR_headers
funnelkind(ce,ip)977 funnelkind(ce, ip) cexp *ce; int *ip;
978 #else
979 funnelkind(cexp *ce, int *ip)
980 #endif
981 {
982 int i, j, k, nzc0, rv;
983 int *vr, *vre;
984
985 ce->vref = 0;
986 if (!(nzc0 = nzc))
987 return 0;
988 for(i = k = rv = 0; i < nzc; i++)
989 if ((j = zci[i]) < nv0) {
990 if (k >= maxfwd)
991 goto done;
992 vrefx[k++] = j;
993 }
994 else {
995 if (!(vr = cexps[j-nv0].vref))
996 goto done;
997 vre = vr + *vr;
998 while(++vr <= vre)
999 if (!zc[*vr]++)
1000 zci[nzc++] = *vr;
1001 }
1002 if (k >= nvref) {
1003 nvref = (maxfwd + 1)*(ncom_togo < vrefGulp
1004 ? ncom_togo : vrefGulp);
1005 vrefnext = (int *)M1alloc(nvref*Sizeof(int));
1006 }
1007 vr = ce->vref = vrefnext;
1008 *vr = k;
1009 vrefnext += i = k + 1;
1010 nvref -= i;
1011 for(i = 0; i < k; i++)
1012 *++vr = vrefx[i];
1013 if (nderp > 3*k && !nocopy) {
1014 *ip = k;
1015 return 2;
1016 }
1017 else {
1018 done:
1019 if (nocopy || nderp > 3*nzc0)
1020 rv = 1;
1021 }
1022 while(nzc > nzc0)
1023 zc[zci[--nzc]] = 0;
1024 return rv;
1025 }
1026
1027 static void
1028 #ifdef KR_headers
cexp_read(R,k,nlin)1029 cexp_read(R, k, nlin) EdRead *R; int k; int nlin;
1030 #else
1031 cexp_read(EdRead *R, int k, int nlin)
1032 #endif
1033 {
1034 int a, fk, i, j, la0, na;
1035 funnel *f, **fp;
1036 linpart *L, *Le;
1037 expr *e;
1038 cplist *cl, *cl0;
1039 cexp *ce;
1040
1041 ce = cexps + k - nv0;
1042 L = ce->L = linpt_read(R, ce->nlin = nlin);
1043 nocopy = 0;
1044 last_d = 0;
1045 last_e = 0;
1046 ce->z.i = la0 = lasta;
1047 nderps += nderp;
1048 nderp = 0;
1049 e = ce->e = eread(R, 1);
1050 if (la0 == lasta) {
1051 a = lasta++;
1052 if (e->op != f_OPNUM)
1053 new_derp(e->a, a, &edagread_one);
1054 }
1055 else
1056 a = e->a;
1057 var_e[k].a = a;
1058 ce->zlen = lasta - la0;
1059 for(Le = L + nlin; L < Le; L++) {
1060 new_derp(i = L->v.i, a, &L->fac);
1061 if (!zc[i]++)
1062 zci[nzc++] = i;
1063 }
1064 if (fk = funnelkind(ce, &i)) {
1065 /* arrange to funnel */
1066 fp = k < nv0b ? &f_b : k < nv0c ? &f_c : &f_o;
1067 ce->funneled = f = (funnel *)mem(sizeof(funnel));
1068 f->next = *fp;
1069 *fp = f;
1070 f->ce = ce;
1071 if (imap_len < lasta)
1072 imap_alloc();
1073 if (fk == 1) {
1074 f->fulld = last_d;
1075 a = lasta00;
1076 for(i = nzc; --i >= 0; )
1077 if ((j = zci[i]) >= nv0)
1078 imap[var_e[j].a] = a++;
1079 if ((na = ce->zlen) || a > lasta00)
1080 na += a - nv1;
1081 f->fcde.zaplen = na*sizeof(real);
1082 i = nzc;
1083 derpadjust(last_d, a, 0);
1084 }
1085 else {
1086 f->fulld = 0;
1087 f->fcde.e = e;
1088 comsubs(ce->zlen, &f->fcde);
1089 memcpy(zci, vrefx, i*sizeof(int));
1090 }
1091 last_d = 0;
1092 cl0 = 0;
1093 while(i > 0)
1094 if ((a = var_e[zci[--i]].a) != nv1) {
1095 new_derp(a, lasta0, 0);
1096 cl = (cplist *)mem(sizeof(cplist));
1097 cl->next = cl0;
1098 cl0 = cl;
1099 cl->ca.i = imap[last_d->a.i];
1100 cl->cfa = last_d->c.rp =
1101 (real *)mem(sizeof(real));
1102 }
1103 f->cl = cl0;
1104 var_e[k].a = lasta0++;
1105 lasta = lasta0;
1106 }
1107 lasta0 = lasta;
1108 ce->d = last_d;
1109 ce->ee = last_e;
1110 ce->cref = crefs();
1111 --ncom_togo;
1112 }
1113
1114 static void
1115 #ifdef KR_headers
cexp1_read(R,j,k,nlin)1116 cexp1_read(R, j, k, nlin) EdRead *R; int j; int k; int nlin;
1117 #else
1118 cexp1_read(EdRead *R, int j, int k, int nlin)
1119 #endif
1120 {
1121 linpart *L, *Le;
1122 cexp1 *ce = cexps1 + (k - nv01);
1123 expr *e;
1124 expr_v *v;
1125 int la0;
1126
1127 n_com1++;
1128 L = ce->L = linpt_read(R, ce->nlin = nlin);
1129
1130 if (!lastj) {
1131 last_d = 0;
1132 if (amax1 < lasta)
1133 amax1 = lasta;
1134 lasta = lasta0;
1135 lastj = j;
1136 }
1137 last_e = 0;
1138 la0 = lasta;
1139 e = ce->e = eread(R, 1);
1140 ce->ee = last_e;
1141 v = var_e + k;
1142 if (la0 == lasta) {
1143 j = lasta++;
1144 if (e->op != f_OPNUM)
1145 new_derp(e->a, j, &edagread_one);
1146 }
1147 else
1148 j = e->a;
1149 v->a = j;
1150 v->dO.r = 0;
1151 for(Le = L + nlin; L < Le; L++)
1152 new_derp(L->v.i, j, &L->fac);
1153 }
1154
1155 static expr *
1156 #ifdef KR_headers
hvadjust(e)1157 hvadjust(e) expr *e;
1158 #else
1159 hvadjust(expr *e)
1160 #endif
1161 {
1162 expr *e0;
1163
1164 for(e0 = 0; e; e = e->bak) {
1165 e->fwd = e0;
1166 e0 = e;
1167 e->a = e->dO.i;
1168 }
1169 return e0;
1170 }
1171
1172 static void
1173 #ifdef KR_headers
co_adjust(d,n)1174 co_adjust(d, n) cde *d; int n;
1175 #else
1176 co_adjust(cde *d, int n)
1177 #endif
1178 {
1179 cde *de1;
1180
1181 for(de1 = d + n; d < de1; d++)
1182 d->ef = hvadjust(d->ee);
1183 }
1184
1185 static void
1186 #ifdef KR_headers
ifadjust(e)1187 ifadjust(e) expr_if *e;
1188 #else
1189 ifadjust(expr_if *e)
1190 #endif
1191 {
1192 for(; e; e = e->next) {
1193 e->Tv.rp = &adjoints[e->Tv.i];
1194 e->Fv.rp = &adjoints[e->Fv.i];
1195 e->Tf = hvadjust(e->Te);
1196 e->Ff = hvadjust(e->Fe);
1197 }
1198 }
1199
1200 static void
1201 #ifdef KR_headers
vargadjust(e)1202 vargadjust(e) expr_va *e;
1203 #else
1204 vargadjust(expr_va *e)
1205 #endif
1206 {
1207 de *d;
1208
1209 for(; e; e = e->next) {
1210 for(d = e->L.d; d->e; d++) {
1211 d->dv.rp = &adjoints[d->dv.i];
1212 d->ef = hvadjust(d->ee);
1213 }
1214 }
1215 }
1216
1217 static void
1218 #ifdef KR_headers
funneladj1(f)1219 funneladj1(f) funnel *f;
1220 #else
1221 funneladj1(funnel *f)
1222 #endif
1223 {
1224 real *a = adjoints;
1225 derp *d;
1226 cplist *cl;
1227
1228 for(a = adjoints; f; f = f->next) {
1229 if (d = f->fulld) {
1230 f->fcde.d = d;
1231 do {
1232 d->a.rp = &a[d->a.i];
1233 d->b.rp = &a[d->b.i];
1234 }
1235 while(d = d->next);
1236 }
1237 for(cl = f->cl; cl; cl = cl->next)
1238 cl->ca.rp = &a[cl->ca.i];
1239 }
1240 }
1241
1242 static void
funneladjust(VOID)1243 funneladjust(VOID)
1244 {
1245 cexp *c, *ce;
1246 linpart *L, *Le;
1247 c = cexps;
1248 for(ce = c + ncom0; c < ce; c++) {
1249 if (L = c->L)
1250 for(Le = L + c->nlin; L < Le; L++)
1251 L->v.vp = (Char*)&var_e[L->v.i];
1252 c->ef = hvadjust(c->ee);
1253 }
1254
1255 funneladj1(f_b);
1256 funneladj1(f_c);
1257 funneladj1(f_o);
1258 }
1259
1260 static void
com1adjust(VOID)1261 com1adjust(VOID)
1262 {
1263 cexp1 *c, *ce;
1264 linpart *L, *Le;
1265
1266 for(c = cexps1, ce = c + ncom1; c < ce; c++) {
1267 for(L = c->L, Le = L + c->nlin; L < Le; L++)
1268 L->v.vp = (Char*)&var_e[L->v.i];
1269 c->ef = hvadjust(c->ee);
1270 }
1271 }
1272
1273 static void
goff_comp(VOID)1274 goff_comp(VOID)
1275 {
1276 int *ka = A_colstarts + 1;
1277 cgrad **cgx, **cgxe;
1278 cgrad *cg;
1279
1280 cgx = Cgrad;
1281 cgxe = cgx + asl->i.n_con0;
1282 while(cgx < cgxe)
1283 for(cg = *cgx++; cg; cg = cg->next)
1284 cg->goff = ka[cg->varno]++;
1285 }
1286
1287 static void
colstart_inc(VOID)1288 colstart_inc(VOID)
1289 {
1290 int *ka, *kae;
1291 ka = A_colstarts;
1292 kae = ka + asl->i.n_var0;
1293 while(ka <= kae)
1294 ++*ka++;
1295 }
1296
1297 static void
zerograd_chk(VOID)1298 zerograd_chk(VOID)
1299 {
1300 int j, n, nv, *z, **zg;
1301 ograd *og, **ogp, **ogpe;
1302
1303 if (!(nv = asl->i.nlvog))
1304 nv = nv0;
1305 zerograds = 0;
1306 ogp = Ograd;
1307 ogpe = ogp + (j = n_obj);
1308 while(ogp < ogpe) {
1309 og = *ogp++;
1310 n = 0;
1311 while(og) {
1312 j += og->varno - n;
1313 n = og->varno + 1;
1314 if (n >= nv)
1315 break;
1316 og = og->next;
1317 }
1318 if (n < nv)
1319 j += nv - n;
1320 }
1321 if (j == n_obj)
1322 return;
1323 zerograds = zg = (int **)mem(n_obj*sizeof(int*)+j*sizeof(int));
1324 z = (int*)(zg + n_obj);
1325 ogp = Ograd;
1326 while(ogp < ogpe) {
1327 *zg++ = z;
1328 og = *ogp++;
1329 n = 0;
1330 while(og) {
1331 while(n < og->varno)
1332 *z++ = n++;
1333 og = og->next;
1334 if (++n >= nv)
1335 break;
1336 }
1337 while(n < nv)
1338 *z++ = n++;
1339 *z++ = -1;
1340 }
1341 }
1342
1343 static void
adjust_compl_rhs(VOID)1344 adjust_compl_rhs(VOID)
1345 {
1346 cde *C;
1347 expr *e;
1348 int *Cvar, i, j, nc, stride;
1349 real *L, *U, t;
1350
1351 L = LUrhs;
1352 if (U = Urhsx)
1353 stride = 1;
1354 else {
1355 U = L + 1;
1356 stride = 2;
1357 }
1358 C = con_de;
1359 Cvar = cvar;
1360 nc = n_con;
1361 for(i = nlc; i < nc; i++)
1362 if (Cvar[i] && (e = C[i].e) && e->op == f_OPNUM
1363 && (t = ((expr_n*)e)->v) != 0.) {
1364 ((expr_n*)e)->v = 0.;
1365 if (L[j = stride*i] > negInfinity)
1366 L[j] -= t;
1367 if (U[j] < Infinity)
1368 U[j] -= t;
1369 }
1370 }
1371
1372 static void
1373 #ifdef KR_headers
adjust(flags)1374 adjust(flags) int flags;
1375 #else
1376 adjust(int flags)
1377 #endif
1378 {
1379 derp *d, **dp;
1380 real *a = adjoints;
1381 relo *r;
1382
1383 for(r = relolist; r; r = r->next) {
1384 dp = &r->D;
1385 while(d = *dp) {
1386 d->a.rp = &a[d->a.i];
1387 d->b.rp = &a[d->b.i];
1388 dp = &d->next;
1389 }
1390 *dp = r->Dnext;
1391 }
1392 ifadjust(iflist);
1393 vargadjust(varglist);
1394 if (ncom0)
1395 funneladjust();
1396 com1adjust();
1397 co_adjust(con_de, n_con);
1398 co_adjust(obj_de, n_obj);
1399 if (n_obj)
1400 zerograd_chk();
1401 if (k_seen)
1402 if (!A_vals)
1403 goff_comp();
1404 else if (Fortran)
1405 colstart_inc();
1406 if (n_cc > nlcc && nlc < n_con
1407 && !(flags & ASL_no_linear_cc_rhs_adjust))
1408 adjust_compl_rhs();
1409 }
1410
1411 static void
1412 #ifdef KR_headers
br_read(R,nc,nc1,Lp,U,Cvar,nv)1413 br_read(R, nc, nc1, Lp, U, Cvar, nv)
1414 EdRead *R; real **Lp, *U; int nc, nc1, *Cvar, nv;
1415 #else
1416 br_read(EdRead *R, int nc, int nc1, real **Lp, real *U, int *Cvar, int nv)
1417 #endif
1418 {
1419 int i, inc, j, k;
1420 real *L;
1421 ASL *asl = R->asl;
1422
1423 if (!(L = *Lp)) {
1424 if (nc1 < nc)
1425 nc1 = nc;
1426 L = *Lp = (real *)M1alloc(2*sizeof(real)*nc1);
1427 }
1428 if (U)
1429 inc = 1;
1430 else {
1431 U = L + 1;
1432 inc = 2;
1433 }
1434 xscanf(R, ""); /* purge line */
1435 for(i = 0; i < nc; i++, L += inc, U += inc) {
1436 switch(edag_peek(R) - '0') {
1437 case 0:
1438 if (xscanf(R,"%lf %lf",L,U)!= 2)
1439 badline(R);
1440 break;
1441 case 1:
1442 if (xscanf(R, "%lf", U) != 1)
1443 badline(R);
1444 *L = negInfinity;
1445 break;
1446 case 2:
1447 if (xscanf(R, "%lf", L) != 1)
1448 badline(R);
1449 *U = Infinity;
1450 break;
1451 case 3:
1452 *L = negInfinity;
1453 *U = Infinity;
1454 xscanf(R, ""); /* purge line */
1455 break;
1456 case 4:
1457 if (xscanf(R, "%lf", L) != 1)
1458 badline(R);
1459 *U = *L;
1460 break;
1461 case 5:
1462 if (Cvar) {
1463 if (xscanf(R, "%d %d", &k, &j) == 2
1464 && j > 0 && j <= nv) {
1465 Cvar[i] = j;
1466 *L = k & 2 ? negInfinity : 0.;
1467 *U = k & 1 ? Infinity : 0.;
1468 break;
1469 }
1470 }
1471 default:
1472 badline(R);
1473 }
1474 }
1475 }
1476
1477 static expr *
1478 #ifdef KR_headers
aholread(R)1479 aholread(R) EdRead *R;
1480 #else
1481 aholread(EdRead *R)
1482 #endif
1483 {
1484 int i, k;
1485 expr_h *rvh;
1486 char *s, *s1;
1487 FILE *nl = R->nl;
1488
1489 k = getc(nl);
1490 if (k < '1' || k > '9')
1491 badline(R);
1492 i = k - '0';
1493 while((k = getc(nl)) != ':') {
1494 if (k < '0' || k > '9')
1495 badline(R);
1496 i = 10*i + k - '0';
1497 }
1498 rvh = (expr_h *)mem(memadj(sizeof(expr_h) + i));
1499 for(s1 = rvh->sym;;) {
1500 if ((k = getc(nl)) < 0) {
1501 fprintf(Stderr,
1502 "Premature end of file in aholread, line %ld of %s\n",
1503 R->Line, R->asl->i.filename_);
1504 exit_ASL(R,1);
1505 }
1506 if (k == '\n') {
1507 R->Line++;
1508 if (!i)
1509 break;
1510 }
1511 if (--i < 0)
1512 badline(R);
1513 *s1++ = k;
1514 }
1515 *s1 = 0;
1516 rvh->op = f_OPHOL;
1517 rvh->a = nv1;
1518 return (expr *)rvh;
1519 }
1520
1521 static expr *
1522 #ifdef KR_headers
bholread(R)1523 bholread(R) EdRead *R;
1524 #else
1525 bholread(EdRead *R)
1526 #endif
1527 {
1528 int i;
1529 expr_h *rvh;
1530 char *s;
1531
1532 if (xscanf(R, "%d", &i) != 1)
1533 badline(R);
1534 rvh = (expr_h *)mem(memadj(sizeof(expr_h) + i));
1535 s = rvh->sym;
1536 if (fread(s, i, 1, R->nl) != 1)
1537 badline(R);
1538 s[i] = 0;
1539 rvh->op = f_OPHOL;
1540 rvh->a = nv1;
1541 for(;;) switch(*s++) {
1542 case 0: goto break2; /* so we return at end of fcn */
1543 case '\n': R->Line++;
1544 }
1545 break2:
1546 return (expr *)rvh;
1547 }
1548
1549 int
1550 #ifdef KR_headers
fgh_read_ASL(a,nl,flags)1551 fgh_read_ASL(a, nl, flags) ASL *a; FILE *nl; int flags;
1552 #else
1553 fgh_read_ASL(ASL *a, FILE *nl, int flags)
1554 #endif
1555 {
1556 int i, j, k, maxfwd1, nc, nco, ncom, nlcon, nlin, no, nv, nvc, nvo, nz;
1557 int *ka, readall;
1558 unsigned x;
1559 expr_v *e;
1560 cgrad *cg, **cgp;
1561 ograd *og, **ogp;
1562 char fname[128];
1563 func_info *fi;
1564 real t;
1565 EdRead ER, *R;
1566 Jmp_buf JB;
1567
1568 ASL_CHECK(a, ASL_read_fgh, "fgh_read");
1569 asl = (ASL_fgh*)a;
1570 #define asl ((ASL_fgh*)a)
1571
1572 ed_reset(a);
1573 R = EdReadInit_ASL(&ER, a, nl, 0);
1574 if (flags & ASL_return_read_err) {
1575 a->i.err_jmp_ = &JB;
1576 i = setjmp(JB.jb);
1577 if (i) {
1578 a->i.err_jmp_ = 0;
1579 return i;
1580 }
1581 }
1582
1583 nlcon = a->i.n_lcon_;
1584 if (nlcon && !(flags & ASL_allow_CLP)) {
1585 if (a->i.err_jmp_)
1586 return ASL_readerr_CLP;
1587 sorry_CLP(R, "logical constraints");
1588 }
1589
1590 if ((readall = flags & ASL_keep_all_suffixes)
1591 && a->i.nsuffixes)
1592 readall |= 1;
1593 if (nfunc)
1594 func_add(a);
1595 if (binary_nl) {
1596 holread = bholread;
1597 xscanf = bscanf;
1598 }
1599 else {
1600 holread = aholread;
1601 xscanf = ascanf;
1602 }
1603
1604 ncom = comb + comc + como + comc1 + como1;
1605 nc = n_con;
1606 no = n_obj;
1607 nvc = c_vars;
1608 nvo = o_vars;
1609 if (no < 0 || (nco = nc + no + nlcon) <= 0)
1610 scream(R, ASL_readerr_corrupt,
1611 "ed2read: nc = %d, no = %d, nlcon = %d\n",
1612 nc, no, nlcon);
1613 nv1 = nv0 = nvc > nvo ? nvc : nvo;
1614 max_var = nv = nv0 + ncom;
1615 combc = comb + comc;
1616 ncom0 = ncom_togo = combc + como;
1617 nzclim = ncom0 >> 3;
1618 ncom1 = comc1 + como1;
1619 nv0b = nv0 + comb;
1620 nv0c = nv0b + comc;
1621 nv01 = nv0 + ncom0;
1622 amax = lasta = lasta0 = lasta00 = nv1 + 1;
1623 ka = 0;
1624 if ((maxfwd1 = maxfwd + 1) > 1)
1625 nvref = maxfwd1*((ncom0 < vrefGulp ? ncom0 : vrefGulp) + 1);
1626 x = nco*sizeof(cde) + nc*sizeof(cgrad *) + no*sizeof(ograd *)
1627 + nv*(sizeof(expr_v) + 2*sizeof(int))
1628 + ncom0*sizeof(cexp)
1629 + ncom1*sizeof(cexp1)
1630 + nfunc*sizeof(func_info *)
1631 + nvref*sizeof(int)
1632 + no;
1633 nvar0 = a->i.n_var0;
1634 if (!(nvinc = a->i.n_var_ - nvar0))
1635 nvar0 += ncom0 + ncom1;
1636 if (pi0) {
1637 memset(pi0, 0, nc*sizeof(real));
1638 if (havepi0)
1639 memset(havepi0, 0, nc);
1640 }
1641 if (X0)
1642 memset(X0, 0, nv0*sizeof(real));
1643 if (havex0)
1644 memset(havex0, 0, nv0);
1645 e = var_e = (expr_v *)M1zapalloc(x);
1646 var_ex = e + nv0;
1647 var_ex1 = var_ex + ncom0;
1648 con_de = (cde *)(e + nv);
1649 lcon_de = con_de + nc;
1650 for(k = 0; k < nv; e++) {
1651 e->op = f_OPVARVAL;
1652 e->a = k++;
1653 }
1654 obj_de = lcon_de + nlcon;
1655 Cgrad = (cgrad **)(obj_de + no);
1656 Ograd = (ograd **)(Cgrad + nc);
1657 cexps = (cexp *)(Ograd + no);
1658 cexpsc = cexps + comb;
1659 cexpso = cexpsc + comc;
1660 cexps1 = (cexp1 *)(cexpse = cexps + ncom0);
1661 funcs = (func_info **)(cexps1 + ncom1);
1662 zc = (int *)(funcs + nfunc);
1663 zci = zc + nv;
1664 vrefx = zci + nv;
1665 objtype = (char *)(vrefx + nvref);
1666 if (nvref) {
1667 vrefnext = vrefx + maxfwd1;
1668 nvref -= maxfwd1;
1669 }
1670 if (n_cc && !cvar)
1671 cvar = (int*)M1alloc(nc*sizeof(int));
1672 if (cvar)
1673 memset(cvar, 0, nc*sizeof(int));
1674 last_d = 0;
1675 nz = 0;
1676 for(;;) {
1677 ER.can_end = 1;
1678 i = edag_peek(R);
1679 if (i == EOF) {
1680 free(imap);
1681 adjoints = (real *)M1zapalloc(amax*Sizeof(real));
1682 adjoints_nv1 = &adjoints[nv1];
1683 adjust(flags);
1684 nzjac = nz;
1685 if (!Lastx)
1686 Lastx = (real *)M1alloc(nv0*sizeof(real));
1687 fclose(nl);
1688 nderps += nderp;
1689 a->p.Objval = obj2val_ASL;
1690 a->p.Objgrd = obj2grd_ASL;
1691 a->p.Conval = con2val_ASL;
1692 a->p.Jacval = jac2val_ASL;
1693 a->p.Hvcomp = hv2comp_ASL;
1694 a->p.Conival = con2ival_ASL;
1695 a->p.Congrd = con2grd_ASL;
1696 a->p.Lconval = lcon2val_ASL;
1697 a->p.Xknown = x2known_ASL;
1698 a->i.err_jmp_ = 0;
1699 return 0;
1700 }
1701 ER.can_end = 0;
1702 k = -1;
1703 switch(i) {
1704 case 'C':
1705 xscanf(R, "%d", &k);
1706 if (k < 0 || k >= nc)
1707 badline(R);
1708 co_read(R, con_de + k, 1);
1709 break;
1710 case 'F':
1711 if (xscanf(R, "%d %d %d %127s",
1712 &i, &j, &k, fname) != 4
1713 || i < 0 || i >= nfunc)
1714 badline(R);
1715 if (fi = func_lookup(a, fname,0)) {
1716 if (fi->nargs != k && fi->nargs >= 0
1717 && (k >= 0 || fi->nargs < -(k+1)))
1718 scream(R, ASL_readerr_argerr,
1719 "function %s: disagreement of nargs: %d and %d\n",
1720 fname,fi->nargs, k);
1721 }
1722 else {
1723 fi = (func_info *)mem(sizeof(func_info));
1724 fi->ftype = j;
1725 fi->nargs = k;
1726 fi->funcp = 0;
1727 fi->name = (Const char *)strcpy((char*)
1728 mem(memadj(strlen(fname)+1)),
1729 fname);
1730 }
1731 if (!fi->funcp && !(fi->funcp = dynlink(fname)))
1732 scream(R, ASL_readerr_unavail,
1733 "function %s not available\n",
1734 fname);
1735 funcs[i] = fi;
1736 break;
1737 case 'G':
1738 if (xscanf(R, "%d %d", &j, &k) != 2
1739 || j < 0 || j >= no || k <= 0 || k > nvo)
1740 badline(R);
1741 ogp = Ograd + j;
1742 while(k--) {
1743 *ogp = og = (ograd *)mem(sizeof(ograd));
1744 ogp = &og->next;
1745 if (xscanf(R, "%d %lf", &og->varno,
1746 &og->coef) != 2)
1747 badline(R);
1748 }
1749 *ogp = 0;
1750 break;
1751 case 'J':
1752 if (xscanf(R, "%d %d", &j, &k) != 2
1753 || j < 0 || j >= nc || k <= 0 || k > nvc)
1754 badline(R);
1755 nz += k;
1756 if (ka) {
1757 if (!A_vals)
1758 goto cg_read;
1759 j += Fortran;
1760 while(k--) {
1761 if (xscanf(R, "%d %lf",
1762 &i, &t) != 2)
1763 badline(R);
1764 i = ka[i]++;
1765 A_vals[i] = t;
1766 A_rownos[i] = j;
1767 }
1768 break;
1769 }
1770 cg_read:
1771 cgp = Cgrad + j;
1772 j = 0;
1773 while(k--) {
1774 *cgp = cg = (cgrad *)mem(sizeof(cgrad));
1775 cgp = &cg->next;
1776 if (ka) {
1777 if (xscanf(R, "%d %lf",
1778 &cg->varno,
1779 &cg->coef) != 2)
1780 badline(R);
1781 }
1782 else
1783 if (xscanf(R, "%d %d %lf",
1784 &cg->varno, &j,
1785 &cg->coef) != 3)
1786 badline(R);
1787 cg->goff = j;
1788 }
1789 *cgp = 0;
1790 break;
1791 case 'L':
1792 xscanf(R, "%d", &k);
1793 if (k < 0 || k >= nlcon)
1794 badline(R);
1795 co_read(R, lcon_de + k, 0);
1796 break;
1797 case 'O':
1798 if (xscanf(R, "%d %d", &k, &j) != 2
1799 || k < 0 || k >= no)
1800 badline(R);
1801 objtype[k] = j;
1802 co_read(R, obj_de + k, 1);
1803 break;
1804 case 'V':
1805 if (xscanf(R, "%d %d %d", &k, &nlin, &j) != 3)
1806 badline(R);
1807 if (k >= nvar0)
1808 k += nvinc;
1809 if (k < nv0 || k >= nv)
1810 badline(R);
1811 if (j)
1812 cexp1_read(R, j, k, nlin);
1813 else
1814 cexp_read(R, k, nlin);
1815 break;
1816 case 'S':
1817 Suf_read_ASL(R, readall);
1818 break;
1819 case 'r':
1820 br_read(R, asl->i.n_con0, nc, &LUrhs,
1821 Urhsx, cvar, nv0);
1822 break;
1823 case 'b':
1824 br_read(R, asl->i.n_var0, nv0, &LUv,
1825 Uvx, 0, 0);
1826 break;
1827 case 'k':
1828 k_seen++;
1829 k = asl->i.n_var0;
1830 if (!xscanf(R,"%d",&j) || j != k - 1)
1831 badline(R);
1832 if (!(ka = A_colstarts)) {
1833 if ((i = k) < n_var)
1834 i = n_var;
1835 ka = A_colstarts = (int *)
1836 M1alloc((i+1)*Sizeof(int));
1837 }
1838 *ka++ = 0;
1839 *ka++ = 0; /* sic */
1840 while(--k > 0)
1841 if (!xscanf(R, "%d", ka++))
1842 badline(R);
1843 ka = A_colstarts + 1;
1844 break;
1845 case 'x':
1846 if (!xscanf(R,"%d",&k)
1847 || k < 0 || k > nv)
1848 badline(R);
1849 if (!X0 && want_xpi0 & 1) {
1850 x = nv0*sizeof(real);
1851 if (want_xpi0 & 4)
1852 x += nv0;
1853 X0 = (real *)M1zapalloc(x);
1854 if (want_xpi0 & 4)
1855 havex0 = (char*)(X0 + nv0);
1856 }
1857 while(k--) {
1858 if (xscanf(R, "%d %lf", &j, &t) != 2)
1859 badline(R);
1860 if (X0) {
1861 X0[j] = t;
1862 if (havex0)
1863 havex0[j] = 1;
1864 }
1865 }
1866 break;
1867 case 'd':
1868 if (!xscanf(R,"%d",&k)
1869 || k < 0 || k > nc)
1870 badline(R);
1871 if (!pi0 && want_xpi0 & 2) {
1872 x = nc*sizeof(real);
1873 if (want_xpi0 & 4)
1874 x += nc;
1875 pi0 = (real *)M1zapalloc(x);
1876 if (want_xpi0 & 4)
1877 havepi0 = (char*)(pi0 + nc);
1878 }
1879 while(k--) {
1880 if (xscanf(R, "%d %lf", &j, &t) != 2
1881 || j < 0 || j >= nc)
1882 badline(R);
1883 if (pi0) {
1884 pi0[j] = t;
1885 if (havepi0)
1886 havepi0[j] = 1;
1887 }
1888 }
1889 break;
1890 default:
1891 badline(R);
1892 }
1893 }
1894 }
1895 #ifdef __cplusplus
1896 }
1897 #endif
1898