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

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

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

revision 1098 by johnpye, Wed Jan 10 05:41:38 2007 UTC revision 1099 by johnpye, Wed Jan 10 06:26:02 2007 UTC
# Line 282  void integrator_lsode_create(IntegratorS Line 282  void integrator_lsode_create(IntegratorS
282      d->ydot_vars=NULL;      d->ydot_vars=NULL;
283      d->rlist=NULL;      d->rlist=NULL;
284      d->dydx_dx=NULL;      d->dydx_dx=NULL;
     CONSOLE_DEBUG("Created LSODE enginedata at %p",d);  
285      blsys->enginedata=(void*)d;      blsys->enginedata=(void*)d;
286      integrator_lsode_params_default(blsys);      integrator_lsode_params_default(blsys);
287    
# Line 313  void integrator_lsode_free(void *engined Line 312  void integrator_lsode_free(void *engined
312      if(d.dydx_dx!=NULL){      if(d.dydx_dx!=NULL){
313          lsode_densematrix_destroy(d.dydx_dx, d.n_eqns);          lsode_densematrix_destroy(d.dydx_dx, d.n_eqns);
314          d.dydx_dx =  NULL;          d.dydx_dx =  NULL;
         CONSOLE_DEBUG("Cleared dydx_dx");  
315      }      }
316    
317      d.n_eqns = 0L;      d.n_eqns = 0L;
# Line 348  int integrator_lsode_params_default(Inte Line 346  int integrator_lsode_params_default(Inte
346      slv_destroy_parms(p);      slv_destroy_parms(p);
347    
348      if(p->parms==NULL){      if(p->parms==NULL){
         CONSOLE_DEBUG("params NULL");  
349          p->parms = ASC_NEW_ARRAY(struct slv_parameter, LSODE_PARAMS_SIZE);          p->parms = ASC_NEW_ARRAY(struct slv_parameter, LSODE_PARAMS_SIZE);
350          if(p->parms==NULL)return -1;          if(p->parms==NULL)return -1;
351          p->dynamic_parms = 1;          p->dynamic_parms = 1;
352      }else{      }else{
353          asc_assert(p->num_parms == LSODE_PARAMS_SIZE);          asc_assert(p->num_parms == LSODE_PARAMS_SIZE);
         CONSOLE_DEBUG("reusing parm memory");  
354      }      }
355    
356      /* reset the number of parameters to zero so that we can check it at the end */      /* reset the number of parameters to zero so that we can check it at the end */
# Line 403  int integrator_lsode_params_default(Inte Line 399  int integrator_lsode_params_default(Inte
399      );      );
400    
401      asc_assert(p->num_parms == LSODE_PARAMS_SIZE);      asc_assert(p->num_parms == LSODE_PARAMS_SIZE);
   
     CONSOLE_DEBUG("Created %d params", p->num_parms);  
   
402      return 0;      return 0;
403  }    }  
404    
# Line 669  int integrator_lsode_derivatives(Integra Line 662  int integrator_lsode_derivatives(Integra
662    
663    asc_assert(blsys!=NULL);    asc_assert(blsys!=NULL);
664    enginedata = (IntegratorLsodeData *)blsys->enginedata;    enginedata = (IntegratorLsodeData *)blsys->enginedata;
   CONSOLE_DEBUG("blsys at %p",blsys);  
   CONSOLE_DEBUG("Enginedata at %p",enginedata);  
665    asc_assert(enginedata!=NULL);    asc_assert(enginedata!=NULL);
666    asc_assert(enginedata->dydx_dx!=NULL);    asc_assert(enginedata->dydx_dx!=NULL);
667    asc_assert(enginedata->input_indices!=NULL);    asc_assert(enginedata->input_indices!=NULL);
# Line 810  static void LSODE_FEX( int *n_eq ,double Line 801  static void LSODE_FEX( int *n_eq ,double
801  #endif  #endif
802      lsodedata->stop = 1;      lsodedata->stop = 1;
803      lsodedata->status = lsode_nok;      lsodedata->status = lsode_nok;
     ERROR_REPORTER_HERE(ASC_PROG_NOTE,"lsodedata->status = %d",lsodedata->status);  
804    }else{    }else{
805      lsodedata->status = lsode_ok;      lsodedata->status = lsode_ok;
806      /* ERROR_REPORTER_HERE(ASC_PROG_NOTE,"lsodedata->status = %d",lsodedata->status); */      /* ERROR_REPORTER_HERE(ASC_PROG_NOTE,"lsodedata->status = %d",lsodedata->status); */
# Line 851  static void LSODE_JEX(int *neq ,double * Line 841  static void LSODE_JEX(int *neq ,double *
841     * Make the real call.     * Make the real call.
842     */     */
843    
   CONSOLE_DEBUG("blsys at %p",l_lsode_blsys);  
   CONSOLE_DEBUG("Enginedata at %p",lsodedata);  
   
844    nok = integrator_lsode_derivatives(l_lsode_blsys    nok = integrator_lsode_derivatives(l_lsode_blsys
845          , *neq          , *neq
846          , *nrpd          , *nrpd
# Line 867  static void LSODE_JEX(int *neq ,double * Line 854  static void LSODE_JEX(int *neq ,double *
854      return;      return;
855    }else{    }else{
856      lsodedata->status = lsode_ok;      lsodedata->status = lsode_ok;
     /* ERROR_REPORTER_HERE(ASC_PROG_NOTE,"lsodedata->status = %d",lsodedata->status); */  
857      lsodedata->lastcall = lsode_derivative;      lsodedata->lastcall = lsode_derivative;
858    }    }
859    /*    /*
# Line 923  int integrator_lsode_solve(IntegratorSys Line 909  int integrator_lsode_solve(IntegratorSys
909      d->input_indices = ASC_NEW_ARRAY_CLEAR(int, d->n_eqns);      d->input_indices = ASC_NEW_ARRAY_CLEAR(int, d->n_eqns);
910      d->output_indices = ASC_NEW_ARRAY_CLEAR(int, d->n_eqns);      d->output_indices = ASC_NEW_ARRAY_CLEAR(int, d->n_eqns);
911      d->dydx_dx = lsode_densematrix_create(d->n_eqns,d->n_eqns);      d->dydx_dx = lsode_densematrix_create(d->n_eqns,d->n_eqns);
     CONSOLE_DEBUG("Allocated dydx_dx into enginedata at %p",d);  
912    
913      d->y_vars = ASC_NEW_ARRAY(struct var_variable *,d->n_eqns+1);      d->y_vars = ASC_NEW_ARRAY(struct var_variable *,d->n_eqns+1);
914      d->ydot_vars = ASC_NEW_ARRAY(struct var_variable *, d->n_eqns+1);      d->ydot_vars = ASC_NEW_ARRAY(struct var_variable *, d->n_eqns+1);
# Line 1063  int integrator_lsode_solve(IntegratorSys Line 1048  int integrator_lsode_solve(IntegratorSys
1048    
1049  # ifndef NO_SIGNAL_TRAPS  # ifndef NO_SIGNAL_TRAPS
1050      }else{      }else{
1051        FPRINTF(stderr,        ERROR_REPORTER_HERE(ASC_PROG_ERR,"Integration terminated due to float error in LSODE call.");
        "Integration terminated due to float error in LSODE call.\n");  
1052        lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);        lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);
1053        d->status = lsode_ok;     /* clean up before we go */        d->status = lsode_ok;     /* clean up before we go */
       ERROR_REPORTER_HERE(ASC_PROG_NOTE,"lsodedata->status = %d",d->status);  
1054        d->lastcall = lsode_none;        d->lastcall = lsode_none;
1055        return 6;        return 6;
1056      }      }
# Line 1101  int integrator_lsode_solve(IntegratorSys Line 1084  int integrator_lsode_solve(IntegratorSys
1084        ERROR_REPORTER_HERE(ASC_PROG_ERR,"Integration terminated due to an error in derivative computations.");        ERROR_REPORTER_HERE(ASC_PROG_ERR,"Integration terminated due to an error in derivative computations.");
1085        lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);        lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);
1086        d->status = lsode_ok;     /* clean up before we go */        d->status = lsode_ok;     /* clean up before we go */
       ERROR_REPORTER_HERE(ASC_PROG_NOTE,"d->status = %d",d->status);  
1087        d->lastcall = lsode_none;        d->lastcall = lsode_none;
1088        integrator_output_close(blsys);        integrator_output_close(blsys);
1089        return 8;        return 8;
# Line 1119  int integrator_lsode_solve(IntegratorSys Line 1101  int integrator_lsode_solve(IntegratorSys
1101          ERROR_REPORTER_HERE(ASC_USER_ERROR,"Integration cancelled");          ERROR_REPORTER_HERE(ASC_USER_ERROR,"Integration cancelled");
1102          lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);          lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);
1103          d->status = lsode_ok;          d->status = lsode_ok;
         ERROR_REPORTER_HERE(ASC_PROG_NOTE,"d->status = %d",d->status);  
1104          d->lastcall = lsode_none;          d->lastcall = lsode_none;
1105          integrator_output_close(blsys);          integrator_output_close(blsys);
1106          return 9;          return 9;
# Line 1150  int integrator_lsode_solve(IntegratorSys Line 1131  int integrator_lsode_solve(IntegratorSys
1131          ERROR_REPORTER_HERE(ASC_PROG_ERR,"Integration terminated due to float error in LSODE FEX call.");          ERROR_REPORTER_HERE(ASC_PROG_ERR,"Integration terminated due to float error in LSODE FEX call.");
1132          lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);          lsode_free_mem(y,reltol,abtol,rwork,iwork,obs,dydx);
1133          d->status = lsode_ok;               /* clean up before we go */          d->status = lsode_ok;               /* clean up before we go */
         ERROR_REPORTER_HERE(ASC_PROG_NOTE,"d->status = %d",d->status);  
1134          d->lastcall = lsode_none;          d->lastcall = lsode_none;
1135          integrator_output_close(blsys);          integrator_output_close(blsys);
1136          return 10;          return 10;
# Line 1174  int integrator_lsode_solve(IntegratorSys Line 1154  int integrator_lsode_solve(IntegratorSys
1154     */     */
1155    
1156    d->status = lsode_ok;    d->status = lsode_ok;
   ERROR_REPORTER_HERE(ASC_PROG_NOTE,"d->status = %d",d->status);  
1157    d->lastcall = lsode_none;    d->lastcall = lsode_none;
1158    
1159    integrator_output_close(blsys);    integrator_output_close(blsys);
Legend:
Removed from v.1098  
changed lines
  Added in v.1099
john.pye@anu.edu.au
 ViewVC Help
Powered by ViewVC 1.1.22