WENO3

  • WENO3

\[\begin{align} \phi^-_{x,i}&=\frac{1}{2h}(\Delta^+\phi_{i-1}+\Delta^+ \phi_i)-\frac{1}{2h}\omega_{-}(\Delta^+ \phi_{i-2}-2\Delta^+\phi_{i-1}+\Delta^+ \phi_i)\\ \omega^- &=\frac{1}{1+2r^2_{-}}\\ r_{-} &=\frac{\epsilon+(\Delta^-\Delta^+ \phi_{i-1})^2}{\epsilon+(\Delta^-\Delta^+ \phi_i)^2}\\ \phi^+_{x,i}&=\frac{1}{2h}(\Delta^+ \phi_{i-1}+\Delta^+\phi_i)-\frac{1}{2h}\omega^+(\Delta^+\phi_{i+1} -2\Delta^+\phi_i+\Delta^+\phi_{i-1})\\ \omega_{+} &=\frac{1}{1+2r^2_{+}}\\ r_{+} &=\frac{\epsilon+(\Delta^-\Delta^+ \phi_{i+1})^2}{\epsilon+(\Delta^-\Delta^+ \phi_{i})^2} \end{align} \]

  • 1D WENO Type Extrapolation
    Assume that we have a stencil of three points \(x_0=0,x_1=h,x_2=2h\) with point values \(u_j,j=0,1,2\), We aim to obtain a \(3-k\)th order approximation of
    \(\displaystyle \frac{d^k u}{dx^k}|_{x=-h/2},k=0,1,2\) denoted by $ u^{*k}$. We have three candidate substencils given by

\[S_r=\{x_0,..x_r\},r=0,1,2 \]

On each stencil \(S_r\) we have a \(r\) degree \(p_r(x)\)

\[ \begin{align} p_0(x) &=u_0\\ p_1(x) &=u_0+\frac{1}{h} x(u_1-u_0)\\ p_2(x) &=u_0+\frac{1}{2h}x(-3u_0+4u_1-u_2)+\frac{1}{2h^2}x^2(u_0-2u_1+u_2) \end{align} \]

Suppose \(u(x)\) is smooth on \(S_r\), the \(u^{*k}\) can be extrapolated by

\[ \begin{align} u^{*k} &=\sum^2_{r=0}d_rp^{(k)}(x)|_{x=\frac{-h}{2}}\\ d_0 &=h^2 \\ d_1 &=h\\ d_2 &=1-h-h^2 \end{align} \]

We now look for WENO type extrapolation in the form

\[ \begin{align} u^{*k} &=\sum^2_{r=0}\omega_rp_r^{(k)}(x)|_{x=\frac{-h}{2}}\\ \beta_0 &=h^2\\ \beta_1 &=\sum_{i=1}^2 \int_{-h}^0 h^{(2i-1)}(p^{(i)}_1(x))^2 dx=(u_1-u_0)^2 \\ \beta_2 &=\sum_{i=1}^2 \int_{-h}^0 h^{(2i-1)}(p^{(i)}_2(x))^2 dx\\ &=\frac{1}{12}(61u_0^2+160u_1^2+74u_0u_2+25u_2^2-196u_0u_1+124u_1u_2)\\ \alpha_r&=\frac{d_r}{(\epsilon+\beta_r)^2}\\ \omega_r&=\frac{\alpha_r}{\sum^2_{s=0} \alpha_s}\\ d_0 &=h^2 \\ d_1 &=h\\ d_2 &=1-h-h^2 \end{align} \]

充分证明这一技术是可靠的,成功的。但是,我的边界是光滑的,使用Lagrange , Hermit 插值是可以的。

function T81
N=500;
eps=1/10;
a=60 ;
x=linspace(0,2*pi,N)';
r0=(x(1)+x(end))/2;
x0=r0;
h=x(2)-x(1);
df=@(r)1./(power(eps,2) + power(-r0 + r,2));
mean_df=@(r)df(r)./sqrt(df(r).^2+1);
d1u=df(x);
ss=@(x)((-2*x + 2*x0)./(sqrt(1 + power(power(eps,2) + power(x - x0,2),-2)).*...
     (1+power(eps,4)+power(x,4) + 2*power(eps,2).*power(x - x0,2) - 4*power(x,3)*x0 +... 
       6*power(x,2).*power(x0,2) - 4*x.*power(x0,3) + power(x0,4))));
An=-(atan((r0 - x)/eps)/eps);
S=-ss(x);
D_BC=[mean_df(x(1:2)-2*h),mean_df(x(end-1:end)+2*h)];
N_BC=[df(x(1)),df(x(end))];
u=An*0.99;
t=0;
dt=0.5*h^2;
t_end=2000;

%============= Runge-Kutta ================
while t<t_end
    u1=u+dt*L(u);
    u2=3/4*u+1/4*u1+1/4*dt*L(u1);
    u=1/3*u+2/3*u2+2/3*dt*L(u2);
    t=t+dt;
    subplot(1,2,1)
    plot(x(1:5:end),u(1:5:end),'b.',x,An,'r')
    subplot(1,2,2)
    plot(x,u-An)
    title(['t=',num2str(t)])
    drawnow
end

function y=L(u)
        %-----1.内部使用中心差分-------
         du=D3(x,u,N_BC,'N');
        % du=WENO3(x,u,1,'no');
         p=du./sqrt(1+d1u.^2);
        %-----2. 外部使用WENO5---------
         d2u_p=WENO3(x,p,1,'no');
         d2u_m=WENO3(x,p,-1,'no');
           d2u=(d2u_p+d2u_m)/2;
        %  d2u=d2u_m;
         y=d2u+S;
    end
end
posted @ 2019-10-06 23:16  yuewen_chen  阅读(214)  评论(0编辑  收藏  举报