Revert "*Bug #8099 fixed - the cumprod function did not simplify the result when...
[scilab.git] / scilab / modules / polynomials / sci_gateway / fortran / sci_f_psimp.f
1 c Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
2 c Copyright (C) ????-2008 - INRIA
3 c
4 c This file must be used under the terms of the CeCILL.
5 c This source file is licensed as described in the file COPYING, which
6 c you should have received as part of this distribution.  The terms
7 c are also available at
8 c http://www.cecill.info/licences/Licence_CeCILL_V2.1-en.txt
9       subroutine intpsimp(id)
10       INCLUDE 'stack.h'
11       integer iadr, sadr
12       integer id(nsiz)
13       integer vola,volb
14       logical chkvar
15 c
16       iadr(l)=l+l-1
17       sadr(l)=(l/2)+1
18 c
19
20 c     simp(num,den)
21 c
22       if(lhs.ne.2) then
23          call error(41)
24          return
25       endif
26 c
27       ila=iadr(lstk(top+1-rhs))
28
29       il=ila
30       if(istk(il).lt.0) il=iadr(istk(il+1))
31       if(istk(il).gt.2) then
32          fun=-1
33          call funnam(ids(1,pt+1),'simp',il)
34          return
35       endif
36 c
37       ilb=iadr(lstk(top))
38
39       il=ilb
40       if(istk(il).lt.0) il=iadr(istk(il+1))
41       if(istk(il).gt.2) then
42          fun=-1
43          call funnam(ids(1,pt+1),'simp',il)
44          return
45       endif
46 c
47       if(istk(ila).lt.0) then
48 c     move b and copy value of a in place of its reference
49          k=istk(ila+2)
50          vola=lstk(k+1)-lstk(k)
51          volb=lstk(top+1)-lstk(top)
52          err=lstk(top)+vola+volb-lstk(bot)
53          if(err.gt.0) then
54             call error(17)
55             return
56          endif
57          call unsfdcopy(volb,stk(lstk(top)),-1,stk(lstk(top-1)+vola),-1)
58          call unsfdcopy(vola,stk(lstk(k)),1,stk(lstk(top-1)),1)
59          lstk(top)=lstk(top-1)+vola
60          ilb=iadr(lstk(top))
61          lstk(top+1)=lstk(top)+volb
62       endif
63
64       if(istk(ilb).lt.0) then
65 c      copy value of b in place of its reference
66          k=istk(ilb+2)
67          volb=lstk(k+1)-lstk(k)
68          err=lstk(top)+volb-lstk(bot)
69          if(err.gt.0) then
70             call error(17)
71             return
72          endif
73          call unsfdcopy(volb,stk(lstk(k)),1,stk(lstk(top)),1)
74          lstk(top+1)=lstk(top)+volb
75       endif
76       lw=lstk(top+1)
77
78 c
79       if(istk(ila+3).ne.0.or.istk(ilb+3).ne.0) then
80          fun=-1
81          call funnam(ids(1,pt+1),'simp',ilb)
82          return
83       endif
84 c
85       mna=istk(ila+1)*istk(ila+2)
86       id(1)=0
87       if(istk(ila).eq.2) then
88          ida=ila+8
89          la=sadr(ida+mna+1)
90          call icopy(4,istk(ila+4),1,id,1)
91       else
92          la=sadr(ila+4)
93          ida=iadr(lw)
94          lw=sadr(ida+mna+1)
95          err=lw-lstk(bot)
96          if(err.gt.0) then
97             call error(17)
98             return
99          endif
100          do 74 i=1,mna+1
101             istk(ida+i-1)=i
102  74      continue
103       endif
104 c
105       mnb=istk(ilb+1)*istk(ilb+2)
106       if(istk(ilb).eq.2) then
107          idb=ilb+8
108          lb=sadr(idb+mnb+1)
109          if(id(1).eq.0) then
110             call icopy(4,istk(ilb+4),1,id,1)
111          else
112             if(.not.chkvar(id,istk(ilb+4))) then
113                call error(43)
114                return
115             endif
116          endif
117       else
118          lb=sadr(ilb+4)
119          idb=iadr(lw)
120          lw=sadr(idb+mna+1)
121          err=lw-lstk(bot)
122          if(err.gt.0) then
123             call error(17)
124             return
125          endif
126          do 75 i=1,mna+1
127             istk(idb+i-1)=i
128  75      continue
129       endif
130
131       if(mnb.ne.mna)then
132          call error(60)
133          return
134       endif
135 c
136 c     determine max of the degrees
137       na=0
138       nb=0
139       do 76 i=1,mna
140          na=max(na,istk(ida+i)-istk(ida-1+i))
141          nb=max(nb,istk(idb+i)-istk(idb-1+i))
142  76   continue
143 c     preserve adress of the beginning of a and b coefficients
144
145       lar=la
146       lbr=lb
147 c     allocate memory for intermediate results
148       law=lw
149       lbw=law+na+1
150       lw=lbw+nb+1
151
152 c     simplify
153       la1=la
154       lb1=lb
155       do 79 i=1,mna
156          na=istk(ida+i)-istk(ida-1+i)-1
157          nb=istk(idb+i)-istk(idb-1+i)-1
158          ierr=lstk(bot)-lw
159          call  dpsimp(stk(la),na,stk(lb),nb,stk(law),nnum,
160      $        stk(lbw),nden,stk(lw),ierr)
161          if(ierr.eq.1) then
162             call error(27)
163             return
164          elseif(ierr.eq.2) then
165             call msgs(43,i)
166          endif
167 c     .  copy overwrite initial polynomials with simplified ones
168          call dcopy(nnum,stk(law),1,stk(la1),1)
169          call dcopy(nden,stk(lbw),1,stk(lb1),1)
170
171          la=la+na+1
172          lb=lb+nb+1
173          la1=la1+nnum
174          lb1=lb1+nden
175          istk(ida-1+i)=nnum
176          istk(idb-1+i)=nden
177  79   continue
178 c
179 c     form vector of pointers from vector of degrees+1
180       ma=1
181       mb=1
182       do 80 i=1,mna+1
183          na=istk(ida-1+i)
184          nb=istk(idb-1+i)
185          istk(ida-1+i)=ma
186          istk(idb-1+i)=mb
187          ma=ma+na
188          mb=mb+nb
189  80   continue
190 c
191 c     compute position of the a and b simplified in the result
192       lstk(top)=lar+istk(ida+mna)-1
193       il=iadr(lstk(top))
194 c
195 c     put new b variable in place
196       if(istk(ilb).eq.2) then
197 c     b matrice de polynome
198          l=sadr(il+9+mna)
199 c     .  move b data up
200          call icopy(9+mna,istk(ilb),1,istk(il),1)
201          call unsfdcopy(istk(il+8+mna),stk(lbr),1,stk(l),1)
202          l=l+istk(il+8+mna)-1
203       else
204 c     b matrice de scalaires
205          call icopy(4,istk(ilb),1,istk(il),1)
206          l=sadr(il+4)
207          call unsfdcopy(mna,stk(lbr),1,stk(l),1)
208          l=l+mna
209       endif
210       lstk(top+1)=l
211       end
212
213 c                       =======================================