{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Chapter 28: Numerical Solution Of Partial Differential Equations"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.1, page no. 725"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"D=Bˆ2−4AC\n",
"if D<0 then elliptic if D=0 then parabolic if D>0 then hyperboic\n",
"(i) A=xˆ2, B1−yˆ2 D=4ˆ2−4∗1∗4=0 so the equation is PARABOLIC \n",
"(ii) D=4xˆ2(yˆ2−1)\n",
"for −inf<x<inf and −1<y<1 D<0 \n",
"So the equation is ELLIPTIC \n",
"(iii) A=1+xˆ2, B=5+2xˆ2, C=4+xˆ2\n",
"D=9>0\n",
"So the equation is HYPERBOLIC \n"
]
}
],
"source": [
"print \"D=Bˆ2−4AC\"\n",
"print \"if D<0 then elliptic if D=0 then parabolic if D>0 then hyperboic\"\n",
"print \"(i) A=xˆ2, B1−yˆ2 D=4ˆ2−4∗1∗4=0 so the equation is PARABOLIC \"\n",
"print \"(ii) D=4xˆ2(yˆ2−1)\"\n",
"print \"for −inf<x<inf and −1<y<1 D<0 \"\n",
"print \"So the equation is ELLIPTIC \"\n",
"print \"(iii) A=1+xˆ2, B=5+2xˆ2, C=4+xˆ2\"\n",
"print \"D=9>0\"\n",
"print \"So the equation is HYPERBOLIC \""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.2, page no. 726"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"See figure in question\n",
"From symmetry u7=u1, u8=u2, u9=u3, u3=u1, u6=u4, u9=u7\n",
"u5=1/4∗(2000+2000+1000+1000)=1500\n",
"u1=1/4(0=1500+1000+2000)=1125\n",
"u2=1/4∗(1125+1125+1000+1500)=1188\n",
"u4=1/4(2000+1500+1125+1125)=1438\n",
"1125 1188 1438 1500\n",
"Iterations :\n",
"\n",
"1301 1414 1164 1031\n",
"\n",
"1250 1337 1087 1019\n",
"\n",
"1199 1312 1062 981\n",
"\n",
"1174 1287 1037 968\n",
"\n",
"1155 1274 1024 956\n",
"\n",
"1144 1265 1015 949\n"
]
}
],
"source": [
"print \"See figure in question\"\n",
"print \"From symmetry u7=u1, u8=u2, u9=u3, u3=u1, u6=u4, u9=u7\"\n",
"print \"u5=1/4∗(2000+2000+1000+1000)=1500\"\n",
"u5 = 1500\n",
"print \"u1=1/4(0=1500+1000+2000)=1125\"\n",
"u1 = 1125\n",
"print \"u2=1/4∗(1125+1125+1000+1500)=1188\"\n",
"u2 = 1188\n",
"print \"u4=1/4(2000+1500+1125+1125)=1438\"\n",
"u4 = 1438\n",
"print u1,u2,u4,u5 \n",
"print \"Iterations :\"\n",
"for i in range(1,7):\n",
" u11 = (1000+u2+500+u4)/4\n",
" u22 = (u11+u1+1000+u5)/4 \n",
" u44 = (2000+u5+u11+u1)/4\n",
" u55 = (u44+u4+u22+u2)/4\n",
" print \"\"\n",
" print u55,u44,u22,u11\n",
" u1 = u11 \n",
" u2 = u22 \n",
" u4 = u44 \n",
" u5 = u55"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.3, page no. 727"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"See figure in question\n",
"To find the initial values of u1 u2 u3 u4 we assume u4=0 \n",
"u1=1/4∗(1000+0+1000+2000)=1000\n",
"u2=1/4(1000+500+1000+500)=625 \n",
"u3=1/4∗(2000+0+1000+500)=875\n",
"u4=1/4(875+0+625+0)=375\n",
"1000 625 875 375\n",
"Iterations:\n",
"\n",
"437 1000 750 1125\n",
"\n",
"453 1031 781 1187\n",
"\n",
"457 1039 789 1203\n",
"\n",
"458 1041 791 1207\n",
"\n",
"458 1041 791 1208\n",
"\n",
"458 1041 791 1208\n"
]
}
],
"source": [
"print \"See figure in question\"\n",
"print \"To find the initial values of u1 u2 u3 u4 we assume u4=0 \"\n",
"print \"u1=1/4∗(1000+0+1000+2000)=1000\"\n",
"u1 = 1000\n",
"print \"u2=1/4(1000+500+1000+500)=625 \"\n",
"u2 = 625\n",
"print \"u3=1/4∗(2000+0+1000+500)=875\"\n",
"u3 = 875\n",
"print \"u4=1/4(875+0+625+0)=375\"\n",
"u4 = 375\n",
"print u1,u2,u3,u4 \n",
"print \"Iterations:\"\n",
"for i in range(1,7):\n",
" u11 = (2000+u2+1000+u3)/4 \n",
" u22 = (u11+500+1000+u4)/4\n",
" u33 = (2000+u4+u11+500)/4\n",
" u44 = (u33+0+u22+0)/4\n",
" print \"\"\n",
" print u44,u33,u22,u11\n",
" u1 = u11 \n",
" u2 = u22 \n",
" u4 = u44 \n",
" u3 = u33"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.5, page no. 729"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Here cˆ2=4, h=1, k=1/8, therefore alpha=(cˆ2)∗k/(hˆ2)\n",
"Using bendre−schmidits recurrence relation ie u(i)(j +1)=t∗u(i−1)(j)+t∗u(i+1)(j)+(1−2t)∗u(i,j)\n",
"Now since u(0,t)=0=u(8,t) therefore u(0,i)=0 and u(8,j)=0 and u(x,0)=4x−1/2xˆ2 \n",
"0.0\n",
"-4.0\n",
"0.0\n",
"4.0\n",
"8.0\n",
"12.0\n",
"16.0\n"
]
}
],
"source": [
"import numpy\n",
"print \"Here cˆ2=4, h=1, k=1/8, therefore alpha=(cˆ2)∗k/(hˆ2)\"\n",
"print \"Using bendre−schmidits recurrence relation ie u(i)(j +1)=t∗u(i−1)(j)+t∗u(i+1)(j)+(1−2t)∗u(i,j)\"\n",
"print \"Now since u(0,t)=0=u(8,t) therefore u(0,i)=0 and u(8,j)=0 and u(x,0)=4x−1/2xˆ2 \"\n",
"c = 2\n",
"h = 1\n",
"k = 1/8\n",
"t = (c**2)*k/(h**2)\n",
"A = numpy.ones((9,9))\n",
"for i in range(0,9):\n",
" for j in range (0,9):\n",
" A[0,i] = 0\n",
" A[8,i] = 0\n",
" A[i,0] = 4*(i-1)-1/2*(i-1)**2\n",
"for i in range(1,8):\n",
" for j in range(1,7):\n",
" A[i,j] = t*A[i-1,j-1]+t*A[i+1,j-1]+(1-2*t)*A[i-1,j-1]\n",
"for i in range(1,8):\n",
" j = 2\n",
" print A[i,j]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.6, page no. 730"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Here cˆ2=1, h=1/3, k=1/36, therefore t=(cˆ2)∗k/(hˆ2)=1/4 \n",
"So bendre−schmidits recurrence relation ie u(i)(j+1)=1/4(u(i−1)(j)+u(i+1)(j)+2u(i,j)\n",
"Now since u(0,t)=0=u(1,t) therefore u(0,i)=0 and u(1,j)=0 and u(x,0)=sin(%pi)x\n",
"-0.433012701892\n",
"0.216506350946\n",
"0.649519052838\n",
"0.433012701892\n",
"-0.216506350946\n",
"-0.649519052838\n",
"-0.433012701892\n"
]
}
],
"source": [
"import numpy,math\n",
"\n",
"print \"Here cˆ2=1, h=1/3, k=1/36, therefore t=(cˆ2)∗k/(hˆ2)=1/4 \"\n",
"print \"So bendre−schmidits recurrence relation ie u(i)(j+1)=1/4(u(i−1)(j)+u(i+1)(j)+2u(i,j)\"\n",
"print \"Now since u(0,t)=0=u(1,t) therefore u(0,i)=0 and u(1,j)=0 and u(x,0)=sin(%pi)x\"\n",
"c = 1.\n",
"h = 1./3\n",
"k = 1./36\n",
"t = (c**2)*k/(h**2)\n",
"A = numpy.ones((9,9))\n",
"for i in range(0,9):\n",
" for j in range(0,9):\n",
" A[0,i] = 0\n",
" A[1,i] = 0\n",
" A[i,0] = math.sin(math.pi/3*(i-1)) \n",
"for i in range(1,8):\n",
" for j in range(1,8):\n",
" A[i,j] = t*A[i-1,j-1]+t*A[i+1,j-1]+(1-2*t)*A[i-1,j-1]\n",
"for i in range(1,8):\n",
" j = 1\n",
" print A[i,j]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.7, page no. 732"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Here cˆ2=16, taking h=1, finding k such that cˆ2tˆ2=1 \n",
"So bendre−schmidits recurrence relation ie u(i)(j+1)=(16tˆ2(u(i−1)(j)+u(i+1)(j))+2(1−16∗tˆ2u(i,j)−u(i)(j−1)\n",
"Now since u(0,t)=0=u(5,t) therefore u(0,i)=0 and u(5,j)=0 and u(x,0)=xˆ2(5−x)\n",
"Also from 1st derivative (u(i)(j+1)−u(i,j−1))/2k=g(x) and g(x)=0 in this case\n",
"So if j=0 this gives u(i)(1)=1/2∗(u(i−1)(0)+u(i+1)(0))\n",
" 0.0 0.0 0.0 0.0 0.0 \n",
" 0.0 4.0 8.0 12.0 16.0 \n",
" 0.0 12.0 24.0 36.0 48.0 \n",
" 0.0 18.0 36.0 54.0 72.0 \n",
" 0.0 16.0 0.0 0.0 0.0 \n"
]
}
],
"source": [
"import numpy\n",
"\n",
"print \"Here cˆ2=16, taking h=1, finding k such that cˆ2tˆ2=1 \"\n",
"print \"So bendre−schmidits recurrence relation ie u(i)(j+1)=(16tˆ2(u(i−1)(j)+u(i+1)(j))+2(1−16∗tˆ2u(i,j)−u(i)(j−1)\"\n",
"print \"Now since u(0,t)=0=u(5,t) therefore u(0,i)=0 and u(5,j)=0 and u(x,0)=xˆ2(5−x)\"\n",
"c = 4\n",
"h =1 \n",
"k = (h/c)\n",
"t = k/h\n",
"A = numpy.zeros((6,6))\n",
"print \"Also from 1st derivative (u(i)(j+1)−u(i,j−1))/2k=g(x) and g(x)=0 in this case\"\n",
"print \"So if j=0 this gives u(i)(1)=1/2∗(u(i−1)(0)+u(i+1)(0))\"\n",
"for i in range(0,6):\n",
" for j in range(1,9):\n",
" A[0,i] = 0\n",
" A[5,i] = 0;\n",
" A[i,1] = (i)**2*(5-i)\n",
"for i in range(0,4):\n",
" A[i+1,2] =1/2*(A[i,1]+A[i+2,1])\n",
"for i in range(2,5):\n",
" for j in range(2,5):\n",
" A[i-1,j] = (c*t)**2*(A[i-2,j-1]+A[i,j-1])+2*(1-(c*t)**2)*A[i-1,j-1]-A[i-1,j-2]\n",
"for i in range(0,5):\n",
" for j in range(0,5):\n",
" print \" \",A[i,j],\n",
" print \"\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Example 28.8, page no. 734"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Here cˆ2=4, taking h=1, finding k such that cˆ2tˆ2=1 \n",
"So bendre−schmidits recurrence relation ie u(i)(j+1)=(16tˆ2(u(i−1)(j)+u(i+1)(j))+2(1−16∗tˆ2u(i,j)−u(i)(j−1)\n",
"Now since u(0,t)=0=u(4,t) therefore u(0,i)=0 and u(4,j)=0 and u(x,0)=x(4−x) \n",
"Also from 1st derivative (u(i)(j+1)−u(i,j−1))/2 k=g(x) and g(x)=0 in this case\n",
"So if j=0 this gives u(i)(1)=1/2∗(u(i−1)(0)+u(i+1)(0))\n",
" 0.0 0.0 0.0 0.0 0.0 \n",
" 3.0 0.0 -3.0 -6.0 -9.0 \n",
" 4.0 0.0 -4.0 -8.0 -12.0 \n",
" 3.0 0.0 -3.0 -6.0 -9.0 \n",
" 0.0 0.0 0.0 0.0 0.0 \n"
]
}
],
"source": [
"import numpy\n",
"\n",
"print \"Here cˆ2=4, taking h=1, finding k such that cˆ2tˆ2=1 \"\n",
"print \"So bendre−schmidits recurrence relation ie u(i)(j+1)=(16tˆ2(u(i−1)(j)+u(i+1)(j))+2(1−16∗tˆ2u(i,j)−u(i)(j−1)\"\n",
"print \"Now since u(0,t)=0=u(4,t) therefore u(0,i)=0 and u(4,j)=0 and u(x,0)=x(4−x) \"\n",
"c = 2\n",
"h = 1\n",
"k = (h/c)\n",
"t = k/h\n",
"A = numpy.zeros((6,6))\n",
"print \"Also from 1st derivative (u(i)(j+1)−u(i,j−1))/2 k=g(x) and g(x)=0 in this case\"\n",
"print \"So if j=0 this gives u(i)(1)=1/2∗(u(i−1)(0)+u(i+1)(0))\"\n",
"for i in range(0,6):\n",
" for j in range(1,9):\n",
" A[0,i] = 0\n",
" A[4,i] = 0\n",
" A[i,0] = (i)*(4-i)\n",
"for i in range(0,4):\n",
" A[i+1,2] = 1/2*(A[i,1]+A[i+2,1])\n",
"for i in range(2,5):\n",
" for j in range(2,5):\n",
" A[i-1,j] = (c*t)**2*(A[i-2,j-1]+A[i,j-1])+2*(1-(c*t)**2)*A[i-1,j-1]-A[i-1,j-2]\n",
"for i in range (0,5):\n",
" for j in range(0,5):\n",
" print \" \",A[i,j],\n",
" print \"\""
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.11+"
}
},
"nbformat": 4,
"nbformat_minor": 0
}