4 elliptic

import%20marimo%0A%0A__generated_with%20%3D%20%220.9.27%22%0Aapp%20%3D%20marimo.App(width%3D%22medium%22)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%20%C2%A77%20%E6%A4%AD%E5%9C%86%E6%96%B9%E7%A8%8B%E7%9A%84%E5%B7%AE%E5%88%86%E6%A0%BC%E5%BC%8F%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%20%E5%87%86%E5%A4%87%E5%B7%A5%E4%BD%9C%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%20%20%20%20return%20(mo%2C)%0A%0A%0A%40app.cell%0Adef%20__()%3A%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20from%20numpy%20import%20linalg%2C%20newaxis%2C%20pi%2C%20sin%0A%0A%20%20%20%20np.set_printoptions(precision%3D3%2C%20suppress%3DTrue)%0A%20%20%20%20return%20linalg%2C%20newaxis%2C%20np%2C%20pi%2C%20sin%0A%0A%0A%40app.cell%0Adef%20__()%3A%0A%20%20%20%20from%20matplotlib%20import%20pyplot%20as%20plt%0A%20%20%20%20return%20(plt%2C)%0A%0A%0A%40app.cell%0Adef%20__()%3A%0A%20%20%20%20from%20util%20import%20show_files%2C%20typst%0A%20%20%20%20return%20show_files%2C%20typst%0A%0A%0A%40app.cell%0Adef%20__(Axes3D%2C%20Figure%2C%20np%2C%20plt)%3A%0A%20%20%20%20from%20matplotlib%20import%20cm%0A%0A%0A%20%20%20%20def%20plot_surface(%0A%20%20%20%20%20%20%20%20x%3A%20np.ndarray%2C%0A%20%20%20%20%20%20%20%20y%3A%20np.ndarray%2C%0A%20%20%20%20%20%20%20%20u%3A%20np.ndarray%2C%0A%20%20%20%20%20%20%20%20title%3A%20str%20%7C%20None%20%3D%20None%2C%0A%20%20%20%20%20%20%20%20invert_x_axis%3DTrue%2C%0A%20%20%20%20%20%20%20%20**kwargs%2C%0A%20%20%20%20)%20-%3E%20tuple%5BFigure%2C%20Axes3D%5D%3A%0A%20%20%20%20%20%20%20%20%22%22%22%E4%B8%AA%E4%BA%BA%E4%B9%A0%E6%83%AF%E7%89%88%60Axes3D.plot_surface%60%0A%0A%20%20%20%20%20%20%20%20Params%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20x%5B%23x%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B%23y%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20u%5B%23x%2C%20%23y%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20title%0A%20%20%20%20%20%20%20%20%20%20%20%20invert_x_axis%3A%20%E6%98%AF%E5%90%A6%E5%8F%8D%E8%BD%ACx%E8%BD%B4%E3%80%82%E6%9C%89%E6%97%B6%E5%8F%8D%E8%BD%AC%E5%8F%AF%E9%81%BF%E5%85%8D%E9%81%AE%E6%8C%A1%EF%BC%8C%E6%9B%B4%E6%B8%85%E6%99%B0%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20assert%20x.ndim%20%3D%3D%201%0A%20%20%20%20%20%20%20%20assert%20y.ndim%20%3D%3D%201%0A%20%20%20%20%20%20%20%20assert%20u.shape%20%3D%3D%20(x.size%2C%20y.size)%0A%0A%20%20%20%20%20%20%20%20ax%3A%20Axes3D%0A%20%20%20%20%20%20%20%20fig%2C%20ax%20%3D%20plt.subplots(%0A%20%20%20%20%20%20%20%20%20%20%20%20layout%3D%22constrained%22%2C%20subplot_kw%3D%7B%22projection%22%3A%20%223d%22%7D%0A%20%20%20%20%20%20%20%20)%20%20%23%20type%3A%20ignore%0A%20%20%20%20%20%20%20%20ax.plot_surface(%0A%20%20%20%20%20%20%20%20%20%20%20%20x%5B%3A%2C%20np.newaxis%5D%2C%20y%5Bnp.newaxis%2C%20%3A%5D%2C%20u%2C%20cmap%3Dcm.coolwarm%2C%20**kwargs%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20ax.xaxis.set_inverted(invert_x_axis)%0A%20%20%20%20%20%20%20%20ax.set_xlabel(%22%24x%24%22)%0A%20%20%20%20%20%20%20%20ax.set_ylabel(%22%24y%24%22)%0A%20%20%20%20%20%20%20%20ax.set_zlabel(%22%24u%24%22)%0A%0A%20%20%20%20%20%20%20%20if%20title%20is%20not%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20ax.set_title(title)%0A%0A%20%20%20%20%20%20%20%20return%20fig%2C%20ax%0A%20%20%20%20return%20cm%2C%20plot_surface%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%20%E9%97%AE%E9%A2%98%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(typst)%3A%0A%20%20%20%20typst(r%22%22%22%0A%20%20%20%20%23import%20%22%40preview%2Fphysica%3A0.9.3%22%3A%20laplacian%2C%20eval%0A%0A%20%20%20%20%24%0A%20%20%20%20-%20laplacian%20u%20%26%3D%20(pi%5E2%20-%201)%20e%5Ex%20sin(pi%20y)%2C%0A%20%20%20%20%20%20quad%20(x%2Cy)%20in%20%5B1%2C2%5D%20times%20%5B0%2C1%5D.%20%5C%0A%20%20%20%20eval(u)_(x%3D1)%20%26%3D%20e%20sin(pi%20y).%20%5C%0A%20%20%20%20eval(u)_(x%3D2)%20%26%3D%20e%5E2%20sin(pi%20y).%20%5C%0A%20%20%20%20eval(u)_(y%3D0)%20%26%3D%20eval(u)_(y%3D1)%20%3D%200.%0A%20%20%20%20%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__()%3A%0A%20%20%20%20dx%20%3D%200.1%0A%20%20%20%20dy%20%3D%200.05%0A%20%20%20%20return%20dx%2C%20dy%0A%0A%0A%40app.cell%0Adef%20__(dx%2C%20dy%2C%20np)%3A%0A%20%20%20%20x%20%3D%20np.arange(1%2C%202%20%2B%20dx%2C%20dx)%0A%20%20%20%20y%20%3D%20np.arange(0%2C%201%20%2B%20dy%2C%20dy)%0A%20%20%20%20return%20x%2C%20y%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%E4%B8%8D%E8%83%BD%E7%BC%93%E5%AD%98%E8%BF%99%E4%BA%9B%E5%87%BD%E6%95%B0%EF%BC%81%E5%9B%A0%E4%B8%BA%E5%8F%AF%E8%83%BD%E6%9C%89%E5%8E%9F%E5%9C%B0%E6%9B%B4%E6%94%B9%E7%9A%84%E9%9C%80%E6%B1%82%E3%80%82%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(np%2C%20pi%2C%20sin)%3A%0A%20%20%20%20def%20ref(x%3A%20np.ndarray%2C%20y%3A%20np.ndarray)%20-%3E%20np.ndarray%3A%0A%20%20%20%20%20%20%20%20%22%22%22%E7%9C%9F%E8%A7%A3%0A%0A%20%20%20%20%20%20%20%20Params%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20x%5B%23x%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B%23y%5D%0A%20%20%20%20%20%20%20%20Returns%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20u%5B%23x%2C%20%23y%5D%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20assert%20x.ndim%20%3D%3D%201%0A%20%20%20%20%20%20%20%20assert%20y.ndim%20%3D%3D%201%0A%20%20%20%20%20%20%20%20return%20np.exp(x%5B%3A%2C%20np.newaxis%5D)%20*%20sin(pi%20*%20y%5Bnp.newaxis%2C%20%3A%5D)%0A%20%20%20%20return%20(ref%2C)%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20ref%2C%20x%2C%20y)%3A%0A%20%20%20%20plot_surface(x%2C%20y%2C%20ref(x%2C%20y)%2C%20title%3D%22%E7%9C%9F%E8%A7%A3%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%200%20%E5%9F%BA%E7%A1%80%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(np%2C%20pi%2C%20ref)%3A%0A%20%20%20%20def%20ref_rhs(x%3A%20np.ndarray%2C%20y%3A%20np.ndarray)%20-%3E%20np.ndarray%3A%0A%20%20%20%20%20%20%20%20%22%22%22%E5%8F%B3%E7%AB%AF%E9%A1%B9%0A%0A%20%20%20%20%20%20%20%20Params%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20x%5B%23x%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B%23y%5D%0A%20%20%20%20%20%20%20%20Returns%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20f%5B%23x%2C%20%23y%5D%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20return%20-(pi**2%20-%201)%20*%20ref(x%2C%20y)%0A%20%20%20%20return%20(ref_rhs%2C)%0A%0A%0A%40app.cell%0Adef%20__(np%2C%20pi%2C%20sin)%3A%0A%20%20%20%20def%20setup_boundary(x%3A%20np.ndarray%2C%20y%3A%20np.ndarray)%20-%3E%20np.ndarray%3A%0A%20%20%20%20%20%20%20%20%22%22%22%E6%A0%B9%E6%8D%AE%E8%BE%B9%E7%95%8C%E6%9D%A1%E4%BB%B6%E5%87%86%E5%A4%87%E9%A2%84%E5%A4%87%E8%A7%A3%0A%0A%20%20%20%20%20%20%20%20Params%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20x%5B%23x%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20y%5B%23y%5D%0A%20%20%20%20%20%20%20%20Returns%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20u%5B%23x%2C%20%23y%5D%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20assert%20x.ndim%20%3D%3D%201%0A%20%20%20%20%20%20%20%20assert%20y.ndim%20%3D%3D%201%0A%0A%20%20%20%20%20%20%20%20u%20%3D%20np.zeros((x.size%2C%20y.size))%0A%20%20%20%20%20%20%20%20u%5B0%2C%20%3A%5D%20%3D%20np.exp(1)%20*%20sin(pi%20*%20y)%0A%20%20%20%20%20%20%20%20u%5B-1%2C%20%3A%5D%20%3D%20np.exp(2)%20*%20sin(pi%20*%20y)%0A%0A%20%20%20%20%20%20%20%20return%20u%0A%20%20%20%20return%20(setup_boundary%2C)%0A%0A%0A%40app.cell%0Adef%20__(np%2C%20ref%2C%20ref_rhs%2C%20setup_boundary)%3A%0A%20%20%20%20class%20Solver%3A%0A%20%20%20%20%20%20%20%20u%3A%20np.ndarray%0A%20%20%20%20%20%20%20%20x%3A%20np.ndarray%0A%20%20%20%20%20%20%20%20y%3A%20np.ndarray%0A%0A%20%20%20%20%20%20%20%20rhs%3A%20np.ndarray%0A%0A%20%20%20%20%20%20%20%20dx%3A%20float%0A%20%20%20%20%20%20%20%20dy%3A%20float%0A%0A%20%20%20%20%20%20%20%20def%20__init__(self%2C%20*%2C%20x%3A%20np.ndarray%2C%20y%3A%20np.ndarray)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20x.ndim%20%3D%3D%201%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20y.ndim%20%3D%3D%201%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20self.dx%20%3D%20np.diff(x).mean()%0A%20%20%20%20%20%20%20%20%20%20%20%20self.dy%20%3D%20np.diff(y).mean()%0A%20%20%20%20%20%20%20%20%20%20%20%20self.x%20%3D%20x%0A%20%20%20%20%20%20%20%20%20%20%20%20self.y%20%3D%20y%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20self.u%20%3D%20setup_boundary(x%2C%20y)%0A%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%20%3D%20ref_rhs(x%2C%20y)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20self.post_init()%0A%0A%20%20%20%20%20%20%20%20def%20post_init(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%22%22%22Prepare%20after%20%60__init__%60%22%22%22%0A%20%20%20%20%20%20%20%20%20%20%20%20pass%0A%0A%20%20%20%20%20%20%20%20def%20solve(self)%20-%3E%20None%3A%20...%0A%0A%20%20%20%20%20%20%20%20def%20validate(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20dv_x_2%20%3D%20(self.u%5B2%3A%2C%20%3A%5D%20%2B%20self.u%5B%3A-2%2C%20%3A%5D%20-%202%20*%20self.u%5B1%3A-1%2C%20%3A%5D)%5B%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%3A%2C%201%3A-1%0A%20%20%20%20%20%20%20%20%20%20%20%20%5D%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20dv_y_2%20%3D%20(self.u%5B%3A%2C%202%3A%5D%20%2B%20self.u%5B%3A%2C%20%3A-2%5D%20-%202%20*%20self.u%5B%3A%2C%201%3A-1%5D)%5B%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%201%3A-1%2C%20%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%5D%20%2F%20self.dy**2%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20dv_x_2%20%2B%20dv_y_2%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%23%20%E8%BF%99%E9%87%8C%E5%AE%9A%E5%BE%97%E5%BE%88%E6%9D%BE%EF%BC%8C%E8%BF%99%E6%A0%B7%E5%85%B6%E5%AE%83%E6%96%B9%E6%B3%95%E7%9A%84%E8%A7%A3%E4%B9%9F%E8%83%BD%E9%80%9A%E8%BF%87%E9%AA%8C%E8%AF%81%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20rtol%3D1e-2%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20def%20error(self)%20-%3E%20np.ndarray%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20return%20self.u%20-%20ref(self.x%2C%20self.y)%0A%0A%20%20%20%20%20%20%20%20def%20max_error(self)%20-%3E%20float%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20return%20abs(self.error()).max()%0A%20%20%20%20return%20(Solver%2C)%0A%0A%0A%40app.cell%0Adef%20__(Solver%2C%20ref%2C%20x%2C%20y)%3A%0A%20%20%20%20solver%20%3D%20Solver(x%3Dx%2C%20y%3Dy)%0A%20%20%20%20solver.u%20%3D%20ref(solver.x%2C%20solver.y)%0A%20%20%20%20solver.validate()%0A%20%20%20%20return%20(solver%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(Axes%2C%20Collection%2C%20Figure%2C%20Solver%2C%20mo%2C%20np%2C%20plt)%3A%0A%20%20%20%20from%20collections%20import%20deque%0A%0A%20%20%20%20from%20pandas%20import%20DataFrame%0A%20%20%20%20from%20seaborn%20import%20lineplot%0A%0A%0A%20%20%20%20def%20benchmark(%0A%20%20%20%20%20%20%20%20solver_cls%3A%20type%5BSolver%5D%2C%0A%20%20%20%20%20%20%20%20*%2C%0A%20%20%20%20%20%20%20%20dx_dy_list%3A%20Collection%5Btuple%5Bfloat%2C%20float%5D%5D%2C%0A%20%20%20%20)%20-%3E%20DataFrame%3A%0A%20%20%20%20%20%20%20%20%22%22%22Benchmark%0A%0A%20%20%20%20%20%20%20%20Returns%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%E5%88%97%E4%B8%BAdx%E3%80%81dy%E3%80%81%E6%9C%80%E5%A4%A7%E8%AF%AF%E5%B7%AE%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20assert%20issubclass(solver_cls%2C%20Solver)%0A%0A%20%20%20%20%20%20%20%20%23%20(dx%2C%20dy%2C%20max_error)%5B%5D%0A%20%20%20%20%20%20%20%20stat%3A%20deque%5Btuple%5Bfloat%2C%20float%2C%20float%5D%5D%20%3D%20deque()%0A%0A%20%20%20%20%20%20%20%20for%20dx%2C%20dy%20in%20mo.status.progress_bar(dx_dy_list)%3A%20%20%23%20type%3A%20ignore%0A%20%20%20%20%20%20%20%20%20%20%20%20dx%3A%20float%0A%20%20%20%20%20%20%20%20%20%20%20%20dy%3A%20float%0A%20%20%20%20%20%20%20%20%20%20%20%20x%20%3D%20np.arange(1%2C%202%20%2B%20dx%2C%20dx)%0A%20%20%20%20%20%20%20%20%20%20%20%20y%20%3D%20np.arange(0%2C%201%20%2B%20dy%2C%20dy)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20solver%20%3D%20solver_cls(x%3Dx%2C%20y%3Dy)%0A%20%20%20%20%20%20%20%20%20%20%20%20solver.solve()%0A%20%20%20%20%20%20%20%20%20%20%20%20stat.append((dx%2C%20dy%2C%20solver.max_error()))%0A%0A%20%20%20%20%20%20%20%20return%20DataFrame(%0A%20%20%20%20%20%20%20%20%20%20%20%20list(stat)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20columns%3D%5B%22dx%22%2C%20%22dy%22%2C%20%22%E6%9C%80%E5%A4%A7%E8%AF%AF%E5%B7%AE%22%5D%2C%0A%20%20%20%20%20%20%20%20)%0A%0A%0A%20%20%20%20def%20plot_benchmark(%0A%20%20%20%20%20%20%20%20data%3A%20DataFrame%2C%20title%3A%20str%20%7C%20None%20%3D%20None%0A%20%20%20%20)%20-%3E%20tuple%5BFigure%2C%20Axes%5D%3A%0A%20%20%20%20%20%20%20%20%22%22%22Plot%20the%20benchmark%20result%0A%0A%20%20%20%20%20%20%20%20Params%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%60df%60%3A%20Output%20of%20%60benchmark()%60%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20fig%2C%20axs%20%3D%20plt.subplots(nrows%3D2%2C%20layout%3D%22constrained%22)%0A%0A%20%20%20%20%20%20%20%20lineplot(ax%3Daxs%5B0%5D%2C%20data%3Ddata%2C%20x%3D%22dx%22%2C%20y%3D%22%E6%9C%80%E5%A4%A7%E8%AF%AF%E5%B7%AE%22%2C%20markers%3DTrue)%0A%20%20%20%20%20%20%20%20axs%5B0%5D.set_xlabel(r%22%24%5Cmathrm%7Bd%7Dx%24%22)%0A%20%20%20%20%20%20%20%20lineplot(ax%3Daxs%5B1%5D%2C%20data%3Ddata%2C%20x%3D%22dy%22%2C%20y%3D%22%E6%9C%80%E5%A4%A7%E8%AF%AF%E5%B7%AE%22%2C%20markers%3DTrue)%0A%20%20%20%20%20%20%20%20axs%5B1%5D.set_xlabel(r%22%24%5Cmathrm%7Bd%7Dy%24%22)%0A%0A%20%20%20%20%20%20%20%20for%20ax%20in%20axs%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20ax.set_xscale(%22log%22)%0A%20%20%20%20%20%20%20%20%20%20%20%20ax.set_yscale(%22log%22)%0A%20%20%20%20%20%20%20%20%20%20%20%20ax.grid(True)%0A%0A%20%20%20%20%20%20%20%20if%20title%20is%20not%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20fig.suptitle(title)%0A%0A%20%20%20%20%20%20%20%20return%20fig%2C%20axs%0A%20%20%20%20return%20DataFrame%2C%20benchmark%2C%20deque%2C%20lineplot%2C%20plot_benchmark%0A%0A%0A%40app.cell%0Adef%20__(np)%3A%0A%20%20%20%20benchmark_kwargs%20%3D%20dict(%0A%20%20%20%20%20%20%20%20dx_dy_list%3D%5B(_dx%2C%20_dx)%20for%20_dx%20in%202.0%20**%20np.arange(-4%2C%20-1%2C%201)%5D%2C%0A%20%20%20%20)%0A%20%20%20%20return%20(benchmark_kwargs%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%201%20%E4%BA%94%E7%82%B9%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%23%20%E5%BC%80%E5%8F%91%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver%2C%20linalg%2C%20np)%3A%0A%20%20%20%20class%20Solver_5(Solver)%3A%0A%20%20%20%20%20%20%20%20def%20post_init(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20x%2Cy%20without%20boundaries%20%E2%86%90%20x%2Cy%20with%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20t%20%3D%20np.zeros((self.x.size%2C%20self.y.size%2C%20self.x.size%2C%20self.y.size))%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20i_x%20in%20range(1%2C%20self.x.size%20-%201)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20all_y%20%3D%20np.arange(1%2C%20self.y.size%20-%201)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20-%201%2C%20all_y%5D%20%3D%201%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%2C%20all_y%5D%20%3D%20-2%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20%2B%201%2C%20all_y%5D%20%3D%201%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20i_y%20in%20range(1%2C%20self.y.size%20-%201)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20all_x%20%3D%20np.arange(1%2C%20self.x.size%20-%201)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Ball_x%2C%20i_y%2C%20all_x%2C%20i_y%20-%201%5D%20%2B%3D%201%20%2F%20self.dy**2%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Ball_x%2C%20i_y%2C%20all_x%2C%20i_y%5D%20%2B%3D%20-2%20%2F%20self.dy**2%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Ball_x%2C%20i_y%2C%20all_x%2C%20i_y%20%2B%201%5D%20%2B%3D%201%20%2F%20self.dy**2%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20x%2Cy%20without%20boundaries%20%E2%86%90%20x%2Cy%20with%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%20%3D%20t%5B1%3A-1%2C%201%3A-1%2C%20...%5D%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20boundary%20terms%20affecting%20x%2Cy%20without%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20self.boundary%20%3D%20np.einsum(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22xyuv%2Cuv-%3Exy%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%23%20Select%20x_min%20and%20x_max%2C%20including%20the%20y%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%5B...%2C%20%3A%3A%20self.x.size%20-%201%2C%20%3A%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B%3A%3A%20self.x.size%20-%201%2C%20%3A%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%20%2B%20np.einsum(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22xyuv%2Cuv-%3Exy%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%23%20Select%20y_min%20and%20y_max%2C%20excluding%20the%20x%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%5B...%2C%201%3A-1%2C%20%3A%3A%20self.y.size%20-%201%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B1%3A-1%2C%20%3A%3A%20self.y.size%20-%201%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20def%20solve(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20All%20to%20RHS%20%2B%20flatten%0A%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B1%3A-1%2C%201%3A-1%5D.flat%20%3D%20linalg.solve(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%5B...%2C%201%3A-1%2C%201%3A-1%5D.reshape(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(-1%2C%20(self.x.size%20-%202)%20*%20(self.y.size%20-%202))%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(self.rhs%5B1%3A-1%2C%201%3A-1%5D%20-%20self.boundary).flat%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20def%20validate(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Original%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20np.einsum(%22xyuv%2Cuv-%3Exy%22%2C%20self.laplacian%2C%20self.u)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Flatten%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian.reshape((-1%2C%20self.x.size%20*%20self.y.size))%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%40%20self.u.flat%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D.flat%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Apart%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20np.einsum(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22xyuv%2Cuv-%3Exy%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%5B...%2C%201%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.boundary%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20All%20to%20RHS%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20np.einsum(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%22xyuv%2Cuv-%3Exy%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%5B...%2C%201%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%20-%20self.boundary%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20All%20to%20RHS%20%2B%20flatten%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%5B...%2C%201%3A-1%2C%201%3A-1%5D.reshape(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(-1%2C%20(self.x.size%20-%202)%20*%20(self.y.size%20-%202))%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%40%20self.u%5B1%3A-1%2C%201%3A-1%5D.flat%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(self.rhs%5B1%3A-1%2C%201%3A-1%5D%20-%20self.boundary).flat%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20return%20(Solver_5%2C)%0A%0A%0A%40app.cell%0Adef%20__(Solver_5%2C%20np)%3A%0A%20%20%20%20_x%20%3D%20np.arange(7)%0A%20%20%20%20_from%20%3D%20Solver_5(x%3D_x%2C%20y%3D_x).laplacian%5B2%2C%202%2C%20%3A%5D%0A%20%20%20%20print(_from)%0A%0A%20%20%20%20assert%20np.count_nonzero(_from)%20%3D%3D%205%0A%20%20%20%20assert%20np.isclose(_from.sum()%2C%200)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%23%20%E5%8D%95%E6%AC%A1%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver_5%2C%20x%2C%20y)%3A%0A%20%20%20%20solver_5%20%3D%20Solver_5(x%3Dx%2C%20y%3Dy)%0A%20%20%20%20solver_5.solve()%0A%20%20%20%20solver_5.validate()%0A%20%20%20%20return%20(solver_5%2C)%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20solver_5%2C%20x%2C%20y)%3A%0A%20%20%20%20plot_surface(x%2C%20y%2C%20solver_5.u%2C%20title%3D%22%E8%BF%91%E4%BC%BC%E8%A7%A3%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20solver_5%2C%20x%2C%20y)%3A%0A%20%20%20%20plot_surface(x%2C%20y%2C%20solver_5.error()%2C%20title%3D%22%E8%AF%AF%E5%B7%AE%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(solver_5)%3A%0A%20%20%20%20solver_5.max_error()%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%23%20%E7%BB%9F%E8%AE%A1%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver_5%2C%20benchmark%2C%20benchmark_kwargs%2C%20plot_benchmark)%3A%0A%20%20%20%20_b%20%3D%20benchmark(Solver_5%2C%20**benchmark_kwargs)%0A%20%20%20%20plot_benchmark(_b%2C%20title%3D%22%E4%BA%94%E7%82%B9%22)%5B0%5D%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20r%22%22%22%0A%20%20%20%20%20%20%20%20%23%23%202%20%E4%B9%9D%E7%82%B9%0A%0A%20%20%20%20%20%20%20%20%E8%BF%99%E5%9B%9E%E7%94%A8%E7%A8%80%E7%96%8F%E7%9F%A9%E9%98%B5%E8%AF%95%E8%AF%95%E3%80%82%0A%0A%20%20%20%20%20%20%20%20%3E%20%5BSparse%20arrays%20currently%20must%20be%20two-dimensional.%5D(https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fscipy%2Freference%2Fsparse.html)%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%23%20%E5%BC%80%E5%8F%91%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__()%3A%0A%20%20%20%20from%20scipy.sparse%20import%20diags_array%2C%20eye_array%2C%20kron%0A%20%20%20%20from%20scipy.sparse.linalg%20import%20spsolve%0A%20%20%20%20return%20diags_array%2C%20eye_array%2C%20kron%2C%20spsolve%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(typst)%3A%0A%20%20%20%20typst(r%22%22%22%0A%20%20%20%20%23import%20%22%40preview%2Fphysica%3A0.9.3%22%3A%20pdv%2C%20laplacian%0A%0A%20%20%20%20%E9%A6%96%E5%85%88%E8%A6%81%E5%A4%84%E7%90%86%E4%B8%80%E4%B8%8B%20%24h_x%20!%3D%20h_y%24%20%E7%9A%84%E9%97%AE%E9%A2%98%E3%80%82%0A%0A%20%20%20%20%23let%20xx%20%3D%20%24cal(X)%24%0A%20%20%20%20%23let%20yy%20%3D%20%24cal(Y)%24%0A%0A%20%20%20%20%E8%AE%BE%20%24xx%20%3A%3D%20h_x%20pdv(%2Cx)%2C%20thick%20yy%20%3A%3D%20h_y%20pdv(%2Cy)%24%EF%BC%8C%E8%A6%81%E5%87%91%E5%87%BA%0A%20%20%20%20%24%20laplacian%20%3D%20(xx%20%2F%20h_x)%5E2%20%2B%20(yy%20%2F%20h_y)%5E2.%20%24%0A%0A%20%20%20%20%E6%B3%A8%E6%84%8F%0A%20%20%20%20%24%0A%20%20%20%20%20%20e%5Exx%20-%20e%5E(-xx)%20%26%3D%202%20cosh%20xx%20%26~%202%20(1%20%2B%20xx%5E2%2F2).%20%5C%0A%20%20%20%20%20%20e%5Eyy%20-%20e%5E(-yy)%20%26%3D%202%20cosh%20yy%20%26~%202%20(1%20%2B%20yy%5E2%2F2).%20%5C%0A%0A%20%20%20%20%20%20e%5E(xx%2Byy)%20%2B%20e%5E(xx-yy)%20%2B%20e%5E(-xx%2Byy)%20%2B%20e%5E(-xx-yy)%0A%20%20%20%20%20%20%26%3D%20(e%5Exx%20%2B%20e%5E(-xx))%20(e%5Eyy%20%2B%20e%5E(-yy))%0A%20%20%20%20%20%20%3D%204%20cosh%20xx%20cosh%20yy%0A%20%20%20%20%20%20%26~%204%20(1%20%2B%20xx%5E2%20%2F%202%20%2B%20yy%5E2%20%2F2).%20%5C%0A%20%20%20%20%24%0A%0A%20%20%20%20%E7%AD%89%E4%B8%80%E4%B8%8B%EF%BC%8C%E7%9C%9F%E8%83%BD%E5%87%91%E5%87%BA%E6%9D%A5%E5%90%97%EF%BC%9F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver%2C%20diags_array%2C%20eye_array%2C%20kron%2C%20np%2C%20pi%2C%20spsolve)%3A%0A%20%20%20%20class%20Solver_9(Solver)%3A%0A%20%20%20%20%20%20%20%20def%20post_init(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.isclose(self.dx%2C%20self.dy)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20without%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20n_x%20%3D%20self.x.size%20-%202%0A%20%20%20%20%20%20%20%20%20%20%20%20n_y%20%3D%20self.y.size%20-%202%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20center%20%3D%20-10%20%2F%203%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20edge%20%3D%202%20%2F%203%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20corner%20%3D%201%20%2F%206%20%2F%20self.dx**2%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20xy%20without%20boundaries%20%E2%86%90%20xy%20without%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%20%3D%20(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20kron(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%23%20x_out%20%3D%20x_in%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20eye_array(n_x)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20diags_array(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5Bcenter%2C%20edge%2C%20edge%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20offsets%3D%5B0%2C%20-1%2C%201%5D%2C%20%20%23%20y_out%20-%20y_in%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20shape%3D(n_y%2C%20n_y)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20kron(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%23%20%23%20x_out%20%3D%20x_in%20%C2%B1%201%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20diags_array(%5B1%2C%201%5D%2C%20offsets%3D%5B-1%2C%201%5D%2C%20shape%3D(n_x%2C)%20*%202)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20diags_array(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%5Bedge%2C%20corner%2C%20corner%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20offsets%3D%5B0%2C%20-1%2C%201%5D%2C%20%20%23%20y_out%20-%20y_in%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20shape%3D(n_y%2C%20n_y)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20).tocsc()%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20To%20perform%20inversion%2C%20first%20convert%20to%20either%20CSC%20or%20CSR%20format.%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20boundary%20terms%20affecting%20x%2Cy%20without%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20self.boundary%20%3D%20np.zeros((n_x%2C%20n_y))%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Select%20x_min%20and%20x_max%0A%20%20%20%20%20%20%20%20%20%20%20%20self.boundary%5B%3A%3A%20n_x%20-%201%2C%20%3A%5D%20%2B%3D%20(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20edge%20*%20self.u%5B%3A%3A%20self.x.size%20-%201%2C%201%3A-1%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20corner%20*%20self.u%5B%3A%3A%20self.x.size%20-%201%2C%20%3A-2%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20corner%20*%20self.u%5B%3A%3A%20self.x.size%20-%201%2C%202%3A%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Select%20y_min%20and%20y_max%0A%20%20%20%20%20%20%20%20%20%20%20%20self.boundary%5B%3A%2C%20%3A%3A%20n_y%20-%201%5D%20%2B%3D%20(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20edge%20*%20self.u%5B1%3A-1%2C%20%3A%3A%20self.y.size%20-%201%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20corner%20*%20self.u%5B%3A-2%2C%20%3A%3A%20self.y.size%20-%201%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20corner%20*%20self.u%5B2%3A%2C%20%3A%3A%20self.y.size%20-%201%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20RHS%20%E8%A6%81%E6%94%B9%E6%88%90%20f%20%2B%201%2F12%20(h%20%E2%88%87)%C2%B2%20f%20%3D%20(1%20%2B%20h%C2%B2%2F12%20(1-%CF%80%C2%B2))%20f%0A%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%20*%3D%201%20%2B%20self.dx**2%20%2F%2012%20*%20(1%20-%20pi**2)%0A%0A%20%20%20%20%20%20%20%20def%20solve(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B1%3A-1%2C%201%3A-1%5D.flat%20%3D%20spsolve(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.laplacian%2C%20(self.rhs%5B1%3A-1%2C%201%3A-1%5D%20-%20self.boundary).flat%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20def%20ref_laplacian(self)%20-%3E%20np.ndarray%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20center%20%3D%20-10%20%2F%203%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20edges%20%3D%202%20%2F%203%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20corners%20%3D%201%20%2F%206%20%2F%20self.dx**2%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20x%2Cy%20without%20boundaries%20%E2%86%90%20x%2Cy%20with%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20t%20%3D%20np.zeros((self.x.size%2C%20self.y.size%2C%20self.x.size%2C%20self.y.size))%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20i_x%20in%20range(1%2C%20self.x.size%20-%201)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20all_y%20%3D%20np.arange(1%2C%20self.y.size%20-%201)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%2C%20all_y%5D%20%3D%20center%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20-%201%2C%20all_y%5D%20%3D%20edges%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20%2B%201%2C%20all_y%5D%20%3D%20edges%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%2C%20all_y%20-%201%5D%20%3D%20edges%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%2C%20all_y%20%2B%201%5D%20%3D%20edges%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20-%201%2C%20all_y%20%2B%201%5D%20%3D%20corners%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20-%201%2C%20all_y%20-%201%5D%20%3D%20corners%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20%2B%201%2C%20all_y%20%2B%201%5D%20%3D%20corners%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20t%5Bi_x%2C%20all_y%2C%20i_x%20%2B%201%2C%20all_y%20-%201%5D%20%3D%20corners%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20x%2Cy%20without%20boundaries%20%E2%86%90%20x%2Cy%20with%20boundaries%0A%20%20%20%20%20%20%20%20%20%20%20%20return%20t%5B1%3A-1%2C%201%3A-1%2C%20...%5D%0A%0A%20%20%20%20%20%20%20%20def%20validate(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20laplacian%20%3D%20self.laplacian.toarray().reshape(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.x.size%20-%202%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.y.size%20-%202%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.x.size%20-%202%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.y.size%20-%202%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20ref_laplacian%20%3D%20self.ref_laplacian()%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(laplacian%2C%20ref_laplacian%5B...%2C%201%3A-1%2C%201%3A-1%5D)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Apart%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20np.einsum(%22xyuv%2Cuv-%3Exy%22%2C%20laplacian%2C%20self.u%5B1%3A-1%2C%201%3A-1%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.boundary%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Original%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20np.einsum(%22xyuv%2Cuv-%3Exy%22%2C%20ref_laplacian%2C%20self.u)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20return%20(Solver_9%2C)%0A%0A%0A%40app.cell%0Adef%20__(Solver_9%2C%20np)%3A%0A%20%20%20%20_x%20%3D%20np.arange(7)%0A%20%20%20%20_l%20%3D%20Solver_9(x%3D_x%2C%20y%3D_x).laplacian.toarray().reshape((_x.size%20-%202%2C)%20*%204)%0A%20%20%20%20print(_l%5B1%2C%200%2C%20...%5D)%0A%0A%20%20%20%20_from%20%3D%20_l%5B1%2C%202%2C%20...%5D%0A%20%20%20%20print(_from)%0A%20%20%20%20assert%20np.count_nonzero(_from)%20%3D%3D%209%0A%20%20%20%20assert%20np.isclose(_from.sum()%2C%200)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%23%20%E5%8D%95%E6%AC%A1%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver_9%2C%20x%2C%20y)%3A%0A%20%20%20%20solver_9%20%3D%20Solver_9(x%3Dy%20-%20y%5B0%5D%20%2B%20x%5B0%5D%2C%20y%3Dy)%0A%20%20%20%20solver_9.solve()%0A%20%20%20%20solver_9.validate()%0A%20%20%20%20return%20(solver_9%2C)%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20solver_9)%3A%0A%20%20%20%20plot_surface(solver_9.x%2C%20solver_9.y%2C%20solver_9.u%2C%20title%3D%22%E8%BF%91%E4%BC%BC%E8%A7%A3%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20solver_9)%3A%0A%20%20%20%20plot_surface(solver_9.x%2C%20solver_9.y%2C%20solver_9.error()%2C%20title%3D%22%E8%AF%AF%E5%B7%AE%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(solver_9)%3A%0A%20%20%20%20solver_9.max_error()%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%23%20%E7%BB%9F%E8%AE%A1%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver_9%2C%20benchmark%2C%20benchmark_kwargs%2C%20plot_benchmark)%3A%0A%20%20%20%20_b%20%3D%20benchmark(Solver_9%2C%20**benchmark_kwargs)%0A%20%20%20%20plot_benchmark(_b%2C%20title%3D%22%E4%B9%9D%E7%82%B9%22)%5B0%5D%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20r%22%22%22%0A%20%20%20%20%20%20%20%20%23%23%20%E7%95%AA%E5%A4%96%0A%0A%20%20%20%20%20%20%20%20%E5%88%98%E5%8A%9B%E8%B1%AA%E5%90%8C%E5%AD%A6%E5%91%8A%E8%AF%89%E6%88%91%EF%BC%9A%E4%B8%8D%E5%86%99%E7%B3%BB%E6%95%B0%E7%9F%A9%E9%98%B5%E4%B9%9F%E8%83%BD%E8%A7%A3%E4%BA%94%E7%82%B9%E6%88%96%E4%B9%9D%E7%82%B9%E6%A0%BC%E5%BC%8F%E2%80%94%E2%80%94%E8%BF%AD%E4%BB%A3%E6%B1%82%E8%A7%A3%EF%BC%81%E8%BF%99%E5%A4%A7%E7%BA%A6%E5%BD%92%E5%8A%9F%E4%BA%8E%E6%A4%AD%E5%9C%86%E6%96%B9%E7%A8%8B%E7%9A%84%E4%BC%98%E8%89%AF%E6%80%A7%E8%B4%A8%E3%80%82%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(Solver_9%2C%20deque%2C%20mo%2C%20np%2C%20pi)%3A%0A%20%20%20%20class%20Solver_9_NoMatrix(Solver_9)%3A%0A%20%20%20%20%20%20%20%20def%20post_init(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.isclose(self.dx%2C%20self.dy)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20RHS%20%E8%A6%81%E6%94%B9%E6%88%90%20f%20%2B%201%2F12%20(h%20%E2%88%87)%C2%B2%20f%20%3D%20(1%20%2B%20h%C2%B2%2F12%20(1-%CF%80%C2%B2))%20f%0A%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%20*%3D%201%20%2B%20self.dx**2%20%2F%2012%20*%20(1%20-%20pi**2)%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20(difference%2C%20error)%5B%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20self.steps%20%3D%20deque()%0A%0A%20%20%20%20%20%20%20%20def%20solve(self%2C%20*%2C%20n_steps%3A%20int)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20center%20%3D%20-10%20%2F%203%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20edge%20%3D%202%20%2F%203%20%2F%20self.dx**2%0A%20%20%20%20%20%20%20%20%20%20%20%20corner%20%3D%201%20%2F%206%20%2F%20self.dx**2%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20center%20%2B%20edges%20%2B%20corners%20%3D%20rhs%0A%0A%20%20%20%20%20%20%20%20%20%20%20%20for%20_%20in%20mo.status.progress_bar(range(n_steps))%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20next_u%20%3D%20(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20-%20edge%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20*%20(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B%3A-2%2C%201%3A-1%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.u%5B2%3A%2C%201%3A-1%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.u%5B1%3A-1%2C%20%3A-2%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.u%5B1%3A-1%2C%202%3A%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20-%20corner%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20*%20(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B%3A-2%2C%20%3A-2%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.u%5B2%3A%2C%20%3A-2%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.u%5B%3A-2%2C%202%3A%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%2B%20self.u%5B2%3A%2C%202%3A%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%20%2F%20center%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.steps.append(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(abs(self.u%5B1%3A-1%2C%201%3A-1%5D%20-%20next_u).max()%2C%20self.max_error())%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.u%5B1%3A-1%2C%201%3A-1%5D%20%3D%20next_u%0A%0A%20%20%20%20%20%20%20%20def%20validate(self)%20-%3E%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20Original%0A%20%20%20%20%20%20%20%20%20%20%20%20assert%20np.allclose(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20np.einsum(%22xyuv%2Cuv-%3Exy%22%2C%20self.ref_laplacian()%2C%20self.u)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20self.rhs%5B1%3A-1%2C%201%3A-1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20return%20(Solver_9_NoMatrix%2C)%0A%0A%0A%40app.cell%0Adef%20__(Solver_9_NoMatrix%2C%20x%2C%20y)%3A%0A%20%20%20%20solver_liu%20%3D%20Solver_9_NoMatrix(x%3Dy%20-%20y%5B0%5D%20%2B%20x%5B0%5D%2C%20y%3Dy)%0A%20%20%20%20solver_liu.solve(n_steps%3D1234)%0A%20%20%20%20solver_liu.validate()%0A%20%20%20%20return%20(solver_liu%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(plt%2C%20solver_5%2C%20solver_9%2C%20solver_liu)%3A%0A%20%20%20%20_fig%2C%20_ax%20%3D%20plt.subplots()%0A%20%20%20%20_ax.semilogy(%0A%20%20%20%20%20%20%20%20list(solver_liu.steps)%2C%0A%20%20%20%20%20%20%20%20label%3D%5B%22%E8%BF%AD%E4%BB%A3%E5%89%8D%E5%90%8E%E6%9C%80%E5%A4%A7%E6%94%B9%E5%8F%98%E9%87%8F%22%2C%20%22%E6%9C%80%E5%A4%A7%E8%AF%AF%E5%B7%AE%EF%BC%88%E4%B8%8E%E7%9C%9F%E8%A7%A3%E7%9B%B8%E6%AF%94%EF%BC%8C%E4%B8%8B%E5%90%8C%EF%BC%89%22%5D%2C%0A%20%20%20%20)%0A%20%20%20%20_ax.set(xlabel%3D%22%E8%BF%AD%E4%BB%A3%E6%AC%A1%E6%95%B0%22)%0A%20%20%20%20_ax.grid(True)%0A%20%20%20%20_ax.axhline(%0A%20%20%20%20%20%20%20%20y%3Dsolver_liu.max_error()%2C%0A%20%20%20%20%20%20%20%20linestyle%3D%22--%22%2C%0A%20%20%20%20%20%20%20%20color%3D%22orange%22%2C%0A%20%20%20%20%20%20%20%20label%3D%22%E6%9C%80%E7%BB%88%E6%9C%80%E5%A4%A7%E8%AF%AF%E5%B7%AE%22%2C%0A%20%20%20%20)%0A%20%20%20%20_ax.axhline(%0A%20%20%20%20%20%20%20%20y%3Dsolver_9.max_error()%2C%0A%20%20%20%20%20%20%20%20linestyle%3D%22%3A%22%2C%0A%20%20%20%20%20%20%20%20color%3D%22green%22%2C%0A%20%20%20%20%20%20%20%20label%3D%22%E8%A7%A3%E6%96%B9%E7%A8%8B%E6%B3%95%E7%9A%84%E8%AF%AF%E5%B7%AE%22%2C%0A%20%20%20%20)%0A%0A%20%20%20%20_ax.axhline(%0A%20%20%20%20%20%20%20%20y%3Dsolver_5.max_error()%2C%0A%20%20%20%20%20%20%20%20linestyle%3D%22%3A%22%2C%0A%20%20%20%20%20%20%20%20color%3D%22purple%22%2C%0A%20%20%20%20%20%20%20%20label%3D%22%E4%BA%94%E7%82%B9%E6%A0%BC%E5%BC%8F%E7%9A%84%E8%AF%AF%E5%B7%AE%22%2C%0A%20%20%20%20)%0A%20%20%20%20_ax.legend()%0A%20%20%20%20_fig%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20solver_liu)%3A%0A%20%20%20%20plot_surface(solver_liu.x%2C%20solver_liu.y%2C%20solver_liu.u%2C%20title%3D%22%E8%BF%91%E4%BC%BC%E8%A7%A3%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(plot_surface%2C%20solver_liu)%3A%0A%20%20%20%20plot_surface(solver_liu.x%2C%20solver_liu.y%2C%20solver_liu.error()%2C%20title%3D%22%E8%AF%AF%E5%B7%AE%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20__(solver_liu)%3A%0A%20%20%20%20solver_liu.max_error()%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%E8%AF%9A%E4%B8%8D%E6%88%91%E6%AC%BA%EF%BC%81%E5%B0%B1%E6%98%AF%E6%AF%8F%E6%AC%A1%E8%BF%AD%E4%BB%A3%E5%8F%AA%E8%83%BD%E4%BC%A0%E9%80%92%E4%B8%80%E6%A0%BC%EF%BC%8C%E6%94%B6%E6%95%9B%E5%A4%AA%E6%85%A2%E4%BA%86%EF%BC%9B%E8%80%8C%E4%B8%94%E9%9A%BE%E4%BB%A5%E4%BB%8E%E8%BF%AD%E4%BB%A3%E6%94%B9%E5%8F%98%E9%87%8F%E4%BC%B0%E8%AE%A1%E8%AF%AF%E5%B7%AE%EF%BC%8C%E4%B8%8D%E5%A5%BD%E7%A1%AE%E5%AE%9A%E8%BF%AD%E4%BB%A3%E5%A4%9A%E5%B0%91%E6%AC%A1%E3%80%82%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%23%23%20%E9%99%84%E5%BD%95%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20__(show_files)%3A%0A%20%20%20%20show_files(%5B%22pyproject.toml%22%2C%20%22util.py%22%5D)%0A%20%20%20%20return%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A