|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
MATHEMATICA讲座第六讲 <BR>线性方程组的表达方式和解法 <BR> <BR>一.向量和矩阵的输入 <BR> <BR>1.Range[n] <BR>Range[n,m] <BR>Range[n,m,h] <BR>执行算例 <BR>Range[5] <BR>Range[1,10,2] <BR>2.Table[(1/2)^n,{n,0,10}] (*由通项构造表*) <BR>Table[{f1(n),F2(n)},{n,n1,n2,h}] <BR>执行算例 <BR>Table[(1/2)^n,{n,0,10}] <BR>A=Table[Normal[Series[Sin[x],{x,0,n}]],{n,1,11,2}] <BR>(*Sin[x]的六个正规化后的幂级数展开式的表*) <BR>执行算例 <BR>Table[{x,x^2},{x,-1.0,1.0,0.2}] <BR>Table[Expand[(x^i+i*x)^2],{i,2,5}] <BR>Table[Mod[n,2],{n,0,17}] <BR>Table[Sum[x^i,{i,0,n}],{n,1,5}] <BR>Table[Random[],{10}] <BR>Table[Random[Real,{1,10}],{10}] <BR>3.Array[函数,n] (*由函数表达式构造表*) <BR>Array[函数,{n1,n2,n3,...}] <BR>执行算例 <BR>Array[Exp,5](*等价于Table[Exp[x],{x,5}]*) <BR>Array[Mod,{10,10}] <BR>(*等价于Table[Mod[n,m]{n,0,10},{m,0,10}]*) <BR>4.NestList[f[n],n0,k] 用递推公式建立表元素 <BR>Clear[n,n0,rnt,fnt,t1,t2] <BR>t1=(Sqrt[5]+1)/2 <BR>t2=(1-Sqrt[5])/2 <BR>fnt=Table[(t1^(n+1)-t2^(n+1))/Sqrt[5],{n,0,40}]//N <BR>rnt=Table[fnt[[n-1]]/fnt[[n]],{n,2,12}] <BR>5.向量与矩阵的标准输入法 <BR>A={x1,x2,x3,...} <BR>A={{a11,a12,a13},{a21,a22,a23},{a31,a32,a33}} <BR>ColumnForm[{a1,a2,...,an}]把向量用列方式输出 <BR>MatrixForm[A] 用矩阵方式显示 <BR>IdentityMatrix[n] 生成n阶单位阵 <BR>DiagonalMatrix[{a11,a22,...,ann}]生成对角阵 <BR>执行算例 <BR>A1={1,2,3} <BR>ColumnForm[A1] <BR>A2=DiagonalMatrix[{1,2,3,4,5}] <BR>MatrixForm[A2] <BR>IdentityMatrix[3] <BR>DiagonalMatrix[{1,2,3,4,5}] <BR>MatrixForm[%] <BR> <BR>二.行列式 <BR>Det[A] <BR>执行算例 <BR>A={{1,2,3},{4,5,6},{7,8,9}} <BR>MatrixForm[A] <BR>Det[A] <BR> <BR>三.矩阵求逆,求特征值,特征向量 <BR>Inverse[A] <BR>Eigen<I>value</I>[A] <BR>Eigenvector[A] <BR>执行算例 <BR>Clear[A] <BR>A={{4,6,0},{-3,-5,0},{-3,-6,1}} <BR>Eigen<I>value</I>s[A] <BR>Eigenvectors[A] <BR> <BR>四.恰定方程求解 <BR>问题1 x1+6x2+36x3=104 <BR>x1+10x2+100x3=160 <BR>x1+20x2+400x3=370 <BR>程序 <BR>Clear[A1,b] <BR>A1={{1,6,36},{1,10,100},{1,20,400}}; <BR>b={104,160,370}; <BR>LinearSolve[A1,b](*求方程组的解*) <BR>X=Inverse[A1].b (*用求逆矩阵方法求解*) <BR> <BR>五.欠定方程求解 <BR>问题2 2x1+x2-x3+x4=1 <BR>x1+2x2+x3-x4=2 <BR>x1+x2+2x3+x4=3 <BR>程序 <BR>A2={{2,1,-1,1},{1,2,1,-1},{1,1,2,1}} <BR>b2={1,2,3}; <BR>X={x1,x2,x3,x4}; <BR>Solve[A2.X==b2,{x1,x2,x3,x4}] <BR> <BR> <BR>MATHEMATICA讲座第七讲 <BR>函数的插值 <BR> <BR>一.拉格朗日插值 <BR>L={List} <BR>InterpolatingPolynomial[L,x] <BR>执行算例1 两点线性插值 <BR>L={{0,0.3},{0.2,0.45}} <BR>I=InterpolatingPolynomial[L,x] <BR>执行算例2 三点抛物插值 <BR>L1={{0,0.3},{0.2,0.45},{0.4,0.15}} <BR>I1=InterpolatingPolynomial[L1,x] <BR>执行算例3 多点拉格朗日插值 <BR>L2={{0,0.3},{0.2,0.45},{0.3,0.47}, <BR>{0.52,0.50},{0.64,0.38},{0.7,0.33},{1.0,0.24}} <BR>I2=InterpolatingPolynomial[L2,x] <BR>Plot[%,{x,-0.25,1.05}] <BR>执行算例4 作正弦在0,P上五点插值函数图形 <BR>g0=Plot[Sin[x],{x,0,Pi}] <BR>L=Line[Table[{x,Sin[x]},{x,0,Pi,Pi/4}]] <BR>g=Graphics[L] <BR>Show[g0,g] <BR>sinAp[n_]:=Graphics[{Line[Table[{x,Sin[x]}, <BR>{x,0,Pi,Pi/(n+1)}]]}] <BR>sinAp[2] <BR>Show[g0,%] <BR> <BR>二.龙格现像演示 <BR>L=Table[{x,1/(1+25*x^2)},{x,-1,1,0.2}] <BR>a=InterpolatingPolynomial[L,x] <BR>b=Plot[1/(1+25*x^2),{x,-1,1}, <BR>PlotStyle->{RGBColor[1,0,0]}] <BR>c=Plot[a,{x,-1,1}] <BR>Show[b,c] <BR> <BR>三. 两点三次Hermite插值 <BR>执行算例5 <BR>Clear[x,y,h0,h1,H0,H1] <BR>x1={0,1};y1={1,2};m={1/2,1/2}; <BR>h0[x_]:=(1+2*x)*(x-1)^2; <BR>h1[x_]=(1-2(x-1))*(x/(x-1))^2 <BR>H0[x_]=x*(x-1)^2 <BR>H1[x_]=(x-1)*x^2 <BR>H[x_]=y1[[1]]*h0[x]+y1[[2]]*h1[x] <BR>+m[[1]]*H0[x]+m[[2]]*H1[x] <BR>%/.{x->0.55} <BR> <BR>四. N+1个节点的2N+1次Hermite插值 <BR>执行算例6 <BR>Clear[x0,y,bb,w,w1,w2,L,h,H,Hm] <BR>x0={0.4,0.5,0.6,0.7,0.8} <BR>y=Table[Log[x],{x,0.40,0.80,0.10}] <BR>m=Table[1/x,{x,0.40,0.80,0.10}] <BR>bb[x_]=InterpolatingPolynomial[y,x] <BR>Simplify[bb[x]] <BR>bb[0.55] <BR>w[x_]=(x-x0[[1]])*(x-x0[[2]])*(x-x0[[3]])*(x-x0[[4]])*(x-x0[[5]]) <BR>w1[x_]=D[w[x],x]; <BR>Simplify[w1[x]] <BR>w2[x_]=D[w[x],{x,2}]; <BR>Simplify[w2[x]] <BR>For[ i=1,i<=5,i++, <BR>L[i_,x_]:=w[x]/((x-x0[])*w1[x0[]])]; <BR>h[i_,x_]:=(1-w2[x0[]]*(x-x0[])/w1[x0[]])*L[i,x]^2; <BR>H[i_,x_]:=L[i,x]^2*(x-x0[]); <BR>] <BR>Hm[x_]=Sum[y[]*h[i,x]+m[]*H[i,x],{i,1,5,1}]; <BR>Simplify[Hm[x]] <BR>Hm[0.55] <BR> <BR>拟合 <BR> <BR>一.一元线性拟合 <BR>执行算例 <BR>b2={{100,45},{110,51},{120,54},{130,61},{140,66}, <BR>{150,70},{160,74},{170,78},{180,85},{190,89}} <BR>fp=ListPlot[b2,PlotStyle->{PointSize[0.03], <BR>RGBColor[0,0,1]}] <BR>ft1=Fit[b2,{1,x},x] <BR>gp=Plot[ft1,{x,100,190},PlotStyle->{RGBColor[1,0,0]}]; <BR>Show[fp,gp] <BR> <BR>二.抛物线拟合 <BR>执行算例 <BR>B=Table[Prime[n],{n,20}] <BR>t1=ListPlot[B,PlotStyle->{RGBColor[0,1,1], <BR>PointSize[0.04]}] <BR>f=Fit[B,{1,x,x^2},x] <BR>t2=Plot[f,{x,0,20}] <BR>Show[t1,t2] <BR> <BR>三.多项式拟合 <BR>执行算例 <BR>data={{0,1.2},{1,1.4},{2,1.3},{3,1.5},{4,1.3},{5,1.3},{6,1.1}}; <BR>t2=ListPlot[data,PlotStyle->{PointSize[0.05], <BR>RGBColor[1,0,0]}] <BR>fx=Fit[data,{1,x,x^2,x^3,x^4,x^5},x] <BR>t1=Plot[fx,{x,0,6},PlotStyle->RGBColor[1,0,1]] <BR>Show[t1,t2] <BR>执行算例 <BR>b3={{1,4},{2,6.4},{3,8.0},{4,8.4},{5,9.28},{6,9.5}, <BR>{7,9.7},{8,9.86},{9,10.0},{10,10.2},{11,10.32}, <BR>{12,10.42},{13,10.5},{14,10.55},{15,10.58}, <BR>{16,10.6}} <BR>gp=ListPlot[b3,PlotStyle->{RGBColor[0,1,0], <BR>PointSize[0.04]}] <BR>ft2=Fit[b3,Table[x^i,{i,0,4}],x] <BR>fp=Plot[ft2,{x,0,17},PlotStyle->{RGBColor[1,0,0]}] <BR>Show[gp,fp] <BR> <BR>四.非线性拟合---指数拟合 <BR>执行算例 求一个经验函数,型如 <BR>x 1 2 3 4 5 6 7 8 <BR>y 15.3 20.5 27.4 36.5 49.1 65.5 87.8 117.6 <BR>程序 <BR>b4={{1,15.3},{2,20.5},{3,27.4},{4,36.6},{5,49.1}, <BR>{6,65.5},{7,87.8},{8,117.6}}; <BR>gb4=ListPlot[b4,PlotStyle->{RGBColor[0,0,1], <BR>PointSize[0.05]}] <BR>y4=Table[b4[[i,2]],{i,1,8}]; <BR>ly4=Log[y4]; <BR>fy4=Fit[ly4,{1,x},x] <BR>s[x_]:=Exp[fy4] <BR>ty=Plot[s[x],{x,1,8},Axes->{2,60}, <BR>AspectRatio->1, <BR>PlotStyle->{RGBColor[1,0,0]}, <BR>PlotRange->{10,120}] <BR>Show[ty,gb4, Axes->{2,60},AxesLabel->{"x","y"}, <BR>AspectRatio->1, <BR>Ticks->{{1,2,3,4,5,6,7,8},{0,30,60,90,120}}] <BR>执行算例 求一个经验函数,型如y=a*exp(-bx)与所给数据拟合 <BR>x 0.4 0.5 0.6 0.7 <BR>y 1.75 1.34 1.00 0.74 <BR>程序 <BR>Clear[fx,fy,biao,nb.ft,ft1,t1] <BR>fx[x_]:=x <BR>fy[y_]:=Log[y] <BR>biao={{0.4,1.75},{0.5,1.34},{0.6,1.00},{0.7,0.74}}; <BR>nb=Table[{fx[biao[[i,1]]],fy[biao[[i,2]]]},{i,1,4}]; <BR>(*拟合方程*) <BR>ft=Fit[nb,{1,x},x] <BR>ft1=Exp[ft] <BR>(*拟合曲线*) <BR>t1=Plot[ft1,{x,0,1.0},AxesLabel->{"x","y"}, <BR>PlotStyle->{RGBColor[1,0,0]}] <BR>t2=ListPlot[biao,PlotStyle->{RGBColor[0,0,1],PointSize[0.04]}] <BR>Show[t1,t2,PlotRange->{0,2}] <BR> <BR>五.用正交多项式作拟合 <BR>[0,1]区间上的勒让得多项式 <BR>(*定义勒让得函数(n=10)*) <BR>Clear[x,t,s] <BR>n=10 <BR>P0[x_]:=1 <BR>P1[x_]:=1-2x/n <BR>P2[x_]:=1-6 x/n+6 x*(x-1)/(n*(n-1)) <BR>P3[x_]:=1-12 x/n+30x*(x-1)/(n*(n-1))-20 x*(x-1)*(x-2)/(n*(n-1)*(n-2)) <BR>P4[x_]:=1-2 x+x*(x-1)-140x*(x-1)*(x-2)/720+70x*(x-1)*(x-2)*(x-3)/(10*9*8*7) <BR>P5[x_]:=1-30x/n+210x*(x-1)/(n*(n-1))-560x*(x-1)*(x-2)/(n*(n-1)*(n-2)) <BR>+630x*(x-1)*(x-2)*(x-3)/(n*(n-1)*(n-2)*(n-3)) <BR>-252 x*(x-1)*(x-2)*(x-3)*(x-4)/(n*(n-1)*(n-2)*(n-3)*(n-4)) <BR> <BR> <BR>(*输入初始数据*) <BR>t={0,5,10,15,20,25,30,35,40,45,50}; <BR>y={0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.60,4.66}; <BR>(*做变量替换*) <BR>x=t/5; <BR>(*计算各多项式在节点处的值*) <BR>A={P0[x],P1[x],P2[x],P3[x],P4[x],P5[x]} <BR>(*计算每一行元素平方的和*) <BR>s=Table[0,{i,1,6}]; <BR>For[k=1,k<=6,k++, <BR>s[[k]]=0; <BR>For [i=1,i<=11,i++, <BR>s[[k]]=s[[k]]+A[[k,i]]^2] <BR>] <BR>N[s,6] <BR>(*计算Pk(xi)*yi*) <BR>r=Table[0,{i,1,6},{j,1,11}] <BR>For[i=1,i<=11,i++, <BR>r[[1,i]]=b0[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[1,i]]=b0[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[1,i]]=b0[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[1,i]]=b0[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[2,i]]=b1[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[3,i]]=b2[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[4,i]]=b3[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[5,i]]=b4[]*y[] <BR>] <BR>For[i=1,i<=11,i++, <BR>r[[6,i]]=b5[]*y[] <BR>] <BR>r <BR> MATHEMATICA讲座第八讲 <BR>线性规划与非线性规划 <BR> <BR>模型 min(or max) f=c'X=c1x1+c2x2+...+cnxn <BR>s.t. MX>=b,X>=0 <BR> <BR>命令格式 <BR>ConstrainedMax[f,{inequalities},{x,y,...}] <BR>ConstrainedMin[f,{inequalities},{x,y,...}] <BR>LinearProgramming[c,M,b] <BR>执行算例1 求解线性规划问题 <BR>max f=2x+3y <BR>s.t. x+2y<=8 <BR>0<=x<=4,0<=y<=3 <BR>ConstrainedMax[2*x+3*y,{x+2 y<=8,x<=4,y<=3},{x,y}] <BR>第一个参数是目标函数的最优值,第二个参数是决策变量的取值。 <BR>执行算例2 <BR>min f=x+2y+3z <BR>s.t. 2x-y =1 <BR>x +z=1 <BR>x,y,z>=0 <BR>(*解法一*) <BR>ConstrainedMin[x+2y+3z,{2x-y==1,x+z==1},{x,y,z}] <BR>(*解法二*) <BR>c={2,-3} <BR>M={{-1,-1},{1,-1},{1,0}} <BR>b={-10,2,1} <BR>q=LinearProgramming[c,M,b] <BR>执行算例3 求解线性规划问题 <BR>min -f =- 0.40x1-0.28x2-0.32x3-0.72x4-0.64x5-0.61x6 <BR>s.t. -0.01x1-0.01x2-0.01x3-0.03x4-0.03x5-0.03x6>=-850 <BR>-0.02x1 -0.05x4 >=-700 <BR>-0.02x2 -0.05x5 >=-100 <BR>-0.03x3 -0.08x6>=-900 <BR>x1,x2,x3,x4,x5,x6>=0 <BR>解: <BR>c={-0.4,-0.28,-0.32,-0.72,-0.64,-0.6}; <BR>A={{-0.01,-0.01,-0.01,-0.03,-0.03,-0.03}, <BR>{-0.02,0,0,-0.05,0,0},{0,-0.02,0,0,-0.05,0}, <BR>{0,0,-0.03,0,0,-0.08}}; <BR>b={-850,-700,-100,-900}; <BR>LinearProgramming[c,A,b] <BR>(* 此命令只给出决策变量。*) <BR> <BR>线性约束条件下的非线性规划问题 <BR>线性逼近法(FW法) <BR>模型: (NLP)minf(x) <BR>S.t.: E={x|AX>=b,X>=0} <BR>执行算例4 求解非线性规划问题 <BR>用线性逼近法求解非线性规划问题 <BR>目标函数 f(x1,x2)=(x1-1)^2+(x2-2)^2 <BR>约束条件 0<=x1<=2,0<=x2<=3 <BR>下面是第一次迭代 <BR>Clear[a,b,c,d,s,e,pf] <BR>f[x1_,x2_]:=(x1-1)^2+(x2-2)^2; <BR>gradf={D[f[x1,x2],x1],D[f[x1,x2],x2]}; <BR>c={0.7,1.25}; (*C即基可行解X0*) <BR>s=gradf/.{x1->Part[c,1],x2->Part[c,2]} (*求X0点梯度*) <BR>{-0.6, -1.5} <BR>x1=.; (*清除X1,X2的值*) <BR>x2=.; <BR>p[u_,v_]:=s.{u,v} <BR>a=ConstrainedMin[p[x1,x2],{x1<=2,x2<=3},{x1,x2}] <BR>(*求出最优解Y0*) <BR>b={x1,x2}/.Part[a,2]; <BR>e=b-c; <BR>pf=s.e; <BR>If[Abs[pf]>0.01, <BR>g=c+d*e; <BR>t[d_]:=f[Part[g,1],Part[g,2]]; <BR>w=FindMinimum[t[d],{d,1,0,1}]; <BR>c=c+(d/.Part[w,2])e; <BR>] <BR>Print["c1=",c] <BR>(* 得到初始点c1,将其替换c,运算后的结果继续代替C,直到w的绝对 <BR>值小于0.01为止。*) *) <BR> |
|