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