# [Scilab-users] Find the points of intersection of two curves

12 messages
Open this post in threaded view
|

## [Scilab-users] Find the points of intersection of two curves

 Hello, How could I improve the accuracy of intersection points? clear; funcprot(0); timer(); function Sys=F(x) x0=x(1);x1=x(2); Sys(1)=x1^2-x0*x1+x0^2-1; Sys(2)=sin(4*x1^2)+sin(5*x0^2); endfunction   function D=Draghilev(x) x0=x(1);x1=x(2);x2=x(3);x00=ics(1);x01=ics(2); D=[8*x01^2*x1*cos(4*x1^2)-8*x00*x01*x1*cos(4*x1^2)+8*x00^2*x1*cos(4*x1^2)-8*x1*cos(4*x1^2)-2*sin(4*x01^2)*x1-2*sin(5*x00^2)*x1+x0*sin(4*x01^2)+x0*sin(5*x00^2); -sin(4*x01^2)*x1-sin(5*x00^2)*x1+2*x0*sin(4*x01^2)-10*x0*cos(5*x0^2)*x01^2+10*x0*cos(5*x0^2)*x00*x01+2*x0*sin(5*x00^2)-10*x0*cos(5*x0^2)*x00^2+10*x0*cos(5*x0^2); -8*x1^2*cos(4*x1^2)+16*x0*x1*cos(4*x1^2)-20*x0*cos(5*x0^2)*x1+10*x0^2*cos(5*x0^2)]; endfunction ics=[0;1] ics=[ics;1]; disp(ics) function [dxdt]=odes(t,x) [dxdt]=Draghilev(x); endfunction N=100; smin=0.0; smax=5; h=0.0001; step=smax/N;     t=smin:step:smax;     t0=0;     atol=h;     LL= ode("fix",ics,t0,t,atol,odes) disp(timer(),"tiempo de ejecucion C"); a=get("current_axes")//get the handle of the newly created axes a.axes_visible="on"; // makes the axes visible a.font_size=3; //set the tics label font size a.x_location="middle"; //set the x axis position plot(t,LL(1,1:\$),"-cya",t,LL(2,1:\$),"-g",t,LL(3,1:\$),"-r") s=LL(3,:)-LL(1,:); ym=s(1:\$-1).*s(2:\$); z=find(ym <= h); yi=interpln([t;LL(1,:)],t(z+1)); plot(t(z+1),yi,"b*") // s=LL(3,:)-LL(2,:); ym=s(1:\$-1).*s(2:\$); z=find(ym <= 0); I1=LL(2,find(ym <= h)); plot(t(z+1),I1,"g*") Gracias -- Sent from: http://mailinglists.scilab.org/Scilab-users-Mailing-Lists-Archives-f2602246.html_______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: {EXT} Find the points of intersection of two curves

 Hello Hermes, If I understand well,, the core of your question : > De : Hermes > Envoyé : vendredi 12 janvier 2018 10:47 > > How could I improve the accuracy of intersection points? Lies here: > s=LL(3,:)-LL(1,:); > ym=s(1:\$-1).*s(2:\$); > z=find(ym <= h); So you have 2 curves which y-values are LL(1,:) and LL(3,:) and you search where the sign of the difference changes. Correct? So you suppose that both curves are continuous. The refinement of the solution is more a mathematical topic than a Scilab one and I'm afraid I don't have the skill and time to analyze your functions. There are some general methods you could use: use a polynomial interpolation instead of a linear interpolation, i.e. choose an interval around the z values found above, perform a polynomial regression - you have an example of such code here : https://fr.wikibooks.org/wiki/D%C3%A9couvrir_Scilab/Calcul_num%C3%A9rique#R%C3%A9gression_polynomiale - and search the local root of the polynomial. You can also refine the step of calculation around the z values. -- Christophe Dang Ngoc Chan Mechanical calculation engineer This e-mail may contain confidential and/or privileged information. If you are not the intended recipient (or have received this e-mail in error), please notify the sender immediately and destroy this e-mail. Any unauthorized copying, disclosure or distribution of the material in this e-mail is strictly forbidden. _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: {EXT} Find the points of intersection of two curves

 Hi, This problem is equivalent to finding n-zeros over a given interval. Matlab function fnzeros uses the following reference, which uses spline interpolation: Mørken K. and Reimers, M. [2007] An unconditionally convergent method for computing zeros of splines and polynomials, Math. Comp., 76, 845--865 Regards, Rafael -----Original Message----- From: users [mailto:[hidden email]] On Behalf Of Dang Ngoc Chan, Christophe Sent: Friday, January 12, 2018 11:03 AM To: Users mailing list for Scilab <[hidden email]> Subject: Re: [Scilab-users] {EXT} Find the points of intersection of two curves Hello Hermes, If I understand well,, the core of your question : > De : Hermes > Envoyé : vendredi 12 janvier 2018 10:47 > > How could I improve the accuracy of intersection points? Lies here: > s=LL(3,:)-LL(1,:); > ym=s(1:\$-1).*s(2:\$); > z=find(ym <= h); So you have 2 curves which y-values are LL(1,:) and LL(3,:) and you search where the sign of the difference changes. Correct? So you suppose that both curves are continuous. The refinement of the solution is more a mathematical topic than a Scilab one and I'm afraid I don't have the skill and time to analyze your functions. There are some general methods you could use: use a polynomial interpolation instead of a linear interpolation, i.e. choose an interval around the z values found above, perform a polynomial regression - you have an example of such code here : https://fr.wikibooks.org/wiki/D%C3%A9couvrir_Scilab/Calcul_num%C3%A9rique#R%C3%A9gression_polynomiale - and search the local root of the polynomial. You can also refine the step of calculation around the z values. -- Christophe Dang Ngoc Chan Mechanical calculation engineer This e-mail may contain confidential and/or privileged information. If you are not the intended recipient (or have received this e-mail in error), please notify the sender immediately and destroy this e-mail. Any unauthorized copying, disclosure or distribution of the material in this e-mail is strictly forbidden. _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users_______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 In reply to this post by Hermes Hi Hermes,   Here below a simple solution using linear interpolation that works fine with your data but would need some further refinement for general arbitrary inputs.   // START OF CODE clear; clf();   function Sys=F(x)     x0=x(1);x1=x(2);     Sys(1)=x1^2-x0*x1+x0^2-1;     Sys(2)=sin(4*x1^2)+sin(5*x0^2); endfunction   function D=Draghilev(x)     x0=x(1);x1=x(2);x2=x(3);x00=ics(1);x01=ics(2);     D=[8*x01^2*x1*cos(4*x1^2)-8*x00*x01*x1*cos(4*x1^2)+8*x00^2*x1*cos(4*x1^2)-8*x1*cos(4*x1^2)-2*sin(4*x01^2)*x1-2*sin(5*x00^2)*x1+x0*sin(4*x01^2)+x0*sin(5*x00^2);     -sin(4*x01^2)*x1-sin(5*x00^2)*x1+2*x0*sin(4*x01^2)-10*x0*cos(5*x0^2)*x01^2+10*x0*cos(5*x0^2)*x00*x01+2*x0*sin(5*x00^2)-10*x0*cos(5*x0^2)*x00^2+10*x0*cos(5*x0^2);     -8*x1^2*cos(4*x1^2)+16*x0*x1*cos(4*x1^2)-20*x0*cos(5*x0^2)*x1+10*x0^2*cos(5*x0^2)]; endfunction   function [dxdt]=odes(t, x)     [dxdt]=Draghilev(x); endfunction   ics=[0;1;1] N=100; smin=0.0; smax=5; h=0.0001;   step=smax/N; t=smin:step:smax; t0=0; atol=h; LL= ode("fix",ics,t0,t,atol,odes)   plot(t,LL(1,1:\$),"-cya",t,LL(2,1:\$),"-g",t,LL(3,1:\$),"-r") a=gca();   s= LL(3,:) - LL(1,:); ym= s(1:\$-1).*s(2:\$); z=find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); y0 = interpln([t;LL(1,:)],t0); plot(t0,y0,"b*")   s= LL(3,:) - LL(2,:); ym= s(1:\$-1).*s(2:\$); z=find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); y0 = interpln([t;LL(2,:)],t0); plot(t0,y0,"black*") // END OF CODE   Regards, Rafael       _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 Hi Hermes, Find here below a simple solution using linear interpolation that works fine with your data but would need some further refinement for general arbitrary inputs. // START OF CODE clear; clf(); function Sys=F(x)     x0=x(1);x1=x(2);     Sys(1)=x1^2-x0*x1+x0^2-1;     Sys(2)=sin(4*x1^2)+sin(5*x0^2); endfunction function D=Draghilev(x)     x0=x(1);x1=x(2);x2=x(3);x00=ics(1);x01=ics(2);     D=[8*x01^2*x1*cos(4*x1^2)-8*x00*x01*x1*cos(4*x1^2)+8*x00^2*x1*cos(4*x1^2)-8*x1*cos(4*x1^2)-2*sin(4*x01^2)*x1-2*sin(5*x00^2)*x1+x0*sin(4*x01^2)+x0*sin(5*x00^2);     -sin(4*x01^2)*x1-sin(5*x00^2)*x1+2*x0*sin(4*x01^2)-10*x0*cos(5*x0^2)*x01^2+10*x0*cos(5*x0^2)*x00*x01+2*x0*sin(5*x00^2)-10*x0*cos(5*x0^2)*x00^2+10*x0*cos(5*x0^2);     -8*x1^2*cos(4*x1^2)+16*x0*x1*cos(4*x1^2)-20*x0*cos(5*x0^2)*x1+10*x0^2*cos(5*x0^2)]; endfunction function [dxdt]=odes(t, x)     [dxdt]=Draghilev(x); endfunction ics=[0;1;1] N=100; smin=0.0; smax=5; h=0.0001; step=smax/N; t=smin:step:smax; t0=0; atol=h; LL= ode("fix",ics,t0,t,atol,odes) plot(t,LL(1,1:\$),"-cya",t,LL(2,1:\$),"-g",t,LL(3,1:\$),"-r") a=gca(); s= LL(3,:) - LL(1,:); ym= s(1:\$-1).*s(2:\$); z=find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); y0 = interpln([t;LL(1,:)],t0); plot(t0,y0,"b*") s= LL(3,:) - LL(2,:); ym= s(1:\$-1).*s(2:\$); z=find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); y0 = interpln([t;LL(2,:)],t0); plot(t0,y0,"black*") // END OF CODE Regards, Rafael _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 In reply to this post by Hermes Hello Rafael, with your proposal the result is as expected. Thank you I have delayed responding because I wanted to arrive to the next step by my means. But without results: I want to achieve the following graph from the data obtained in LL. I have played with several configurations of the following graphical function: plot3d (LL (1, :), LL (2, :), 0 * ones (1,101), alpha = 35, theta = 45), for the circle parallel to the xy plane. plot3d (LL (1, :), LL (2, :), LL (3, :), alpha = 35, theta = 45), for the function in space. then it would be to put the intersection points. -- Sent from: http://mailinglists.scilab.org/Scilab-users-Mailing-Lists-Archives-f2602246.html_______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 The correct graph for the script shown should be the following:   my apologies -- Sent from: http://mailinglists.scilab.org/Scilab-users-Mailing-Lists-Archives-f2602246.html_______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 Hi Hermes, Check this solution out: // START OF CODE clear; clf(); function Sys=F(x)     x0=x(1);x1=x(2);     Sys(1)=x1^2-x0*x1+x0^2-1;     Sys(2)=sin(4*x1^2)+sin(5*x0^2); endfunction function D=Draghilev(x)     x0=x(1);x1=x(2);x2=x(3);x00=ics(1);x01=ics(2);     D=[8*x01^2*x1*cos(4*x1^2)-8*x00*x01*x1*cos(4*x1^2)+8*x00^2*x1*cos(4*x1^2)-8*x1*cos(4*x1^2)-2*sin(4*x01^2)*x1-2*sin(5*x00^2)*x1+x0*sin(4*x01^2)+x0*sin(5*x00^2);     -sin(4*x01^2)*x1-sin(5*x00^2)*x1+2*x0*sin(4*x01^2)-10*x0*cos(5*x0^2)*x01^2+10*x0*cos(5*x0^2)*x00*x01+2*x0*sin(5*x00^2)-10*x0*cos(5*x0^2)*x00^2+10*x0*cos(5*x0^2);     -8*x1^2*cos(4*x1^2)+16*x0*x1*cos(4*x1^2)-20*x0*cos(5*x0^2)*x1+10*x0^2*cos(5*x0^2)]; endfunction function [dxdt]=odes(t,x)     [dxdt]=Draghilev(x); endfunction ics=[0;1;1] N=100; smin=0.0; smax=5; h=0.0001; step=smax/N; t=smin:step:smax; t0=0; atol=h; LL= ode("fix",ics,t0,t,atol,odes) s= LL(3,:); ym= s(1:\$-1).*s(2:\$); z= find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); x0 = interpln([t;LL(1,:)],t0); y0 = interpln([t;LL(2,:)],t0); param3d(LL(1,:), LL(2,:), 0 * LL(1,:), alpha=35, theta=45); // circle in the xy plane e = gce(); e.foreground =17; param3d(LL(1,:), LL(2,:), LL(3,:), alpha=35, theta=45) // parametric function in 3d-space param3d(x0, y0, 0*x0, alpha=35, theta=45) // zero-crossings e = gce(); e.line_mode = "off"; e.mark_mode = "on" e.mark_style = 3; e.mark_foreground =5 ; // END OF CODE Regards, Rafael _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 In previous email the code formatting got messed up, fixed it here below: // START OF CODE clear; clf(); function Sys=F(x)     x0=x(1);x1=x(2);     Sys(1)=x1^2-x0*x1+x0^2-1;     Sys(2)=sin(4*x1^2)+sin(5*x0^2); endfunction function D=Draghilev(x)     x0=x(1);x1=x(2);x2=x(3);x00=ics(1);x01=ics(2);     D=[8*x01^2*x1*cos(4*x1^2)-8*x00*x01*x1*cos(4*x1^2)+8*x00^2*x1*cos(4*x1^2)-8*x1*cos(4*x1^2)-2*sin(4*x01^2)*x1-2*sin(5*x00^2)*x1+x0*sin(4*x01^2)+x0*sin(5*x00^2);     -sin(4*x01^2)*x1-sin(5*x00^2)*x1+2*x0*sin(4*x01^2)-10*x0*cos(5*x0^2)*x01^2+10*x0*cos(5*x0^2)*x00*x01+2*x0*sin(5*x00^2)-10*x0*cos(5*x0^2)*x00^2+10*x0*cos(5*x0^2);     -8*x1^2*cos(4*x1^2)+16*x0*x1*cos(4*x1^2)-20*x0*cos(5*x0^2)*x1+10*x0^2*cos(5*x0^2)]; endfunction function [dxdt]=odes(t,x)     [dxdt]=Draghilev(x); endfunction ics=[0;1;1] N=100; smin=0.0; smax=5; h=0.0001; step=smax/N; t=smin:step:smax; t0=0; atol=h; LL= ode("fix",ics,t0,t,atol,odes) s= LL(3,:); ym= s(1:\$-1).*s(2:\$); z= find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); x0 = interpln([t;LL(1,:)],t0); y0 = interpln([t;LL(2,:)],t0); param3d(LL(1,:), LL(2,:), 0 * LL(1,:), alpha=35, theta=45); // circle in the xy plane e = gce(); e.foreground =17; param3d(LL(1,:), LL(2,:), LL(3,:), alpha=35, theta=45); // parametric function in 3d-space param3d(x0, y0, 0*x0, alpha=35, theta=45); // zero-crossings e = gce(); e.line_mode = "off"; e.mark_mode = "on" e.mark_style = 3; e.mark_foreground =5 ; // END OF CODE Regards, Rafael _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users
Open this post in threaded view
|

## Re: Find the points of intersection of two curves

 (still problems with the cut and pasting of code from SciNotes, removed some of the comments now) // START OF CODE clear; clf(); function Sys=F(x)     x0=x(1);x1=x(2);     Sys(1)=x1^2-x0*x1+x0^2-1;     Sys(2)=sin(4*x1^2)+sin(5*x0^2); endfunction function D=Draghilev(x)     x0=x(1);x1=x(2);x2=x(3);x00=ics(1);x01=ics(2);     D=[8*x01^2*x1*cos(4*x1^2)-8*x00*x01*x1*cos(4*x1^2)+8*x00^2*x1*cos(4*x1^2)-8*x1*cos(4*x1^2)-2*sin(4*x01^2)*x1-2*sin(5*x00^2)*x1+x0*sin(4*x01^2)+x0*sin(5*x00^2);     -sin(4*x01^2)*x1-sin(5*x00^2)*x1+2*x0*sin(4*x01^2)-10*x0*cos(5*x0^2)*x01^2+10*x0*cos(5*x0^2)*x00*x01+2*x0*sin(5*x00^2)-10*x0*cos(5*x0^2)*x00^2+10*x0*cos(5*x0^2);     -8*x1^2*cos(4*x1^2)+16*x0*x1*cos(4*x1^2)-20*x0*cos(5*x0^2)*x1+10*x0^2*cos(5*x0^2)]; endfunction function [dxdt]=odes(t,x)     [dxdt]=Draghilev(x); endfunction ics=[0;1;1] N=100; smin=0.0; smax=5; h=0.0001; step=smax/N; t=smin:step:smax; t0=0; atol=h; LL= ode("fix",ics,t0,t,atol,odes) s= LL(3,:); ym= s(1:\$-1).*s(2:\$); z= find(ym <= 0); t0 = t(z) - s(z).*(t(z+1)-t(z))./(s(z+1)-s(z)); x0 = interpln([t;LL(1,:)],t0); y0 = interpln([t;LL(2,:)],t0); param3d(LL(1,:), LL(2,:), 0 * LL(1,:), alpha=35, theta=45); e = gce(); e.foreground =17; param3d(LL(1,:), LL(2,:), LL(3,:), alpha=35, theta=45); param3d(x0, y0, 0*x0, alpha=35, theta=45); e = gce(); e.line_mode = "off"; e.mark_mode = "on"; e.mark_style = 3; e.mark_foreground =5; // END OF CODE Regards, Rafael _______________________________________________ users mailing list [hidden email] http://lists.scilab.org/mailman/listinfo/users