/[ascend]/trunk/base/generic/solver/slv3.c
ViewVC logotype

Diff of /trunk/base/generic/solver/slv3.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 153 by johnpye, Thu Dec 22 10:58:33 2005 UTC revision 154 by johnpye, Thu Dec 22 15:18:02 2005 UTC
# Line 1  Line 1 
1  /*  /*
2   *  SLV: Ascend Nonlinear Solver      SLV: Ascend Nonlinear Solver
3   *  by Karl Michael Westerberg      by Karl Michael Westerberg
4   *  Created: 2/6/90      Created: 2/6/90
5   *  Version: $Revision: 1.87 $      
6   *  Version control file: $RCSfile: slv3.c,v $      This file is part of the SLV solver.
7   *  Date last modified: $Date: 2000/01/25 02:27:32 $  
8   *  Last modified by: $Author: ballan $      Copyright (C) 1990 Karl Michael Westerberg
9   *      Copyright (C) 1993 Joseph Zaher
10   *  This file is part of the SLV solver.      Copyright (C) 1994 Joseph Zaher, Benjamin Andrew Allan
11   *  
12   *  Copyright (C) 1990 Karl Michael Westerberg      This program is free software; you can redistribute it and/or modify
13   *  Copyright (C) 1993 Joseph Zaher      it under the terms of the GNU General Public License as published by
14   *  Copyright (C) 1994 Joseph Zaher, Benjamin Andrew Allan      the Free Software Foundation; either version 2 of the License, or
15   *      (at your option) any later version.
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      This program is distributed in the hope that it will be useful,
18   *  published by the Free Software Foundation; either version 2 of the      but WITHOUT ANY WARRANTY; without even the implied warranty of
19   *  License, or (at your option) any later version.      MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
20   *      GNU General Public License for more details.
21   *  The SLV solver is distributed in hope that it will be  
22   *  useful, but WITHOUT ANY WARRANTY; without even the implied warranty of      You should have received a copy of the GNU General Public License
23   *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU      along with this program; if not, write to the Free Software
24   *  General Public License for more details.      Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
25   *      This file is part of the SLV solver.
26   *  You should have received a copy of the GNU General Public License  */
27   *  along with the program; if not, write to the Free Software Foundation,  /* Last *CVS* version ballan 2000/01/25 02:27:32 */
  *  Inc., 675 Mass Ave, Cambridge, MA 02139 USA.  Check the file named  
  *  COPYING.  COPYING is found in ../compiler.  
  *  
  */  
28    
29  #include <math.h>  #include <math.h>
30  #include <stdarg.h>  #include <stdarg.h>
# Line 719  static boolean calc_J( slv3_system_t sys Line 715  static boolean calc_J( slv3_system_t sys
715    double time0;    double time0;
716    real64 resid;    real64 resid;
717    
   CONSOLE_DEBUG("jacobian");  
   
718    if( sys->J.accurate )    if( sys->J.accurate )
719      return TRUE;      return TRUE;
720    
# Line 754  static void calc_nominals( slv3_system_t Line 748  static void calc_nominals( slv3_system_t
748    int32 col;    int32 col;
749    FILE *fp = MIF(sys);    FILE *fp = MIF(sys);
750    
   CONSOLE_DEBUG("nominals...");  
   
751    if( sys->nominals.accurate ) return;    if( sys->nominals.accurate ) return;
752    fp = MIF(sys);    fp = MIF(sys);
753    col = sys->nominals.rng->low;    col = sys->nominals.rng->low;
# Line 815  static void calc_weights( slv3_system_t Line 807  static void calc_weights( slv3_system_t
807    mtx_coord_t nz;    mtx_coord_t nz;
808    real64 sum;    real64 sum;
809    
   CONSOLE_DEBUG("...");  
   
810    if( sys->weights.accurate )    if( sys->weights.accurate )
811      return;      return;
812    
# Line 891  static void jacobian_scaled(slv3_system_ Line 881  static void jacobian_scaled(slv3_system_
881  static void scale_variables( slv3_system_t sys)  static void scale_variables( slv3_system_t sys)
882  {  {
883    int32 col;    int32 col;
   CONSOLE_DEBUG("...");  
884    
885    if( sys->variables.accurate ) return;    if( sys->variables.accurate ) return;
886    
# Line 914  static void scale_residuals( slv3_system Line 903  static void scale_residuals( slv3_system
903   **/   **/
904  {  {
905    int32 row;    int32 row;
   CONSOLE_DEBUG("...");  
906    
907    if( sys->residuals.accurate ) return;    if( sys->residuals.accurate ) return;
908    
# Line 942  static void calc_relnoms(slv3_system_t s Line 930  static void calc_relnoms(slv3_system_t s
930    struct rel_relation *rel;    struct rel_relation *rel;
931    real64 *var_list;    real64 *var_list;
932    
   CONSOLE_DEBUG("...");  
933    
934    var_list = create_array(sys->cap,real64);    var_list = create_array(sys->cap,real64);
935    col = 0;    col = 0;
# Line 988  static real64 col_max_ratio(mtx_matrix_t Line 975  static real64 col_max_ratio(mtx_matrix_t
975    real64 num, denom, dummy;    real64 num, denom, dummy;
976    mtx_coord_t coord;    mtx_coord_t coord;
977    
   CONSOLE_DEBUG("...");  
   
978    max_ratio = 0;    max_ratio = 0;
979    for(coord.col = reg->col.low;    for(coord.col = reg->col.low;
980      coord.col <= reg->col.high; coord.col++) {      coord.col <= reg->col.high; coord.col++) {
# Line 1025  static real64 row_max_ratio(mtx_matrix_t Line 1010  static real64 row_max_ratio(mtx_matrix_t
1010    mtx_coord_t coord;    mtx_coord_t coord;
1011    max_ratio = 0;    max_ratio = 0;
1012    
   CONSOLE_DEBUG("...");  
   
1013    for(coord.row = reg->row.low;    for(coord.row = reg->row.low;
1014      coord.row <= reg->row.high; coord.row++) {      coord.row <= reg->row.high; coord.row++) {
1015      ratio = 0;      ratio = 0;
# Line 1064  static real64 calc_fourer_scale(mtx_matr Line 1047  static real64 calc_fourer_scale(mtx_matr
1047    real64 min, max, dummy;    real64 min, max, dummy;
1048    real64 scale;    real64 scale;
1049    
   CONSOLE_DEBUG("...");  
   
1050    if(option == 0){    if(option == 0){
1051      if((loc < reg.row.low) || (loc > reg.row.high)){      if((loc < reg.row.low) || (loc > reg.row.high)){
1052        return 1;        return 1;
# Line 1119  static void scale_J_iterative(slv3_syste Line 1100  static void scale_J_iterative(slv3_syste
1100    real64 *rowvec = sys->weights.vec;    real64 *rowvec = sys->weights.vec;
1101    real64 rowscale, colscale;    real64 rowscale, colscale;
1102    
   CONSOLE_DEBUG("...");  
   
1103    rho_col_old = col_max_ratio(&(sys->J.mtx),&(sys->J.reg));    rho_col_old = col_max_ratio(&(sys->J.mtx),&(sys->J.reg));
1104    rho_row_old = row_max_ratio(&(sys->J.mtx),&(sys->J.reg));    rho_row_old = row_max_ratio(&(sys->J.mtx),&(sys->J.reg));
1105    k = 0;    k = 0;
# Line 1166  static void scale_system( slv3_system_t Line 1145  static void scale_system( slv3_system_t
1145   *** Scale system dependent on interface parameters   *** Scale system dependent on interface parameters
1146   **/   **/
1147  {  {
   CONSOLE_DEBUG("...");  
   
1148    if(strcmp(SCALEOPT,"NONE") == 0){    if(strcmp(SCALEOPT,"NONE") == 0){
1149      if(sys->J.accurate == FALSE){      if(sys->J.accurate == FALSE){
1150        calc_nominals(sys);        calc_nominals(sys);
# Line 1240  static boolean calc_gradient(slv3_system Line 1217  static boolean calc_gradient(slv3_system
1217   * This entire function needs to be reimplemented with relman_diffs.   * This entire function needs to be reimplemented with relman_diffs.
1218   *   *
1219   */   */
   CONSOLE_DEBUG("...");  
   
1220    if( sys->gradient.accurate ) return TRUE;    if( sys->gradient.accurate ) return TRUE;
1221    
1222    calc_ok = TRUE;    calc_ok = TRUE;
# Line 1286  static void create_update( slv3_system_t Line 1261  static void create_update( slv3_system_t
1261  {  {
1262    struct hessian_data *update;    struct hessian_data *update;
1263    
   CONSOLE_DEBUG("...");  
   
1264    if( !OPTIMIZING(sys) )    if( !OPTIMIZING(sys) )
1265       return;       return;
1266    
# Line 1405  static int calc_pivots(slv3_system_t sys Line 1378  static int calc_pivots(slv3_system_t sys
1378    linsolqr_system_t lsys = sys->J.sys;    linsolqr_system_t lsys = sys->J.sys;
1379    FILE *fp = LIF(sys);    FILE *fp = LIF(sys);
1380    
   CONSOLE_DEBUG("...");  
   
1381    oldtiming = g_linsolqr_timing;    oldtiming = g_linsolqr_timing;
1382    g_linsolqr_timing =LINTIME;    g_linsolqr_timing =LINTIME;
1383    linsolqr_factor(lsys,sys->J.fm);          /* factor */    linsolqr_factor(lsys,sys->J.fm);          /* factor */
# Line 1507  static void calc_ZBZ(slv3_system_t sys) Line 1478  static void calc_ZBZ(slv3_system_t sys)
1478  {  {
1479     mtx_coord_t nz;     mtx_coord_t nz;
1480    
   CONSOLE_DEBUG("...");  
   
1481     if( sys->ZBZ.accurate ) return;     if( sys->ZBZ.accurate ) return;
1482    
1483     for( nz.row = 0; nz.row < sys->ZBZ.order; nz.row++ ) {     for( nz.row = 0; nz.row < sys->ZBZ.order; nz.row++ ) {
# Line 1582  static void calc_rhs(slv3_system_t sys, Line 1551  static void calc_rhs(slv3_system_t sys,
1551   ***  rhs.  Caller is responsible for initially zeroing the rhs vector.   ***  rhs.  Caller is responsible for initially zeroing the rhs vector.
1552   **/   **/
1553  {  {
   CONSOLE_DEBUG("...");  
   
1554    if( transpose ) {     /* vec is indexed by col */    if( transpose ) {     /* vec is indexed by col */
1555      int32 col;      int32 col;
1556      for( col=vec->rng->low; col<=vec->rng->high; col++ ) {      for( col=vec->rng->low; col<=vec->rng->high; col++ ) {
# Line 1604  static void calc_multipliers(slv3_system Line 1571  static void calc_multipliers(slv3_system
1571   ***  Computes the lagrange multipliers for the equality constraints.   ***  Computes the lagrange multipliers for the equality constraints.
1572   **/   **/
1573  {  {
   CONSOLE_DEBUG("...");  
   
1574         if( sys->multipliers.accurate )         if( sys->multipliers.accurate )
1575        return;        return;
1576    
# Line 3260  static SlvClientToken slv3_create(slv_sy Line 3225  static SlvClientToken slv3_create(slv_sy
3225  {  {
3226    slv3_system_t sys;    slv3_system_t sys;
3227    
   CONSOLE_DEBUG("...");  
   
3228    sys = (slv3_system_t)asccalloc(1, sizeof(struct slv3_system_structure) );    sys = (slv3_system_t)asccalloc(1, sizeof(struct slv3_system_structure) );
3229    if (sys==NULL) {    if (sys==NULL) {
3230      *statusindex = 1;      *statusindex = 1;
# Line 3300  static SlvClientToken slv3_create(slv_sy Line 3263  static SlvClientToken slv3_create(slv_sy
3263    /* we don't give a damn about the objective list or the pars or    /* we don't give a damn about the objective list or the pars or
3264     * bounds or extrels or any of the other crap.     * bounds or extrels or any of the other crap.
3265     */     */
   CONSOLE_DEBUG("...");  
3266    slv_check_var_initialization(server);    slv_check_var_initialization(server);
3267    *statusindex = 0;    *statusindex = 0;
3268    
# Line 3893  static void slv3_iterate(slv_system_t se Line 3855  static void slv3_iterate(slv_system_t se
3855    int               minor = 0,ds_status=0, rank_defect=0;    int               minor = 0,ds_status=0, rank_defect=0;
3856    double            time0;    double            time0;
3857    
   CONSOLE_DEBUG("QRSlv start");  
   
3858    sys = SLV3(asys);    sys = SLV3(asys);
3859    mif = MIF(sys);    mif = MIF(sys);
3860    lif = LIF(sys);    lif = LIF(sys);
# Line 4271  static void slv3_iterate(slv_system_t se Line 4231  static void slv3_iterate(slv_system_t se
4231  static void slv3_solve(slv_system_t server, SlvClientToken asys)  static void slv3_solve(slv_system_t server, SlvClientToken asys)
4232  {  {
4233    slv3_system_t sys;    slv3_system_t sys;
   CONSOLE_DEBUG("starting");  
4234    sys = SLV3(asys);    sys = SLV3(asys);
4235    if (server == NULL || sys==NULL) return;    if (server == NULL || sys==NULL) return;
4236    if (check_system(sys)) return;    if (check_system(sys)) return;
# Line 4302  static int slv3_destroy(slv_system_t ser Line 4261  static int slv3_destroy(slv_system_t ser
4261    
4262  int slv3_register(SlvFunctionsT *sft)  int slv3_register(SlvFunctionsT *sft)
4263  {  {
   CONSOLE_DEBUG("registering");  
   
4264    if (sft==NULL)  {    if (sft==NULL)  {
4265      FPRINTF(stderr,"slv3_register called with NULL pointer\n");      FPRINTF(stderr,"slv3_register called with NULL pointer\n");
4266      return 1;      return 1;

Legend:
Removed from v.153  
changed lines
  Added in v.154

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