利用matlab的PDE工具箱求解Neumann边界的Poisson方程之3
针对于 利用matlab的PDE工具箱求解Neumann边界的Poisson方程之2 出现的问题,后面察觉是由于 对法向导数理解的误解。
首先明白:Poisson方程的第二边值问题:
n是外法线方向向量.
例如:对于下图所示
下边界时条件为 ![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAM4AAABhCAIAAAARNwmoAAAD0ElEQVR4nO2cy5bjIAwF/f8/7Vn0JhNAlgQSKK5a9Ykf3HDLTkJO57oBUrh2B4C3gGqQBKpBEqgGSaAaJIFqkASqQRKoBkmgGiSBapCEWbXrfyIyQRwb67MNdg0ICgdr2VufYZhRUGwrwfb6tGN0M6FaFU6oTzWGkAbbzueQ+hYMgGqlqaTaLV43cD6oBkn8jmqjNwS8zzsEx5s5R3ez62pC1nZnVNuLvr6I7syLHZ+nFqJoMqFaGr761nZnWOzwDam5LFAtFHd9a7t73k8vk3ysEFeZFRy461ve3cOumpNqnob8TPRxwcRMfcu726OacMXAQiJUc3cXqJqQCc9ycNcX0Z20t/KkaZcFmJipL6K7WdWErd3DrwZTXNAzU19Ed6j2s/yUag7PNKeFJbjrC+rOr5rmcvnc5/OPdhPaLcddX1B32k+gyse7cdvr4/NBPIvDV19Qd57FjkfP2rjC08CzOHz1BXVn/g5UOYCw8xLJcFSPtb6g7qq2hWrlqNoWqpWjaluoVo6qbaFaOaq2hWrlqNoWqpWjaluoVo6qbaFaOaq2hWrlqNpWjmrtOruS6GAVqTopqFaOqpNCo+V4/veTXTzkVuxzOKEzs5d+sDOz/vE4ofpeDyRiWmZOu5ZOsJhpDGf0fOBYqraVo9rCaxqqTgqqlaPqpHw1Sv3nU7WArj0jq1DtBKoWgGrlqFqASTU4gaqtoFo5qraCauWo2gqqlcP+j6NnrCCgmptdDRp/t2hAULjHJKN4+XmqsLFBy28uiMRFFMKMEiaHqcLeBtU/I9MLVFG1jZn3sr1Bw8/DCJvym5u8q71NtRMaXHD2c24SqOajjGr3Ac1dPR73T4t3PqgWRbnA0fyCasKdZuMrL6p9MZqQ7uO+2ZtdVxuN3X1F+9pf+WIXwWtVUzaouUFYJ9C82KEMJBwoKJjGC1XTNyi3E6uacHafapOh53mbao4GNYfYMsykfBy4vTi6nmW2fjWkDb0Ld4OPttliTKaURx15du9T7W3MNCirZk4SF/QWfcKzHN6uGre0NCJU81W24JtpTVDhloZqcSxs8OtxT5jJoMJW4XA8yyGiQXdlG1TjlpbG8gZn+ppSTZlSUM2dGzTMNNg9fI9q1pTCJrQLYqbB9vDJjrSfQJWPj1IKm/AsFHeD9+q3Op7FDpNnsmp4Fo2vwfbYcNW6Q5oiWrdCBL7Le21TlA1D1t4RUA2GrH3lQTUYgmoQQivW2jfTqAb3vXq1tj/E2tNBUXwfUW1DLD8jFCV6BQrVIAlUgyRQDZJANUgC1SAJVIMkUA2S+AdR+70DservOgAAAABJRU5ErkJggg==)
对于左边边界为 ![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFcAAAAzCAIAAABtzwN3AAABjUlEQVRoge2Y2xLDIAhE/f+f3j6E6aQCDiJekrJPTQLremqd1IIUUHYHOEJJAUgKl2IolPJsmkkBSAqXhPTlVyaXMyg4klOjaFSZWoa3DzlJvuRUb7dulxnHW6akAIxTKLq0AQ6h4EhO90WjylS8X3XFzschX3K6z43udUac2ym4k9Mj7sUvn0KBX06hwJF/KwOmMiB3cqrXvNpGld1RFLqSU4tox00PXwvwJqfKqAQhPruUFICkcCkpABspFF0hkfrCxLhMi74GlvzKade3NzBTr0aSk0NUjsDoIZH6wsS45O6IpEAuSQGvpODYq86h4Nto33MSj4HD+NeeQaPnGDopAI6TeHGlnUOh/Samcek7ib9/rrrC5+OQGJ5/i0Ijd7nXaWuBDx87H4e08DD8NF5yEg89PH8k9FqM/pdCG0GR9oXCNDixLvGc96hacqoXvTSQ2vTa7JapCm9JTk+jhr8PFuK5UhMpPIjI3LXwXxTIS9pTAv3naW7KpND333avHhBxgZICkBQuJQUA+AB3VqTi/njz7wAAAABJRU5ErkJggg==)
在左边边界设定的时候,前面添加一负号。
原方程在PDE工具箱中的设置
上边界与下边界 sin(3*pi*x+pi/4)*sin(pi/4)
(注意:此处不能为'x',后同)
左边界 -3*pi*sin(2*pi*y+pi/4)*cos(pi/4)
右边界 -3*pi*sin(2*pi*y+pi/4)*cos(pi/4)
方程设置:13*pi*pi.*sin(3*pi*x+pi/4).*sin(2*pi*y+pi/4)
最后的数值结果为:
X-Y平面图为:
与精切解吻合比较好。
如果为一导数边界,而其他边界为第二类,具体代码实现如下:
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
function pdemodel [pde_fig,ax]=pdeinit; pdetool('appl_cb',1); set(ax,'DataAspectRatio',[1 2.583333333333333 1.6666666666666665]); set(ax,'PlotBoxAspectRatio',[1 1 1]); set(ax,'XLim',[-0.10000000000000001 1.1000000000000001]); set(ax,'YLim',[-0.10000000000000001 3]); set(ax,'XTickMode','auto'); set(ax,'YTickMode','auto'); pdetool('gridon','on'); % Geometry description: pderect([0 1 1 0],'R1'); set(findobj(get(pde_fig,'Children'),'Tag','PDEEval'),'String','R1') % Boundary conditions: pdetool('changemode',0) pdesetbd(4,... 'dir',... 1,... '1',... 'sin(3*pi*x+pi/4)*sin(pi/4)') pdesetbd(3,... 'neu',... 1,... '0',... '-2*pi.*sin(3*pi*x+pi/4)*sin(pi/4)') pdesetbd(2,... 'neu',... 1,... '0',... '-3*pi.*sin(2*pi.*y+pi/4)*cos(pi/4)') pdesetbd(1,... 'neu',... 1,... '0',... '2*pi.*sin(3*pi*x+pi/4)*sin(pi/4)') % Mesh generation: setappdata(pde_fig,'Hgrad',1.3); setappdata(pde_fig,'refinemethod','regular'); pdetool('initmesh') % PDE coefficients: pdeseteq(1,... '1.0',... '0.00',... '13*pi*pi.*sin(3*pi*x+pi/4).*sin(2*pi*y+pi/4)',... '1.0',... '0:10',... '0.0',... '0.0',... '[0 100]') setappdata(pde_fig,'currparam',... ['1.0 ';... '0.00 ';... '13*pi*pi.*sin(3*pi*x+pi/4).*sin(2*pi*y+pi/4)';... '1.0 ']) % Solve parameters: setappdata(pde_fig,'solveparam',... str2mat('0','1000','10','pdeadworst',... '0.5','longest','0','1E-4','','fixed','Inf')) % Plotflags and user data strings: setappdata(pde_fig,'plotflags',[1 1 1 1 1 1 1 1 0 0 0 1 1 1 1 0 0 1]); setappdata(pde_fig,'colstring',''); setappdata(pde_fig,'arrowstring',''); setappdata(pde_fig,'deformstring',''); setappdata(pde_fig,'heightstring','');