# HP Forums

You're currently viewing a stripped down version of our content. View the full version with proper formatting.
Hello,

First of all sorry for my bad English .
I'm making a program to calculate integrals by gaussian quadrature and I have the following problem , I've done two different ways.
Program 1 , due to the unique character of this calculator , sometimes it works, sometimes not , so i try to make the program 2
Program 2 , do not put the correct result in a0 .

Program 1
Code:
```#pragma mode( separator(.,;) integer(h64) )  EXPORT Cuadratura_Gaussiana() // Calculo de las raices de los polinomios de Legendre BEGIN     local i,aux;     local nr:=3;     M1:=makemat(0,nr,2); // 1ª col factores de ponderación ,2ª col raices     // Cálculo de las raices  de los polinomios de Legendre     M2:=CAS.zeros(legendre(nr));           // Cálculo de los coeficientes      F1:=CAS.diff(legendre(nr));     for i from 1 to nr do          M1(i,2):=M2(i);          aux:=F1(M2(i));          M1(i,2):=2/((1-M2(i)^2)*(aux^2)); // coeficiente     end;       end;```

Programa 2
Code:
```#pragma mode( separator(.,;) integer(h64) ) coeficientes(); EXPORT Cuadratura_Gaussiana() BEGIN       local n:=3;     M1:=makemat(0,n,4);           coeficientes(n);         msgbox("Orden de las columnas: c_i, r_i , F(r_i), c_i*F(ri)");     editmat(M1); end; #cas coeficientes(nr):= BEGIN     // Cálculo de las raices y los coeficientes para la cuadratura gaussiana de de nr términos     local i,a0;     M2:=zeros(legendre(nr)); //raices del polinomio     // Cálculo de los coeficientes de ponderacion     for i from 1 to nr do          M1(i,2):=M2(i);//  raiz del polinomio;           a0:=subst(diff(legendre(nr)),x=M2(i));          L1(i):=a0;           M1(i,1):=2/((1-M2(i)^2)*(a0^2));//coeficiente;     end;      END; #end```

Let's see if someone can lend me a hand, thanks .
(08-23-2015 07:01 PM)Antonio Wrote: [ -> ]Program 1 , due to the unique character of this calculator , sometimes it works, sometimes not , so i try to make the program 2
Can you explain what happens when it doesn’t work? On my side it seems to works fine.

(08-23-2015 07:01 PM)Antonio Wrote: [ -> ]Program 2 , do not put the correct result in a0 .
I’ve done some tests with your program 2 and if I’m correct it seems that there is an issue when the CAS function legendre() is called from a CAS program:
• When called from CAS it returns a polynom of x, for ex. legendre(3) returns 5/2*x^3-3/2*x
• However when called from a CAS program it returns a polynom of X instead of x: 5/2*X^3-3/2*X

For example if you add the following line in your code: f:=legendre(nr); just before a0:=subst(diff(legendre(nr)),x=M2(i)); then after running your program, if you type CAS.f you will get 5/2*X^3-3/2*X.
Same for diff(legendre(nr)) which returns: 15/2*X^2-3/2
As X is a global predefined variable, then a0 is always evaluated to the same value and just depends on the value of X (-3/2 if X=0, its default value).
Hello Didier,

Regarding the program 1 , what I mean is that sometimes it runs and calculates the roots and coefficients and sometimes fails to execute
M2:=CAS.zeros(legendre(nr)); or
F1:=CAS.diff(legendre(nr));
does not give any error message unless it does not.

Regarding the program 2, I tried to change x by X , as you suggested, I had to add one more line to assign the value to the X ,as it does not properly
a0:=subst(diff(legendre(nr)),X=M2(i));

they work as

X:=M2(i);
a0:=subst(diff(legendre(nr)),X=0);

Another problem is that not running the loop, and I 've been manually deactivated by assigning values ​​to the variable i, work when i = 2 or i = 3 , but not when i = 1, the error is in

M1(i,1):=2/((1-M2(i)^2)*(a0^2));
and I dont see the failure more I look.
(08-25-2015 09:23 AM)Antonio Wrote: [ -> ]Regarding the program 1 , what I mean is that sometimes it runs and calculates the roots and coefficients and sometimes fails to execute
M2:=CAS.zeros(legendre(nr)); or
F1:=CAS.diff(legendre(nr));
does not give any error message unless it does not.
On my side every time I've tried the program 1 it has worked correctly so I can't reproduce your issue.

(08-25-2015 09:23 AM)Antonio Wrote: [ -> ]Another problem is that not running the loop, and I 've been manually deactivated by assigning values ​​to the variable i, work when i = 2 or i = 3 , but not when i = 1, the error is in

M1(i,1):=2/((1-M2(i)^2)*(a0^2));
and I dont see the failure more I look.
I think this is because M1 is predefined as a matrix of reals and you're trying to assign to it an exact number when a0 is -3/2. Try with adding a dot to the first 2 to get an approximate number (or create your own matrix variable, different from the predefined ones M0-M9):

M1(i,1):=2./((1-M2(i)^2)*(a0^2));

As you said I put the point to take the number as Real rather than exact and have solved two problems , which did not work when i = 1 and the problem of the loop. Although I do not see the relationship between the loop and the latter problem that the variable a0 can be real or exact.

The corrections program that works is
Code:
```#pragma mode( separator(.,;) integer(h64) ) coeficientes(); EXPORT Cuadratura_Gaussiana() BEGIN       local n:=3;     M1:=makemat(0,n,4);           coeficientes(n);         msgbox("Orden de las columnas: c_i, r_i , F(r_i), c_i*F(ri)");     editmat(M1); end; #cas coeficientes(nr):= BEGIN     // Cálculo de las raices y los coeficientes para la cuadratura gaussiana de de nr términos     local i,a0;     M2:=zeros(legendre(nr)); //raices del polinomio     // Cálculo de los coeficientes de ponderacion     for i from 1 to nr do          M1(i,2):=M2(i);//  raiz del polinomio;          X:=M2(i);          a0:=subst(diff(legendre(nr)),X=0);            M1(i,1):=2.[/b]/((1-M2(i)^2)*(a0^2));//coeficiente;     end;      END; #end```

In relation to program 1 , I have not noticed any particular pattern in failure . It happens also in other programs which use CAS.funtion( ) sometimes makes the call to CAS successfully and some not.

Every time I like least this calculator , since to do something that comes out of the most basic , even as simple as this program is a headache..

In any case thank you very much for your help Didier
Reference URL's
• HP Forums: https://www.hpmuseum.org/forum/index.php
• :