哇叽菜 2021-09-22 18:43 采纳率: 100%
浏览 134
已结题

一维抛物型方程的差分方法

急求一维抛物型方程的差分方法,要有数值解,精确解和误差分析,急求!

  • 写回答

1条回答 默认 最新

  • 技术专家团-Joel 2021-09-22 18:54
    关注

    你好!比如一维热传导方程就是抛物线型的

    % Number of Control Volumes
    N = 5;
    % Domain length
    L = 0.02; %[m]
    % Grid Spacing
    h = L/(N);
    % Diffusivity (Thermal conductivity)
    k = 0.5;  %[W/m-K]
    % Uniform Heat Generation, A Source Term
    q = 1000e3;     %  [W/m^(3)]
    % The center of the first cell is at dx/2 & the last cell is at L-dx/2.
    x = h/2 : h : L-(h/2);
    % Cross sectional area of the 1D domain
    A = 1;          %[m^2]
    % Left Boundary Condition
    T_a = 100;    %[\circ c]
    % Right Boundary Condition
    T_b = 200;    %[\circ c]
    % Initialising Temperature
    Nh = 2*N + 1; % the  number of node with temperature
    T=zeros(Nh,1);
    % Left Boundary Condition
    T(1)= T_a;    %[\circ c]
    % Right Boundary Condition
    T(Nh) = T_b;    %[\circ c]
    % Interior of Domain
    index_x = 2:Nh-1;
    %% Formation of Coefficient Matrices, Numerical Solution
    % Tolerence
    eps = 1.e-6;
    % Iteration
    It = 0;
    T_old = T;
    Error = 2*eps;
    % Calculation
    while   (Error > eps)
        It = It +1;
        T(index_x)= 0.5*(T(index_x-1) + T(index_x+1)+ (q*h*h/4)/(k)); % here divided by 4, cause h/2*h/2
        Error = max(abs(T(:)-T_old(:)));
        if any(isnan(T(:))) || any(isinf(T(:)))
            fprintf(' Iteration Diverge \n');
            return;
        end
        T_old = T;
    end
    T_Numeric = T(2:2:end);
    %% Analytical Results
    X_G = x;
    T_Analytical=((T_b-T_a)/L + (q/(2*k))*( L - X_G )).*(X_G) + T_a;%解析解
    %% Error Analysis
    M_Error(N)=0;
    for i=1:N
        % This x is only use to compare results on Versteeg book
        k=find(abs(X_G - x(i))<=1e-3);
        % No need to find k (You can but no Need)
        % Absolute Error
        M_Error(i) = abs(T_Numeric(i)-T_Analytical(i));
        % Relative Error
        Relative_Error(i) = abs(T_Numeric(i)-T_Analytical(i))/abs(T_Numeric(i));
        % Percentage Error
        Percentage_Error(i) = Relative_Error(i).*100;
    end
    %% Post Processing
    % Data Table Comparison of Different Results
    Nodes = 1 : N;
    fprintf('Node\tGrid_{Loc.}\t\tT_{Num.}\t\tT_{An.}\t\tError\n\n');
    fprintf('===\t\t=======\t\t=========\t\t=========\t\t=======\n');
    for i = 1: N;
        fprintf('%3.0f\t\t%3.5f\t\t%3.5f\t\t%3.5f\t\t%3.5f\n',[Nodes(i); X_G(i); T_Numeric(i); T_Analytical(i); Percentage_Error(i)]);
    end
    % Graphical Representation
    plot(100*x,T_Numeric,'-o','MarkerSize',5,'MarkerFaceColor','b');
    % hold on
    % plot([0, 100.*L], [T_a, T_b], 'd', 'LineWidth',3,'MarkerFaceColor','b');
    hold on
    plot(100*X_G,T_Analytical,'r--','LineWidth',2);
    % hold on
    % plot([0, 100.*L], [T_a, T_b], 'd', 'LineWidth',3,'MarkerFaceColor','r');
    % xlim([0  2]);
    % ylim([50 300]);
    %Boundary Points Inclusion
    xlabel('Distance (m)');
    ylabel('Temperature [\circ C]');
    legend('T_{Numeric}','T_{Analytical}');
    title ('Numerical & Analytical Results Comparison');
    grid on;
    %通过不断迭代误差控制在1e-6内
    % Number of Control Volumes
    N = 5;
    % Domain length
    L = 0.02; %[m]
    % Grid Spacing
    h = L/(N);
    % Diffusivity (Thermal conductivity)
    k = 0.5;  %[W/m-K]
    % Uniform Heat Generation, A Source Term
    q = 1000e3;     %  [W/m^(3)]
    % The center of the first cell is at dx/2 & the last cell is at L-dx/2.
    x = h/2 : h : L-(h/2);
    % Cross sectional area of the 1D domain
    A = 1;          %[m^2]
    % Left Boundary Condition
    T_a = 100;    %[\circ c]
    % Right Boundary Condition
    T_b = 200;    %[\circ c]
    % Initialising Temperature
    Nh = 2*N + 1; % the  number of node with temperature
    T=zeros(Nh,1);
    % Left Boundary Condition
    T(1)= T_a;    %[\circ c]
    % Right Boundary Condition
    T(Nh) = T_b;    %[\circ c]
    % Interior of Domain
    index_x = 2:Nh-1;
    %% Formation of Coefficient Matrices, Numerical Solution
    % Tolerence
    eps = 1.e-6;
    % Iteration
    It = 0;
    T_old = T;
    Error = 2*eps;
    % Calculation
    while   (Error > eps)
        It = It +1;
        T(index_x)= 0.5*(T(index_x-1) + T(index_x+1)+ (q*h*h/4)/(k)); % here divided by 4, cause h/2*h/2
        Error = max(abs(T(:)-T_old(:)));
        if any(isnan(T(:))) || any(isinf(T(:)))
            fprintf(' Iteration Diverge \n');
            return;
        end
        T_old = T;
    end
    T_Numeric = T(2:2:end);
    %% Analytical Results
    X_G = x;
    T_Analytical=((T_b-T_a)/L + (q/(2*k))*( L - X_G )).*(X_G) + T_a;%解析解
    %% 误差分析
    M_Error(N)=0;
    for i=1:N
        % This x is only use to compare results on Versteeg book
        k=find(abs(X_G - x(i))<=1e-3);
        % No need to find k (You can but no Need)
        % Absolute Error
        M_Error(i) = abs(T_Numeric(i)-T_Analytical(i));
        % Relative Error
        Relative_Error(i) = abs(T_Numeric(i)-T_Analytical(i))/abs(T_Numeric(i));
        % Percentage Error
        Percentage_Error(i) = Relative_Error(i).*100;
    end
    %% Post Processing
    % Data Table Comparison of Different Results
    Nodes = 1 : N;
    fprintf('Node\tGrid_{Loc.}\t\tT_{Num.}\t\tT_{An.}\t\tError\n\n');
    fprintf('===\t\t=======\t\t=========\t\t=========\t\t=======\n');
    for i = 1: N;
        fprintf('%3.0f\t\t%3.5f\t\t%3.5f\t\t%3.5f\t\t%3.5f\n',[Nodes(i); X_G(i); T_Numeric(i); T_Analytical(i); Percentage_Error(i)]);
    end
    % Graphical Representation
    plot(100*x,T_Numeric,'-o','MarkerSize',5,'MarkerFaceColor','b');
    % hold on
    % plot([0, 100.*L], [T_a, T_b], 'd', 'LineWidth',3,'MarkerFaceColor','b');
    hold on
    plot(100*X_G,T_Analytical,'r--','LineWidth',2);
    % hold on
    % plot([0, 100.*L], [T_a, T_b], 'd', 'LineWidth',3,'MarkerFaceColor','r');
    % xlim([0  2]);
    % ylim([50 300]);
    %Boundary Points Inclusion
    xlabel('Distance (m)');
    ylabel('Temperature [\circ C]');
    legend('T_{Numeric}','T_{Analytical}');
    title ('Numerical & Analytical Results Comparison');
    grid on;
    

    误差分析

    Node    Grid_{Loc.}        T_{Num.}        T_{An.}        Error
    
    ===        =======        =========        =========        =======
      1        0.00200        146.00000        146.00000        0.00000
      2        0.00600        213.99999        214.00000        0.00000
      3        0.01000        249.99999        250.00000        0.00000
      4        0.01400        253.99999        254.00000        0.00000
      5        0.01800        226.00000        226.00000        0.00000
    
    

    img

    本回答被题主选为最佳回答 , 对您是否有帮助呢?
    评论

报告相同问题?

问题事件

  • 系统已结题 10月1日
  • 已采纳回答 9月23日
  • 修改了问题 9月22日
  • 创建了问题 9月22日

悬赏问题

  • ¥15 fpga自动售货机数码管(相关搜索:数字时钟)
  • ¥20 Python安装cvxpy库出问题
  • ¥15 用前端向数据库插入数据,通过debug发现数据能走到后端,但是放行之后就会提示错误
  • ¥15 python天天向上类似问题,但没有清零
  • ¥30 3天&7天&&15天&销量如何统计同一行
  • ¥30 帮我写一段可以读取LD2450数据并计算距离的Arduino代码
  • ¥15 C#调用python代码(python带有库)
  • ¥15 活动选择题。最多可以参加几个项目?
  • ¥15 飞机曲面部件如机翼,壁板等具体的孔位模型
  • ¥15 vs2019中数据导出问题