|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
<P><FONT size=4>原来的程序如下<BR><BR>/*e.g:y'=y-2x/y,x∈[0,0.6]</FONT></P>
<P><FONT size=4>/* y(0)=1<BR>/*使用经典四阶龙格-库塔算法进行高精度求解 <BR>/* y(i+1)=yi+( K1+ 2*K2 +2*K3+ K4)/6<BR>/* K1=h*f(xi,yi)<BR>/* K2=h*f(xi+h/2,yi+K1/2)<BR>/* K3=h*f(xi+h/2,yi+K2/2)<BR>/* K4=h*f(xi+h,yi+K3)<BR>*/<BR>#include <stdio.h><BR>#include <math.h></FONT></P>
<P><FONT size=4>float f(float den,float p0) //要求解的微分方程的右部的函数 e.g: y-2x/y<BR>{<BR> float rus;<BR> // den=w/(W0+sl);<BR> // rus=k*A*f/Wp*sqrt(RTp)-(k-1)*.....<BR> rus=p0-2*den/p0;<BR> return(rus);<BR>}</FONT></P>
<P><FONT size=4>//float fden()<BR>//{</FONT></P>
<P><FONT size=4>//}</FONT></P>
<P><BR><FONT size=4>void main()<BR>{<BR> float x0; //范围上限<BR> int x1; //范围下限<BR> float h; //步长<BR> int n; //计算出的点的个数<BR> float k1,k2,k3,k4;<BR> float y[3]; //用于存放计算出的常微分方程数值解<BR> int i=0;<BR> int j;</FONT></P>
<P><FONT size=4>//以下为函数的接口<BR> printf("intput x0:");<BR> scanf("%f",&x0);</FONT></P>
<P><FONT size=4> printf("input x1:");<BR> scanf("%f",&x1);</FONT></P>
<P><FONT size=4> printf("input y[0]:");<BR> scanf("%f",&y[0]);</FONT></P>
<P><FONT size=4> printf("input h:");<BR> scanf("%f",&h);</FONT></P>
<P><FONT size=4>// 以下为核心程序<BR> n=((x1-x0)/h);<BR> n=3;<BR> <BR> for(j=0;j<n;j++)<BR> { </FONT></P>
<P><FONT size=4> k1=h*f(x0,y); //求K1<BR> k2=h*f((x0+h/2),(y+k1/2)); //求K2<BR> k3=h*f((x0+h/2),(y+k2/2)); //求K3<BR> k4=h*f((x0+h),(y+k3)); //求K4<BR> <BR> y[i+1]=y+((k1+2*k2+2*k3+k4)/6); //求yi+1<BR> x0+=float(0.2);<BR> printf("y[%f]=%f\n",x0,y[i+1]); <BR> ++i;<BR> <BR> }<BR> <BR>}</FONT></P>
<P><BR><FONT size=4> y'是微分的意思,要把原来的方程改成y'=-(400+4000*cos(x))+(3200000-8000000*sin(x))/y,,,x∈[-180°,180°],,,中间要取100个点,也就是h=3.6°,,,y(0)=16000,度数也可以用实数表示。<BR><BR>那个大师会帮帮我,感激不尽!加我QQ也行157544051</FONT></P> |
|