/[ascend]/trunk/base/generic/integrator/idaanalyse.c
ViewVC logotype

Diff of /trunk/base/generic/integrator/idaanalyse.c

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

revision 1249 by johnpye, Sat Jan 27 04:14:58 2007 UTC revision 1250 by johnpye, Sat Jan 27 04:50:35 2007 UTC
# Line 222  static int integrator_ida_create_lists(I Line 222  static int integrator_ida_create_lists(I
222            
223      /* allocate space (more than enough) */      /* allocate space (more than enough) */
224      sys->y = ASC_NEW_ARRAY(struct var_variable *,sys->n_y);      sys->y = ASC_NEW_ARRAY(struct var_variable *,sys->n_y);
     sys->y_id = ASC_NEW_ARRAY_CLEAR(int,sys->n_y);  
225      sys->ydot = ASC_NEW_ARRAY_CLEAR(struct var_variable *,sys->n_y);      sys->ydot = ASC_NEW_ARRAY_CLEAR(struct var_variable *,sys->n_y);
226      j = 0;      sys->n_ydot = 0;
227    
228      CONSOLE_DEBUG("Passing through chains...");      CONSOLE_DEBUG("Passing through chains...");
229    
230      /* add the variables from the derivative chains */      /* create the lists y and ydot, ignoring 'bad' vars */
231      for(i=0; i<diffvars->nseqs; ++i){      for(i=0; i<diffvars->nseqs; ++i){
   
232          CONSOLE_DEBUG("i = %d",i);          CONSOLE_DEBUG("i = %d",i);
233                            
234          seq = diffvars->seqs[i];          seq = diffvars->seqs[i];
235          asc_assert(seq.n >= 1);          asc_assert(seq.n >= 1);
236          v = seq.vars[0];          v = seq.vars[0];
237            j = var_sindex(v);
         varname = var_make_name(sys->system, v);  
         CONSOLE_DEBUG("alg '%s'",varname);  
         ASC_FREE(varname);  
   
238    
239          if(!var_apply_filter(v,&integrator_ida_nonderiv)){          if(!var_apply_filter(v,&integrator_ida_nonderiv)){
240              continue;              continue;
241          }          }
242    
243          varname = var_make_name(sys->system, v);          sys->y[j] = v;
244          CONSOLE_DEBUG("alg '%s' is GOOD",varname);          VARMSG("'%s' is good non-deriv");
         ASC_FREE(varname);  
245    
246          if(seq.n > 1 && var_apply_filter(seq.vars[1],&integrator_ida_deriv)){          if(seq.n > 1 && var_apply_filter(seq.vars[1],&integrator_ida_deriv)){
247              asc_assert(var_sindex(seq.vars[1]) >= sys->n_y);              v = seq.vars[1];
248              asc_assert(var_sindex(seq.vars[1])-sys->n_y < sys->n_y);              asc_assert(var_sindex(v) >= sys->n_y);
249                asc_assert(var_sindex(v)-sys->n_y < sys->n_y);
250              varname = var_make_name(sys->system, seq.vars[1]);  
251              CONSOLE_DEBUG("diff '%s' IS GOOD",varname);              sys->ydot[j] = v;
252              ASC_FREE(varname);              sys->n_ydot++;
253                VARMSG("'%s' is good deriv");
             sys->y_id[var_sindex(seq.vars[1]) - sys->n_y] = j;  
             sys->ydot[j] = seq.vars[1];  
254          }else{          }else{
255              asc_assert(sys->ydot[j]==NULL);              asc_assert(sys->ydot[j]==NULL);
256          }          }
257        }
258    
259          sys->y[j] = v;      CONSOLE_DEBUG("Found %d good non-derivs",j);
260    
261        /* create the list y_id by looking at non-NULLs from ydot */
262        sys->y_id = ASC_NEW_ARRAY(int,sys->n_ydot);
263        for(i=0,j=0; i <  sys->n_y; ++i){
264            if(sys->ydot[i]==NULL)continue;
265            v = sys->ydot[i]; VARMSG("deriv '%s'...");
266            v = sys->y[i]; VARMSG("diff '%s'...");
267            sys->y_id[var_sindex(sys->ydot[i]) - sys->n_y] = i;
268          j++;          j++;
269      }      }
270        
271        asc_assert(j==sys->n_ydot);
272    
273      return 0;      return 0;
274  }  }
# Line 596  static int integrator_ida_check_diffinde Line 598  static int integrator_ida_check_diffinde
598          v = list[i];          v = list[i];
599          if(!var_apply_filter(v, &integrator_ida_nonderiv)){          if(!var_apply_filter(v, &integrator_ida_nonderiv)){
600              msg = "'%s' (in first n_y vars) fails non-deriv filter"; goto finish;              msg = "'%s' (in first n_y vars) fails non-deriv filter"; goto finish;
         }else if(v != sys->y[i]){  
             msg = "'%s' not matched in y vector"; goto finish;  
601          }else if(var_sindex(v) != i){          }else if(var_sindex(v) != i){
602              msg = "'%s' has wrong var_sindex"; goto finish;              msg = "'%s' has wrong var_sindex"; goto finish;
603            }else if(v != sys->y[i]){
604                msg = "'%s' not matched in y vector"; goto finish;
605          }          }
606          /* meets filter, matches in y vector, has correct var_sindex. */          /* meets filter, matches in y vector, has correct var_sindex. */
607      }      }
# Line 659  finish: Line 661  finish:
661  }  }
662    
663  /**  /**
664      @TODO provide a macro implementation of this      Given a derivative variable, return the index of its corresponding differential
665        variable in the y vector (and equivalently the var_sindex of the diff var)
666  */  */
667  int integrator_ida_diffindex(const IntegratorSystem *sys, const struct var_variable *deriv){  int integrator_ida_diffindex(const IntegratorSystem *sys, const struct var_variable *deriv){
     int diffindex;  
 #ifdef DIFFINDEX_DEBUG  
     asc_assert( var_deriv    (deriv));  
     asc_assert( var_active   (deriv));  
     asc_assert( var_incident (deriv));  
     asc_assert( var_svar     (deriv));  
   
     asc_assert(!var_fixed    (deriv));  
   
668      asc_assert(var_sindex(deriv) >= sys->n_y);      asc_assert(var_sindex(deriv) >= sys->n_y);
669      asc_assert(diffindex == var_sindex(sys->y[diffindex]));      asc_assert(var_sindex(deriv) < sys->n_y + sys->n_ydot);
670  #endif      return sys->y_id[var_sindex(deriv) - sys->n_y];
     asc_assert(var_sindex(deriv) >= sys->n_y);  
     diffindex = sys->y_id[var_sindex(deriv) - sys->n_y];  
   
     return diffindex;  
671  }  }
672    
673    
Legend:
Removed from v.1249  
changed lines
  Added in v.1250
john.pye@anu.edu.au
 ViewVC Help
Powered by ViewVC 1.1.22