[Scilab-users] root calculation

[Paul Carrico]

Hi Calixte,

 

Your suggestion looks encouraging 
 I’ll have a deeper look on that specific
issue  when all the aspects of my new project will be bracketed ! Many
thanks 


 

And Yes:

-          The exponents are always negative

-          The coefficients are always positive

 

To go ahead, in one aspect of the project, only the last constant change (it
won(t be the case in a second stage), I mean :

0.403*X^(-0.121) + 60.5*X^(-0.73) è always the same == constant !

– 0.1839 è it changes 
 may  be either positive or negative in a first
understanding

 

Imagine I’ve a matrix (n x1) where reach line is independent from the
others, and where it’s necessary to solve such equation 
 I’ve been thinking
using a loop:

For i=1:n

   Solve my non-linear equation

End

(n from hundred thousand’s to million’s)

 

It may be interesting to slip that loop onto several processor (solving by
“blocks”), isn’t it ?

 

Paul

 

De : users-bounces at lists.scilab.org [mailto:users-bounces at lists.scilab.org]
De la part de Calixte Denizet
Envoyé : lundi 5 novembre 2012 16:10
À : users at lists.scilab.org
Objet : Re: [Scilab-users] root calculation

 

Hi Paul,

If the exponent are always negative and the coefficients are always
positive, the function is strictly decreasing (so the derivative is
non-null), you can use a Newton method.

function y=newton(fun, dfun, x0, eps)
    y=fun(x0)
    while abs(y) > eps
        x0 = x0 - y / dfun(x0);
        y = fun(x0)
    end
    y = x0;
endfunction

deff('y=foo(x)','y=0.403*x.^(-0.121)+60.5*x.^(-0.73)-0.1839')
deff('y=dfoo(x)','y=0.403*-0.121*x.^(-1.121)+60.5*-0.73*x.^(-1.73)')

To find a "good" starting point:
If x is a solution then 0.403*x^(-0.121)<=0.1839 and 60.5*x^(-0.73)<=0.1839,
so x >=max((0.1839/0.403)^(-1/0.121), (0.1839/60.5)^(-1/0.73))

newton(foo,dfoo,2806,1e-10)

A faster way:
function y=newton2(coeffs, expo, cste, x0, eps)
    y = coeffs * (x0 .^ expo)' - cste;
    dcoeffs = coeffs .* expo;
    dexpo = expo - 1;
    while abs(y) > eps
        x0 = x0 - y / (dcoeffs * (x0 .^ dexpo)');
        y = coeffs * (x0 .^ expo)' - cste;
    end
    y = x0;
endfunction

newton2([0.403 60.5], [-0.121 -0.73],0.1839,2806,1e-10)

Calixte

On 05/11/2012 13:13, Paul Carrico wrote: 

Dear all,

 

This a stupid question, but how can I solve directly in, Scilab  an equation
such as : 

0.403*X^(-0.121) + 60.5*X^(-0.73) – 0.1839 = 0

?

Is-it necessary to code a function ? from memory : dichotomy method, secant
method, Brent one etc. 


 

Regards

 

Paul






_______________________________________________
users mailing list
users at lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users






-- 
Calixte Denizet
Software Development Engineer
-----------------------------------------------------------
Scilab Enterprises
143bis rue Yves Le Coz - 78000 Versailles, France
http://www.scilab-enterprises.com
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.scilab.org/pipermail/users/attachments/20121105/79d0d003/attachment.htm>