{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Strong (Dirichlet) boundary conditions" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "Let's move on from the Helmholtz problem to Poisson:\n", "$$-\\nabla^2u = f$$\n", "for some known function $f$. To have a well-posed problem, we must impose Dirichlet conditions over at least part of the domain boundary:\n", "$$u(x) = g(x) \\quad \\forall x \\in \\Gamma_D,\\\\\n", "\\nabla u(x)\\cdot \\vec{n} = h(x) \\quad \\forall x \\in \\Gamma_N.$$\n", "\n", "As before, the Neumann condition is imposed weakly by setting the boundary integral over the relevant part of the boundary. The Dirichlet condition is imposed strongly by modifying the function space from which $u$ is drawn.\n", "\n", "Suppose $V\\ni u$ is a finite element space. Then we write:\n", "$$V_D = \\{u\\in V \\,|\\,u(x) = g \\quad \\forall x \\in \\Gamma_D\\}$$\n", "$$V_0 = \\{u\\in V \\,|\\,u(x) = 0 \\quad \\forall x \\in \\Gamma\\setminus\\Gamma_D\\}$$\n", "\n", "The weak form of the problem then becomes, find $u\\in V_D$ such that:\n", "$$\\int_\\Omega \\nabla u \\cdot \\nabla v\\, \\mathrm{d}x - \\int_{\\Gamma_N} hv\\, \\mathrm{d}s = \\int_\\Omega fv\\mathrm{d}x \\quad\\forall v\\in V_0.$$\n", "\n", "In Firedrake, the weak form is always formulated in the unmodified space $V$. The modification to the function spaces implied by the Dirichlet conditions is then imposed by supplying a [`DirichletBC` object](http://firedrakeproject.org/firedrake.html#firedrake.bcs.DirichletBC) to the solve." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As usual, we start by setting up the notebook matplotlib interaction, and importing Firedrake. We then choose a triangulation of our domain and select an appropriate $H^1$-conforming function space." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": true }, "outputs": [], "source": [ "%matplotlib notebook\n", "from firedrake import *\n", "mesh = UnitSquareMesh(10, 10)\n", "V = FunctionSpace(mesh, \"CG\", 1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let us specify the details of our problem. We will solve\n", "\n", "$$ -\\nabla^2 u = f \\text{ on } \\Omega = [0, 1]\\times[0, 1] \\\\\n", " u = 0 \\text{ on } \\Gamma$$.\n", " \n", "We'll choose an $f$ such that the exact solution $u^* = \\sin(\\pi x)\\sin(\\pi y)$, which satisfies the zero Dirichlet conditions. With a little algebra, we determine that we need\n", "\n", "$$ f = 2\\pi^2\\sin(\\pi x) \\sin(\\pi y) $$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "x, y = SpatialCoordinate(mesh)\n", "u_exact = sin(pi*x)*sin(pi*y)\n", "f = 2*pi**2*u_exact" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As before, we'll require test and trial functions to define our problem" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": true }, "outputs": [], "source": [ "u = TrialFunction(V)\n", "v = TestFunction(V)\n", "a = dot(grad(u), grad(v))*dx\n", "L = f*v*dx" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we need to tell Firedrake how to modify the function space to take into account the Dirichlet condition. This is done by constructing a [DirichletBC object](http://firedrakeproject.org/firedrake.html#firedrake.bcs.DirichletBC). As usual, we can ask for help using Python's builtin documentation facilities" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "?DirichletBC" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We want to apply the zero Dirichlet condition to all exterior nodes of the space, so we need to provide the list of boundary ids that corresponds to the exterior boundary of the domain. For the `UnitSquareMesh`, this is the set `(1, 2, 3, 4)`." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": true }, "outputs": [], "source": [ "boundary_ids = (1, 2, 3, 4)\n", "bcs = DirichletBC(V, 0, boundary_ids)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We supply the boundary condition as an argument to `solve`." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "uh = Function(V)\n", "solve(a == L, uh, bcs=bcs, solver_parameters={\"ksp_type\": \"cg\"})" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support.' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " if (mpl.ratio != 1) {\n", " fig.send_message(\"set_dpi_ratio\", {'dpi_ratio': mpl.ratio});\n", " }\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " this.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var backingStore = this.context.backingStorePixelRatio ||\n", "\tthis.context.webkitBackingStorePixelRatio ||\n", "\tthis.context.mozBackingStorePixelRatio ||\n", "\tthis.context.msBackingStorePixelRatio ||\n", "\tthis.context.oBackingStorePixelRatio ||\n", "\tthis.context.backingStorePixelRatio || 1;\n", "\n", " mpl.ratio = (window.devicePixelRatio || 1) / backingStore;\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width * mpl.ratio);\n", " canvas.attr('height', height * mpl.ratio);\n", " canvas.attr('style', 'width: ' + width + 'px; height: ' + height + 'px;');\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
')\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('