Opened 14 years ago

Closed 14 years ago

#211 closed bug (fixed)

error in proc actionIsProper of ainvar.lib

Reported by: Imad Jaradat (imad_jar@… Owned by: pfister@…
Priority: minor Milestone: 3-1-1
Component: singular-libs Version: 3-1-1
Keywords: ainvar.lib Cc:

Description

Using the procedure actionIsProper of the library ainvar.lib of Singular 3.1.0 we identified a problem: during the computation of the exponential, the formula for the calculation of the denominator works correctly only for the first three terms (iterations). Here it is the code: / proc actionIsProper(matrix m) . . . {

poly inv,delta,tee,j; ideal id=imap(bsr,id); matrix @m[size(id)+1][1]; @m=imap(bsr,m),0;

computes the exp(@t*m)(var(i)) for all i for(i=1;i<=nvars(basering)-1;i++) {

inv=var(i); delta=derivate(@m,inv); j=1; tee=@t; while(delta!=0) {

inv=inv+1/j*delta*tee; j=j*(j+1); Problem Line (denominator) tee=tee*@t; delta=derivate(@m,delta);

} id=id+ideal(inv);

}

. . . /

One way to fix it is by introducing a new variable (auxv) as seen in the following modified version of your code: / proc actionIsProper(matrix m) . . . {

. . . @m=imap(bsr,m),0; int auxv; /Added this line (1/3) to replace the line j=j+1

computes the exp(@t*m)(var(i)) for all i

for(i=1;i<=nvars(basering)-1;i++) {

inv=var(i); delta=derivate(@m,inv); j=1; auxv=0; /Added this line (2/3) to replace j=j+1 tee=@t; while(delta!=0) {

cc=cc+1; inv=inv+1/j*delta*tee; j=j*(auxv+1); /Added this line (3 of 3) to replace j=j+1

below

/ j=j+1 /Commented this line (1/1).

tee=tee*@t; delta=derivate(@m,delta);

} id=id+ideal(inv);

} i=inSubring(@t,id)[1]; setring(bsr); return(i);

}

/

Change History (1)

comment:1 Changed 14 years ago by seelisch

Resolution: fixed
Status: newclosed

fixed like this:

... int auxv;

computes the exp(@t*m)(var(i)) for all i for(i=1;i<=nvars(basering)-1;i++) {

inv=var(i); delta=derivate(@m,inv); j=1; auxv=1; tee=@t; while(delta!=0) {

inv=inv+1/j*delta*tee; auxv=auxv+1; j=j*auxv; tee=tee*@t; delta=derivate(@m,delta);

} id=id+ideal(inv);

}

...

Note: See TracTickets for help on using tickets.