diff --git a/five_point_stencil.ipynb b/five_point_stencil.ipynb index 099e014c0383328e9e37f5d7789f18358158f40b..ef94a3a96985b902cb13f7f4b6009ef5516735df 100644 --- a/five_point_stencil.ipynb +++ b/five_point_stencil.ipynb @@ -53,7 +53,7 @@ }, { "cell_type": "code", - "execution_count": 77, + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ @@ -70,7 +70,7 @@ }, { "cell_type": "code", - "execution_count": 78, + "execution_count": 5, "metadata": {}, "outputs": [], "source": [ @@ -85,82 +85,67 @@ ] }, { - "cell_type": "raw", + "cell_type": "code", + "execution_count": 6, "metadata": {}, + "outputs": [], "source": [ - "# just the initialisation of matrix B\n", "def matrix_B(h):\n", " N = int(1/h)\n", " m = pow(N-1,2)\n", " l = 4*N\n", " B = np.zeros((m,l))\n", + " for i in range(N-1):\n", + " B[(i+1)*(N-1)-1][2*(N+i)]= pow(h,-2)\n", + " B[(N-2)*(N-1)+i][N+2+2*i]=pow(h,-2)\n", " return B" ] }, { - "cell_type": "raw", + "cell_type": "code", + "execution_count": 7, "metadata": {}, + "outputs": [], "source": [ "def vector_g(h):\n", " N = int(1/h)\n", " l = 4*N\n", " g = np.zeros(l)\n", - " g[-1] = u(1,1)\n", + " for i in range(N):\n", + " g[3*N+i] = u((i+1)*h,1)\n", + " for i in range(N-1): \n", + " g[N+2+2*i] = u(1,(i+1)*h)\n", " return g " ] }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Since the exact solution is zero at the boundary, the product B*g is zero." - ] - }, { "cell_type": "code", - "execution_count": 115, + "execution_count": 64, "metadata": {}, "outputs": [], "source": [ - "n = 6\n", + "n = 2\n", "h = pow(2,-n)\n", "N = pow(2,n)" ] }, { "cell_type": "code", - "execution_count": 116, + "execution_count": 65, "metadata": {}, "outputs": [], "source": [ "A=matrix_A(h)\n", "f=vector_f(h)\n", - "appr_u = np.linalg.solve(A,f)" - ] - }, - { - "cell_type": "code", - "execution_count": 129, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "4096.0" - ] - }, - "execution_count": 129, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "A[3968,3905]" + "B=matrix_B(h)\n", + "g=vector_g(h)\n", + "RHS = f+np.matmul(B,g)\n", + "appr_u = np.linalg.solve(A,RHS)" ] }, { "cell_type": "code", - "execution_count": 130, + "execution_count": 66, "metadata": {}, "outputs": [], "source": [ @@ -176,7 +161,7 @@ }, { "cell_type": "code", - "execution_count": 131, + "execution_count": 67, "metadata": {}, "outputs": [], "source": [ @@ -185,7 +170,7 @@ }, { "cell_type": "code", - "execution_count": 132, + "execution_count": 68, "metadata": {}, "outputs": [], "source": [ @@ -194,16 +179,16 @@ }, { "cell_type": "code", - "execution_count": 133, + "execution_count": 69, "metadata": {}, "outputs": [ { "data": { "text/plain": [ - "0.015864236841443172" + "7.389346886114219" ] }, - "execution_count": 133, + "execution_count": 69, "metadata": {}, "output_type": "execute_result" } @@ -214,16 +199,16 @@ }, { "cell_type": "code", - "execution_count": 134, + "execution_count": 70, "metadata": {}, "outputs": [ { "data": { "text/plain": [ - "13.161396848144852" + "16.380407649947266" ] }, - "execution_count": 134, + "execution_count": 70, "metadata": {}, "output_type": "execute_result" }