Re: [Scilab-users] lack of simplification for rational

[Tim Wescott]

I think this is a case of "you can't get there from here", certainly if
you're using Scilab, and absolutely if you're getting your polynomials
by measuring the coefficients, or by, at some point, going from symbolic
representations to some sort of floating-point representation.

The underlying problem is that reducing the fractions the way that you
want to requires finding the roots of the polynomials involved, but the
task of finding the roots of a polynomial gets more and more numerically
unstable as the polynomial order goes up.

You might be able to do this reliably if you use a symbolic solver as,
for instance, Maxima, and if the order of your polynomials doesn't
exceed whatever that solver is capable of.  I'm not sure what Maxima can
do, but fifth-order is the maximum for which an algorithm exists to
factor any polynomial symbolically, and my understanding is that the
algorithm is horrid.  I would be pleasantly surprised if Maxima could
reliable do it with 3rd-order polynomials, and astonished if it could do
fourth- and fifth-order ones.

You may be reduced to doing this job by hand: factor the numerator and
denominator polynomials into their roots, then look for the roots that
come close to matching, then use that and knowledge of whatever it is
that you're modeling to decide which ones "really" cancel.

Just as a note: I'm not sure what physical processes you're trying to
represent with ratios of polynomials, but if you're working with
transfer functions of linear systems, having identical roots in
numerator and denominator can either arise from a mathematical artifact
in your computations, or it can be indicative of a pole-zero
cancellation in your physical system.  The former doesn't matter (other
than as an irritant), but a pole-zero cancellation can cause problems if
the pole in question is unstable or if it decays more slowly than it
receives excitation in the real world.  I don't know if similar cautions
apply in other realms with canceling roots in the numerator and
denominator -- but given my experience with control systems, I'd at
least check the physics behind any such cancellations.

On Wed, 2016-11-23 at 19:37 +0100, philippe wrote:
> Hi,
> 
> I  have a problem  to simplify automatically rational with an imaginary 
> part. In my case  I only manipulate polynomials with real coefficients 
> but some imaginary part arise due to rounding errors. See this example :
> 
> -->X=poly(0,'x');
> 
> -->A=(X-1)^2;
> 
> -->B=(X-1)*(X-2);
> 
> -->R=A/B// = (X-1)/(X-2)
>   R  =
> 
>- 1 + x
>  -
>- 2 + x
> 
> -->A=clean(A+%eps*(1+%i))//  add a rounding error
>   A  =
> 
> Partie réelle
> 
> 
>2
>  1 - 2x + x
> Partie imaginaire
> 
> 
>  0
> 
> -->R=A/B// = no simplifications
>   R  =
> 
>  2
>  1 - 2x +  1x
> --
>2
>  2 - 3x + x
> 
> 
> I can't simplify R, even using simp(R), I've found some work around 
> applying clean  to numerator/denominator coefficients  TWICE (I don't 
> known why twice?!?!?)  :
> 
> -->R=poly(clean(coeff(numer(R))),'x')/poly(clean(coeff(denom(R))),'x')
>   R  =
> 
>   2
>- 2 + x + x
>  -
>   2
>- 6 + x + x
> 
> -->R=poly(clean(coeff(numer(R))),'x')/poly(clean(coeff(denom(R))),'x')
>   R  =
> 
>  2 + x
>  -
>  6 + x
> 
> but of course this doesn't work  for polynomials with complex 
> coefficients like :
> 
> 
> -->(X+%i)/(X^2+1) // = 1/(X-%i)
>   ans  =
> 
> 
> i +  1x
> 
>   2
>  1 + x
> 
> I would like to find a solution for all polynomials with real/complex 
> coefficients , any idea ?
> 
> Philippe
> 
> ___
> users mailing list
> users@lists.scilab.org
> http://lists.scilab.org/mailman/listinfo/users

-- 

Tim Wescott
www.wescottdesign.com
Control & Communications systems, circuit & software design.
Phone: 503.631.7815
Cell:  503.349.8432


___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[Samuel Gougeon]

Hello Philippe,

Le 23/11/2016 19:37, philippe a écrit :

Hi,

I  have a problem  to simplify automatically rational with an 
imaginary part. In my case  I only manipulate polynomials with real 
coefficients but some imaginary part arise due to rounding errors. See 
this example :

.../...
I would like to find a solution for all polynomials with real/complex 
coefficients , any idea ?

.
For all rationals, no (see Tim's post). But for some of them, yes: you 
should use Scilab 6.0-b.

Some issues of this kind were reported and fixed quite recently:
http://bugzilla.scilab.org/8493 https://codereview.scilab.org/#/c/18420
http://bugzilla.scilab.org/13893

X = poly(0,'x');
A = (X-1)^2;
B = (X-1)*(X-2);
A = A+%i*0
A/B

-->A/B   // Scilab 5.5.2
 ans  =
2
1 - 2x +  1x
   --
  2
2 - 3x + x

--> A/B  // Scilab 6
 ans  =
   -1 + x
   ---
   -2 + x

HTH
Samuel

___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[Serge Steer]

First of all if you compute the rational fraction numerator or 
denominator from their roots (and if they are complex conjugate) the 
result is a complex polynomial even if the imaginary part is not shown.


Second: the simplification of numerical rational fraction is a very 
difficult problem. As Tim said  the roots computed by numerical 
algorithm cannot be use for that purpose.

Just try:
p=poly(ones(1,n),"x","c");
r=roots(p)
plot(real(r)-1,imag(r),"+")
You will see the roots on a circle which radius is of order %eps^(1/n) !
It is not a bug just a consequence of the numerical computations.

The simplication algorithm present in Scilab tryes to find reasonable 
simplification based on the bezout equation error. But it has only be 
written for the real case.


The bug fixes Samuel pointed to you  allow simplification for complex 
polynomials having zero  imaginary parts.

Serge
Le 23/11/2016 à 19:37, philippe a écrit :

Hi,

I  have a problem  to simplify automatically rational with an 
imaginary part. In my case  I only manipulate polynomials with real 
coefficients but some imaginary part arise due to rounding errors. See 
this example :


-->X=poly(0,'x');

-->A=(X-1)^2;

-->B=(X-1)*(X-2);

-->R=A/B// = (X-1)/(X-2)
 R  =

  - 1 + x
-
  - 2 + x

-->A=clean(A+%eps*(1+%i))//  add a rounding error
 A  =

Partie réelle


  2
1 - 2x + x
Partie imaginaire


0

-->R=A/B// = no simplifications
 R  =

2
1 - 2x +  1x
   --
  2
2 - 3x + x


I can't simplify R, even using simp(R), I've found some work around 
applying clean  to numerator/denominator coefficients TWICE (I don't 
known why twice?!?!?)  :


-->R=poly(clean(coeff(numer(R))),'x')/poly(clean(coeff(denom(R))),'x')
 R  =

 2
  - 2 + x + x
-
 2
  - 6 + x + x

-->R=poly(clean(coeff(numer(R))),'x')/poly(clean(coeff(denom(R))),'x')
 R  =

2 + x
-
6 + x

but of course this doesn't work  for polynomials with complex 
coefficients like :



-->(X+%i)/(X^2+1) // = 1/(X-%i)
 ans  =


   i +  1x
   
 2
1 + x

I would like to find a solution for all polynomials with real/complex 
coefficients , any idea ?


Philippe

___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users



___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[Tim Wescott]

Correction:

p=poly(ones(1,n),"x","r");  // Note 'r', not 'c'
r=roots(p)
plot(real(r)-1,imag(r),"+")

And my yes, it does illustrate the problem quite well.

On Wed, 2016-11-23 at 21:20 +0100, Serge Steer wrote:
> p=poly(ones(1,n),"x","c");
> r=roots(p)
> plot(real(r)-1,imag(r),"+")

-- 

Tim Wescott
www.wescottdesign.com
Control & Communications systems, circuit & software design.
Phone: 503.631.7815
Cell:  503.349.8432


___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[philippe]

Le 23/11/2016 à 20:08, Tim Wescott a écrit :

[...]

You might be able to do this reliably if you use a symbolic solver as,
for instance, Maxima, and if the order of your polynomials doesn't
exceed whatever that solver is capable of.


I try to do this with scilab only because I can't no longer use maxima 
via the scimax toolboxe (since I've switch to ubuntu 16.04 I couldn't 
compile scimax anymore)



[...]
Just as a note: I'm not sure what physical processes you're trying to
represent with ratios of polynomials,


I'm working on partial fraction decomposition  algorithms.

thanks

Philippe


___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[Samuel Gougeon]

Le 23/11/2016 23:19, philippe a écrit :

.../...
Le 23/11/2016 à 20:08, Tim Wescott a écrit :


[...]
Just as a note: I'm not sure what physical processes you're trying to
represent with ratios of polynomials,


I'm working on partial fraction decomposition  algorithms.

.
Like there?:
https://fileexchange.scilab.org/toolboxes/451000

See also 
http://mailinglists.scilab.org/Users-fr-Racines-multiples-d-un-polynome-tt4032750.html


Cheers
Samuel

___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[philippe]

Le 23/11/2016 à 21:20, Serge Steer a écrit :
> First of all if you compute the rational fraction numerator or
> denominator from their roots (and if they are complex conjugate) the
> result is a complex polynomial even if the imaginary part is not shown.
>


Yes I know these problems, but I was surprised that the simplification 
failed just because of a rounding error. In the last months I've got 
several times problems with rounding errors generating complex values.



>
> The simplication algorithm present in Scilab tryes to find reasonable
> simplification based on the bezout equation error. But it has only be
> written for the real case.

so that's the root of the problem!

Thanks,

Philippe

___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[philippe]

Le 23/11/2016 à 20:57, Samuel Gougeon a écrit :
> For all rationals, no (see Tim's post). But for some of them, yes: you
> should use Scilab 6.0-b.
> Some issues of this kind were reported and fixed quite recently:
> http://bugzilla.scilab.org/8493
> https://codereview.scilab.org/#/c/18420
> http://bugzilla.scilab.org/13893

thanks samuel, so the problem is already fixed in scilab 6!


Philippe

___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[philippe]

Le 23/11/2016 à 23:30, Samuel Gougeon a écrit :

Like there?:
https://fileexchange.scilab.org/toolboxes/451000

See also
http://mailinglists.scilab.org/Users-fr-Racines-multiples-d-un-polynome-tt4032750.html


great! It seems you've already done it :-) I've already write functions 
for  the decomposition in C(X) writing the different steps , but your 
code will be usefull.


Philippe.

___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[Samuel Gougeon]

Le 23/11/2016 23:32, philippe a écrit :

Le 23/11/2016 à 20:57, Samuel Gougeon a écrit :
> For all rationals, no (see Tim's post). But for some of them, yes: you
> should use Scilab 6.0-b.
> Some issues of this kind were reported and fixed quite recently:
> http://bugzilla.scilab.org/8493
> https://codereview.scilab.org/#/c/18420
> http://bugzilla.scilab.org/13893

thanks samuel, so the problem is already fixed in scilab 6!



Only for complex-encoded decimal numbers, that is to say for truly null 
imaginary parts.
Otherwise, either clean() must be priorly applied or a complex 
processing be done.


___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users

Re: [Scilab-users] lack of simplification for rational

[Tim Wescott]

Out of curiosity, how does that relate to pfss?

On Wed, 2016-11-23 at 23:30 +0100, Samuel Gougeon wrote:
> Le 23/11/2016 23:19, philippe a écrit :
> > .../...
> > Le 23/11/2016 à 20:08, Tim Wescott a écrit :
> >
> >> [...]
> >> Just as a note: I'm not sure what physical processes you're trying to
> >> represent with ratios of polynomials,
> >
> > I'm working on partial fraction decomposition  algorithms.
> .
> Like there?:
> https://fileexchange.scilab.org/toolboxes/451000
> 
> See also 
> http://mailinglists.scilab.org/Users-fr-Racines-multiples-d-un-polynome-tt4032750.html
> 
> Cheers
> Samuel
> 
> ___
> users mailing list
> users@lists.scilab.org
> http://lists.scilab.org/mailman/listinfo/users

-- 

Tim Wescott
www.wescottdesign.com
Control & Communications systems, circuit & software design.
Phone: 503.631.7815
Cell:  503.349.8432


___
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users