| 1 |
/* |
| 2 |
* SLV: Ascend Nonlinear Solver |
| 3 |
* by Karl Michael Westerberg |
| 4 |
* Created: 2/6/90 |
| 5 |
* Version: $Revision: 1.29 $ |
| 6 |
* Version control file: $RCSfile: relman.h,v $ |
| 7 |
* Date last modified: $Date: 1998/04/23 23:56:24 $ |
| 8 |
* Last modified by: $Author: ballan $ |
| 9 |
* |
| 10 |
* This file is part of the SLV solver. |
| 11 |
* |
| 12 |
* Copyright (C) 1990 Karl Michael Westerberg |
| 13 |
* Copyright (C) 1993 Joseph Zaher |
| 14 |
* Copyright (C) 1994 Joseph Zaher, Benjamin Andrew Allan |
| 15 |
* |
| 16 |
* The SLV solver is free software; you can redistribute |
| 17 |
* it and/or modify it under the terms of the GNU General Public License as |
| 18 |
* published by the Free Software Foundation; either version 2 of the |
| 19 |
* License, or (at your option) any later version. |
| 20 |
* |
| 21 |
* The SLV solver is distributed in hope that it will be |
| 22 |
* useful, but WITHOUT ANY WARRANTY; without even the implied warranty of |
| 23 |
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
| 24 |
* General Public License for more details. |
| 25 |
* |
| 26 |
* You should have received a copy of the GNU General Public License along with |
| 27 |
* the program; if not, write to the Free Software Foundation, Inc., 675 |
| 28 |
* Mass Ave, Cambridge, MA 02139 USA. Check the file named COPYING. |
| 29 |
* COPYING is found in ../compiler. |
| 30 |
*/ |
| 31 |
|
| 32 |
/************************* |
| 33 |
* Contents: Relation manipulator module |
| 34 |
* |
| 35 |
* Authors: Karl Westerberg |
| 36 |
* Joseph Zaher |
| 37 |
* |
| 38 |
* Dates: 06/90 - original version |
| 39 |
* 06/93 - separated out exprman |
| 40 |
* 11/93 - added relman_calc_satisfied |
| 41 |
* |
| 42 |
* Description: This module will provide supplemental operations for |
| 43 |
* relations such as simplification, evaluation, and |
| 44 |
* differentiation. |
| 45 |
*************************/ |
| 46 |
#ifndef relman__already_included |
| 47 |
#define relman__already_included |
| 48 |
|
| 49 |
/* requires #include "mtx.h" */ |
| 50 |
/* requires #include "var.h" */ |
| 51 |
/* requires #include "rel.h" */ |
| 52 |
|
| 53 |
|
| 54 |
#define relman_is_linear(a,b) (FALSE) |
| 55 |
/** |
| 56 |
extern boolean relman_is_linear(struct rel_relation *,var_filter_t *); |
| 57 |
* islinear = relman_is_linear(rel,filter) |
| 58 |
* boolean islinear; |
| 59 |
* struct rel_relation *rel; |
| 60 |
* var_filter_t *filter; |
| 61 |
* |
| 62 |
* Determines whether or not the given relation is linear in |
| 63 |
* all of the variables which pass through the variable filter, treating |
| 64 |
* those variables which fail to pass as constants. |
| 65 |
* |
| 66 |
* Example: |
| 67 |
* x1 + x2 >= 3 is a linear relation. |
| 68 |
**/ |
| 69 |
|
| 70 |
extern real64 relman_linear_coef(struct rel_relation *, |
| 71 |
struct var_variable *,var_filter_t *); |
| 72 |
/** |
| 73 |
* coef = relman_linear_coef(rel,var,filter) |
| 74 |
* real64 coef; |
| 75 |
* struct rel_relation *rel; |
| 76 |
* struct var_variable *var; |
| 77 |
* var_filter_t *filter; |
| 78 |
* |
| 79 |
* Computes the coefficient of the given variable in a linear |
| 80 |
* relation. If var=NULL, then the "RHS" is returned instead. More |
| 81 |
* nprecisely, "a[var]" is returned where: |
| 82 |
* |
| 83 |
* rel : sum[variables v](a[v] * v) COMPARISON a[NULL] |
| 84 |
* |
| 85 |
* It is assumed that the relation is linear, otherwise, |
| 86 |
* something will be returned, but it won't be very meaningful. |
| 87 |
**/ |
| 88 |
|
| 89 |
extern relman_simplify(struct rel_relation *,int); |
| 90 |
/** NOT IMPLEMENTED |
| 91 |
* relman_simplify(rel,opt_level) |
| 92 |
* struct rel_relation *rel; |
| 93 |
* int opt_level; |
| 94 |
* |
| 95 |
* The left and right hand sides of the given relation are simplified |
| 96 |
* to the extent given by opt_level. The effect of varying values of |
| 97 |
* opt_level is given in the description for exprman_simplify(). |
| 98 |
**/ |
| 99 |
|
| 100 |
extern relman_dissolve_vars(struct rel_relation *,var_filter_t *); |
| 101 |
/** NOT IMPLEMENTED |
| 102 |
* relman_dissolve_vars(rel,filter) |
| 103 |
* struct rel_relation *rel; |
| 104 |
* var_filter_t *filter; |
| 105 |
* |
| 106 |
* Variables which pass through the filter are replaced in the |
| 107 |
* relation by their current values. |
| 108 |
**/ |
| 109 |
|
| 110 |
extern void relman_decide_incidence(struct rel_relation *); |
| 111 |
/** |
| 112 |
* relman_decide_incidence(rel) |
| 113 |
* struct rel_relation *rel; |
| 114 |
* |
| 115 |
* Sets the incidence field of each variable which is found to be |
| 116 |
* incident in the relation rel to TRUE. If these variables make |
| 117 |
* up a subset of some larger variable list, it is important to first |
| 118 |
* set the incidence field of all of the variables to FALSE before |
| 119 |
* using this function in order to determine the unattached variables. |
| 120 |
**/ |
| 121 |
|
| 122 |
extern void relman_get_incidence(struct rel_relation *, |
| 123 |
var_filter_t *,mtx_matrix_t); |
| 124 |
/** |
| 125 |
* relman_get_incidence(rel,filter,matrix) |
| 126 |
* struct rel_relation *rel; |
| 127 |
* var_filter_t *filter; |
| 128 |
* mtx_matrix_t matrix; |
| 129 |
* |
| 130 |
* Upon return, coefficient (rel_n,var_n) (using original row and column |
| 131 |
* numbers) is non-zero if and only if the relation rel with index rel_n |
| 132 |
* depends on a variable with index var_n. |
| 133 |
**/ |
| 134 |
|
| 135 |
extern real64 relman_eval(struct rel_relation *, int32 *, int); |
| 136 |
/** |
| 137 |
* residual = relman_eval(rel,status,safe) |
| 138 |
* struct rel_relation *rel; |
| 139 |
* int32 *status; |
| 140 |
* int safe; |
| 141 |
* real64 residual; |
| 142 |
* |
| 143 |
* The residual of the relation is calculated and returned. In addition |
| 144 |
* to returning the residual, the residual field of the relation is |
| 145 |
* updated. Residual := LHS - RHS regardless of comparison. The status |
| 146 |
* value can be monitored to |
| 147 |
* determine if any calculation errors were encountered. It will be set |
| 148 |
* 0 if a calculation results in an error. |
| 149 |
* If the value of safe is nonzero, "safe" functions will be used to |
| 150 |
* calculate the residual. |
| 151 |
* The residual field of the relation is not updated when an error occurs. |
| 152 |
* |
| 153 |
* This function should be surrounded by Asc_SignalHandlerPush/Pop both |
| 154 |
* with arguments (SIGFPE,SIG_IGN). If it is being called in a loop, |
| 155 |
* the push/pop should be _outside_ the loop. |
| 156 |
**/ |
| 157 |
|
| 158 |
extern int32 relman_obj_direction(struct rel_relation *); |
| 159 |
/** |
| 160 |
* direction = relman_eval_obj(rel,status,safe) |
| 161 |
* struct rel_relation *rel; |
| 162 |
* |
| 163 |
* Returns: |
| 164 |
* direction = -1 if objective is minimization |
| 165 |
* direction = 1 if objective is maximization |
| 166 |
* direction = 0 otherwise. (ie. if not an objective) |
| 167 |
**/ |
| 168 |
|
| 169 |
extern real64 relman_scale(struct rel_relation *); |
| 170 |
/** |
| 171 |
* residual = relman_scale(rel) |
| 172 |
* struct rel_relation *rel; |
| 173 |
* real64 residual; |
| 174 |
* |
| 175 |
* Calculates relation nominal scaling factor for |
| 176 |
* current values stored in the relations variables. |
| 177 |
* Fills the relations nominal field and also returns |
| 178 |
* the relations nominal. |
| 179 |
**/ |
| 180 |
|
| 181 |
#define relman_diff(a,b,c,d) (abort(),1) |
| 182 |
/** needs to be reimplemented. |
| 183 |
ExTERn int relman_diff(struct rel_relation *,struct var_variable*,real64*,int); |
| 184 |
*! status = relman_diff(rel,var,pd,safe); |
| 185 |
*! struct rel_relation *rel; |
| 186 |
*! struct var_variable *var; |
| 187 |
*! real64 *pd; |
| 188 |
*! int safe; |
| 189 |
*! int status; |
| 190 |
*! |
| 191 |
*! Calculates the derivative of the relation residual with respect to |
| 192 |
*! the specified variable and stuffs it in pd. if problem with |
| 193 |
*! calculation, returns 1, else 0. |
| 194 |
*! If the value of safe is nonzero, "safe" functions will be used to |
| 195 |
*! calculate the residual. |
| 196 |
*! Needs compiler side work. |
| 197 |
**/ |
| 198 |
|
| 199 |
extern int relman_diff2(struct rel_relation *, var_filter_t *, |
| 200 |
real64 *, int32 *, int32 *, int32); |
| 201 |
/** |
| 202 |
* status = relman_diff2(rel,filter,derivatives,variables,count,safe) |
| 203 |
* struct rel_relation *rel; |
| 204 |
* var_filter_t *filter; |
| 205 |
* real64 *derivatives; |
| 206 |
* int32 *variables; |
| 207 |
* int32 *count; |
| 208 |
* int32 status; |
| 209 |
* |
| 210 |
* Calculates the row of the jacobian matrix (the transpose gradient of |
| 211 |
* the relation residual grad^T(f) ) corresponding to the relation |
| 212 |
* rel. The filter determines which variables actually contribute to the |
| 213 |
* jacobian. |
| 214 |
* If an error is encountered in the calculation, the status returned is |
| 215 |
* 1. Status = 0 is OK. |
| 216 |
* If the value of safe is nonzero, "safe" functions are used to for |
| 217 |
* the calucaltions. |
| 218 |
* The calling function should allocate the output vectors 'derivatives' |
| 219 |
* and 'variables'. 'count' will be set to the number of elements |
| 220 |
* assigned upon exit. |
| 221 |
* derivative(I) will contain the derivative of the relation with |
| 222 |
* respect to the variable whose solver index is stored in |
| 223 |
* variables(I). |
| 224 |
**/ |
| 225 |
|
| 226 |
extern int relman_diff_grad(struct rel_relation *, var_filter_t *, |
| 227 |
real64 *, int32 *, int32 *, int32 *, real64 *, |
| 228 |
int32); |
| 229 |
/** |
| 230 |
*** status = relman_diff_grad(rel,filter,derivatives,variables_master, |
| 231 |
*** variables_solver,count,resid,safe) |
| 232 |
*** struct rel_relation *rel; |
| 233 |
*** var_filter_t *filter; |
| 234 |
*** real64 *derivatives; |
| 235 |
*** real64 *resid; |
| 236 |
*** int32 *variables_master; |
| 237 |
*** int32 *variables_solver; |
| 238 |
*** int32 *count; |
| 239 |
*** int32 status; |
| 240 |
*** |
| 241 |
*** Calculates the row of the jacobian matrix (the transpose gradient of |
| 242 |
*** the relation residual grad^T(f) ) corresponding to the relation |
| 243 |
*** rel. The filter determines which variables actually contribute to the |
| 244 |
*** jacobian. The residual of the relation is also computed and returned. |
| 245 |
*** If an error is encountered in the calculation, the status returned is |
| 246 |
*** 1. Status = 0 is OK. |
| 247 |
*** If the value of safe is nonzero, "safe" functions are used to for |
| 248 |
*** the calculations. |
| 249 |
*** The calling function should allocate the output vectors 'derivatives', |
| 250 |
*** 'variables_master' and 'variables_solver'. 'count' will be set to |
| 251 |
*** the number of elements assigned upon exit. |
| 252 |
*** derivative(i) will contain the derivative of the relation with |
| 253 |
*** respect to the variable whose master index is stored in |
| 254 |
*** variables_master(i). The solver index of each variable is stored in |
| 255 |
*** variables_solver(i). |
| 256 |
*** |
| 257 |
*** Ther are two differences wrt to relman_diff2: |
| 258 |
*** 1) the master index (solver independent) is obtained |
| 259 |
*** 2) the residual is evaluated |
| 260 |
**/ |
| 261 |
|
| 262 |
|
| 263 |
extern int relman_diffs(struct rel_relation *, var_filter_t *, |
| 264 |
mtx_matrix_t, real64 *, int); |
| 265 |
/** |
| 266 |
* status = relman_diffs(rel,filter,mtx,resid,safe) |
| 267 |
* struct rel_relation *rel; |
| 268 |
* var_filter_t *filter; |
| 269 |
* real64 *resid; |
| 270 |
* mtx_matrix_t mtx; |
| 271 |
* int safe; |
| 272 |
* int status; |
| 273 |
* |
| 274 |
* Calculates the row of the jacobian matrix (the transpose gradient of |
| 275 |
* the relation residual grad^T(f) ) corresponding to the relation |
| 276 |
* rel. The filter determines which variables actually contribute to the |
| 277 |
* jacobian. The residual of the relation is also computed and returned. |
| 278 |
* If an error is encountered in the calculation, the status returned is |
| 279 |
* 1 and the residual is set to some number we managed to calculate, |
| 280 |
* while the gradient is discarded. status = 0 is OK. |
| 281 |
* If the value of safe is nonzero, "safe" functions are used to for |
| 282 |
* the calucaltions. |
| 283 |
* It doesn't matter how you have permuted the columns and rows: |
| 284 |
* for the vars which pass the filter you send we |
| 285 |
* fill the org row determined by rel_sindex and the org cols |
| 286 |
* determined by var_sindex. |
| 287 |
* |
| 288 |
* NOTE: The row of the mtx corresponding to rel should be cleared |
| 289 |
* before calling this function, since this FILLS with the gradient. |
| 290 |
* |
| 291 |
* CHANGE: This operator used to just ADD on top of any incidence already |
| 292 |
* in the row. This is not TRUE now. |
| 293 |
* |
| 294 |
* Bugs: none known except this operator really needs to be redesigned |
| 295 |
* so it can deal with harwellian matrices, glassbox rels and blackbox. |
| 296 |
**/ |
| 297 |
|
| 298 |
extern int32 relman_diff_harwell(struct rel_relation **, |
| 299 |
var_filter_t *, rel_filter_t *, |
| 300 |
int32, int32, int32, |
| 301 |
real64 *, int32 *, int32 *); |
| 302 |
/* |
| 303 |
* err = relman_diff_harwell(rlist,vfilter,rfilter,rlen,bias,mors, |
| 304 |
* avec,ivec,jvec); |
| 305 |
* This fills an "a-i-j" sparse matrix in the avec/ivec/jvec given. |
| 306 |
* struct rel_relation **rlist; list of relations rlen long. |
| 307 |
* var_filter_t *vfilter; stuffs gradient for matching variables only. |
| 308 |
* int32 rlen; length of list of relations. |
| 309 |
* int32 bias; 0 = row grouped together, 1 = column grouped together. |
| 310 |
* There is a substantial penalty for bias = 1. we MODEL by row. |
| 311 |
* int32 mors; 0 = master var index of columns, master rel index of rows |
| 312 |
* 1 = solver var index of columns, master rel index of rows |
| 313 |
* 2 = master var index of columns, solver rel index of rows |
| 314 |
* 3 = solver var index of columns, solver rel index of rows |
| 315 |
* Size of avec,ivec,jvec given is assumed big enough. |
| 316 |
* big_enough = relman_jacobian_count(rlist,rlen,vfilter,rfilter,&dummy); |
| 317 |
* If ivec or jvec given is NULL, then neither is stuffed, though avec is. |
| 318 |
* int32 err; |
| 319 |
* err = 1 --> unrecoverable error/bad input. caller should probably punt. |
| 320 |
* err = 0 --> ok; |
| 321 |
* err < 0 --> -(number of floating point errors in evaluation). |
| 322 |
* The matrix will contain an approximation only. |
| 323 |
* |
| 324 |
* Bugs: |
| 325 |
* bias == 1 is not yet implemented. |
| 326 |
*/ |
| 327 |
|
| 328 |
extern int32 relman_jacobian_count(struct rel_relation **, int32, |
| 329 |
var_filter_t *, rel_filter_t *, int32 *); |
| 330 |
/* |
| 331 |
* nnz = relman_jacobian_count(rlist, rlen, vfilter, rfilter, rhomax); |
| 332 |
* Return the number of nonzero gradient entries in the equations |
| 333 |
* given. Only equations passing rfilter and entries passing vfilter |
| 334 |
* are counted. rlen is the length of the relation list. |
| 335 |
* *rhomax is the largest row count on return. |
| 336 |
*/ |
| 337 |
|
| 338 |
extern boolean relman_calc_satisfied(struct rel_relation *,real64); |
| 339 |
extern boolean relman_calc_satisfied_scaled(struct rel_relation *,real64); |
| 340 |
/** |
| 341 |
* satisfied = relman_calc_satisfied(rel,tolerance) |
| 342 |
* satisfied = relman_calc_satisfied_scaled(rel,tolerance) |
| 343 |
* boolean satisfied; |
| 344 |
* real64 tolerance; |
| 345 |
* |
| 346 |
* relman_calc_satisfied: |
| 347 |
* Returns TRUE or FALSE depending on whether the relation whose residual |
| 348 |
* has been previously calculated is satisfied based on the value stored |
| 349 |
* in the residual field. The satisfied field of the relation is also |
| 350 |
* updated. A tolerance specification allows equalities to be declared |
| 351 |
* satisfied as long as their residuals are close to zero. |
| 352 |
* |
| 353 |
* relman_calc_satisfied_scaled: |
| 354 |
* This definition of satisfaction includes the notion |
| 355 |
* of scaling by the relation nominal before comparison. |
| 356 |
**/ |
| 357 |
|
| 358 |
#define relman_directly_solve(r,s,a,n) \ |
| 359 |
relman_directly_solve_new(r,s,a,n,1.0e-8) |
| 360 |
extern real64 *relman_directly_solve_new(struct rel_relation *,struct var_variable *, |
| 361 |
int *,int *,real64); |
| 362 |
/** |
| 363 |
* solution_list = relman_directly_solve(rel,solvefor,able,nsolns) |
| 364 |
* solution_list = relman_directly_solve_new(rel,solvefor,able,nsolns,tol) |
| 365 |
* real64 *solution_list; |
| 366 |
* struct rel_relation *rel; |
| 367 |
* struct var_variable *solvefor; |
| 368 |
* int *able; |
| 369 |
* int *nsolns; |
| 370 |
* real64 tol; |
| 371 |
* |
| 372 |
* Attempts to solve the given equation for the given variable. If this |
| 373 |
* function is able to determine the solution set, then *able is set to |
| 374 |
* 1 and a newly allocated solution list is returned: *nsolns will be |
| 375 |
* set to the length of this array. Otherwise *able is 0 and NULL |
| 376 |
* is returned. NULL *may* also be returned if the solution set is empty. |
| 377 |
* A return of able == 1, solution_list != NULL, and nsolns == 0 is |
| 378 |
* possible for certain classes of floating point exceptions. |
| 379 |
* It is assumed that the relation is a condition of equality. |
| 380 |
* |
| 381 |
* relman_directly_solve_new handles passing in a tolerance for glassbox |
| 382 |
* relations so a rootfinder can do the work rather than leaving it to |
| 383 |
* someone else. The rootfinder is based on Brent's algorithm. Old clients of |
| 384 |
* relman_directly solve are grandfathered at a default tolerance of 1e-8. |
| 385 |
* Once all these clients are gone, go back to the old name space for |
| 386 |
* the new function. |
| 387 |
* |
| 388 |
* Do NOT free or keep a persistent pointer to the solution_list this |
| 389 |
* function returns. |
| 390 |
**/ |
| 391 |
|
| 392 |
#define relman_make_string_infix(s,r) \ |
| 393 |
relman_make_vstring_infix((s),(r),TRUE) |
| 394 |
#define relman_make_string_postfix(s,r) \ |
| 395 |
relman_make_vstring_postfix((s),(r),TRUE) |
| 396 |
#define relman_make_xstring_infix(s,r) \ |
| 397 |
relman_make_vstring_infix((s),(r),FALSE) |
| 398 |
#if 0 /* needs compiler side work */ |
| 399 |
#define relman_make_xstring_postfix(s,r) \ |
| 400 |
relman_make_vstring_postfix((s),(r),FALSE) |
| 401 |
#else |
| 402 |
#define relman_make_xstring_postfix(s,r) \ |
| 403 |
dummy_rel_string(s,r,0) |
| 404 |
#endif |
| 405 |
/** |
| 406 |
* string = relman_make_string_infix(sys,rel) |
| 407 |
* string = relman_make_string_postfix(sys,rel) |
| 408 |
* string = relman_make_xstring_infix(sys,rel) |
| 409 |
* string = relman_make_xstring_postfix(sys,rel) // not working |
| 410 |
* char *string; |
| 411 |
* struct rel_relation *rel; |
| 412 |
* |
| 413 |
* Creates a sufficiently large string (you must free it when you're |
| 414 |
* done with it) into which it converts a relation. The string will be |
| 415 |
* terminated with a '\0' character. |
| 416 |
* |
| 417 |
* The xstring versions of this call make strings where all the variables |
| 418 |
* are written as x<varindex> (e.g. x23) rather than as object names. |
| 419 |
* The MASTER index (var_mindex) is used, not the solver's sindex. |
| 420 |
* Currently the compiler does not support xstring postfix format, |
| 421 |
* but could easily do so if needed. |
| 422 |
* |
| 423 |
* The name of a var is context dependent, so you have to provide the |
| 424 |
* slv_system_t from which you got the relation. |
| 425 |
**/ |
| 426 |
extern char *relman_make_vstring_infix(slv_system_t, |
| 427 |
struct rel_relation *,int); |
| 428 |
extern char *relman_make_vstring_postfix(slv_system_t, |
| 429 |
struct rel_relation *,int); |
| 430 |
/* |
| 431 |
* tmp function to placehold unimplemented io functions. |
| 432 |
*/ |
| 433 |
extern char *dummyrelstring(slv_system_t, struct rel_relation *,int); |
| 434 |
|
| 435 |
/* |
| 436 |
* relman_free_reused_mem(void); |
| 437 |
* Call when desired to free memory cached internally. |
| 438 |
*/ |
| 439 |
extern void relman_free_reused_mem(void); |
| 440 |
|
| 441 |
#endif |
Parent Directory
| 
Revision Log