Post new topic Reply to topic  [ 3 posts ] 
Author Message
 Post subject: Singular and Maxima gives different solutions
PostPosted: Thu Jun 30, 2011 5:38 am 
Hi,

I've a system of linear-quadratic equations to solve. I tried Maxima and it gives me the correct solution. Singular gives me a complex solution, which is not accurate. Below is my code and outputs in Maxima and Singular.

Singular:
Code:
number rax, ray, raz = 8.076651, -5.766802, 3.363392;
number rbx, rby, rbz = 6.532685, -6.184719, 4.475940;
number rcx, rcy, rcz = 7.850585, -4.163334, 4.243473;
number dax, day, daz = 10.079651, -2.603802, 0.338392;
number dbx, dby, dbz = 8.535685, -3.021719, 1.450940;
number dcx, dcy, dcz = 9.853585, -1.000334, 1.218473;
number la = 5.189281;
number lb = 2.145810;
number lc = 2.236297;
number ta = 4.788316;
number tb = 0.935814;
number tc = 4.252196;
number fab = 0.093240;
number fac = 0.512108;
number fbc = 1.061882;

poly f1 = raz*maz+ray*may+rax*max - la;
poly f2 = daz*maz+day*may+dax*max - ta;
poly f3 = rbz*mbz+rby*mby+rbx*mbx - lb;
poly f4 = dbz*mbz+dby*mby+dbx*mbx - tb;
poly f5 = rcz*mcz+rcy*mcy+rcx*mcx - lc;
poly f6 = dcz*mcz+dcy*mcy+dcx*mcx - tc;
poly f7 = -raz*rbz*maz*mbz-raz*rby*may*mbz-raz*rbx*max*mbz+la*raz*mbz
                         -ray*rbz*maz*mby-ray*rby*may*mby-ray*rbx*max*mby
                         +la*ray*mby-rax*rbz*maz*mbx-rax*rby*may*mbx
                         -rax*rbx*max*mbx+la*rax*mbx+lb*rbz*maz+lb*rby*may
                         +lb*rbx*max-la*lb+fab
           ;
poly f8 = -raz*rcz*maz*mcz-raz*rcy*may*mcz-raz*rcx*max*mcz+la*raz*mcz
                         -ray*rcz*maz*mcy-ray*rcy*may*mcy-ray*rcx*max*mcy
                         +la*ray*mcy-rax*rcz*maz*mcx-rax*rcy*may*mcx
                         -rax*rcx*max*mcx+la*rax*mcx+lc*rcz*maz+lc*rcy*may
                         +lc*rcx*max-la*lc+fac
           ;
poly f9 = -rbz*rcz*mbz*mcz-rbz*rcy*mby*mcz-rbz*rcx*mbx*mcz+lb*rbz*mcz
                         -rby*rcz*mbz*mcy-rby*rcy*mby*mcy-rby*rcx*mbx*mcy
                         +lb*rby*mcy-rbx*rcz*mbz*mcx-rbx*rcy*mby*mcx
                         -rbx*rcx*mbx*mcx+lb*rbx*mcx+lc*rcz*mbz+lc*rcy*mby
                         +lc*rcx*mbx-lb*lc+fbc
           ;

ideal I = f1, f2, f3, f4, f5, f6, f7, f8, f9;
ideal G = groebner(I);
dim(G);
vdim(G);
G;

def s = solve(I, 12);


Maxima:
Code:
ra : [8.076651, -5.766802, 3.363392];
rb : [6.532685, -6.184719, 4.475940];
rc : [7.850585, -4.163334, 4.243473];
da : [10.079651, -2.603802, 0.338392];
db : [8.535685, -3.021719, 1.450940];
dc : [9.853585, -1.000334, 1.218473];
la : 5.189281;
lb : 2.145810;
lc : 2.236297;
ta : 4.788316;
tb : 0.935814;
tc : 4.252196;
fab : 0.093240;
fac : 0.512108;
fbc : 1.061882;


ma : [max, may, maz];
mb : [mbx, mby, mbz];
mc : [mcx, mcy, mcz];

Hab : expand((ma . rb - la) * (lb - mb . ra));
Hac : expand((ma . rc - la) * (lc - mc . ra));
Hbc : expand((mb . rc - lb) * (lc - mc . rb));

eqns : [ma . ra = la, ma . da = ta, mb . rb = lb, mb . db = tb, mc . rc = lc, mc . dc = tc, Hab + fab = 0, Hac + fac = 0, Hbc + fbc = 0];
eqns; 
sols : solve(eqns, [max, may, maz, mbx, mby, mbz, mcx, mcy, mcz]);

solsf : bfloat(sols);


The ground truth solution is
Code:
ma: 0.397328 -0.288674, 0.093798
mb: -0.007418 -0.283056, 0.099119
mc: 0.459667 0.599035, 0.264318


Output of Singular:
Code:
[1]:
   [1]:
      (0.396501063225-i*0.0585142482085)
   [2]:
      (-0.292459889853-i*0.267965150227)
   [3]:
      (0.0892914674573-i*0.318934797814)
   [4]:
      (-0.0228603746306-i*0.0180139060626)
   [5]:
      (-0.380524257884-i*0.113697415565)
   [6]:
      (-0.0130221380534-i*0.130812174086)
   [7]:
      (0.467285747848+i*0.00508868957766)
   [8]:
      (0.302313482452-i*0.198182075022)
   [9]:
      (-0.0408951640129-i*0.20385360867)
[2]:
   [1]:
      (0.396501063225+i*0.0585142482085)
   [2]:
      (-0.292459889853+i*0.267965150227)
   [3]:
      (0.0892914674573+i*0.318934797814)
   [4]:
      (-0.0228603746306+i*0.0180139060626)
   [5]:
      (-0.380524257884+i*0.113697415565)
   [6]:
      (-0.0130221380534+i*0.130812174086)
   [7]:
      (0.467285747848-i*0.00508868957766)
   [8]:
      (0.302313482452+i*0.198182075022)
   [9]:
      (-0.0408951640129+i*0.20385360867)


Output of Maxima:
Code:
[[max = 3.463461b-1,may = -5.2214398b-1,maz = -1.8408099b-1,
          mbx = -2.3113602b-2,mby = -3.8212235b-1,mbz = -1.4860691b-2,
          mcx = 4.6593653b-1,mcy = 3.5486107b-1,mcz = 1.3156148b-2],

         [max = 3.973295b-1,may = -2.8866626b-1,maz = 9.3806532b-2,
          mbx = -7.4168349b-3,mby = -2.8304991b-1,mbz = 9.9125014b-2,
          mcx = 4.5966645b-1,mcy = 5.9905295b-1,mcz = 2.6433629b-1]]


Did I do something wrong in Singular? Thanks for your help.


Report this post
Top
  
Reply with quote  
 Post subject: Re: Singular and Maxima gives different solutions
PostPosted: Fri Jul 01, 2011 9:23 am 

Joined: Thu Jun 30, 2011 6:23 am
Posts: 1
Previously, the declared ring is
Code:
ring r1 = real, (max, may, maz, mbx, mby, mbz, mcx, mcy, mcz), lp;


It seems like there are something wrong when the ring is declared with real type. I switch to the following ring
Code:
ring r1 = 0, (max, may, maz, mbx, mby, mbz, mcx, mcy, mcz), lp;

and convert all floating point numbers to rational form. The solver gives solution similar to Maxima now.

The full working code is as follows:
Code:
ring r1 = 0, (max, may, maz, mbx, mby, mbz, mcx, mcy, mcz), lp;

number rax, ray, raz = 2318/287, -1459/253, 2499/743;
number rbx, rby, rbz = 2698/413, -971/157, 837/187;
number rcx, rcy, rcz = 8722/1111, -1249/300, 6501/1532;
number dax, day, daz = 1139/113, -1781/684, 223/659;
number dbx, dby, dbz = 10166/1191, -3200/1059, 695/479;
number dcx, dcy, dcz = 3163/321, -2994/2993, 686/563;
number la = 1453/280;
number lb = 1972/919;
number lc = 1836/821;
number ta = 6311/1318;
number tb = 1283/1371;
number tc = 1450/341;
number fab = 491/5266;
number fac = 571/1115;
number fbc = 2042/1923;

poly f1 = raz*maz+ray*may+rax*max - la;
poly f2 = daz*maz+day*may+dax*max - ta;
poly f3 = rbz*mbz+rby*mby+rbx*mbx - lb;
poly f4 = dbz*mbz+dby*mby+dbx*mbx - tb;
poly f5 = rcz*mcz+rcy*mcy+rcx*mcx - lc;
poly f6 = dcz*mcz+dcy*mcy+dcx*mcx - tc;
poly f7 = -raz*rbz*maz*mbz-raz*rby*may*mbz-raz*rbx*max*mbz+la*raz*mbz
                         -ray*rbz*maz*mby-ray*rby*may*mby-ray*rbx*max*mby
                         +la*ray*mby-rax*rbz*maz*mbx-rax*rby*may*mbx
                         -rax*rbx*max*mbx+la*rax*mbx+lb*rbz*maz+lb*rby*may
                         +lb*rbx*max-la*lb+fab
           ;
poly f8 = -raz*rcz*maz*mcz-raz*rcy*may*mcz-raz*rcx*max*mcz+la*raz*mcz
                         -ray*rcz*maz*mcy-ray*rcy*may*mcy-ray*rcx*max*mcy
                         +la*ray*mcy-rax*rcz*maz*mcx-rax*rcy*may*mcx
                         -rax*rcx*max*mcx+la*rax*mcx+lc*rcz*maz+lc*rcy*may
                         +lc*rcx*max-la*lc+fac
           ;
poly f9 = -rbz*rcz*mbz*mcz-rbz*rcy*mby*mcz-rbz*rcx*mbx*mcz+lb*rbz*mcz
                         -rby*rcz*mbz*mcy-rby*rcy*mby*mcy-rby*rcx*mbx*mcy
                         +lb*rby*mcy-rbx*rcz*mbz*mcx-rbx*rcy*mby*mcx
                         -rbx*rcx*mbx*mcx+lb*rbx*mcx+lc*rcz*mbz+lc*rcy*mby
                         +lc*rcx*mbx-lb*lc+fbc
           ;

ideal I = f1, f2, f3, f4, f5, f6, f7, f8, f9;
def s = solve(I, 12);


Could anybody help explain why I cannot define a ring with real data type to use with the solve function? Thanks.


Report this post
Top
 Profile  
Reply with quote  
 Post subject: Re: Singular and Maxima gives different solutions
PostPosted: Mon Jul 04, 2011 11:14 am 

Joined: Wed May 25, 2005 4:16 pm
Posts: 275
As soon as one calls groebner or std in a ring with inexact coefficients
(like real or complex),
most probably the result will be unusable.


Report this post
Top
 Profile  
Reply with quote  
Display posts from previous:  Sort by  
Post new topic Reply to topic  [ 3 posts ] 

You can post new topics in this forum
You can reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot post attachments in this forum

It is currently Fri May 13, 2022 10:56 am
cron
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group