/[ascend]/trunk/tcltk98/generic/interface/SlvProc.c
ViewVC logotype

Contents of /trunk/tcltk98/generic/interface/SlvProc.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 411 - (show annotations) (download) (as text)
Sun Apr 2 07:05:54 2006 UTC (16 years, 1 month ago) by ben.allan
File MIME type: text/x-csrc
File size: 32008 byte(s)
Restored autotools to working, parsers to typ_ and zz_,
Fixed many missing initializations, many casting insanities
that have been creeping in, many missing forward declarations
in preparation for fixing external relations.

1 /*
2 * SlvProc.c
3 * by Ben Allan
4 * Created: 1/94
5 * Version: $Revision: 1.6 $
6 * Version control file: $RCSfile: SlvProc.c,v $
7 * Date last modified: $Date: 2003/08/23 18:43:08 $
8 * Last modified by: $Author: ballan $
9 *
10 * This file is part of the ASCEND Tcl/Tk interface
11 *
12 * Copyright 1997, Carnegie Mellon University
13 *
14 * The ASCEND Tcl/Tk interface is free software; you can redistribute
15 * it and/or modify it under the terms of the GNU General Public License as
16 * published by the Free Software Foundation; either version 2 of the
17 * License, or (at your option) any later version.
18 *
19 * The ASCEND Tcl/Tk interface is distributed in hope that it will be
20 * useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
22 * General Public License for more details.
23 *
24 * You should have received a copy of the GNU General Public License
25 * along with the program; if not, write to the Free Software Foundation,
26 * Inc., 675 Mass Ave, Cambridge, MA 02139 USA. Check the file named
27 * COPYING. COPYING is found in ../compiler.
28 */
29
30 #include <math.h>
31 #include <tcl.h>
32 #include <utilities/ascConfig.h>
33 #include <utilities/ascMalloc.h>
34 #include <general/list.h> /* needed? */
35 #include <compiler/instance_enum.h>
36 #include <solver/slv_types.h>
37 #include <solver/mtx.h>
38 #include <solver/var.h>
39 #include <solver/rel.h>
40 #include <solver/discrete.h>
41 #include <solver/conditional.h>
42 #include <solver/logrel.h>
43 #include <solver/bnd.h>
44 #include <solver/slv_common.h>
45 #include <solver/linsol.h>
46 #include <solver/linsolqr.h>
47 #include <solver/slv_client.h>
48 #include "HelpProc.h"
49 #include "SolverGlobals.h"
50 #include "Commands.h" /* for registration function */
51 #include "SlvProc.h"
52
53 #ifndef lint
54 static CONST char SlvProcID[] = "$Id: SlvProc.c,v 1.6 2003/08/23 18:43:08 ballan Exp $";
55 #endif
56
57
58 #ifndef MAXIMUM_STRING_LENGTH
59 #define MAXIMUM_STRING_LENGTH 256
60 #endif
61
62
63 /* hellaciously bad assumption */
64 #define MAXIMUM_ID_LENGTH 80
65
66
67 /** monitor stuff *******************************************/
68
69 /*
70 * cast a pointer to monitor.
71 */
72 #define SMC(cdata) ((struct SlvMonitor *)cdata)
73
74 /* the usual free and reassign operator */
75 #undef free_unless_null
76 #define free_unless_null(p) if ((p)!=NULL) { ascfree(p); } (p)=NULL
77
78 /*
79 * a bunch of bit values for the updated field of the struct below.
80 */
81 #define MON_ALLCLEAR 0x0
82 #define MON_VARCHANGE 0x1
83 #define MON_RELCHANGE 0x2
84 #define MON_VARSPEED 0x4
85 #define MON_RELLOG 0x8
86
87 #define MONALLSET(u) (u) = \
88 (MON_VARCHANGE | MON_RELCHANGE | MON_VARSPEED | MON_RELLOG)
89
90 #define MONCLEAR(u,bit) (u) &= (~(bit))
91
92 /*
93 * our basic data package.
94 * Functions using the package should pay attention to the
95 * update flag and if it is set, clear it while doing the appropriate
96 * action.
97 */
98 struct SlvMonitor {
99 char *interface_id; /* a constant symbolic handle for this object */
100 real64 *lastrelres; /* last scaled real relation residual */
101 real64 *lastvarval; /* last scaled variable value */
102 real64 *lastvardel; /* last scaled variable delta */
103 unsigned int sys_id; /* slv system unique identifier */
104 int32 nrels; /* size of lastrelres */
105 int32 nvars; /* size of lastvarval,lastvardel */
106 int ulx, uly; /* upper left corner of plot region */
107 int w, h; /* WxH plot region size, pixels. */
108 int updated; /* geometry change or wholesale data change */
109 real64 varmax; /* largest plottable value of any scaled variable */
110 real64 relmax; /* largest plottable value of any scaled residual */
111 real64 relmin; /* smallest plottable value of any scaled residual */
112 };
113
114
115
116
117 /*********************************************************************/
118
119 /*
120 * Analyze Stuff. Everything below here is independent
121 * of everything above and should be in a separate file.
122 * Maybe a ascSolverQuery.[ch]
123 */
124
125 /*
126 * This is the beginning of some analysis routines.
127 * For upper bounds, lower bounds, nominals, scaling
128 * langrange multipliers, residuals etc.
129 */
130
131 #define SMALL_NUMBER 1e-12
132
133 enum Bounds_Enum {
134 b_close, b_far, b_equal,
135 b_invalid,
136 b_lower, b_upper, b_nominal, b_othervalue,
137 b_residual
138 };
139
140 static enum Bounds_Enum CloseToBound(real64 value, real64 bound,
141 real64 tolerance,
142 enum Bounds_Enum bounds_type,
143 int relative_check)
144 {
145 double error;
146
147 switch(bounds_type) {
148 case b_lower:
149 if (value <= bound-SMALL_NUMBER) {
150 return b_invalid;
151 }
152 break;
153 case b_upper:
154 if (value >= bound+SMALL_NUMBER) {
155 return b_invalid;
156 }
157 break;
158 case b_nominal:
159 case b_residual:
160 case b_othervalue:
161 break;
162 default:
163 return b_invalid;
164 }
165 /*
166 * The rest of this code does not care about
167 * the type of bound that we are checking.
168 */
169 if ((fabs(value) < SMALL_NUMBER) || (relative_check==0)) {
170 error = fabs(value - bound);
171 } else { /* do relative_check all other cases */
172 error = fabs((value - bound)/value);
173 }
174 if (error <= tolerance) {
175 return b_close;
176 } else {
177 return b_far;
178 }
179 }
180
181
182 /*
183 *
184 * This is the analyze routine.
185 * We also look for relations.
186 * We extract the necessary values, and pass them on to the
187 * analyze routine.
188 * If b_type==b_othervalue we should do the
189 * comparison. For example, one could be checking the value of a
190 * lower_bound against some other value. HOWEVER, if we have
191 * b_type!=b_othervalue, we just simply ignore it.
192 *
193 */
194
195
196 static void DoVarAnalyze(Tcl_Interp *interp,
197 struct var_variable **v,
198 unsigned long low, unsigned long high,
199 enum Bounds_Enum b_type,
200 real64 tolerance,
201 real64 othervalue, /* an arbitrary value */
202 int relative_check)
203 {
204 enum Bounds_Enum b_result; /* result of the query */
205 real64 value =0,checkvalue =0;
206 unsigned long c;
207 char tmp[80];
208
209 for (c=low; c<=high; c++) {
210 switch(b_type) {
211 case b_nominal:
212 value = var_value(v[c]);
213 checkvalue = var_nominal(v[c]);
214 break;
215 case b_lower:
216 value = var_value(v[c]);
217 checkvalue = var_lower_bound(v[c]);
218 break;
219 case b_upper:
220 value = var_value(v[c]);
221 checkvalue = var_upper_bound(v[c]);
222 break;
223 case b_othervalue:
224 value = var_value(v[c]);
225 checkvalue = othervalue;
226 break;
227 default:
228 b_result = b_invalid; /* should not be here */
229 break;
230 }
231 b_result = CloseToBound(value,checkvalue,tolerance,
232 b_type,relative_check);
233 if (b_result==b_close) {
234 sprintf(tmp,"%lu b_close",c);
235 Tcl_AppendElement(interp,tmp);
236 }
237 }
238 }
239
240
241
242 static void DoRelAnalyze(Tcl_Interp *interp,
243 struct rel_relation **r,
244 unsigned long low, unsigned long high,
245 enum Bounds_Enum b_type,
246 real64 tolerance,
247 real64 othervalue,
248 int relative_check)
249 {
250 enum Bounds_Enum b_result;
251 real64 value,checkvalue;
252 unsigned long c;
253 char tmp[80];
254
255 (void)othervalue; /* stop gcc whine about unused parameter */
256
257 for (c=low; c<=high; c++) {
258 if (b_type==b_residual) {
259 value = rel_residual(r[c]);
260 checkvalue = 0.0; /* see if the residual is zero */
261 b_result = CloseToBound(value,checkvalue,tolerance,
262 b_type,relative_check);
263 if (b_result==b_close) {
264 sprintf(tmp,"%lu b_close",c);
265 Tcl_AppendElement(interp,tmp);
266 }
267 } else {
268 b_result = b_invalid;
269 }
270 }
271 }
272
273
274
275 /*
276 * Here is the version of the implementation.
277 * Usage: __var_analyze low high \
278 * scaling?lower?upper?other tolerance rel?abs <otherval>.
279 * Usage: __rel_analyze low high \
280 * residual?other tolerance rel?abs <otherval>.
281 */
282
283
284
285 int Asc_VarAnalyzeCmd(ClientData cdata, Tcl_Interp *interp,
286 int argc, CONST84 char *argv[])
287 {
288 struct var_variable **vp;
289 unsigned long maxvar;
290
291 unsigned long low,high; /* range checking */
292 enum Bounds_Enum b_type; /* query type */
293 real64 tolerance = 0; /* tolerance value */
294 real64 othervalue = 0; /* a check an arbitrary value
295 * which is not a child attribute.*/
296 int relative_check = 1; /* what type of check;
297 * relative is the default */
298
299 (void)cdata; /* stop gcc whine about unused parameter */
300
301 if ( g_solvsys_cur == NULL ) {
302 FPRINTF(stderr, "Asc_VarAnalyzeCmd called with NULL pointer\n");
303 Tcl_SetResult(interp, "Asc_VarAnalyzeCmd called without slv_system",
304 TCL_STATIC);
305 return TCL_ERROR;
306 }
307
308 vp=slv_get_solvers_var_list(g_solvsys_cur);
309 maxvar = (unsigned long) slv_get_num_solvers_vars(g_solvsys_cur);
310
311 if ( argc < 6 ) {
312 Tcl_AppendResult(interp,"wrong # args: Usage :",
313 "\" __var_analyze\" low high \n",
314 "scaling?lower?upper?other tolerance rel?abs <othervalue>",
315 (char *)NULL);
316 return TCL_ERROR;
317 }
318
319 low = atol(argv[1]);
320 high = atol(argv[2]);
321 if (!((low>0)&&(high>0)&&(high<=maxvar))) {
322 Tcl_SetResult(interp, "Invalid index ranges in __var_analyze", TCL_STATIC);
323 return TCL_ERROR;
324 }
325
326 /*
327 * We should now have a valid solver_var description.
328 * Now go for the analysis type. Lagrange multipliers
329 * could be added here, for relations as well as say
330 * relation residuals.
331 */
332 if (strncmp(argv[3],"scaling",3)==0) {
333 b_type = b_nominal;
334 } else if (strncmp(argv[3],"lower",3)==0) {
335 b_type = b_lower;
336 } else if (strncmp(argv[3],"upper",3)==0) {
337 b_type = b_upper;
338 } else if (strncmp(argv[3],"other",3)==0) {
339 b_type = b_othervalue;
340 if ( argc != 7 ) {
341 Tcl_AppendResult(interp,"A \"other value\" analysis requires an ",
342 " additional arg which is the comparison value",
343 (char *)NULL);
344 return TCL_ERROR;
345 } else {
346 othervalue = atof(argv[6]);
347 }
348 } else {
349 Tcl_SetResult(interp, "Invalid analyze type requested\n", TCL_STATIC);
350 return TCL_ERROR;
351 }
352
353 /*
354 * We should now have a valid analyze query. The last 2 things
355 * that we need are a relative or absolute check and a tolerance.
356 */
357
358 tolerance = atof(argv[4]);
359 if (strncmp(argv[5],"relative",3)==0) {
360 relative_check = 1;
361 } else {
362 relative_check = 0;
363 }
364
365 DoVarAnalyze(interp,vp,
366 low,high,
367 b_type,
368 tolerance,
369 othervalue,
370 relative_check);
371 return TCL_OK;
372 }
373
374
375 int Asc_RelAnalyzeCmd(ClientData cdata, Tcl_Interp *interp,
376 int argc, CONST84 char *argv[])
377 {
378 struct rel_relation **rp;
379 unsigned long maxrel;
380
381 unsigned long low,high; /* range checking */
382 enum Bounds_Enum b_type; /* query type */
383 double tolerance = 0; /* tolerance value */
384 double othervalue = 0; /* a check an arbitrary value
385 * which is not a child attribute.*/
386 int relative_check = 1; /* what type of check;
387 * relative is the default */
388
389 (void)cdata; /* stop gcc whine about unused parameter */
390
391 if ( g_solvsys_cur == NULL ) {
392 FPRINTF(stderr, "Asc_RelAnalyzeCmd called with NULL pointer\n");
393 Tcl_SetResult(interp, "Asc_RelAnalyzeCmd called without slv_system",
394 TCL_STATIC);
395 return TCL_ERROR;
396 }
397
398 rp=slv_get_solvers_rel_list(g_solvsys_cur);
399 maxrel = (unsigned long) slv_get_num_solvers_rels(g_solvsys_cur);
400
401 if ( argc < 6 ) {
402 Tcl_AppendResult(interp,"wrong # args: Usage :",
403 "\" __rel_analyze\" low high\n",
404 "residual?other tolerance rel?abs <othervalue>",
405 (char *)NULL);
406 return TCL_ERROR;
407 }
408
409 low = atol(argv[1]);
410 high = atol(argv[2]);
411 if (!((low>0)&&(high>0)&&(high<=maxrel))) {
412 Tcl_SetResult(interp, "Invalid index ranges in __rel_analyze", TCL_STATIC);
413 return TCL_ERROR;
414 }
415
416 if (strncmp(argv[3],"residual",3)==0) {
417 b_type = b_residual;
418 } else if (strncmp(argv[3],"other",3)==0) {
419 b_type = b_othervalue;
420 if ( argc != 7 ) {
421 Tcl_AppendResult(interp,"A \"other value\" analysis requires an ",
422 " additional arg which is the comparison value",
423 (char *)NULL);
424 return TCL_ERROR;
425 } else {
426 othervalue = atof(argv[6]);
427 }
428 } else {
429 Tcl_SetResult(interp, "Invalid rel analyze type requested\n", TCL_STATIC);
430 return TCL_ERROR;
431 }
432
433 tolerance = atof(argv[4]);
434 if (strncmp(argv[5],"relative",3)==0) {
435 relative_check = 1;
436 } else {
437 relative_check = 0;
438 }
439
440 DoRelAnalyze(interp,rp,low,high,b_type,tolerance,othervalue,relative_check);
441 return TCL_OK;
442 }
443
444
445 /*
446 * start of a solver monitor.
447 * These functions are to be generic to every NLP client, so that
448 * we DO NOT NEED to depend on block or mtx information.
449 */
450
451 static
452 void MonDestroy(struct SlvMonitor *m)
453 {
454 if (m==NULL) {
455 return;
456 }
457 free_unless_null(m->interface_id);
458 free_unless_null(m->lastrelres);
459 free_unless_null(m->lastvarval);
460 free_unless_null(m->lastvardel);
461 m->sys_id = m->nrels = m->nvars = 0;
462 m->updated = 103050301;
463 ascfree(m);
464 }
465
466 /*
467 * Updates the values of vars and rels, if the sys
468 * id or sizes may have changed, possibly becoming non-NULL.
469 * Returns 0 if ok, or 1 if insufficient memory.
470 * If the system identifier is the same and its master list
471 * sizes have not changed, this function does nothing.
472 */
473 static
474 int MonUpdateData(struct SlvMonitor *m, slv_system_t sys)
475 {
476 struct rel_relation **rp;
477 struct var_variable **vp;
478 int len,i;
479 if (m==NULL) {
480 return 1;
481 }
482 if (sys==NULL) {
483 free_unless_null(m->lastrelres);
484 free_unless_null(m->lastvarval);
485 free_unless_null(m->lastvardel);
486 m->sys_id = m->nrels = m->nvars = 0;
487 return 0;
488 }
489 /* handle monitoring new system */
490 if (m->sys_id != slv_serial_id(sys)) {
491 m->sys_id = slv_serial_id(sys);
492
493 free_unless_null(m->lastrelres);
494 len = m->nrels = slv_get_num_master_rels(sys);
495 if (len > 0) {
496 rp = slv_get_master_rel_list(sys);
497 m->lastrelres = (real64 *)ascmalloc(sizeof(real64)*len);
498 if (m->lastrelres == NULL) {
499 return 1;
500 }
501 for (i = 0; i < len; i++) {
502 m->lastrelres[i] = rel_residual(rp[i]);
503 }
504 MONALLSET(m->updated);
505 }
506
507 free_unless_null(m->lastvarval);
508 free_unless_null(m->lastvardel);
509 len = m->nvars = slv_get_num_master_vars(sys);
510 if (len > 0) {
511 vp = slv_get_master_var_list(sys);
512 m->lastvarval = (real64 *)ascmalloc(sizeof(real64)*len);
513 if (m->lastvarval == NULL) {
514 return 1;
515 }
516 m->lastvardel = (real64 *)ascmalloc(sizeof(real64)*len);
517 if (m->lastvardel == NULL) {
518 return 1;
519 }
520 for (i = 0; i < len; i++) {
521 m->lastvarval[i] = var_value(vp[i]);
522 m->lastvardel[i] = 0.0;
523 }
524 MONALLSET(m->updated);
525 }
526
527 return 0;
528 }
529
530 /* if, oddly, the rel list changed size, fix up */
531 len = slv_get_num_master_rels(sys);
532 if (len > m->nrels) {
533 free_unless_null(m->lastrelres);
534 m->nrels = len;
535 m->lastrelres = (real64 *)ascmalloc(sizeof(real64)*len);
536 if (m->lastrelres == NULL) {
537 return 1;
538 }
539 /* update the residual data */
540 rp = slv_get_master_rel_list(sys);
541 for (i = 0; i < len; i++) {
542 m->lastrelres[i] = rel_residual(rp[i]);
543 }
544 MONALLSET(m->updated);
545 }
546
547 /* if, oddly, the var list changed size, fix up */
548 len = slv_get_num_master_vars(sys);
549 if (len > m->nvars) {
550 free_unless_null(m->lastvarval);
551 free_unless_null(m->lastvardel);
552 m->nvars = len;
553 m->lastvarval = (real64 *)ascmalloc(sizeof(real64)*len);
554 if (m->lastvarval == NULL) {
555 return 1;
556 }
557 m->lastvardel = (real64 *)ascmalloc(sizeof(real64)*len);
558 if (m->lastvardel == NULL) {
559 return 1;
560 }
561 /* update the value data */
562 vp = slv_get_master_var_list(sys);
563 for (i = 0; i < len; i++) {
564 m->lastvarval[i] = var_value(vp[i]);
565 m->lastvardel[i] = 0.0;
566 }
567 MONALLSET(m->updated);
568 }
569
570 return 0;
571 }
572
573 /*
574 * if geometry has changed or newness OTHERWISE detected,
575 * updates values. Returns 0 if ok, 1 if check cannot be
576 * completed properly.
577 * If failure, interp has a message appended.
578 */
579 static
580 int MonUpdateCheck(struct SlvMonitor *m, Tcl_Interp *interp,
581 int argc, CONST84 char *argv[], slv_system_t sys)
582 {
583 (void)argc;
584 if (MonUpdateData(m,sys)) {
585 Tcl_AppendResult(interp,argv[0],": malloc failure",(char *)NULL);
586 return 1;
587 }
588 return 0;
589 }
590
591 /*
592 * s change var return scaled values that changed
593 * s change rel return scaled residuals that changed
594 * returns a tcl list whose elements are the changed values
595 * with master list indices {index value}
596 * values are scaled by nominals, or 1 if nominal is 0.
597 *
598 * This function may raise a floating point exception. The
599 * caller should check the exception global and discard the
600 * result of the function if exception occured.
601 * The cause of such exceptions is division of a large value
602 * by a small one.
603 */
604 static
605 int MonChange(struct SlvMonitor *m, Tcl_Interp *interp,
606 int argc, CONST84 char *argv[], slv_system_t sys)
607 {
608 int i,len,do_all=0;
609 struct var_variable **vp;
610 struct rel_relation **rp;
611 real64 tmpd,limit1,limit2,sign;
612 real64 *list;
613 char buf[40];
614 assert(m!=NULL && interp!=NULL && sys != NULL && argc >= 3);
615
616 if (argc != 4) {
617 Tcl_AppendResult(interp,"need option to ",argv[0]," ",argv[1]," ",
618 argv[2],(char *)NULL);
619 return TCL_ERROR;
620 }
621
622 if (MonUpdateCheck(m,interp,argc,argv,sys)) {
623 return TCL_ERROR;
624 }
625
626 switch (argv[3][0]) {
627 case 'v':
628 if (m->updated & MON_VARCHANGE) {
629 MONCLEAR(m->updated,MON_VARCHANGE);
630 do_all = 1;
631 }
632 vp = slv_get_master_var_list(sys);
633 list = m->lastvarval;
634 limit2 = fabs(m->varmax);
635 len = m->nvars;
636 for (i = 0; i < len; i++) {
637 tmpd = var_nominal(vp[i]);
638 if (tmpd == 0.0) {
639 tmpd = 1;
640 }
641 tmpd = var_value(vp[i])/tmpd;
642 if (tmpd > limit2) {
643 tmpd = limit2;
644 } else {
645 if (tmpd < -limit2) {
646 tmpd = (-limit2);
647 }
648 }
649 if (do_all || tmpd != list[i]) {
650 list[i] = tmpd;
651 sprintf(buf,"%d %.18g",i,tmpd); /* do faster with tcl objects */
652 Tcl_AppendElement(interp,buf);
653 }
654 }
655 return TCL_OK;
656 case 'r':
657 if (m->updated & MON_RELCHANGE) {
658 MONCLEAR(m->updated,MON_RELCHANGE);
659 do_all = 1;
660 }
661 rp = slv_get_master_rel_list(sys);
662 list = m->lastrelres;
663 limit1 = m->relmin;
664 limit2 = m->relmax;
665 len = m->nrels;
666 for (i = 0; i < len; i++) {
667 tmpd = rel_nominal(rp[i]);
668 if (tmpd == 0.0) {
669 tmpd = 1;
670 }
671 /* record sign of residual and do work on abs value */
672 tmpd = rel_residual(rp[i])/tmpd;
673 sign = (tmpd <0) ? ((tmpd = -tmpd),-1.0) : 1.0;
674 if (tmpd > limit2) {
675 tmpd = limit2;
676 } else {
677 if (tmpd < limit1) {
678 tmpd = limit1;
679 }
680 }
681 tmpd *= sign;
682 if (do_all || tmpd != list[i]) {
683 list[i] = tmpd;
684 sprintf(buf,"%d %.18g",i,tmpd); /* do faster with tcl objects */
685 Tcl_AppendElement(interp,buf);
686 }
687 }
688 return TCL_OK;
689 default:
690 /* whine about bad arg to change */
691 Tcl_AppendResult(interp,"unknown option to ",argv[0]," ",argv[1]," ",
692 argv[2],": ",argv[3],(char *)NULL);
693 return TCL_ERROR;
694 }
695 }
696
697 /*
698 * s geometry w h x y rmin rmax vmax; sets conversion parameters for plotdata
699 * turns on all need update, presumably because this is called with a window
700 * size or range change that will require new drawing.
701 */
702 static
703 int MonGeometry(struct SlvMonitor *m, Tcl_Interp *interp,
704 int argc, CONST84 char *argv[])
705 {
706 /* parses an int arg n that should look like s and assigns v */
707 #define CSTI(n,s,v) \
708 status = Tcl_GetInt(interp,argv[(n)],&i); \
709 if (status != TCL_OK) { \
710 Tcl_ResetResult(interp); \
711 Tcl_AppendResult(interp,"error parsing ",argv[(n)]," as ",(s), \
712 (char *)NULL); \
713 return TCL_ERROR; \
714 } (v) = i;
715
716 /* parses an real64 arg n that should look like s and assigns v */
717 #define CSTD(n,s,v) \
718 status = Tcl_GetDouble(interp,argv[(n)],&x); \
719 if (status != TCL_OK) { \
720 Tcl_ResetResult(interp); \
721 Tcl_AppendResult(interp,"error parsing ",argv[(n)]," as ",(s), \
722 (char *)NULL); \
723 return TCL_ERROR; \
724 } (v) = fabs(x);
725
726 int status;
727 int i;
728 real64 x;
729
730 assert(m!=NULL);
731 if (argc !=10) {
732 Tcl_AppendResult(interp,argv[0]," ",argv[1]," ", argv[2], " requires ",
733 "4 ints and 3 reals: width height x y minresidual maxresidual maxvar",
734 (char *)NULL);
735 return TCL_ERROR;
736 }
737 CSTI(3,"width",m->w);
738 CSTI(4,"height",m->h);
739 CSTI(5,"x",m->ulx);
740 CSTI(6,"y",m->uly);
741 m->w = abs(m->w);
742 m->h = abs(m->h);
743
744 /* parse positive real limits */
745 CSTD(7,"minresidual",m->relmin);
746 CSTD(8,"maxresidual",m->relmax);
747 CSTD(9,"maxvar",m->varmax);
748 /* swap rel range limits if user is a git. */
749 if (m->relmax < m->relmin) {
750 x = m->relmin;
751 m->relmin = m->relmax;
752 m->relmax = x;
753 }
754 /* enforce display of at least 1 order of magnitude */
755 if (m->relmax < 10*m->relmin) {
756 m->relmax = 10*m->relmin;
757 }
758 MONALLSET(m->updated);
759 return TCL_OK;
760 #undef CSTI
761 #undef CSTD
762 }
763
764 /*
765 * s plotdata value return plot info for scaled values that changed
766 * s plotdata speed return plot info for scaled rates of value change
767 * s plotdata residual return plot info for scaled residuals that changed
768 *
769 * Each option returns a list of {x y index} for changed values of the
770 * variables or relations. The x,y are coordinates at which a point
771 * should be plotted based on a transformation derived from whxy info
772 * last obtained by the Geometry command of the monitor.
773 * The transformation may specify the same coordinate for more than
774 * one relation or variable.
775 * If this function raises a floating point exception, as it may,
776 * the result should be discarded.
777 */
778 static
779 int MonPlotData (struct SlvMonitor *m, Tcl_Interp *interp,
780 int argc, CONST84 char *argv[], slv_system_t sys)
781 {
782 int i,do_all=0;
783 struct var_variable **vp;
784 struct rel_relation **rp;
785 real64 tmpd,limit1,limit2,sign,delta,maxlog,minlog,rlog;
786 int px,py,width,halfheight,x,y,nvars,nrels,center;
787 real64 *list;
788 char buf[40];
789
790 assert(m!=NULL && interp!=NULL && sys != NULL && argc >= 3);
791
792 if (argc != 4) {
793 Tcl_AppendResult(interp,"need option to ",argv[0]," ",argv[1]," ",
794 argv[2],(char *)NULL);
795 return TCL_ERROR;
796 }
797 if (MonUpdateCheck(m,interp,argc,argv,sys)) {
798 return TCL_ERROR;
799 }
800
801 x = m->ulx;
802 y = m->uly;
803 width = m->w;
804 halfheight = m->h/2;
805 center = y + halfheight;
806
807 switch (argv[3][0]) {
808 case 'v':
809 /* plot variables as scaled between bounds +/-varmax.
810 * the variables real bounds do not figure in this picture.
811 */
812 vp = slv_get_master_var_list(sys);
813 if (m->updated & MON_VARCHANGE) {
814 MONCLEAR(m->updated,MON_VARCHANGE);
815 do_all = 1;
816 }
817 nvars = (int)m->nvars;
818 list = m->lastvarval;
819 limit2 = fabs(m->varmax);
820 for (i = 0; i < nvars; i++) {
821 tmpd = var_nominal(vp[i]);
822 if (tmpd == 0.0) {
823 tmpd = 1;
824 }
825 tmpd = var_value(vp[i])/tmpd;
826 if (tmpd > limit2) {
827 tmpd = limit2;
828 } else {
829 if (tmpd < -limit2) {
830 tmpd = (-limit2);
831 }
832 }
833 if (do_all || tmpd != list[i]) {
834 list[i] = tmpd;
835 px = x + (i*width)/nvars;
836 py = center - (int)((tmpd*halfheight)/limit2);
837 sprintf(buf,"%d %d %d",px,py,i); /* do faster with tcl objects */
838 Tcl_AppendElement(interp,buf);
839 }
840 }
841 return TCL_OK;
842 case 's':
843 /* plot variable changes in the scaled space where varmax is the largest
844 * plottable scaled change.
845 */
846 vp = slv_get_master_var_list(sys);
847 list = m->lastvarval;
848 limit2 = fabs(m->varmax);
849 nvars = (int)m->nvars;
850 if (m->updated & MON_VARSPEED) {
851 MONCLEAR(m->updated,MON_VARSPEED);
852 /* all speeds are 0 because we have no good history. update data. */
853 for (i = 0; i < nvars; i++) {
854 tmpd = var_nominal(vp[i]);
855 if (tmpd == 0.0) {
856 tmpd = 1;
857 }
858 tmpd = var_value(vp[i])/tmpd;
859 if (tmpd > limit2) {
860 tmpd = limit2;
861 } else {
862 if (tmpd < -limit2) {
863 tmpd = (-limit2);
864 }
865 }
866 list[i] = tmpd;
867 m->lastvardel[i] = 0.0;
868 px = x + (i*width)/nvars;
869 sprintf(buf,"%d %d %d",px,center,i); /* do faster with tcl objects */
870 Tcl_AppendElement(interp,buf);
871 }
872 return TCL_OK;
873 }
874 /* using the last scaled values, calculate all deltas.
875 * if delta is not the same as last delta, return data for it.
876 * We calculate everything, but we don't return excess stuff for
877 * plotting.
878 */
879 for (i = 0; i < nvars; i++) {
880 tmpd = var_nominal(vp[i]);
881 if (tmpd == 0.0) {
882 tmpd = 1;
883 }
884 tmpd = var_value(vp[i])/tmpd;
885 if (tmpd > limit2) {
886 tmpd = limit2;
887 } else {
888 if (tmpd < -limit2) {
889 tmpd = (-limit2);
890 }
891 }
892 delta = tmpd-list[i];
893 list[i] = tmpd;
894 if (m->lastvardel[i] != delta) {
895 m->lastvardel[i] = delta;
896 px = x + (i*width)/nvars;
897 py = center - (int)((delta*halfheight)/limit2);
898 sprintf(buf,"%d %d %d",px,py,i); /* do faster with tcl objects */
899 Tcl_AppendElement(interp,buf);
900 }
901 }
902 return TCL_OK;
903 case 'r':
904 /* log representation of scaled residuals */
905 rp = slv_get_master_rel_list(sys);
906 if (m->updated & MON_RELLOG) {
907 MONCLEAR(m->updated,MON_RELLOG);
908 do_all = 1;
909 }
910 list = m->lastrelres;
911 limit1 = m->relmin;
912 limit2 = m->relmax;
913 maxlog = log10(limit2);
914 minlog = log10(limit1);
915 delta = maxlog-minlog;
916 nrels = (int)m->nrels;
917 for (i = 0; i < nrels; i++) {
918 tmpd = rel_nominal(rp[i]);
919 if (tmpd == 0.0) {
920 tmpd = 1;
921 }
922 /* record sign of residual and do work on abs value */
923 tmpd = rel_residual(rp[i])/tmpd;
924 sign = (tmpd <0) ? ((tmpd = -tmpd),-1.0) : 1.0;
925 if (tmpd > limit2) {
926 tmpd = limit2;
927 } else {
928 if (tmpd < limit1) {
929 tmpd = limit1;
930 }
931 }
932 rlog = log10(tmpd);
933 tmpd *= sign;
934 if (do_all || tmpd != list[i]) {
935 list[i] = tmpd;
936 px = x + (width*i)/nrels;
937 py = center - (int)(sign*halfheight*((rlog-minlog)/delta));
938 sprintf(buf,"%d %d %d",px,py,i); /* do faster with tcl objects */
939 Tcl_AppendElement(interp,buf);
940 }
941 }
942 return TCL_OK;
943 default:
944 /* whine about bad arg */
945 Tcl_AppendResult(interp,"unknown option to ",argv[0]," ",argv[1]," ",
946 argv[2],": ",argv[3],(char *)NULL);
947 return TCL_ERROR;
948 }
949 }
950
951 /*
952 * Distributor to the options of a monitor command.
953 * destroy destroy monitor
954 * s change var return scaled values that changed
955 * s change rel return scaled residuals that changed
956 * s geometry w h x y rmin rmax vmax; sets conversion parameters for plotdata
957 * s plotdata value return plot info for scaled values that changed
958 * s plotdata speed return plot info for scaled rates of value change
959 * s plotdata residual return plot info for scaled residuals that changed
960 *
961 * The change and plotdata subcommands return information only for those
962 * relations or variables that changed since we last recorded their values
963 * or we last changed geometry in the case of plotdata.
964 *
965 * In the description above, the leading s is the symbolic name of a
966 * slv_system. At present the value of s is ignored and g_solvsys_cur
967 * is assumed. The solver needs to define symbolic handles and a lookup
968 * function for us to interpret s properly.
969 *
970 * Most of the internals of this function and downstream use of
971 * the data it returns in the interpreter could be much faster
972 * if reimplemented in the tcl Object functions when those functions
973 * have stabilized in 8.1.
974 *
975 * This function visits the entire var or rel master list and should not
976 * be called at all on small blocks in large systems.
977 */
978 static
979 int SolveMonitor(ClientData cdata,Tcl_Interp *interp, int argc, CONST84 char *argv[])
980 {
981 slv_system_t sys;
982 char command[80];
983 int i;
984
985 ASCUSE;
986
987 if (cdata==NULL) {
988 /* what the? */
989 return TCL_ERROR;
990 }
991
992 if (argc<2) {
993 /* do usage thing */
994 return TCL_ERROR;
995 }
996 if (argc < 3) {
997 /* must be destroy */
998 if (argv[1][0]=='d') {
999 sprintf(command,"rename %s {}",SMC(cdata)->interface_id);
1000 return Tcl_GlobalEval(interp,command);
1001 }
1002 Tcl_AppendResult(interp,argv[0],": unknown option ",argv[1],(char *)NULL);
1003 return TCL_ERROR;
1004 }
1005
1006 sys = g_solvsys_cur;
1007 /* this needs to be generalized to lookup sys based on argv[1].*/
1008
1009 switch (argv[2][0]) {
1010 case 'c':
1011 if (sys==NULL) {
1012 return TCL_OK;
1013 }
1014 return MonChange(SMC(cdata),interp,argc,argv,sys);
1015 case 'g':
1016 return MonGeometry(SMC(cdata),interp,argc,argv);
1017 case 'p':
1018 if (sys==NULL) {
1019 return TCL_OK;
1020 }
1021 return MonPlotData(SMC(cdata),interp,argc,argv,sys);
1022 default:
1023 break;
1024 }
1025 Tcl_AppendResult(interp,argv[0],": unknown option",(char *)NULL);
1026 for (i=1; i<argc; i++) {
1027 Tcl_AppendResult(interp," ",argv[i],(char *)NULL);
1028 }
1029 return TCL_ERROR;
1030 }
1031
1032 STDHLF(Asc_SolveMonitorCmd,(Asc_SolveMonitorCmdHL1,
1033 Asc_SolveMonitorCmdHL2,
1034 Asc_SolveMonitorCmdHL3,
1035 HLFSTOP));
1036
1037 static
1038 STDHLF(SolveMonitor,(SolveMonitorHL1,
1039 SolveMonitorHL2,
1040 SolveMonitorHL3,
1041 SolveMonitorHL4,
1042 SolveMonitorHL5,
1043 SolveMonitorHL6,
1044 SolveMonitorHL7,
1045 SolveMonitorHL8,
1046 HLFSTOP));
1047 /*
1048 * Creates a monitor and returns its symbolic handle.
1049 * Multiple monitors can exist and are manipulated by their
1050 * symbolic handles.
1051 * A monitor may be used on a series of unrelated slv_system_t.
1052 * Currently, this function gets its slv_system_t from g_solvsys_cur,
1053 * but it should be changed to take a slvsys interface id when
1054 * the solver interface is changed to work by name.
1055 */
1056 int Asc_SolveMonitorCmd(ClientData cdata,Tcl_Interp *interp,
1057 int argc, CONST84 char *argv[])
1058 {
1059 struct SlvMonitor *result;
1060 static unsigned int nextid = 1;
1061
1062 ASCUSE;
1063 if (argc!=1) {
1064 Tcl_AppendResult(interp,argv[0],": no arguments allowed yet",(char *)NULL);
1065 return TCL_ERROR;
1066 }
1067
1068 result = SMC(asccalloc(1,sizeof(struct SlvMonitor)));
1069 if (result==NULL) {
1070 Tcl_AppendResult(interp,argv[0],": insufficient memory",(char *)NULL);
1071 return TCL_ERROR;
1072 }
1073 result->interface_id = (char *)ascmalloc(20+strlen(Asc_SolveMonitorCmdHN));
1074 if (result->interface_id==NULL) {
1075 Tcl_AppendResult(interp,argv[0],": insufficient memory",(char *)NULL);
1076 MonDestroy(result);
1077 return TCL_ERROR;
1078 }
1079 sprintf(result->interface_id,"%s%u",Asc_SolveMonitorCmdHN,nextid++);
1080 if (MonUpdateData(result,g_solvsys_cur)) {
1081 Tcl_AppendResult(interp,argv[0],result->interface_id,
1082 ": insufficient memory",(char *)NULL);
1083 MonDestroy(result);
1084 return TCL_ERROR;
1085 }
1086
1087 result->w = result->h = 1;
1088 MONALLSET(result->updated);
1089
1090 Asc_AddCommand(interp, result->interface_id, SolveMonitor,
1091 (ClientData)result, (Tcl_CmdDeleteProc *)MonDestroy, "solver-monitor",
1092 SolveMonitorHU, SolveMonitorHS, SolveMonitorHLF);
1093
1094 Tcl_AppendResult(interp,result->interface_id,(char *)NULL);
1095 return TCL_OK;
1096 }
1097

john.pye@anu.edu.au
ViewVC Help
Powered by ViewVC 1.1.22