Mathematica代码.zip

  • gzhgzh
    了解作者
  • Mathematica
    开发工具
  • 2KB
    文件大小
  • zip
    文件格式
  • 0
    收藏次数
  • 10 积分
    下载积分
  • 2
    下载次数
  • 2020-06-20 11:52
    上传日期
代码为计算偏微分方程中的有限差分法,软件使用为mathematica
Mathematica代码.zip
  • Mathematica代码.txt
    4KB
内容介绍
(*second order IDE*) (*data*) jan[x_,s_]:={{If[s<=x,(1-x)s,x(1-s)]}}; g[x_,u_]:={Exp[-u[[1]]]}; lamb={1}; (*module*) dane=solveIDE[g,1,0,1,9,bound->{{0,1}}] (*Qh(x) module*) intPol[n_,p_]:=Module[{r,pary,vek}; Clear[x]; pary = Table[Table[{dane[[i,1]],dane[[i,r]]}, {i,1,n+2}], {r,2,p+1}]; vek[r_]:=Part[pary,r]; Table[Expand[InterpolatingPolynomial[vek[k],x]], {k,1,p}] ]; (*command*) intPol[9,1] (*residual error*) resError[g_,p_,n_]:=Module[{p1,pxx,ps,ca1,prawa}; p1[x_]:=intPol[n,p]; pxx=D[D[intPol[n,p],x],x]; ps=intPol[n,p]/.x- >s; ca1=Expand[Integrate[jan[x,s].ps*lamb,{s,0,1}]]; prawa=g[x,p1[x]]; blad=-pxx+ca1-prawa; Table[N[Prepend[blad,x]],{x,0.1,0.9,0.1} ]//TableForm ]; (*a system of IDE*) jan[x_,s_]:={{s*sin[x],2*Exp[-s]}; {s,cos[x]*Exp[-s]}}; lamb = {1,1}; g[x_,u_]:={-1 u[[1]]-Exp[-x],-1 u[[2]]-4}; N[dane=solveIDE[g,2,0,1,9,bound- >{{0,1}; {0,1}}]] (*command*) intPol[9,2] (*module solveIDE*) Options[solveIDE] = { bound -> homogeneous,iters -> twon}; solveIDE[g_, p_, a_, b_, n_, opts_] := Module[ { h, x, boundary,w0,v0, iterNumber, beta, gamma,solution, bcorrections, useboundary,solve3, solve5, oneStep}, h = (b-a)/(n+1); x = Table [ N[a + i h], {i, 1, n} ]; homogeneous = Table [0, {p}, {2} ]; twon = 2n; boundary=bound/.{opts} /.Options[solveIDE]; v0=Table [(boundary[[q,2]]-boundary[[q,1]])* i*h/(b-a)+boundary[[q,1]],{q,1,p},{i,1,n}]; iterNumber = iters/.{opts} /.Options[solveIDE]; beta={0}; Do[AppendTo[beta, 1/(14-Last[beta]) ], {3} ]; Do[AppendTo[beta, N[7 - 4 Sqrt[3]] ], {n-3} ]; gamma=Table [i/(i+1), {i, 1, n} ]; bcorrections=Table [ {{12*h^2* boundary[[q,1]], - h^2*boundary[[q,1]]},{-h^2*boundary[[q,2]], 12*h^2*boundary[[q,2]]} } /h^2, {q,1,p} ]; useboundary[f_, bcorr ]:=Join[Take[f,2] + bcorr[[1]],Take[f, {3,-3} ], Take[f,-2] + bcorr[[2]]]; solve3[f_, d11_, alpha_]:= Module[{ f1, sol }, f1[1]=f[[1]]=d11; f1[i_]:=f1[i]=(f[[i]]+f1[i-1])*alpha[[i]]; sol[n]=If[d11==12, f[[n]]/12, f1[n]]; sol[i_]:=sol[i]=f1[i]+alpha[[i]]*sol[i+1]; Table [sol[i],{i,1,n} ]]; solve5[f_]:=With[{w= solve3[f, 12, beta]}, solve3[w, 2, gamma]]; cab[i_,q_,r_,vv_]:=Module[{ks0,ks1,ks2,ma }, ks0=-Sqrt [3/5]; ks1=0; ks2=Sqrt[3/5]; N[Sum[5/9 jan[a+i h,a+j h+(1+ks0) h/2][[q,r]]* ((1-ks0)*vv[[r,j+1]]+(1+ks0)vv[[r,j+2]])+ 8/9 jan[a+i h,a+j h+(1-ks1) h/2][[q,r]] * ((1-ks1)*vv[[r,j+1]]+(1+ks1)vv[[r,j+2]])+ 5/9 jan[a+i h,a+j h+(1+ks2)h/2][[q,r]]* ((1-ks2)*vv[[r,j+1]]+(1+ks2)vv[[r,j+2]]), {j,0,n} ]*h/4,4] ]; calka[ ii_,qq_,ww_]:=Sum[cab[ii,qq,rr,ww],{rr,1,p} ]; gaussint[i_,q_,r_,vv_]:=Module[{ }, m3=Mod[n+1,3]; last4=If[m3==0, n-1,n-m3-1]; integrandValues[k_, y_]:=Module[ {intpoints, intpol, integrand, c=Sqrt [3/5]}, intpoints = Table [{a+(j-1)*h,y[[j]]}, {j,k,k+3} ]; intpoly = InterpolatingPolynomial[intpoints,s] ; integrand=jan[a+i*h, s][[q,r]]*intpoly; With[ {xj=a+(j-1)*h}, Table [N[{ (5/9)*integrand/.s->xj+(1-c)*h/2, (8/9)*integrand/.s->xj+ h/2, (5/9)*integrand/.s->xj+(1+c)* h/2 } ], {j, k,k+2} ]* h/2 ]]; integral3[k_,y_]:=Apply[Plus, Flatten[integrandValues[k,y]]]; int1[y_]:=Apply[Plus, Map[integral3[#,y]&, Range[1, last4,3]]]; restint[y_]:=If[m3 == 0, 0, Apply[Plus, Flatten[ Take[integrandValues[n-1,y],-m3]]] ]; y=vv[[r]]; int1[y]+restint[y]]; kcalka[i_,q_,ww_]:=Sum[gaussint[i,q,r,ww],{r,1,p} ]; oneStep[v_]:=Module[ {ff, bff,d,q,r,s,y,zz}, ff=12*h^2*N[Transpose[ MapThread[g,{x,Transpose[v]} ]]]; zz=Transpose[boundary]; cc=Table [Append[Prepend[v[[q]], zz[[1,q]]],zz[[2,q]]],{q,1,p} ]; fff=Table [lamb[[q]]*kcalka[i,q,cc],{i,1,n},{q,1,p} ]; ff=ff-12 h^2 Transpose[fff]; bff=MapThread[useboundary,{ff,bcorrections} ]; Map[solve5, bff] ]; solution = Nest[oneStep, v0, iterNumber]; PrependTo[solution, x]; With[ { initiala} = Prepend[Table [boundary[[q,1]], {q,1,p} ],a], endb= Prepend[Table [boundary[[q,2]], {q,1,p} ],b], Append[Prepend[Transpose[solution], initiala],endb]] ]
评论
    相关推荐
    • 求解Bessel方程.rar
      使用差分法解一类Bessel方程在不同边界条件下的离散解
    • vkkarn.zip
      介绍了OPNET无线仿真的实例,主要对接入层相关仿真的介绍
    • polbcy.rar
      介绍了OPNET无线仿真的实例,主要对接入层相关仿真的介绍
    • 数学建模十大算法
      数学建模十大算法
    • matlab图像相嵌代码-daydayup:一些编码练习
      matlab图像相嵌代码 平时练手的一些代码片段 ImageLabel,一个用于图片标注的web系统,用tomcat + java servlet + ajax实现 acm acm竞赛的习题练习 bitmap_picture 图像压缩的编码练习 ...使用AudioTrack和AudioRecord...
    • 视频图matlab代码-DifferentialEquations.jl-0c46a032-eb83-5123-abaf-570
      (随机)偏微分方程((S)PDEs)(同时使用有限差分法有限元方法) 经过良好优化的DifferentialEquations求解器使用一些经典算法以及最近的研究得出的基准,是一些最快的实现方案,这些算法通常优于“标准” C / ...
    • 常微分方程论文求解算法的matlab源码
      有限差分法是一种成熟而有效的求解常微分方程近似解的方法,这种方法是基于差商代替导数(数值微分)或者积分插值(数值积分),然后构造差分格式,通过差分迭代格式求解原微分方程,获得原微分方程的近似解。 这种方法在...
    • hntroduces.zip
      介绍了OPNET无线仿真的实例,主要对接入层相关仿真的介绍
    • nfxyb9.zip
      介绍了OPNET无线仿真的实例,主要对接入层相关仿真的介绍
    • matlabcnhelp.rar
      matlab中文帮助很难找的,快速下载