该楼层疑似违规已被系统折叠 隐藏此楼查看此楼
data 文件
%% data
a=0;b=1;
K=100;% TVB Limiter中的系数M
n=20;h=(b-a)/n;p=pi;
k=1;% 多项式次数
dd=1;
T=0.001;
tao=0.001*h*h;
M=floor(T/tao)+1;
u_in1=[];u_in2=[];u_in3=[];% 存储求出来的系数
q_in1=[];q_in2=[];q_in3=[];
u_11=[];u_12=[];u_21=[];u_22=[];% 存储RK3中间变量求出来的系数
errorL2=0; % 误差
p(1)=-0.9324695142;p(2)=-0.6612093865;p(3)=-0.2386191861;
p(4)=0.2386191861;p(5)=0.6612093865;p(6)=0.9324695142;
aa(1)=0.1713244924;aa(2)=0.3607615730;aa(3)=0.4679139346;
aa(4)=0.4679139346;aa(5)=0.3607615730;aa(6)=0.1713244924;
rhf函数
function rhf=rhf(x,y,t)
rhf=4*pi*cos(2*pi*(x+y+t))-8*pi*pi*(cos(2*pi*(x+y+t)))^2-8*pi*pi*(cos(2*pi*(x+y+t)))^2-2*pi*sin(2*pi*(x+y+t))*cos(2*pi*(x+y+t))+128*pi^3*(cos(2*pi*(x+y+t)))^2*sin(2*pi*(x+y+t))+16*pi*pi*(cos(2*pi*(x+y+t)))^2*4*pi*pi*sin(2*pi*(x+y+t));
end
主程序
for m=1:M
for i=1:n
for j=1:n
bb(1)=0;bb(2)=0;bb(3)=0;
for l=1:6
for ll=1:6
bb(1)=bb(1)+h^2/4*rhf(p(l)*h/2+x(i),p(ll)*h/2+y(j),(m-1)*tao)*aa(l)*aa(ll);
bb(2)=bb(2)+h^2/4*rhf(p(l)*h/2+x(i),p(ll)*h/2+y(j),(m-1)*tao)*p(l)/2*aa(l)*aa(ll);
bb(3)=bb(3)+h^2/4*rhf(p(l)*h/2+x(i),p(ll)*h/2+y(j),(m-1)*tao)*p(ll)/2*aa(l)*aa(ll);
end
end
f1_in(i,j,:,m)=A\bb';
end
end
end