gaussian quadrature - Printable Version +- HP Forums (https://www.hpmuseum.org/forum) +-- Forum: HP Calculators (and very old HP Computers) (/forum-3.html) +--- Forum: HP Prime (/forum-5.html) +--- Thread: gaussian quadrature (/thread-4564.html) gaussian quadrature - Antonio - 08-23-2015 07:01 PM 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 . RE: gaussian quadrature - Didier Lachieze - 08-24-2015 03:18 PM (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). RE: gaussian quadrature - Antonio - 08-25-2015 09:23 AM Hello Didier, Thanks for the reply. 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. RE: gaussian quadrature - Didier Lachieze - 08-25-2015 11:42 AM (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)); RE: gaussian quadrature - Antonio - 08-25-2015 12:26 PM Already solved 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