mirror of
https://github.com/wassname/simpeg.git
synced 2026-10-30 13:50:07 +08:00
Secondary field approach
This commit is contained in:
SEOGI KANG committed
1 parent
cc33d148b6
commit
35d132f52b
1 file changed
+704
@@ -0,0 +1,704 @@
|
||||
{
|
||||
"metadata": {
|
||||
"name": ""
|
||||
},
|
||||
"nbformat": 3,
|
||||
"nbformat_minor": 0,
|
||||
"worksheets": [
|
||||
{
|
||||
"cells": [
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"from SimPEG import *\n",
|
||||
"from test_boundary import spheremodel, MagSphereAnalFun\n",
|
||||
"%pylab inline"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"Warning: Python backend is being used for solver. Run setup.py from the command line.\n",
|
||||
"Warning: mumps solver not available.\n",
|
||||
"Warning: upgrade your scipy to 0.13.0"
|
||||
]
|
||||
},
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"\n",
|
||||
"Populating the interactive namespace from numpy and matplotlib"
|
||||
]
|
||||
},
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 1
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Step:1 Generating mesh and operators"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"hxind = ((10,25,1.3),(41, 12.5),(5,25,1.3))\n",
|
||||
"hyind = ((10,25,1.3),(41, 12.5),(5,25,1.3))\n",
|
||||
"hzind = ((10,25,1.3),(40, 12.5),(1,25,1.3))\n",
|
||||
"hx, hy, hz = Utils.meshTensors(hxind, hyind, hzind)\n",
|
||||
"M3 = Mesh.TensorMesh([hx, hy, hz], [-sum(hx)/2,-sum(hy)/2,-sum(hz)/2])"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 46
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Step2: Compute Boundary indicies and set $\\mathbf{B}_{sBC}$"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"BC = [['neumann', 'neumann'], ['neumann', 'neumann'], ['neumann', 'neumann']]\n",
|
||||
"indxd, indxu, indyd, indyu, indzd, indzu = M3.faceBoundaryInd\n",
|
||||
"ind = np.r_[(indxd | indxu), (indyd | indyu), (indzd | indzu)]"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 47
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"Box = 1 # Primary field in x-direction (background)\n",
|
||||
"Boy = 0 # Primary field in y-direction (background)\n",
|
||||
"Boz = 0 # Primary field in z-direction (background)\n",
|
||||
"B0 = np.r_[Box*np.ones(np.prod(M3.nFx)), Boy*np.ones(np.prod(M3.nFy)), Boz*np.ones(np.prod(M3.nFz))]\n",
|
||||
"Bbcx = np.zeros(np.prod(M3.nFx))\n",
|
||||
"Bbcy = np.zeros(np.prod(M3.nFy))\n",
|
||||
"Bbcz = np.zeros(np.prod(M3.nFz))\n",
|
||||
"\n",
|
||||
"Bbcx[indxd] = 0\n",
|
||||
"Bbcx[indxu] = 0\n",
|
||||
"Bbcy[indyd] = Boy\n",
|
||||
"Bbcy[indyu] = Boy\n",
|
||||
"Bbcz[indzd] = Boz\n",
|
||||
"Bbcz[indzu] = Boz\n",
|
||||
"\n",
|
||||
"Bbc = np.r_[Bbcx, Bbcy, Bbcz]"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 48
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Step3: Generating model $\\mu$\n",
|
||||
"\n",
|
||||
"$\\mu = \\mu_0(1+\\chi)$"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"mu0 = 4*np.pi*1e-7\n",
|
||||
"chibkg = 0.\n",
|
||||
"chiblk = 0.01\n",
|
||||
"chi = np.ones(M3.nC)*chibkg"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 49
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"sph_ind = spheremodel(M3, 0, 0, 0, 100)\n",
|
||||
"chi[sph_ind] = chiblk\n",
|
||||
"mu = (1.+chi)*mu0"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 50
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Additional: Congrous sphere method (for boudary condition)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"V = M3.vol[sph_ind].sum()\n",
|
||||
"gamma = 1/V*(chiblk*M3.vol).sum()\n",
|
||||
"mom = 1/mu0*gamma*V/(1+gamma/3)"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 72
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"4109375.0\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 70
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"# figsize(10,10)\n",
|
||||
"# M3.plotGrid()"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAAjwAAAI8CAYAAAD1D3GaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXWYVNUbx7/3Tu6yLLGwxNINUlLSLQjSIIoIiISkhBIK\nAgIiYCCiSBqgSAsiAtIh+bOQEJBSYlli2emb7++Pw70zs4G7hMh6Ps9znt2ZW2fmxvnOW0cgIgKH\nw+FwOBxOJkZ80B3gcDgcDofDud9wwcPhcDgcDifTwwUPh8PhcDicTA8XPBwOh8PhcDI9XPBwOBwO\nh8PJ9HDBw+FwOBwOJ9PDBQ+Hw+FwOJxMDxc8HA6Hw+FwMj1c8HA4HA6Hw8n0cMHD4XA4HA4n08MF\nD4fD4XA4nEwPFzwcDofD4XAyPVzwcDgcDofDyfRwwcPhcDgcDifTwwUPh8PhcDicTA8XPBwOh8Ph\ncDI9XPBwOBwOh8PJ9HDBw+FwOBwOJ9PDBQ+Hw+FwOJxMDxc8HA6Hw+FwMj1c8HA4HA6Hw8n0cMHD\n4XA4HA4n08MFD4fD4XA4nEwPFzwcDofD4XAyPVzwcDgcDofDyfRwwcPhcDgcDifTwwUPh8PhcDic\nTA8XPBwOh8PhcDI9XPBwOBwOh8PJ9HDBw+FwOBwOJ9PDBQ+Hw+FwOJxMDxc8HA6Hw+FwMj1c8HA4\nHA6Hw8n0cMHD4XA4HA4n08MFD4fD4XA4nEwPFzwcDofD4XAyPVzwcDgcDofDyfRwwcPhcDgcDifT\nwwUPh8PhcDicTA8XPBwOh8PhcDI9XPBwOBwOh8PJ9HDBw+E8ZOi6DlVVQUQPuiscDofz0GB90B3g\ncDjph4ggSRL8fj8sFgtsNhusVissFgtEUYQgCA+6ixwOh/OvRCD+M5HDeShQVRWKogCA+VfXdRAR\nBEGAIAimALJareZ7HA6Hw+GCh8P510NEUFUVqqqaAkaW5RRihohMAQQAoijCZrPBZrPBYrFwAcTh\ncP7TcMHD4fyL0XUdiqJA13VTsBBRqoInFOO2DrUAiaIIq9XKBRCHw/lPwgUPh/MvhIigaZrpuhJF\nMWzZ3wme1PYHMAFkwAUQh8P5L8EFD4fzL4OIoCgKNE1LVYTcieBJ7RjA7QVQqMjicDichx0ueDic\nfxG6rkOW5bBA5OTcC8GT2j6N4xuIohiWCcYFEIfDeZjhaekczr+A5IHJ/7S4MMSTxWIx+wMAXq/X\ntPxwCxCHw3mY4U8sDucBY1hsDLGTXsuNJAHLlt2f3yxGPwzxJYqi2U+v1wuXywWXywW/328GVXM4\nHM6/Ge7S4nAeIIZV53YurOQQEW7ckFGjRhQuXhRRqJCO+vU11K+vol49DXFx9+6WDgQCplsreR+M\nZmCkwYcWQuRwOJx/C1zwcDgPgOQurIzE4/z8M9CtmxX58xOOHhWxYYMfu3ZZsHs3a9mzwxQ/9epp\nyJv3zm/xtARPap+HCyAOh/NvhgseDucfxkg3z4hVBwCIgI8/FjFlihVTp/pRq5aO5s0jceyY11xH\n14Fjx0Ts3m3Brl0W7N1rRe7czAJUty4TQLlzp/+WT6/gSdnXlALI2I/FYjErQXM4HM4/BRc8HM4/\nhFFbx+fzQZIkREdHp3vb69eBF1+04tIlAYsWKShYUMKlSyIaN47EiRPeNLfTNODIEUMAWbFvnwX5\n8zMBVK+ehrp1VeTMmfZx71TwJCe5ACKisABoLoA4HM79hgseDucfILS2jqqqCAQC6RY8e/YI6NnT\nhg4dNEyapMFmY8HDV66IqFcvEqdOpS14kqOqwK+/irfcX1bs329BkSI6fvvNgtatFXz0UQDZswfX\nv1eCJzmhAsiwdCWfDJULIA6Hcy/hgofDuc8kr62jqir8fv/fCh5NA6ZNs2DuXAvmzFHRogXLhDKy\npa5eFVGzZiTOnEm/4EmOJAHdu0dgwwaW7RUVRShRQke9eiwIukoVH7JnF++54EnO7VxgXABxOJx7\nAa/Dw+HcJ+6mts6lS8DzzzORsXevjLi4lOtYLCxm505JSBDw4otO+HxAp04KKlbUMGCAgh9/ZPE/\ns2bZ8b//RaBMGQ0NGjAXWM2aGqKi7vyYaZE8lsmYCDUQCJjvcQHE4XDuBp42weHcB+60to6uAxs2\niKhVy46GDXVs2KCkKnYAQBAIun5ng/7OnRbUqxeJKlU0rF/vR4ECOogE2O1ArVoaRo2S8e23fvz+\n+zWMH++FwwG8844dJUtGoWnTSEycaMf27Rb4fHd0+L/FEIgWi8XM8DIEkMfjQVJSEtxuNwKBgJnW\nz+FwOLeDW3g4nHvMnWZhyTIQHe0AAMTFEfLmJZw7BxQrBqS2izux8GgaMH26HZ9+asOcOQE0bqwB\nAEQx9X1FRAioU0dFo0YyXnsN8PmAgwdZ+vuUKQ4cOMAqM48eLaFePQ3Vq2twOjPWp/TALUAcDudu\n4YKHw7lHpNeFJQhCCovE6dNAt242REQQnE5g3DgVO3aImDzZClEEGjbU0aCBjvr1dRQqxLYRRSZg\n0svlywJ693ZCFIFdu3xh9XkEIX3iKTISaNhQQ506GhQFpuCRJGDcOAd+/11ElSrarRggDdWqabDb\n09/H9JKWAPL7/WHTZHABxOFwDLjg4XDuAbqum1MsZLSQ4LJlIl5+2YoxY1SULk2YNs2KHj109Oih\ngwj44w8BO3YI2LRJxGuvWREVRahXT0T16jq83vQdZ8sWC/r3d6J3bwWvvCLj1pRZJqLI6vykhwsX\nBPTsGYFs2QgbNvgwdKgDb7whA5DhcgH797MU+FdfdeDUKRHVq2tmEHS5cvp9jQEyRCYXQBwOJzlc\n8HA4d4FRW0dRFADIUGCy1wsMG2bFvn0C1q9XUKkSYdcuIcxqIwhAyZKEkiUJffowAXT0KLBtG+Gb\nb9jtW7VqpBlUXK+ehly5gspFUYA337Rj6VIbPv00gLp1UzcJpeXSSs6mTRYMHOjE4MEyBg9WcPKk\nGCaUoqOBZs00NGvGjnPzJrB3LxNATZtmAQA0bKiiUSMN9eqpqFxZh/U+PIVuJ4AURYHNZjPFDxdA\nHM5/Ay54OJw7JLS2TkatOkePWtG/vw01ahD27VNMq4fFcns3lSAAjzxCKFVKQa9eCnLlisInnwSw\na5cFS5bYMHiwEwULsrRym41NLlqpko49e3xhQii1/d7OwqMowKRJdqxcacMXX/hRs6aeru3YNBca\nvv7ahnz5dFy+LKJ3bwW7d1swaJATFy+KqFkzOA9YhQp6CuvTvSBUAEmSBEEQTKFqnLfQQohcAHE4\nmQ8ueDicOyB5bZ2MTA8xf74NkydHYPp0DV27hptV/k7whGIYkypV0lGpko7BgxWoKvDLLyJefdVp\nxtckJBDef9+O+vVV1KqlIWvWlPu6nXC5eFFE376RyJaNsHu3DzExlGy7tD/78eMiunVzokYNHdu3\n+/DYY1nQpo2KNm1UAMDVqwL27GFB0AsX2nHmDPtQ06Yxa9Qjj+i4H1Nwhc4CD4Rb6pILIKvVClEU\nuQDicB5yuODhcDJA8oExIy6sxEQgXz6WhbV0aSLatYtMsU5GBQ+RAKJgFpeuA6tW2XD5soCSJVld\nnXLldOzaZcHMmXb06GFBuXK6aVGpWVNDZCTb1y2vXBhbttgwZEhWDBok46WXlBTiQxAoTaG0ZIkV\nY8Y4MGmShOeeU3HzZsp1cucmtG+vokYNDUeOiDhzRkSOHITjx0XMm2dHYiLMOcCyZCE8/bSK+1ED\nMTUXmFEV24ALIA7n4YYLHg4nnRi1de4kMHnvXsEsJFi3robBg6MxapSIBg10MwMrLi5jgkcQjFo8\nbLuzZ1kwcd68Onbt8mLcOAesVlZXh9XWAfz+YFr51Kl2/PabBZUqaThwwILq1TVIEuBwBF1YK1ZY\n8emnbtSrl1bGWUrLkM8HjBjhwP79Fqxf70e5ckErVmriaOtWC1580Yl+/RQ0aqRB14ExY2QAEi5e\nZBagPn0iAAADBgAdOihwuwV07KigSxc11ZT925Gemj2pZYGFCiBBEMz4Hy6AOJyHAy54OJx0cKe1\ndTQNeOcdCz76yILZs1X072/FokUSIiLcuHIlB3bsEPDNNyJeecWKXLkIOXMyl9TVq0Du3GnvNzj/\nFLPqrFtnxfDhDowYIaNfPwWCkHraekQE0KABq5wMAB4Py6rq0CES+/dbUbRo1K33BWTPTtiz5wby\n5gWA1HPLk38Np04J6NYtAuXL69i50xeWkZV8XTZ1hh2ff27DJ58EUL++hqlTw4+TPTthyxYrypXT\ncOyYBZs2+TBjhh2bN1uxebMVEyboiIwEzpwR8fPPHhQrRukSQBkVJ6kJIFVVzWB1LoA4nH8/XPBw\nOLfhbqaHuHwZeOEFGxSFTQ9RoAAweDBuVUcmlCnDWr9+OnQd+O03AR9+aMH+/cAjj9hRuDDdsgAR\n6tbVzUk9jclHRVGEqmbF0KEO7NljxcqVflSpErSmpKcwYVQU0LSphvHjJbhcQPnyOnr1YtaUAgV0\n1KyZA489pt6qAaShUqWUQcWGwWT5citGjXJg3DgZzz+vpCo8jHWvXmU1gVSV1QTKkyel1eXECRb/\nU6WKjm3bfKhYMQs+/tiG48dZqvtTTymIjAQGDWKVDlu2jITHI8DtFjBxYgDt26soXPj+VGDmAojD\nefjgU0twOGlgBCZndHoIAPj+ewG1atlRp46OjRsVFCjA3mcuq5T7EUWgUiXC8OEaSpfWcemSjI8+\nUhEbC8yebUGJEnbUqWPDqFGEdetkBAIWnD7NlMelSzo2bbqGRx7xh02zkN5Uc4BZW2bMcGD8eAfe\nfjuAChU07N3rw48/JqJbNwkXL4ro18+JokWj8PTTEZg924bffmOPj0AAGDbMgSlTHFi71oLine truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x231e8588>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 51
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"M3.plotImage(np.log10(mu), imageType='CC')"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "pyout",
|
||||
"prompt_number": 52,
|
||||
"text": [
|
||||
"<matplotlib.collections.QuadMesh at 0x1f145d68>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAJWCAYAAACEbPmsAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XtcVHXi//H3DLe85xVZA10vC2KpoIgtbli5Lbl5jdR2\ntfJ++bbihmlpJv3cNS1b++rXC9UD3ZQibNdVM1fznpZ4Lw0kTc0rIl6HvAXM7w9ykgAF/cwwyOvp\ng0dx5szh83l7hnnPOWdGi91utwsAAADGWMt6AAAAAHcbChYAAIBhFCwAAADDKFgAAACGUbAAAAAM\no2ABAAAYRsECcNdYsWKFIiMjVbduXdWvX1+dO3fW5s2bJUknTpxQv3791LRpU9WoUUOtWrXS9OnT\nlZeXV8ajBnA3omABuCu8++67+vOf/6z27dtry5Yt2r9/vwYMGKDk5GRlZGSoZcuWstlsmjlzprKy\nsvT+++9r3bp1ys7OLuuhA7gLWfigUQDlnc1mU8OGDfXiiy/q5ZdfLnT74MGDtXHjRqWnp5fB6ABU\nRBzBAlDu7d27V+fPn1fXrl2LvH3jxo3F3gYAzkDBAlDuHT16VNWqVVOLFi2KvP3YsWP67W9/6+JR\nAajIKFgAyj1/f3/ZbDZ98803xd7+xRdfuHhUACoyChaAcu/+++9XzZo1tXTp0iJvf+ihh7Rs2TIX\njwpARUbBAlDuVatWTVOnTtVbb72lcePG6bvvvtPFixe1ePFixcTEaNKkSTpz5oyefPJJrVy5Uteu\nXdOePXvUo0cPXbhwoayHD+AuRMECcFcYNGiQFi5cqM2bNys8PFy/+c1vlJCQoKefflq+vr76+uuv\nVaVKFY0YMUJ16tRRv3799PDDD6tatWplPXQAdyE+pgEAAMAwjmABAAAYRsECAAAwjIIFAABgmGdZ\nD+BGFoulrIcAAABQYsVdyu5WBStfXFkPwK3Y7RMlSRbLa5LWS+pYhqNxH9dzuVFcXJzi4uJcP5hy\ngnxujYxujnxujnxu7W7L6GYHhjhFCAAAYBgFCwAAwDAKVrnSqKwH4NY6duxY1kNwa+Rza2R0c+Rz\nc+RzaxUpI7f6oNH8c5lxZT0Mt1LwGixcV9Q1WAAAuJLFYin2IneOYAEAABhGwQIAADCMggUAAGAY\nBQsAAMAwChYAAIBhFCwAAADDKFgAAACGUbAAAAAMo2ABAAAYRsECAAAwjIIFAABgGAULAADAMAoW\nAACAYRQsAAAAwyhYAAAAhlGwAAAADKNgAQAAGEbBAgAAMIyCBQAAYBgFCwAAwDAKFgAAgGEULAAA\nAMMoWAAAAIZRsAAAAAyjYAEAABhGwQIAADCMggUAAGAYBQsAAMAwChYAAIBhFCwAAADDKFgAAACG\nUbAAAAAMo2ABAAAYRsECAAAwjIIFAABgGAULAADAMAoWAACAYRQsAAAAwyhYAAAAhlGwAAAADKNg\nAQAAGEbBAgAAMIyCBQAAYBgFCwAAwDAKFgAAgGEULAAAAMMoWAAAAIZRsAAAAAyjYAEAABhGwQIA\nADCMggUAAGAYBQsAAMAwChYAAIBhtyxYP/zwg5599ln95je/UXBwsFJSUmSz2dStWzcFBASoe/fu\nys7Odqw/Y8YMNWvWTMHBwdq0aZNjeVpamkJDQ9W4cWONHz/eObMBAABwA7csWBMnTlRAQIC+/vpr\nff311woKCtKcOXMUEBCg/fv367777tPcuXMlSZmZmZo9e7bWrFmjOXPmaOTIkY7txMbGauzYsdq2\nbZs2bNig7du3O29WAAAAZeiWBWv16tUaN26c7rnnHnl6eqpGjRraunWrBg4cKB8fHw0YMEApKSmS\npJSUFEVFRSkgIECRkZGy2+2Oo1vp6enq3bu3ateurZ49ezruAwAAcLfxvNmNx44d05UrVzR8+HCl\npaWpZ8+eGjlypLZt26agoCBJUlBQkLZu3Sopv2A1b97ccf/AwEClpKSoYcOGqlevnmN5cHCwEhMT\n9T//8z9F/NT1N/x/o5++AAAAytb69eu1fv36Eq1704J15coVffvtt3rzzTfVqVMnDR06VMnJybLb\n7SUejMViKbTs5vfvWOJtAwAAuErHjh3VsWNHx/evvfZaseve9BRh06ZNFRgYqC5duqhSpUp6+umn\n9d///ldhYWFKS0uTlH/xelhYmCQpPDxcqampjvvv27dPYWFhatq0qU6dOuVYnpqaqvbt29/W5AAA\nANzdLa/BatasmVJSUpSXl6fly5erU6dOCg8PV0JCgi5fvqyEhARHWWrXrp1WrlypI0eOaP369bJa\nrapWrZqk/FOJSUlJysrK0uLFixUeHu7cmQEAAJSRm54ilKRp06bpmWee0ZUrV9SpUyf16dNHeXl5\n6tu3rwIDAxUaGqqpU6dKknx9fTV8+HA98sgj8vb2Vnx8fIHt9O3bVy+//LL69Omjtm3bOm9WAAAA\nZchiL80FVU6Wf71WXFkPo0QyMmLVufMH2rnzpDZseE5z5mxXUtLeQus1b15H27YNlre3h7y9/1bg\ntqiopvrznx/QI4/8WhkZ2frss+80c+ZWHT9uc6xjt0+UJFksxZ/ndRcmMnn88abq1auFOnduppMn\nbZo3b7f+938Lv+P0ei4AAJQVi8VS7HXlfJL7bWjSpKYqV/bSrl0n5eVlVdu2v9KmTUcKrVepkqeS\nk5/SmjWH9Mv8a9WqpOTkaGVm/qCoqIUaN26NOnVqrKlTO7loFmaZyKRNGz8tWdJHqamn9eij72vu\n3B2aMOEhvfRSBxfNAgAAMyhYtyEiIkApKcdlt0thYQ105swlHTt2sdB6s2Z11saN3+tf/0rTL99M\n2bVroDw8rBoz5jPt2ZOplSu/04wZW/Xkk8Hy9Cx/fy0mMomJCdeKFQf05ptfaO/eTM2du13z5u1W\nbOyD8vHxcNFMAAC4c7e8Bgs/O3durOx2u3x8PGW1WnT27Bh5eXnIx8dDZ8+Okd0u1a79hiSpX7+W\natPmVwoLe1dPP31/oW2tWvWd7Ha7RowI0/z5u1WvXhX169dSy5d/q5ycPFdP7baZzKR6dR9lZ18r\nsMxmu6ZatSopOLiudu3KcMmcAAC4UxSsUmjZco4sFou2bBmoYcOWa/fuDCUlPakPPtirJUv2OdYL\nCqqjadMeU8eO83XtWm6R2zpxwqbQ0Hf0+ef99Y9//EFWq0Uff5yqPn0+dtV0jDCZybvv7tSiRU+p\nZ8/mWrnygB580F+DBoVIkvz9a1CwAADlRvk7F1WGjh69qBo1fOTl5aFly9J17txltW5dX0lJe3X0\n6EUdPXpR3t4eWrToKb3yylqlpWUVu63g4Lpav/5ZJSXt1cMP/1P9+i1W48Y1tWjRUy6c0Z0zmcmn\nn+7XpEkbNWHCQ7pw4SUlJHTVzJn5/0pAeTqqBwAA7yIsob17hysgoIY8Pa3y8vLQ5cs/ymq1qFIl\nL/3wQ/5prebNZ8nT06qDB2OUm/tzIbBYLLJaLcrNzdOECes0depmTZ/+Bz38cCO1bv3zR1mEhzfQ\nF18MVIsWs7VvX34Rced3EZrO5EY1avjowoWr6to1UIsX91Zg4P/pwIGzjtt5FyEAoKzd7F2EnCIs\noaioRHl7eyghoatWrDig5ORvNHFipK5ezdWUKZskSSdPZstike6/f3aB+3bvHqTXXuuoVq3mKjPz\nB0lSbm5eoaMyubn5f0lu1HlvynQmN7pw4aokqX//1tq9O6NAuQIAwN1RsEro2LGLslotatnSV0OG\nfKJDh87rgQd8FRe3XocOnS+w7i9Pg7VrZyu0fM6c7Ro1qr0mT35Eixalyt+/hsaOjdDatYeUnn7G\n+RMywHQmjRvXVIcOAfryy6MKCfHT6NEPKji4rqKiEp0/GQAADKJglUJISH1dvZqrb789o+rVfdSi\nRV1t3Ph9ie77y4NS3313Tt27f6QePYK0bNnTOnkyW5988q3ef/8rJ4zceUxmYrVa9PzzYZo1q7Ns\ntqvatu2EnntuiVJTTzth5AAAOA/XYLk5d74GqyxxDRYAoKzxSe4AAAAu5HZHsNxoOAAAAMXiCBYA\nAIALud1F7lxrVBDXYBWNXAojk6KRS9HIpWjkUhiZFO1W1wJzBAsAAMAwChYAAIBhFCwAAADDKFgA\nAACGUbAAAAAMo2ABAAAYRsECAAAwjIIFAABgGAULAADAMAoWAACAYRQsAAAAwyhYAAAAhlGwAAAA\nDKNgAQAAGEbBAgAAMIyCBQAAYBgFCwAAwDAKFgAAgGEULAAAAMMoWAAAAIZRsAAAAAyjYAEAABhG\nwQIAADCMggUAAGAYBQsAAMAwChYAAIBhFCwAAADDKFgAAACGUbAAAAAMo2ABAAAYRsECAAAwjIIF\nAABgGAULAADAMAoWAACAYRQsAAAAwyhYAAAAhlGwAAAADKNgAQAAGEbBAgAAMIyCBQAAYBgFCwAA\nwDAKFgAAgGEULAAAAMMoWAAAAIZRsAAAAAyjYAEAABhGwQIAADCMggUAAGAYBQsAAMAwChYAAIBh\nFCwAAADDKFgAAACGUbAAAAAMo2ABAAAYRsECAAAwjIIFAABgGAULAADAMAoWAACAYRQsAAAAwyhY\nAAAAht2yYDVq1EgtW7ZUSEiI2rVrJ0my2Wzq1q2bAgIC1L17d2VnZzvWnzFjhpo1a6bg4GBt2rTJ\nsTwtLU2hoaFq3Lixxo8f74SpAAAAuIdbFiyLxaL169dr165d2rp1qyRpzpw5CggI0P79+3Xfffdp\n7ty5kqTMzEzNnj1ba9as0Zw5czRy5EjHdmJjYzV27Fht27ZNGzZs0Pbt2500JQAAgLJVolOEdru9\nwPdbt27VwIED5ePjowEDBiglJUWSlJKSoqioKAUEBCgyMlJ2u91xdCs9PV29e/dW7dq11bLine truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x1eef3e48>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 52
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Step4: get system matrix $\\mathbf{A}$ and right handside\n",
|
||||
"#### 4.1 Backgrounds: Weak formulation for Poisson equation with cell-centered system\n",
|
||||
"\n",
|
||||
"$$\\frac{1}{\\mu}\\mathbf{B}_s = (\\frac{1}{\\mu}_0-\\frac{1}{\\mu})\\mathbf{B}_0+\\nabla\\phi_s$$\n",
|
||||
"\n",
|
||||
"$$\\nabla \\cdot \\mathbf{B}_s = 0 $$\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"$$ (\\mathbf{B}_s\\cdot{\\bar{\\mathbf{n}}})_{\\partial\\Omega} = B_{sBC}$$\n",
|
||||
"\n",
|
||||
"In discretized form we have \n",
|
||||
"\n",
|
||||
"$$\\mathbf{diag(v)}\\mathbf{D}\\mathbf{P}_{in}^T \\mathbf{P}_{in}\\mathbf{B}_s + \\mathbf{diag(v)}\\mathbf{D} \\mathbf{P}_{out}^T B_{sBC}=0$$\n",
|
||||
"\n",
|
||||
"$$\\mathbf{M}^f_{\\frac{1}{\\mu}}\\mathbf{B}_s = \\mathbf{M}^f_{\\frac{1}{\\mu_0}}\\mathbf{B}_0-\\mathbf{M}^f_{\\frac{1}{\\mu}}\\mathbf{B}_0\n",
|
||||
"-\\mathbf{P}_{in}\\mathbf{P}_{in}^T\\mathbf{D}^T\\mathbf{diag(v)}\\phi_s$$\n",
|
||||
"\n",
|
||||
"$$ \\mathbf{Div} = \\mathbf{diag(v)}\\mathbf{D}\\mathbf{P}_{in}^T \\mathbf{P}_{in}$$\n",
|
||||
"\n",
|
||||
"$$ \\mathbf{Div}\\mathbf{B}_s + \\mathbf{diag(v)}\\mathbf{D} \\mathbf{P}_{out}^T B_{sBC}=0$$\n",
|
||||
"\n",
|
||||
"$$\\mathbf{B}_s = (\\mathbf{M}^f_{\\frac{1}{\\mu}})^{-1}\\mathbf{M}^f_{\\frac{1}{\\mu_0}}\\mathbf{B}_0-\\mathbf{B}_0\n",
|
||||
"-(\\mathbf{M}^f_{\\frac{1}{\\mu}})^{-1}\\mathbf{Div}^T\\phi_s$$\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"$$\\mathbf{A} = - \\mathbf{Div}(\\mathbf{M}^f_{\\frac{1}{\\mu}})^{-1}\\mathbf{Div}^{T}$$\n",
|
||||
"\n",
|
||||
"$$\\mathbf{rhs} = - \\mathbf{Div}(\\mathbf{M}^f_{\\frac{1}{\\mu}})^{-1}\\mathbf{M}^f_{\\frac{1}{\\mu_0}}\\mathbf{B}_0 + \\mathbf{Div}\\mathbf{B}_0-\\mathbf{diag}(v)\\mathbf{D} \\mathbf{P}_{out}^T B_{sBC}$$"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"Dface = M3.faceDiv\n",
|
||||
"Pbc,Pin, Pout = M3.getBCProjWF(BC, discretization='CC')\n",
|
||||
"Mc = Utils.sdiag(M3.vol)\n",
|
||||
"Div = Mc*Dface*Pin.T*Pin\n",
|
||||
"MfmuI = Utils.sdiag(1/M3.getFaceInnerProduct(mu = 1/mu).diagonal())\n",
|
||||
"Mfmu0 = 1/mu0*M3.getFaceInnerProduct()\n",
|
||||
"A = -Div*MfmuI*Div.T\n",
|
||||
"rhs = -Div*MfmuI*Mfmu0*B0 + Div*B0 - Mc*Dface*Pout.T*Bbc[ind]"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 53
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"print MfmuI.shape, Div.T.shape"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"(397851, 397851) (397851, 130050)\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 54
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"# import petsc4py\n",
|
||||
"# import sys\n",
|
||||
"# from sys import getrefcount\n",
|
||||
"# petsc4py.init(sys.argv)\n",
|
||||
"# from petsc4py import PETSc\n",
|
||||
"# import PETScIO as IO\n",
|
||||
"# Apetsc = PETSc.Mat().createAIJ(size=A.shape,csr=(A.indptr, A.indices, A.data))\n",
|
||||
"# bpetsc = IO.arrayToVec(rhs)\n",
|
||||
"# xpetsc = IO.arrayToVec(0*rhs)"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 55
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"# u1 = PETSc.Vec().createSeq(M3.nC)\n",
|
||||
"# u1.set(1)\n",
|
||||
"# u1.normalize()\n",
|
||||
"# basis = [u1]\n",
|
||||
"# nullsp = PETSc.NullSpace().create(False, basis, comm=PETSc.COMM_SELF)"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 56
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"# %%time\n",
|
||||
"# ksp = PETSc.KSP().create()\n",
|
||||
"# pc = PETSc.PC().create()\n",
|
||||
"# ksp.setOperators(Apetsc)\n",
|
||||
"# ksp.setType(ksp.Type.BCGS)\n",
|
||||
"# pc = ksp.getPC()\n",
|
||||
"# pc.setType(pc.Type.SOR)\n",
|
||||
"# OptDB = PETSc.Options()\n",
|
||||
"# OptDB[\"ksp_rtol\"] = 1e-8\n",
|
||||
"# OptDB[\"pc_factor_levels\"] = 1\n",
|
||||
"# ksp.setFromOptions()\n",
|
||||
"# ksp.view()\n",
|
||||
"# ksp.solve(bpetsc, xpetsc)\n",
|
||||
"# print ksp.its# print ksp.its\n",
|
||||
"# phi = IO.vecToArray(xpetsc)\n",
|
||||
"# print np.linalg.norm(A*phi-rhs)/norm(rhs)"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 57
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"%%time\n",
|
||||
"m1 = sp.linalg.interface.aslinearoperator(Utils.sdiag(-1/A.diagonal()))\n",
|
||||
"phi, info = sp.linalg.bicgstab(A, rhs, tol=1e-6, maxiter = 1000, M =m1)\n",
|
||||
"# m1 = sp.linalg.interface.aslinearoperator(Utils.sdiag(np.sqrt(1/(Utils.mkvc(-A.diagonal())))))\n",
|
||||
"# m2 = m1\n",
|
||||
"# phi, info = sp.linalg.qmr(A, rhs, tol = 1e-5, M1 = m1, M2=m2)"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"Wall time: 1.31 s\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 58
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"print info\n",
|
||||
"print M3.nC"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"0\n",
|
||||
"130050\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 59
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"B = MfmuI*Mfmu0*B0-B0-MfmuI*Div.T*phi\n",
|
||||
"rhsa = A*phi\n",
|
||||
"print info\n",
|
||||
"print np.linalg.norm(rhs-rhsa)/np.linalg.norm(rhs)"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"0\n",
|
||||
"9.95677949419e-07\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 60
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"#### BICG, qmr or BICGstab with Jacobi preconditioner works well (scipy), minres does not work ... Need to play with some other iterative solvers. Remeber that our system has null-space, which is constant"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"figsize(16,5)\n",
|
||||
"M3.plotImage(B, imageType='F')"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAA6cAAAE/CAYAAAC6trV3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsvXmcFNXV///u6Z4Z9n2dgQFc2ERQZDFuIEFFDQqIIIom\nX9yXqL8HX/qgGNHEuCR54gMRJTHqY9QQjSLuEomIS2QREIEBVBYB2UEYtlma/v1RXd23q2vtrl5m\nOG9evKaX6jqfunXr1j33nHsrEIlEIgiCIAiCIAiCIAhCDinItQBBEARBEARBEARBEOdUEARBEARB\nEARByDninAqCIAiCIAiCIAg5R5xTQRAEQRAEQRAEIeeIcyoIgiAIgiAIgiDkHHFOBUEQBEEQBEEQ\nhJwjzqmQUzp37kyDBg1o3LgxjRs3pkmTJmzbti3XsgRBEHxH2jtBEI4FpK0T0kGcUyGnBAIB3n77\nbSoqKqioqGD//v20a9cu17IEQRB8R9o7QRCOBaStE9JBnFMh7/jHP/7BcccdR0VFBQDvvfce7du3\nZ/fu3TlWJgiC4B+9evXi7bffjr2vrq6mVatWfPXVVzlUJQiC4D+33XZbLJLauHFjCgsLefDBB3Mt\nS8hDxDkVck4kEkl4P3bsWM444wxuv/12du/ezXXXXcdf//pXWrZsmSOFgiAI/qC2dz//+c958cUX\nY+/fffddSktL6dOnTy6kCYIg+Iaxb/enP/0pFkn95JNPaN68OSNGjMiROiGfCUSMtUcQskjnzp3Z\nvXs3oVAIgHPPPZfXX3+dffv20bt3b5o2bcqZZ57JU089lWOlgiAI6WFs7wYMGMDnn3/O1q1badSo\nEaNHj+b000/nrrvuyrFSQRCE1LHq2wHs3LmT/v378/jjjzNmzJhcyhTylFCuBQjHNoFAgNmzZzNk\nyJCEz5s2bcro0aP54x//GGvQBEEQajNm7d2wYcP45z//yYgRI3j//feZNm1aDhUKgiCkj1Xfrrq6\nmtGjRzN+/HhxTAVLJK1XyEuWLVvGc889x5VXXskvf/nLXMsRBEHICHpq76uvvsoZZ5xB+/btcy1J\nEAQhI/zyl7+kWbNm/OY3v8m1FCGPEedUyDuOHDnC+PHjeeSRR3j22WfZsmWLpPUKglAnGTlyJEuW\nLGHq1Klcc801uZYjCIKQEWbMmMH8+fMT5tkLghninAp5x6RJk+jUqRM33ngjRUVFvPjii0yePJnv\nvvsu19IEQRB8pV69eowaNYoNGzYwatSoXMsRBEHICDNnzmT9+vWUlJTEVux99NFHcy1LyENkQSRB\nEARByCEPPfQQ3377LS+88EKupQiCIAhCThHnVBAEQRByxO7du+nduzfvvfcevXv3zrUcQRAEQcgp\nktYrCIIgCDngL3/5C926dWPChAnimAqCIAgCEjkVBEEQBEEQBEEQ8oC8es5pIBDItQRBEPKQujaG\nJm2dIAhmSFsnCMKxgF1bl1fOKQAv2TTMNT7s3+aII1dqfwOXTYHLpuRUA8CUKVOYMsVCRxYRHaIj\nlzrqaufm3sj9Kf82SDgt2w/xMADnTBnE2VMG5VQDHDt1WXSIDjukrTMnnbbGj7bOLw1w7NRl0SE6\n7HBq62TOqSAIgiAIgiAIgpBzxDkVBEEQBEEQBEEQco44p2b0GJxrBQAMHjw41xIA0WFEdCSSLzoE\n75QN7pRrCUD+1CHRkYjoSCRfdAjekbYuEdGRiOhIJNc68mq13kAgkB9zTl+2+X0W55wKgqC1C3nU\nTPlCIBDIizmn9/GrnGsQBEFD2jpz/JjvmU5b55cGQRA0nNo6iZwKgiAIgiAIgiAIOUecU0EQBEEQ\nBEEQBCHniHPqBT9Sev3cjyAIgiAIgiAIQh1BnFNBEARBEARBEAQh59gszVPLUKORdeeoBEEQBEEQ\nBEEQjgnqRuTUmCYrabOCIAiCIAiCIAi1itrvnFo5orlwUGuU/4IgCIIgCIIgCIJrar9zmi9I9FYQ\nBEEQBEEQBCFlavfszHxxAPNFhyAIgiAIQi0hTDDhfZCwq9+p2xn3IQhC7UYip4IgCIIgCEJWMXMq\n3TiaRgfWrUMrCELtoG47p35HNL3Gma3s1+54tSAIgiAIQsqkGu20ckTFQRWEukPddk4FQRAEQRCE\nWoOV4yoOqCAcG4hzKgiCIAiCIAiCIOQccU7TRRZDEgRBEARBcEWmFjAKEpboqiDUAcQ5FQRBEARB\nEARBEHKOLM0jCIKQA9Id4fcj+pAPGgRBEARBEHQkcioIgiAIgiAIgiDkHHFO00Viz4IgCIIgCIIg\nCGkjzqkgCIIgCIKQFTK1aFE4uiSSIAi1G3FOBUEQBEEQhLzAynkVx1MQjg3qtnOa65TbXNsXBEEQ\nBEHIM/yOnorjKgh1h7rtnPqNX880lWejCoIgCIJwDGPmoDo5rWZOqDimglC3cHRODx48yM9//nO6\ndu1Kz549WbBgARUVFVx66aWUlZUxYsQIDhw4ENt+6tSpnHjiifTs2ZNPP/009nl5eTl9+/bluOOO\n47777vNHvUQmBUHwibxu6wRBEHwin9q6YGymaNh1NDVs+JUgCHULR+f0gQceoKysjOXLl7N8+XK6\nd+/OU089RVlZGd988w0dOnTg6aefBmDHjh1Mnz6duXPn8tRTT3H77bfH9jNx4kTuueceFi1axMcf\nf8zixYs9Cd02Evq20F5/fAFc0dl8ux7N4MDVUPWL5O+GlcDffgJbRsCXw+DRU6C0vgcNl1toMDjJ\nMQ3jk/dxYSk8dwZsHwXLLoQ7urm3LwhC5siXtu5uJtGeEgAmcB0n09t0u9a05gEm8xAPJH13Iidy\nOZdxD3dxKzdxAefRhCZZ1dCVE7mMkUzibm7jFs7gdNf2BUHIHPnS1gFM4m5KaA/AdUygNyebbtea\n1kzmAR7goaTvTuRELuNy7uIebuJWzuMC1+2dP/a7MpLLuJtJ3MJtnM4ZrmwLgmCOo3P64Ycfcu+9\n91KvXj1CoRBNmzZl4cKFXHvttRQXFzNhwgQWLFgAwIIFCxg2bBhlZWUMGjSISCQSG31bs2YNY8eO\npWXLlowaNSr2Gzcc3wgahGDpHigsgH4t4dMd0S8Vx7B+EF45F+ZuhUgkcR8tiuCVQbDjCAybB/d+\nBUPbwWOn+qBB0RHTsC1Zw2ktYfa5sGof/HQuPP0N3N8L/run66IQBCFD5ENb14IWFFLINrYSJEgJ\npWxkY9J2hRQyjjF8xzoiJDY09anPOMZwgIM8z9+Yw4ccz/FcwHlZ01BKCeO5kh3s4FmeZyGLOJfB\nnMPZrstCEITMkA9tHcTbmq1sI0iQUkos25oxjGMd35m2d2MYx0EO8Dee50PmcDzHcx4XZMV+CaVc\nyXh2sIPneZZFLGQw53I253gqC0EQ4tgmxm7evJkjR45w8803U15ezqhRo7j99ttZtGgR3bt3B6B7\n9+4sXLgQ0BqxHj16xH7frVs3FixYQKdOnWjTpk3s8549e/LSSy9x6623Jht9bUr8dY/B0HMwZ7aG\nBbshAvRvCbsrYfOh5J8++ROYvw0W7IQLOyR+d0lHCAbg7mUQjsDXP0LbejBjAPwiADWR5P2pnNka\nFuxy0BCCJ8+A+TuiGkoTv76jO7y3BX63Utt2xT44vjFM7A5/XA2VR+01CMKxwLx585g3b15Wbeai\nrft4yiex150Gl9F5cCfK6MRmNhEhQgmlHOYQ+9mX9NuLGc56NrKZzXTlxITvetCdAAE+YA5HOcp2\nttOIRlzKcF5jFkexb2j80PATTmct3/AJnwGwnR20oAVncQaf8x9qZOK9IBwzbd0nUz6OvS4b3IlO\ngzsD0IkyNrGZCBFKKeEQh9nH/qTfX8xwNrKezWzmRLomfNedHgQIMIcPEtq74VzKLF6zLQu39odz\nsaX90/kJ37CWz9Da8x1spwUtOIOz+A+fS1snCHhv62yd0yNHjrB27Vp+97vfMXToUG688UZeeeUV\nIsaQoA2BQCDpM9vfXzYl9nLvaM0ZLC6AggDsGatFLYuD2usI0PIf2lFc3QlOawX934RxJybvds4P\nWiTzlhPh+XXQph5c3QXe2WLvmLrWAFx9nBYd7f8OjOuSvK8mhXCgOvGzimpoUQw9m8LSvdY6BOFY\nYfDgwQwePDj2/sEHH8y4zVy0dYOmxKOIk5gMQJAgAQJMYjIFFBAixCQmEyHCozwMQB9OoYRSpjPD\nNAXtW74FYCD9WcIyGtKQU+jDGtbaOqZ+aiimmCqqEj6rpJL61Kc1rdnKVksdgnCscKy0dWdPGZTw\nfjKTgHhbM5lJsbZmMpOIEOFhHgXgFPpQQikzmG46vUBv7/ozkGUsoSEN6cMprGWNZXvn1X4pJZb2\npa0TBGe8tnW2zukJJ5xAt27dGD58OADjxo3jhRdeoH///pSXl3PqqadSXl5O//79ARg4cCAffvhh\n7PerV6+mf//+NG7cmO3bt8c+X7VqFaef7jz/qPe7EAC+uABuWgjLdsHMc+Dl9TB7U3y77k3h9wNg\n8ByoskhU/uEw9H0HPhkG/9NXczT/+T1c8ZkHDV/Asj02Gk6LarDo//3lG3h1EIwqgw92wE9awXXH\na991bCDOqSDkily3ddOZBgS4gZt4i9lsZStjGMtyvmI15bHtWtGaC7iQ53iGsMXiIfup4Eme4nqu\n5UKGESDASlbxD17NmobFfMk4xnISPfmGb+lIR/pxGgDNaCodNkHIEblu6wCmMZ0AcBM3MJu32MpW\nxjKGr1hOOatj27WmFRdygW1bU8F+nuJJruV6hnEhAQKsYiWv8g/f7D/Dc5b2v2QxYxlHT07iW76h\nIx05jX4ANKWZtHWCkAKOc05PPPFEFixYwNGjR3nnnXcYOnQoAwcO5Nlnn+Xw4cM8++yzsQLine truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x1f138828>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 61
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Step4: Compute analytic function (sphere in whole space)\n",
|
||||
"Outside of the sphere $(r>R)$\n",
|
||||
"\n",
|
||||
"$$\\mathbf{H}_1 = H_0 \\hat{x} + H_0\\frac{R}{r^5}\\frac{\\mu_2-\\mu_1}{\\mu_2+2\\mu_1}[(2x^2-y^2-z^2)\\hat{x}+(3xy)\\hat{y}+(3xz)\\hat{z}]$$\n",
|
||||
"\n",
|
||||
"$$H_{x1} = H_0 + H_0\\frac{R}{r^5}\\frac{\\mu_2-\\mu_1}{\\mu_2+2\\mu_1}(2x^2-y^2-z^2)$$\n",
|
||||
"\n",
|
||||
"$$H_{y1} = H_0\\frac{R}{r^5}\\frac{\\mu_2-\\mu_1}{\\mu_2+2\\mu_1}(3xy)$$\n",
|
||||
"\n",
|
||||
"$$H_{z1} = H_0\\frac{R}{r^5}\\frac{\\mu_2-\\mu_1}{\\mu_2+2\\mu_1}(3xz)$$\n",
|
||||
"\n",
|
||||
"Inside of the sphere $(r\\le R)$\n",
|
||||
"\n",
|
||||
"$$\\mathbf{H}_2 = H_0\\frac{3\\mu_1}{\\mu_2+2\\mu_1}\\hat{x}$$\n",
|
||||
"\n",
|
||||
"$$H_{x2} = H_0\\frac{3\\mu_1}{\\mu_2+2\\mu_1}$$"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Step5: Projection to receiver plane"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"xr = np.linspace(-400, 400, 41)\n",
|
||||
"yr = np.linspace(-400, 400, 41)\n",
|
||||
"X, Y = np.meshgrid(xr, yr)\n",
|
||||
"Z = np.ones((size(xr), size(yr)))*150"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 62
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"rxLoc = np.c_[Utils.mkvc(X), Utils.mkvc(Y), Utils.mkvc(Z)]\n",
|
||||
"Qfx = M3.getInterpolationMat(rxLoc,'Fx')\n",
|
||||
"Qfy = M3.getInterpolationMat(rxLoc,'Fy')\n",
|
||||
"Qfz = M3.getInterpolationMat(rxLoc,'Fz')"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 63
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"Bxr = np.reshape(Qfx*B, (size(xr), size(yr)), order='F')\n",
|
||||
"Byr = np.reshape(Qfy*B, (size(xr), size(yr)), order='F')\n",
|
||||
"Bzr = np.reshape(Qfz*B, (size(xr), size(yr)), order='F')\n",
|
||||
"H0 = Box/mu0\n",
|
||||
"\n",
|
||||
"flag = 'secondary'\n",
|
||||
"\n",
|
||||
"Bxra, Byra, Bzra = MagSphereAnalFun(X, Y, Z, 100, 0., 0., 0., mu0, mu0*(1+chiblk), H0, flag)\n",
|
||||
"Bxra = np.reshape(Bxra, (size(xr), size(yr)), order='F')\n",
|
||||
"Byra = np.reshape(Byra, (size(xr), size(yr)), order='F')\n",
|
||||
"Bzra = np.reshape(Bzra, (size(xr), size(yr)), order='F')"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"prompt_number": 64
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"figsize(10, 4)\n",
|
||||
"plot(Utils.mkvc(Bxra))\n",
|
||||
"plot(Utils.mkvc(Bxr), 'k:')\n",
|
||||
"plot(Utils.mkvc(Byra))\n",
|
||||
"plot(Utils.mkvc(Byr), 'k:')\n",
|
||||
"plot(Utils.mkvc(Bzra))\n",
|
||||
"plot(Utils.mkvc(Bzr), 'k:')"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "pyout",
|
||||
"prompt_number": 65,
|
||||
"text": [
|
||||
"[<matplotlib.lines.Line2D at 0x2b71c908>]"
|
||||
]
|
||||
},
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAAnUAAAD9CAYAAADai4HIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd4VFX6wPHvTJJJ7z1AAgkpBEKH0KtSRRSxrtgVu+uK\nuurqsvpz1bU3FOu6Ioqi2JAqQighoYYQ0kjvPZNMpmQmc39/3EkMksxcjA08n+fxefZmz5u5EzKT\nd855z3tUkiRJCIIgCIIgCGc19e99A4IgCIIgCELfiaROEARBEAThHCCSOkEQBEEQhHOASOoEQRAE\nQRDOASKpEwRBEARBOAeIpE4QBEEQBOEc0OekLiUlhSFDhhAbG8urr77a45iHHnqI6OhoxowZQ05O\njsPYzz77jKFDh+Lk5MThw4dP+V6vvPIKsbGxJCYmsmfPnr7eviAIgiAIwjmhz0ndPffcw+rVq9m+\nfTuvv/469fX1p/z/6enp7N69m4MHD7JixQpWrFjhMDYpKYkNGzYwbdq0U75XbW0tq1at4vvvv+eN\nN97g7rvv7uvtC4IgCIIgnBP6lNRptVoApk2bRlRUFHPmzCEtLe2UMWlpaSxdupSAgACuvPJKsrOz\nHcYmJCQQFxd32uOlpaUxb948IiMjmT59OpIk0dra2penIAiCIAiCcE5w7kvwgQMHSEhI6LpOTExk\n//79LFy4sOtr6enpLFu2rOs6ODiYgoICioqKHMb+VHp6OkOGDOm6jo+PJz09ndmzZ3d9TaVS9eUp\nCYIgCIIg/KZ+qcO9fvWNEpIknXazPzfx6ulJ9/S9Oh/zz/TfP//5z9/9HsTzFs9bPG/xvMXzFs9b\nPO8z+++X1Kekbty4cadsfMjKymLChAmnjElOTubEiRNd13V1dURHRzN27FiHsT/10++Vk5PDuHHj\n+vIUBEEQBEEQzgl9Sup8fX0BeRdrcXEx27ZtIzk5+ZQxycnJfP755zQ0NLB27dqu5VM/Pz+HsXDq\n7Nz48ePZsmULpaWl7Ny5E7Vajbe3d1+egiAIgiAIwjmhTzV1AC+99BLLly/HbDZz9913ExQUxOrV\nqwFYvnw548ePZ8qUKYwdO5aAgADWrFljNxZgw4YN3H333dTX17Nw4UJGjRrFpk2bCA0N5bbbbmPW\nrFloNJquxxFgxowZv/ct/C7E8/5zEc/7z0U87z+XP+vz/iWppF96Qfd3plKpfvE1akEQBEEQhF/D\nL5m3iBMlBEEQBEEQzgEiqRMEQRAEQTgHiKROEARBEAThHCCSOkEQBEEQhHOASOoEQRAEQRDOASKp\nEwRBEARBOAeIpE4QBEEQBOEcIJI6QRAEQRCEc4BI6gRBEM4G774L9933e9+FIAh/YOJECUEQhLNA\nXUwMrYWFRIv3N0E4p4gTJQRBEM5mVivs23dGIZ/rdDz0K92OIAjnBjFTJwiC8FtLSYHp0+EM3quk\n+HhUeXlnFCMIwh+fmKkTBEE4i7VWV5MNYDAojlGpVL/a/QiCcG4QSZ0gCMJv7GRGBi8AlJcrjrmr\npoY0OLOZuocfhqNHz/DuBEE4W4mkThAE4Tc2yseHBYC2rExxzGUaDVPgjGb3eOopWL36TG9PEISz\nlEjqBEEQfmtGIweBlro6xSFTOzowA7S1KQuwWjkInGhs/Bk3KAjC2UgkdYIgCL+x0poabgMGeHkp\nD9LpICREeVKn05EDlNbX/5xbFAThLCSSOkEQhN/Ylrw8toDyBM1iYaHJRL2vLx1araIQU10du4F5\n3t4/9zYFQTjLiKROEAThN3bz4MEMBgqLi5UFGAzc7+bGzPJyigsKFIWodDqmArS2/sy7FAThbCOS\nOkEQhN+awcAJjYbamhpl441GZnh5kTl2LDGBgYpCNEYjs4FNZ7AZQxCEs5tI6gRBEH5j+XV1LA0O\nZkJIiLIAoxHc3MDdXfnu19ZWWkJCyGhq+vk3KgjCWUUkdYIgCL+xNQUFHHJzA71e0fi6ykqWNjVh\ncHHBrHA59WRuLp86O/N3F5e+3KogCGcRkdQJgiD8xv4VFUWYvz+HFdbUeTs58bfQUK7OyGDHoUPK\nYoDRERGKE0cAzGbw8gKLRXmMIAh/GCKpEwRB+K0ZDJQ4O1OssN2ImyQxKSCAz887j7lxcYpiQjUa\nzo+L4386neLbkrKz+WdbG1J2tuIYQRD+OERSJwiC0FdJSVBUpHj4scZGZsXGsqR/f2UBRiO4up5Z\nTZ3BAL6+pHZ0yDNwCnRkZ5MGWE6cUPYYgiD8oYikThAEoS8kiZuPH6c5JUVxyKtVVRSrVHKypkBm\ndjZ3lpRgdHHBpLCmbkdmJmuLi3nDy0txIuis1bIZcGlpUTReEIQ/FpHUCYIg9EVdHfMAT6WNhIG3\nAwLwCAxkV2mpovEDvL25YeBAHtm/nw9SUxXF9HN1JT4sDDw8lNfVde6UFUeLCcJZSSR1giAIfVFX\nxxyg9Uz6wRmN1KrV5ChMnvycnRkdGsrzixZxy7BhimLiPT2ZGBfHWquVtoYGRTFfpqVxgbMzX6en\nKxovCMIfi0jqBEEQ+kKv503gwwMHFIccaG1lzPDhLA8NVRZgMp15nzqDAdzdOWK1YlJ4tNgQZ2cu\n69ePGLX40yAIZyPxyhUEQeiD3JwcGoF7EhMVx6zU6Whzc1NcU7ft8GH+LyeHdhcXDApr6j7IyGB7\ncTHPDhpEgMJedfFOTlyTlMRQJydF4wVB+GMRSZ0gCEIfBGk0zAFQWlMnSWyUJNo9Pfm2qkpRyPCg\nIJbExbE6LY0n0tIUxSS4udH/TGvq2togLAzOoA2KIAh/HM6/9w0IgiCczQKdnRkP1DU2EqwkwGwG\nZ2darFYyWlq4QEFIqEZDaEQEiWPGKG4MnOzpCTExbGxrY2RZGf0UxDyVlUWgxUJTTg4PKnoUQRD+\nSMRMnSAIQl/o9Wxzc+OJjAxFwzva2kh1cSE2IYFHfHyUPUb3s18VLtl21tTltLejVbhRYpqbG2OG\nDWPCmRwt9tFHEKwonRUE4VcmkjpBEIQ++P7wYco8PXklOlrReFNzM4+azXKSpjBB+196Ou8cP45F\no6FN4dLoU7m5ZNXUcN/w4SQGBiqKmaxWM2bcOKYrGm2zYwfU14ujxQThD0AkdYIgCH0w2MuLSWFh\nimvqPNRqtoeHo7da+aS5WVHM9LAwZg0ZwqaMDG7OzFQUM0ajITAk5Mxq6nQ6CA09o5o6qbGRtYC1\npkZxjCAIv44+J3UpKSkMGTKE2NhYXn311R7HPPTQQ0RHRzNmzBhycnIcxra2trJ48WIiIyO56KKL\n0NneYIqLi3F3d2fUqFGMGjWK22+/va+3LwiC0CdRbm4kRUZSoTBBw2AANzfMTk7sVzhTF6XREB0R\nwaLzzmNtZKSimDkaDWH9+7O3uZmsggJFMbdWV5Oj13NrXZ2i8QDPHD7MfwCDwkbKgiD8evqc1N1z\nzz2sXr2a7du38/rrr1P/kwOq09PT2b17NwcPHmTFihWsWLGi19gGW93HG2+8QWRkJPn5+fTv3583\n33yzK2bw4MEcOXKEI0eOsGrVqr7eviAIQt/o9ZzUaLippETR8LamJg5LEr5hYbwkScoeo7Omzs3t\nzPrUeXhQYjJRo7CmbrEkEZKQwFyLBRTe21SLhe/698dTHC0mCL+7PiV1WltDy2nTphEVFcWcOXNI\n+8l2+7S0NJYuXUpAQABXXnkl2dnZvcbu378fkBPBG2+8EVdXV2644YbTvqcgCMIfxQeHDpEjSWzy\n9VU0vqq8nJebmsDZ1nxAQS3afw4f5pvsbDo0GloVLqXeU1VFg9HIVWPHMis83HGA1cr89nYCIiO5\n2NVVcfI4Wa8nYtgwua5OEITfVZ9amhw4cICEhISu68TERPbv38/ChQu7vpaens6yZcu6roODgyko\nKKCoqKjX2O7fNyEhgfRuR9YUFRUxcuRIkpOTuf322xkxYsRp97Vy5cqu/z1jxgxmzJjRl6cpCILQ\nqzE+PmgGDICdOxWNHxwUxAdDhwLwvpMT1+h0OPn52Y25KDgY74QEssvLWV5RwV4FjzNNpcLdz0+u\nqVPSD8+2LIxaDV5ecl2dh4f9GEmC1lb2urgQkp9PrIL7EoQ/u507d7JT4fvFmfrV+9RJkoT0k2l8\nlUrV49jOr/90fKeIiAjKysrw9/dn06ZNLFu2jGPHjp02rntSJwiC8Gsa5uaGNSGBYp2OgZIEvby/\ndTEa5dYkQLpKxVUKkro4Z2eIiCB81Cj2enoquq9LOjogMJDjzc1oCwuZ7GC8saGB6zo6+AS42WTi\n2YoK/EJC7MZYW1tZCiwwm7FUVoqkThAU+Olk07/+9a9f7Hv3afl13Lhxp2x8yMrKYsKECaeMSU5O\n5sSJE13XdXV1REdHM3bs2NNik5OTu75v5zJtdnY248aNA0Cj0eDv7w/A/PnzcXZ25uTJk315CoIg\nCH2j19Ph7c15VqvcWNiB+poaTrS3A/BGcDCuSmrXuvepU7IsKkldyWOt2Uypgpo6J6ORv9iSy0UB\nAWhMJscPo9VyrY8PN02dynSFbVMEQfj19Cmp87XVkKSkpFBcXMy2bdu6ErNOycnJfP755zQ0NLB2\n7VqGDBkCgJ/tzaOn2OTkZN577z0MBgPvvfdeV6JYX19PR0cHAIcPH8ZgMDB48OC+PAVBEIQ++Xdm\nJoerqjjp66uorcmJwkI+rK6WLxT2qvvriRMcKilB0mjQmkxgtdodb9bpuE6lArWaWaNHc2VYmMPH\ncDGZWGRrInxhRAQetvdae5x0OhaHhICvL9jqpBVZswZSUpSPFwRBkT4vv7700kssX74cs9nM3Xff\nTVBQEKtXrwZg+fLljB8/nilTpjB27FgCAgJYs2aN3ViA2267jauvvpr4+HhGjx7NM888A8Line truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x1eed9438>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 65
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"fig, ax = subplots(3,2, figsize = (10,15))\n",
|
||||
"dat1 = ax[0,0].imshow(Bxr); fig.colorbar(dat1, ax=ax[0,0])\n",
|
||||
"dat2 = ax[0,1].imshow(Bxra); fig.colorbar(dat2, ax=ax[0,1])\n",
|
||||
"dat3 = ax[1,0].imshow(Byr); fig.colorbar(dat3, ax=ax[1,0])\n",
|
||||
"dat4 = ax[1,1].imshow(Byra); fig.colorbar(dat4, ax=ax[1,1])\n",
|
||||
"dat5 = ax[2,0].imshow(Bzr); fig.colorbar(dat5, ax=ax[2,0])\n",
|
||||
"dat6 = ax[2,1].imshow(Bzra); fig.colorbar(dat6, ax=ax[2,1])"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "pyout",
|
||||
"prompt_number": 66,
|
||||
"text": [
|
||||
"<matplotlib.colorbar.Colorbar instance at 0x00000000236DAD88>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAAmYAAANVCAYAAADStTI9AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsvX+sHkd5P/qZ3ffHsWMTYhLC5YaTNqqvf1AKcWI7EimY\nFIhLRB1aaBtqqW2iYhw1TiWSVhGpZCgIRUIqJAQnfCWr3IZcqECBK1UiSuAakyvFJ2nNvZfE6Q1c\naCooAWzHObaPz/vu7tw/Zmb32eedmZ19f51z7P1Ic3bO7szu7O7sZz/PM8/OK6SUEg0aNGjQoEGD\nBg2WHNFSN6BBgwYNGjRo0KCBQiPMGjRo0KBBgwYNlgkaYdagQYMGDRo0aLBM0AizBg0aNGjQoEGD\nZYJGmDVo0KBBgwYNGiwTNMKsQYMGDRo0aNBgmaARZg0aXKA4fPgwNm3ahPXr1+OBBx6wlrnnnntw\n1VVX4ZprrsELL7xQWXd+fh67du3C7Owsbr75Zpw+fRoAcOLECbzrXe/C2rVrcccdd5SOsWPHDmzc\nuBFXX301rr76avz617+ewNk2aNCgwcpAI8waNFhhWLduHYQQtdK6desG9nPnnXfi4YcfxpNPPokH\nH3xwQBDNzc3h+9//Pp599lncdddduOuuu5x1jx8/DgA4cOAAZmdn8eKLL+KKK67AQw89BACYmZnB\npz71KXz2s58daIcQAo8++iiOHj2Ko0eP4tJLLx3n5WrQoMEyw+ox8Nf5jNZSN6BBgwb1cPLkSeyv\nWWf/yZOl/0+dOgUAeMc73gEAeO9734sjR47gpptuysscOXIEH/zgB7Fu3TrccsstuPfee511n376\nadx0002Ym5vDvffei263i1tvvRWf+cxnAACrV6/G29/+drz44ovW9jXzXDdocOFgAajFYZy/znc0\nHrMGDVYgWhXpvwAcJonjmWeewcaNG/P/N2/ejKeffrpUZm5uDps3b87/v+yyy/DjH//YW5du27hx\nI+bm5kr7FEJYz+fP//zP8Z73vAdf/vKXq0++QYMGKx5VHEbThYYL8ZwbNFjxaFds/190MvjuEMeQ\nUg54slzCyqwfxvP1la98BW984xvxn//5n/jQhz6EN7/5zbj22mvrN7hBgwYrBlUcdiGj8Zg1aLAC\nUcfatFlfW7duLQXzP/fcc7juuutKZbZv347nn38+//9Xv/oVrrrqKlx77bUDdbdv357v99ixYwCA\nY8eOYevWrZXn8sY3vhEAcOWVV2L37t147LHHKus0aNBgZaPxmLnRCLMGDVYg2jUTx8UXXwxAfV35\n05/+FE888UQurgy2b9+Ob3zjGzh+/DgeffRRbNq0CQDw2te+1ll3+/btOHjwIBYWFnDw4MEBscc9\namma5h8dvPrqq3jsscfwvve9b7iL0qBBgxWDUfjrfMeFKEYbNFjxGMeD+7nPfQ579uxBv9/Hvn37\ncOmll+Lhhx8GAOzZswfbtm3D9ddfj2uvvRbr1q3DI4884q0LAHv37sXu3buxYcMGbNmyBffdd19e\n5zd+4zcwPz+PXq+Hb37zm3jiiScwOzuLnTt3ot/v4zWveQ3e//734+1vf/sYzq5BgwbLGY34cEPI\n5nOoBg1WFIQQ+B816/wVmi8fGzRosDxQl8Nc/HX48GHs2bMHSZJg3759A3MkAmouxq997Wu45JJL\n8JWvfCX/OMlVd35+Hrt378bRo0exZcsWPPLII1izZg1OnDiBP/qjP8Kzzz6Lv/iLvyjN37hjxw78\n4he/wKpVqwAATzzxxEjT/jRDmQ0arECMGmPWoEGDBkuJcfDX+ToXYyPMGjRYgRg1xqxBgwYNlhKj\n8hedT/HKK6/M52Kk4HMxmg+TbHXNlD9zc3O47bbb8rkYzT7NXIzdbtfannGOSDTCrEGDFYjGY9ag\nQYOVDB9f/b8AvkmSDefzXIwNZzdosALReMEaNGiwkuHjsLfqZPD1IY+xUudibDxmDRqsQDRDmQ0a\nNFjJGJW/zue5GBth1qDBCkQzlNmgQYOVjFH563yei7Hh7AYNViAaL1iDBg1WMsbBYefrXIzNPGYN\nGqwwCCHwf9Ss8y4085g1aNBgeaAuh11o/NV4zBo0WIFoPGYNGjRYyWg4zI1GmDVosALRPLgNGjRY\nyWg4zI3m2jRosAKxaqkb0KBBgwYjoOEwNxph1qDBCkQzDNCgQYOVjIbD3Bh6uozDhw9j06ZNWL9+\nfenHPBs0aDB5NNNljI6Gwxo0WDo0/OXG0F9lXn311fj85z+PK6+8EjfeeCOeeuqp0g93umbXbdCg\nwSDqPIZCCPyk5v5/s+YxLgT4OKzhrwYNwlGXW+py2IXGX0OJUfoDoADyHw+96aabSuXktcD+nwP7\n3+jYkXAsbeUEK+da58qPWHf/C8D+TTX2F7GlLfFz9+2z6roEvEf2Pw3sv666XAmSpTrlqvI0ZWzp\n2ef+54H9Gx3l6hzXt863X9f+bdfEtkT5uRDPOup70AwDjIYQDpM3qOX+nwD7f9OyExdvhMDFA7Zy\nIfuv4KX9zwP73zzicX3cGdK+uudhKbd/Dti/zVG3LlfBUy6UM4bY3/4fAvs31zyGr31Vxw3hqiHq\nmudCfLfGPgkaDnNjqKHMkB8PbdCgweTQDGWOhobDGjRYWjT85cZEz3n/z4BD88B+CexYC+x4zSSP\ntkwQYo2Oun9+rNDj1bXqpaf8uL3Kgh3vPPZaH5rXz8XPht9Hu+6Tmwx/rAsV+/VYy6GTwKHXAjsu\nsRQK7bO250iSbdJTzlZ+UrDxC837RiN8basaHXFtt5VzXUtax+VJ4ver6nrz7Xyda3+T5K+QfdvO\ne0z7OnQSOPSKWu4PqO5CLQ67wPhrKGG2detW3H333fn/zz33HHbu3DlQbv//BBxao0RZMPksU+y4\ntLpMjroCqM5++TFqDKXsuGKIdtmIyBxvEuKMHsOBWvdiGWLHGvVs7FgLQAKf+O/6+2g1wmwkhHCY\nGb50ijKDquegivt8Ak2S9fz5qPEs77gsoI38/zqizMc/rmFR27ErhNmO/9lzHNexq4Y1XfXM9rri\njNdnqLwXLsgh87ZzdNVxlSNldlyiknkuPvETT30PanHYBcZfQw1lhvx4qMGOtTqzwj0gXjEQGmcx\nCkJIsgI73jTEMauOOwlUXMuhiW0ZwYiyYdGO66UGZdTiMJcoq/tys60LiU+idUKPyRD8zPgEWIgo\n4+WrytBYXFtcLj+PKyraXbU/XzvrnJfvOlQZyK/3b7ei7n239akqgVZnGyqMlQA0/OXG0EOZrh8A\nLcF2Y7n1d75hkiKtijAmCZvLfpxi2+aVM+tXuKgvYdhAXIbaHrMGAwjisCqEeJAp5/mGvnzl4Kg3\nDri8VDZuqRIvvv3byld4yIaC6164rmnIe6qO58y2HBdCvGBVHwbU4aAJcm/DYW5M7EfMhRCQb6Mr\n4Lc06HJgZ576PrEyrbpVX2Da6rpIsKqNVevGDcnSMOtced86iVpfao7tuHBsCzlHjsD9iR+g9nQZ\n8nXBxVWd4/WOcaFDCJF/lTn6zmB/zl3lQoXOsDw3Du7z8VfVuS1H/vLVdeUnzUG+fdRp77DlAiG+\nW59b6nLYhcZfy0OzcgujTj2Qur687X9bedqeSWNcoozvb1zg18Ws48epWsfzk4Rk+VFFmW2/ddpR\nRabD4gJ0708dVfesLlfxZ8BW38Vjyxk+MWkTcLROiBAdFlW8RK8xv9fT5KxRYWv7KOVCjzMqGg5z\nYuiZ/8cC24uvroLnL9cqK8VWh+9r0qiyROuKMte2qGYK3W9Iu1z7mAaGsUjhKT9Mf6T7mwTGMF9G\nyMz399xzD6666ipcc801eOGFFyrrzs/PY9euXZidncXNN9+M06dP59vuv/9+rF+/Hps3b8ZTTz2V\nrz927Bi2bNmCq666Ch//+MeHuRpLg1G4qqr+mL0aUwF//n1zOIZ45Ebhr1CO9bV1uaOKz+qWCznO\nONHMl+HE0gozwO2VqOoIVS9UOPLjePkOC5uHaxTxw/c5DCmGEuYw7Qtp/6QRIsqqtsGSdx2rKj8u\njEGY3XnnnXj44Yfx5JNP4sEHH8Svf/3r0va5uTl8//vfx7PPPou77roLd911l7Pu8ePHAQAHDhzA\n7OwsXnzxRVxxxRV46KGHAAC//OUv8cUvfhHf+c53cODAAezbty/f18c+9jH83d/9HZ555hl873vf\nw7PPDjHj7iRg4wifoPfth+ddL1B63JB9LzWqhI4vCL9KlNXlqzriLLR9K0Wg0SXP1y1nq2Py43xv\nNsLMiaUXZsDwRGSrNw7P2TRQR6RViR6+P05AodZmVLG/YdrpOt9Jw0Ugw4iyYSzMYerVwYjCjM58\nf+WVV+Yz31McOXIEH/zgB7Fu3TrccsstOHbsmLOumZx1bm4Ot912G7rdLm699dZ8n0eOHMHOnTsx\nOzuLd77znZBS5t60//iP/8Cf/Mmf4HWvex3+8A//cKAdyxbj4Cq+n9Byyxk2ERQiynx163rL6uy7\nSogtR3Hm6yN1nRZ1jjdONMLMieUhzAyqiMjVOVxWbB1x5tv/qLCJFL6trigLWVfH6qwzlDls2/i5\nTwpV/cDXL3i90D4x7Zdo7E+HFoH9rxSJI2Tm+7m5OWzevDn//7LLLsOPf/xjb126bePGjZibmwOg\nhNmmTZvyOhs2bMCRI0fwox/9CK9//eut+1pyhHjM6t7v0H61VAajDa5n1sdjPuFTxSkhP2c3TJ1h\n2+m7BssBVX3FxYch++XLcXrMKjislC4wrBwtKlE8HKYzTOJh4fse57FcFpmtTB1RxJe27VWg19fW\nJkmWIXV9DyzdX1XZKrhelJNAqKEwDVQ8uTteU/6ljU+cqH8IKeXAl1CuH/c26+t+XWo75oLine truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x1f133518>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 66
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"fig, axes = subplots(1,3, figsize=(14,5))\n",
|
||||
"epsx = np.linalg.norm(Utils.mkvc(Bxr))*1e-6\n",
|
||||
"epsy = np.linalg.norm(Utils.mkvc(Byr))*1e-6\n",
|
||||
"epsz = np.linalg.norm(Utils.mkvc(Bzr))*1e-6\n",
|
||||
"id = 21\n",
|
||||
"axes[0].plot(Y[:,id], Bxra[:,id], 'b', Y[:,id], Bxr[:,id], 'r.')\n",
|
||||
"axes[1].plot(Y[:,id], Byra[:,id], 'b', Y[:,id], Byr[:,id], 'r.')\n",
|
||||
"axes[2].plot(Y[:,id], Bzra[:,id], 'b', Y[:,id], Bzr[:,id], 'r.')\n",
|
||||
"print X[1,id]"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"output_type": "stream",
|
||||
"stream": "stdout",
|
||||
"text": [
|
||||
"20.0\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAA1IAAAE1CAYAAADtSNHMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XtcVHX++PEX3kjSnFAzTcEw4tJWgoxoKlKWqZS6pZmb\n5nd1G7OLruVmbm3pb7c2t9o01xB3v/TdLdHcNdvaUkEL0FpmMCmVi6CBppUXSoXyiuf3xwdG7nKZ\nmTPnzPv5ePAohjln3jNwjud9Pu/P++OnaZqGEEIIIYQQQogma6N3AEIIIYQQQghhNJJICSGEEEII\nIUQzSSIlhBBCCCGEEM0kiZQQQgghhBBCNJMkUkIIIYQQQgjRTJJICSGEEEIIIUQztTqRyszMJCIi\ngtDQUJYtW1bvcxYsWEBISAgDBgygoKDgktuWlZUxbtw4goKCGD9+POXl5c6fvf7664SGhhIZGcm2\nbdtaG74QwsM8ec74/vvvufXWW+ncuTOPP/54jdeIj48nPDycqKgooqKiOHbsmBverRDCnTx9DQJw\n4MABOnXqxKuvvup8LD8/n+joaEJCQnjmmWdc/C6FEF5La6X+/ftrGRkZWklJiRYWFqYdPXq0xs/t\ndrs2ZMgQrbS0VEtJSdESEhIa3PbYsWOapmna4sWLtccee0w7ffq09uijj2ovv/yypmmadvjwYS0s\nLEzbv3+/lp6erkVFRbU2fCGEh3nynPHjjz9q27Zt01asWKE99thjNV4nPj5e+/zzz938boUQ7uTJ\n80mVe++9V7vvvvu0V155xfnY6NGjtTVr1mjHjh3ThgwZomVnZ7vxXQshvEWrRqROnDgBQFxcHMHB\nwYwcORK73V7jOXa7nQkTJhAYGMjkyZPJz89vcNusrCwAHA4HM2bMwN/fn+nTpzv3abfbGTVqFEFB\nQQwfPhxN0ygrK2vNWxBCeJCnzxkBAQEMGTIEf3//euPRZD1yIQzL0+cTgPfee4+QkBAiIyNrvM6e\nPXuYNGkSXbt25Z577qkThxDCnFqVSGVnZxMeHu78PjIy0nkiquJwOGqccLp3786+ffsa3bb6z8LD\nw3E4HIA6IUZERDi3CQsLc/6sip+fn3zJl3x54Zce54zq54X6TJs2jTvuuIO///3v9f5c789MvuRL\nvup+VfH0+aS8vJw//elPLFy4sMZr7N27l6uuuqrROPT+zORLvuSr/q/WcnuzCU3T6tz1bSjwqseb\nc5e4vn1Vvaanvp5//nlTv568pnleT6/XbI76tnHlOaPKqlWr2L17N3/7299Yvnw527dvbzQes/5u\n5DXN8Xq+8prNVd92LT2fLFy4kLlz5xIQEFDjObWf39D2Zv/d+Mpr+sJ79JXP1RValUhZrdYaEzdz\nc3MZNGhQjefExsaSl5fn/P7o0aOEhIQQExNTZ9vY2FjnfquG3/Pz87FarfXuq6CgwPkzIYT38/Q5\nozG9evUCIDg4mClTprB+/fqWvzEhhMd5+nzicDh46qmnuPbaa1m6dCkvvvgib7zxBqGhoRw+fNi5\nr7y8vDpxCCHMqVWJVJcuXQDV+aakpIS0tDTniahKbGws69ato7S0lJSUFGdpnsViaXDb2NhYkpOT\nOXXqFMnJyc4T0sCBA9m0aRMHDhwgPT2dNm3a0Llz59a8BSGEB3n6nFGl9p2niooKZ5e+kydPsn79\nesaMGeP6NyyEcBtPn08yMzMpLi6muLiYX//61zzzzDM88sgjgCoBXLNmDceOHWP9+vV14hBCmJTW\nSunp6Vp4eLjWr18/benSpZqmadqKFSu0FStWOJ8zf/58rW/fvlp0dLSWl5fX6LaapmknT57Uxo4d\nq/Xp00cbN26cVlZW5vzZkiVLtH79+mkRERFaZmZmnXhc8Jaa7ZNPPjH168lrmuf19HrN6selp88Z\nwcHBWmBgoNapUyetd+/eWn5+vvbjjz9qAwYM0G666SZt6NCh2uLFiy8Ztyf4yt+DL7ymL7xHPV6z\n9jHp6fNJlYULF2qvvvqq8/vc3FwtKipK69u3r/b0009fMm5P8IW/Bz1e0xfeox6vqfe1SUv5Ve7I\nNPz8/FxW9yiEcA2jHpdGjVsIszLqMWnUuIUwM1ccl25vNiGEEEIIIYQQZiOJlBBCCCGEEEI0kyRS\nQgghhBBCCNFMkkgJIYQQQgghRDNJIiWEEEIIIYQQzSSJlBBCCCGEEEI0Uzu9A3AHqxXGjIGEBIiJ\ngTaSLgrhcT/+CB9/DB99pHckQgghhBCuZ8oU4+WX4aef4H/+B3r2hGnTYO1aKC/XOzIhzO3AAVi2\nDEaPVsfea6/BddfpHZUQQgghhOuZfkHe4mJ1R/w//4HPPoNRo+D++9WF3mWX6RioECZx+DD885+w\nejXs2QN3361Gg++4A7p0Uc8x6mKURo1bCLMy6jFp1LiFMDNXHJemT6SqO3YM1q2DNWvgyy9h3Dj4\nxS9gxAgp/xOiOcrKVPK0Zg04HHDXXTB5skqeOnSo+3yjXkQYNW4hzMqox6RR4xbCzCSRqkdTP5RD\nh1S53z/+ASdPwqxZ8MtfQteuHghSCIPKzYU33lCjT8OGwQMPqCQqIKDx7Yx6EWHUuIUwK6Mek0aN\nW9Ris3E+r5DS0wH02JwCFoveEYlWcMVx6bPjMNdcA3Pnwo4dsGoV7NwJ/fqpZCo7W+/ohPAe586p\nmw7x8WrEqVs3dbz8+99w332XTqKEEEIIUygspN2nGfT4fAPlD9j0jkZ4AZ9NpKr4+cGgQWpkqqgI\nwsPVxeHgwbBxI8gNJOGrTp2CJUugb181CvXII1BSAosWQe/eekcnhBBCeNbpNurOYXE3K09ZVuoc\njfAGPp9IVde9O8yfD3v3qtGqJ5+UhEr4nqoEql8/yMhQjVrS09UNhvrmPwkhhBC+4LnrUsgJnUjg\n9lT+mWphzx69IxJ689k5Uk1x4QL861/qDnznzrBwIdx5pxrFEsJsTp2CpCT4058gNhaeew6iolyz\nb6PODzBq3EKYlVGPSaPGLS7av1/9m1hQAFddBX/8oypzX71a78hES0mziXq442RVPaHq3h2WL4cb\nbnDpSwihG01Tf99PPKEWsHZlAlXFqBcRRo1bCLMy6jFp1LjFRQ89pK4BX3xRfV9ertZJTE2Fm27S\nNzbRMpJI1cOdJ6uKCkhMVAnV//wPPP88dOrklpcSwiMKC+Gxx+Dbb9UNgrg497yOUS8ijBq3EGZl\n1GPSqHELZe9eNZ++qAiuvPLi46+9BpmZsH69frGJlpOufQ0ZMwaOH3f5btu2VRedu3erRUgjI9Wd\nfDk3CqP56Sd49lm45Ra1SPWOHe5LorBJZyMhhBDGtWgRzJlTM4kCePhh1el5+3Z94hL6M+eIFMDE\niapnM6gLucJC1ac5xXV9/zMzVSeza66BlSshONgluxXCrTZtUif/2Fh49VX19+sSDR1n8fH4ZWQY\n8m6s3EUWwrsY9Zg0atwC8vLU8h9798IVV9T9eWIivP8+bNjg8dBEK8mIVEOsVpXZVCksVO3HNmxw\n6d3xuDjIyVEHmNUqQ7vCu507B089Bb/6lTo81qxxYRIFDR9nstCUEEIIg3r+efjNb+pPogBmzFAN\nKLZt82xcwjuYM5FKTa056lR1IVc9wbLZVAbUyjLA9u1hwQJ1N+KJJ1Tp3+nTLQ9dCHcoLoZhwyA3\nV5Xx3XFHK3bW0LFT33EGanRKCCGEMJicHPj0U3j00Yaf06GDatL07LMy1cMXmTORql26l5KiSv2q\nJ1guHqUaNEgdcN99p9aeKixs9S6FcIl//UuV8d13H3zwgeo61CoNHTv1HWfgslJaIYQQwpOee07d\nLL9UYcXUqfDNN7Bli2fiEt7DnHOkmvKWxoxRF4JWa80Lv1bOp9I0WLFCHXyvvQZTprTgTQjhAqdP\nq4WlU1NVGZ/V2swdNHQsNHTsNMKo8wOMGrcQZmXUY9KocfuykyehZ08oLYXLLrv085ctgy+/hL/9\nzf2xCdeQ9uf1aPKHcvy4ulBcubLmhWB8vLrbDjUbVjTTzp1wzz0qkXr+eVnEV3jW8eMwfjwEBsKb\nb0KXLi3YSUPHQkPHTiOMehFh1LiFMCujHpNGjduXfTvWxrfphUQPbcKNdZuNk9sL2bkvgKH7XdfU\nTLiXNJtoDYtFXRjW/mNvaJ5HM910k6qr/eAD1SGtoqIVsQrRDN98A8OHw403qrK+FiVR0PCx0NCx\nI4QQQphERX4h0WVNnAJSWMgVORkMPbmBil/Jkh++xHcTqYY0NM+jBc0pevSA9HT46iuYMAFOnXJL\nxEI47dkDQ4bApEnw+uvQpqlHeH1/3w0dC0IIIYTJlZ5qxo31yhuPuy6zkjen5TfhhfFIIlVbQ3fb\nW9iconNn+PBDVV97551uWSdYCAAcDjUS9bvfwW9/28xy0vr+vmXkSQjh5TIzM4mIiCA0NJRly5bV\n+5wFCxYQEhLCgAEDKCgouOS2ZWVljBs3jqCgIMaPH095eTkADoeDqKgo+vfvz4gRI9i4caNzm/j4\neMLDw4mKiiIqKopjx4656R0LT3mkSwrf397Em4mVNx6X3Z1KdpH8m+lLJJFqqlaU/HXoAKtWwYAB\nqgX1oUNuiE/4tA0bICFBTXKdPr0FO3BRSasQQnjSnDlzSEpKYvPmzSxfvrxOAuNwONi6dSvbt29n\n3rx5zJs3r8FtS0tLAUhMTCQoKIiioiJ69+7NihUrALjxxhv5/PPP+eKLL0hKSmL27NnO+RV+fn6k\npKSQk5NDTk4O3bp189AnINzh9GnIKbYQ8EETbyZW3ngMH2QhJ8f98QnvIYlUU7Wy5K9NG/jzn1Xz\nieHD4cgR94csfENqKkybBv/+N9x11yWe3NDfq5TxCSEM5sSJEwDExcURHBzMyJEjsdvtNZLine truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x22f6a828>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 67
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [
|
||||
"fig, axes = subplots(1,3, figsize=(14,5))\n",
|
||||
"epsx = np.linalg.norm(Utils.mkvc(Bxr))*1e-6\n",
|
||||
"epsy = np.linalg.norm(Utils.mkvc(Byr))*1e-6\n",
|
||||
"epsz = np.linalg.norm(Utils.mkvc(Bzr))*1e-6\n",
|
||||
"axes[0].plot(Y[:,id], abs((Bxr[:,id]-Bxra[:,id])/(Bxra[:,id]+epsx)), 'r.')\n",
|
||||
"axes[1].plot(Y[:,id], abs((Byr[:,id]-Byra[:,id])/(Byra[:,id]+epsy)), 'r.')\n",
|
||||
"axes[2].plot(Y[:,id], abs((Bzr[:,id]-Bzra[:,id])/(Bzra[:,id]+epsz)), 'r.')"
|
||||
],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "pyout",
|
||||
"prompt_number": 68,
|
||||
"text": [
|
||||
"[<matplotlib.lines.Line2D at 0x1eeea898>]"
|
||||
]
|
||||
},
|
||||
{
|
||||
"metadata": {},
|
||||
"output_type": "display_data",
|
||||
"png": "iVBORw0KGgoAAAANSUhEUgAAAzcAAAE1CAYAAADahh74AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X90VPWd//HX5EAbI56EkCjsUQORaCZaQg6ZTFgxHVkM\nlDS45cfRuAWP6DLYdYMK3bOoPSXuqd2tspLNum7SPZ6qkPVwllbxB8VkPZOAmkmwVtr8aFYloLvF\nTeBLQJE9Ue/3j5QpQzLJTDJ3Zu6d5+McDpm5d+bzvnfu5z3zvj8+12EYhiEAAAAAsLiUeAcAAAAA\nANFAcQMAAADAFihuAAAAANgCxQ0AAAAAW6C4AQAAAGALFDcAAAAAbGHc4mb9+vW64oor9I1vfGPU\n6bt27VJhYaEKCwt1xx13qLe3N+pBArCO1tZWOZ1O5eXlqa6ubsT0l156SYWFhZo/f74qKirU0dER\n9msBJJfxckJPT48WLlyo1NRUbd++PWjaZ599pjvvvFPXXnutCgoK1NbWFquwAcSRY7z73Bw4cEDT\npk3TunXr9Jvf/GbE9LffflsFBQVKT0/Xs88+q+bmZj3//POmBQwgsRUVFam2tlY5OTlaunSpDh48\nqKysrMD0zz77TJdeeqkkqaWlRT/4wQ/U2toa1msBJJfxckJ/f7+OHj2qF198UdOnT9fmzZsD07Zs\n2aJLLrlEDz/8sKZMmaLPPvtM6enp8VgMADE07pGbm266SdOnTw85feHChYFkUVFRoZaWluhFB8BS\nBgcHJUllZWXKyclReXm5/H5/0DznC5vz86empob9WgDJI5yckJ2dreLiYk2dOnXE65ubm/XQQw8p\nNTVVU6ZMobABksSUaL5ZQ0ODKisrRzzvcDii2QyAKBnnwG3EOjo6lJ+fH3h8/lSQioqKoPl+8Ytf\n6IEHHtCnn36qd955J6LXkk+AxBPtXCKFnxNG8/HHH+vcuXO699571d3drZUrV2rTpk2BnSkSuQRI\nVJPNJ1EbUKC5uVk7d+7Uj370o1GnG4YR038//OEPbd0ebdqnvXi1GU/f+c531NfXp6eeekp//ud/\nHvHr7f7Z0KY92kuWNhPRuXPn1Nvbq1WrVsnn86mzs1O7d+8eMZ/dP5tkaTMZljFZ1ms0RKW4OXz4\nsDZu3Ki9e/cqIyMjGm8JwIJcLpd6enoCjzs7O1VaWhpy/ttuu03/8z//o88//1zFxcURvRaAvUWa\nTy40d+5cXXfddaqsrNQll1yiqqoq7du3z6xQASSQSRc3x44d06pVq7Rr1y7NnTs3GjEBsKjz57S3\ntraqr69PTU1NcrvdQfN88MEHgb0zr732mhYsWKBLLrkksGNkrNcCSB7h5JPzRtvjm5eXJ7/fr6++\n+kqvvvqqlixZYmq8ABLDuNfcVFVVqaWlRQMDA7rqqqtUU1OjoaEhSZLX69Wjjz6qkydPauPGjZKk\nqVOnqr293dyow+DxeGzdHm3ap714tWmWHTt2yOv1amhoSNXV1crKylJ9fb2k4ZyxZ88ePffcc5o6\ndaqKior0k5/8ZMzXxluybA/J0GYyLGO82jTLePnk+PHjcrlcOn36tFJSUlRbW6uuri5NmzZNTzzx\nhNatW6dz585pyZIluv322+O8NMmzPdC37dGmVXPJuENBR6URhyNq59EBiA6r9kurxg3YlVX7pFXj\nBuwsGv0yagMKAAAAAEA8UdwAAAAAsAWKGwAAAAC2QHEDAAAAwBYobgAAAADYAsUNAAAAAFuguAEA\nAABgCxQ3AAAAAGyB4gYAAACALVDcALCe5culU6fiHQWADRviHQEABKG4AWA9+/bxowpIBL298Y4A\nAIJQ3ACwHpdLamiIdxQA0tLiHQEABHEYhmGY3ojDoRg0AyACVu2XDodDxv/7f1JGRrxDAXDqlBzT\np1s3l1gwbsDOotEvKW6AJGXVfmnVuAG7smqftGrcgJ1Fo19yWhoAAAAAW6C4AQAAAGALFDcAAAAA\nbIHiBgAAJCfumQUkjijd4oHiBgAAJCfumQUkjijdN4viBgAAJCfumQUkjijdN4uhoIEkZdV+adW4\nAbuyap/knllAgonSfbMoboAkZdV+adW4Abuyap+0atyAnXGfGwAAAAD4A4obAAAAALZAcQMAAADA\nFihuAAAAANgCxQ0AAEhIra2tcjqdysvLU11d3YjpPT09WrhwoVJTU7V9+/YR07/88ksVFRWpsrIy\nFuECSABT4h0AAADAaDZt2qT6+nrl5ORo6dKlqqqqUlZWVmD6jBkzVFdXpxdffHHU19fW1qqgoEBn\nzpyJVcgA4ix2R26WL5dOnYpZcwDGwB25ASS4wcFBSVJZWZlycnJUXl4uv98fNE92draKi4s1derU\nEa//+OOP9dprr+mee+5hyGcgicTuyM2+fcM/qHbvjlmTAELo7Y13BAAwpo6ODuXn5wceFxQUqK2t\nTRUVFWG9/oEHHtDjjz+u06dPh5xn27Ztgb89Ho88Hs9EwwUwAT6fTz6fL6rvGbvixuWSGhpi1hyA\nMaSlxTsCADDNK6+8ossvv1xFRUVj/nC6sLgBEHsX71SoqamZ9HvG7rS011+XMjJi1hyAMTQ2xjsC\nABiTy+VST09P4HFnZ6dKS0vDeu1bb72lvXv3as6cOaqqqtIbb7yhdevWmRUqgAQSu+KGwgZIHPRH\nAAkuPT1d0vCIaX19fWpqapLb7R513ouvqXnsscf00Ucf6ciRI3rhhRe0ePFiPffcc6bHDCD+GC0N\nAAAkpB07dsjr9WpoaEjV1dXKyspSfX29JMnr9er48eNyuVw6ffq0UlJSVFtbq66uLk2bNi3ofRwO\nRzzCBxAHDiMGQ4g4HA5GKgESjFX7pVXjBuzKqn3SqnEDdhaNfslNPAEAAADYAsUNAAAAAFsYs7hZ\nv369rrjiCn3jG98IOc/WrVuVm5urBQsWBI1qAiA5tba2yul0Ki8vT3V1dSOm79q1S4WFhSosLNQd\nd9yh3gvuuTN79mzNmzdPRUVFKikpiWXYAADABsYsbu666y798pe/DDm9vb1dBw4c0KFDh7RlyxZt\n2bIl6gECsJZNmzapvr5ezc3NeuqppzQwMBA0PTc3V62trXrvvfe0dOlS/d3f/V1gmsPhkM/n07vv\nvqv29vZYhw4AACxuzOLmpptu0vTp00NO9/v9Wr16tTIzM1VVVaXu7u6oBwjAOgYHByVJZWVlysnJ\nUXl5ufx+f9A8CxcuDAzxWlFRoZaWlqDpXOALAAAmalJDQbe3t2vt2rWBx9nZ2frggw90zTXXjJj3\nwrsAX3w3UgDm8/l8Y96pOxo6OjqUn58feFxQUKC2tjZVVFSMOn9DQ4MqKysDjx0OhxYvXqw5c+Zo\n/fr1WrFixaivI58A8ROLXAIAEzWp4sYwjBF7WUONJX/hjxEAsXdxEVBTUxO/YCQ1Nzdr586deuut\ntwLPvfnmm5o1a5a6u7tVWVmpkpISzZw5c8RrySdA/CRaLgGAC01qtDS3262urq7A4/7+fuXm5k46\nKADW5HK5ggYW6ezsVGlp6Yj5Dh8+rI0bN2rv3r3KyMgIPD9r1ixJktPp1IoVK/Tyyy+bHzQAALCN\nSRc3e/bs0YkTJ9TY2Cin0xmtuABY0PlraVpbW9XX16empia53e6geY4dO6ZVq1Zp165dmjt3buD5\ns2fP6syZM5KGd5Ts379fy5Yti13wAADA8sY8La2qqkotLS0aGBjQVVddpZqaGg0NDUmSvF6vSkpK\ntGjRIhUXFyszM1M7d+6MSdAAEteOHTvk9Xo1NDSk6upqZWVlqb6+XtJw3nj00Ud18uRJbdy4UZI0\ndepUtbe36/jx41q5cqUkacaMGdq8ebOuuuqquC0HAACwHocRg6GJHA4HIyABCcaq/dKqcQN2ZdU+\nadW4ATuLRr+c1GlpAAAAAJAoKG4AAAAA2ALFDQAAAABboLgBAAAAYAsUNwAAAABsgeIGAAAAgC1Q\n3AAAAACwBYobAAAAALZAcQMAAADAFihuAAAAANgCxQ0AAAAAW6C4AQAAAGALFDcAAAAAbIHiBgAA\nAIAtUNwAAICE1NraKqfTqby8PNXV1Y2Y3tPTo4ULFyo1NVXbt28PPP/RRx/p5ptv1vXXXy+Px6PG\nxsZYhg0gjhyGYRimN+JwKAbNAIiAVfulVeMG7MrMPllUVKTa2lrl5ORo6dKlOnjwoLKysgLT+/v7\ndfToUb344ouaPn26Nm/eLEk6fvy4jh8/rvnz52tgYEAlJSV67733dNlll8UkbgATE41+yZEbAACQ\ncAYHByVJZWVlysnJUXl5ufx+f9A82dnZKi4u1tSpU4OenzlzpubPny9JysrK0vXXX69Dhw7FJnAA\ncUVxAwAAEk5HR4fy8/MDjwsKCtTW1hbx+7z//vvq7OxUSUlJNMMDkKCmxDsAAAAAM5w5c0a33Xab\nnnzySV166aUjpm/bti3wt8fjkcfjiV1wAOTz+eTz+aL6nlxzAyQpq/ZLq8YN2JVZfXJwcFAej0fv\nvvuuJOmv//qvtWzZMlVUVIyYt6amRtOmTQtccyNJQ0NDqqio0PLly3X//ffHLG4AE8c1NwAAwJbS\n09MlDY+Y1tfXp6amJrnd7lHnvfjHkGEYuvvuu3XDDTeMWtgAsC+O3ABJyqr90qpxA3ZlZp9saWnR\nxo0bNTQ0pOrqalVXV6u+vl6S5PV6dfz4cblcLp0+fVopKSm67LLL1NXVpV//+tcqKyvTvHnz5HA4\nJEk//vGPtWzZspjEDWBiotEvKW6AJGXVfmnVuAG7smqftGrcgJ1Z+7S0DRskj0davlw6dSpuYQC2\nRj8DAABJJH7FTW+v1NIi7ds3/AMMQPTRzwAAQBKJX3GTljb8v8slNTTELQzA1uhnAAAgicLine truncated
|
||||
"text": [
|
||||
"<matplotlib.figure.Figure at 0x111e0390>"
|
||||
]
|
||||
}
|
||||
],
|
||||
"prompt_number": 68
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Thoughts\n",
|
||||
"\n",
|
||||
"- It works well with non-uniform mesh!!\n",
|
||||
"- Actual accuray is ~10% relative error in secondary fields. As we pad more we can get better accuracy, since we did not consider secondary field at the boudnary ($\\partial\\Omega$). \n",
|
||||
"- Here, we can try primary secondary field approach to get better accuracy. \n",
|
||||
"- In addition, we can use the congruous sphere method to handle this secondary fields at boundaries"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"collapsed": false,
|
||||
"input": [],
|
||||
"language": "python",
|
||||
"metadata": {},
|
||||
"outputs": []
|
||||
}
|
||||
],
|
||||
"metadata": {}
|
||||
}
|
||||
]
|
||||
}
|
||||
Reference in new issue
Block a user