偏微分方程讲义.ppt

《偏微分方程讲义.ppt》由会员分享，可在线阅读，更多相关《偏微分方程讲义.ppt（76页珍藏版）》请在三一办公上搜索。

1、偏微分方程讲义建模、数值解和Matlab工具箱,温罗生2013.7,提 纲,几个偏微分方程的典型问题建模偏微分方程的基本概念偏微分方程的建模偏微分方程的数值方法偏微分方程求解的数值方法偏微分方程的Matlab求解工具Pdetool的使用方法练习题,1.1 偏微分方程的基本概念,含有未知多元函数及其偏导数的方程，称之为偏微分方程。方程中出现的未知函数偏导数的最高阶数称为偏微分方程的阶。如果方程中对于未知函数和它的所有偏导数都是线性的，这样的方程称为线性偏微分方程，否则称它为非线性偏微分方程。一些常见的偏微分方程：,拉普拉斯方程：uxx+uyy+uzz=0适用于重力场或静电场。,波动方程式：未知函

2、数 u(x,y,z,t):,热传导方程式：其中 k 代表该材料的热导率。,泊松方程：适用于所有物质或电荷的重力场或静电场。,初始条件和边界条件称为定解条件，未附加定解条件的偏微分方程称为泛定方程。对于一个具体的问题，定解条件与泛定方程总是同时提出。定解条件与泛定方程作为一个整体，称为定解问题。,1.2偏微分方程的建模,例子1，弦的振动方程,一根长为L的均匀细线，线密度为(x)，放于x轴上拉紧后，让它离开平衡位置做微小振动。求任意时刻t弦线的位移u(x,t)。使用元素法的思想，在弦上取x,x+x一段进行分析。因为弦的质量相对于张力非常小，因此可以看成弦仅仅受到两个张力，T(x)和T(x+x)，利

3、用力学基本定律，可以得到水平和竖直方向两个方程,这里，a分别是两个力和水平方向的夹角，以及弦线在竖直方向的加速度。,注意到弦仅仅在接近水平位置振动，所以和都是很小的量，于是前一个方程可以近似为,这代表张力和位置无关。直观的理解这是显然的。同时因为和都是很小，也有如下的近似关系,于是竖直方向的方程可以改写为,注意到弧长ds近似等于dx，因此上述方程又可以写成,当弦受竖直方向的外力F(x,t)作用时，竖直方向的方程可以写成,因此方程可以写成,仅仅通过上面的方程并不能完全确定弦的状态，我们必须先给出在初始时刻t=0时的状态，也就是所谓的初始时刻，即,也就是给出初始位置和初始速度。根据方程和上面的初始

4、条件，仍然不能决定，还需要给出边界条件。更加具体的情况，一般有三类边界条件，第一类是两个端点固定，也就是,第二类边界条件描述端点受到外力作用的情况，有,第三类边界条件描述端点受到弹性支撑的情况，结合虎克定理，能得到条件,例子2，热传导方程,热传导定理的推导依据两个定理，一个是傅里叶定理和能量守恒定理。傅里叶定理说热量从高温区域向低温区域传播，热流强度满足,所谓热流强度是单位时间通过单位横截面积的热量。上面dt是时间段，k是热传导系数。是温度场梯度方向面积的投影，u(x,y,z,t)是t时刻在空间点(x,y,z)的温度。取点(x,y,z)的一个小邻域，从时间t1到t2流入该区域的热量为,因为,该

5、区域的问题升高需要的热量满足,再根据能量守恒定理，传入的热量等于温度升高所需要的热量，于是有,后一方程进一步使用高斯公式,所以有,由t和V的任意性知道，被积函数恒为0，特别的，如果k是一个常量时，得到方程,当物体有内部热源的时候，方程为,因为,左边是升高温度需要的热量，右边第一项是传入的热量，第二项是内部热源的效果。,类似与弦振动方程，要确定温度，也需要初始条件和边界条件。,初始条件,第二边界条件，表示表面各点热流已知。,第一边界条件，表示表面各点温度已知。,第三边界条件，表示外界温度为u1,表面的热量和温度差成正比。,2.1 一些常见的偏微分方程,Poisson 方程,带有稳定热源或内部无热

6、源的稳定温度场的温度分布，不可压缩流体的稳定无旋流动及静电场的电势等均满足这类方程。下面的方程是Poisson 方程的第一边值问题。,第二第三边界条件,在研究热传导过程，气体扩散现象及电磁场的传播等随时间变化的非定常物理问题时，常常会遇到抛物型方程。下面是初边值问题,Cauchy问题,常见的一维振动与波动问题可用该方程表示。下面是初边值问题,双曲型方程,2.2 偏微分方程数值解,差分方法又称为有限差分方法或网格法，是求偏微分方程定解问题的数值解中应用最广泛的方法之一。它的基本思想是：先对求解区域作网格剖分，将自变量的连续变化区域用有限离散点（网格点）集代替；将问题中出现的连续变量的函数用定义在

7、网格点上离散变量的函数代替；通过用网格点上函数的差商代替导数，将含连续变量的偏微分方程定解问题化成只含有限个未知数的代数方程组（称为差分格式）。如果差分格式有解，且当网格无限变小时其解收敛于原微分方程定解问题的解，则差分格式的解就作为原问题的近似解（数值解）。因此，用差分方法求偏微分方程定解问题一般需要解决以下问题：（i）选取网格；（ii）对微分方程及定解条件选择差分近似，列出差分格式；（iii）求解差分格式；,椭圆型方程第一边值问题的差分解法,首先将区域通过一系列横纵平行线划分，考察方程在交点处的值。注意到在边界的交点值是已知的，而在内部的交点的值是要求的。为计算u在内部交点处的值，需要列出

8、方程。使用二阶中心差商代替方程中的二阶导数，可以得到方程组,这里h和是横纵方向上的步长，右端是已知函数f在交点kj处的值，左边分子是待求的未知量。,一些内点的四个邻居都属于集合，这种点称为正则点，另一些内点的部分邻居不属于集合，这种点属于非正则点。对于每个正则点，可以得到一个方程，但是非正则点不能得到方程，这样最后得到的方程个数少于未知量个数。因此，对于非正则点要使用一些方法得到方程，常用的方法是直接转移和线性插值的方法。前者就是直接用内临近的点的值代替，后者是用内附近点的线性组合代替。下面以例子讲解具体的求解过程。,例3，用五点菱形格式求解 Laplace 方程第一边值问题,首先画出网格，如

9、右图。可以看到，内点一共4个，都是正则点，其他的是边界点，其值已知，使用五点菱形格式，得到方程组如下：,代人边界条件的值，得到方程组,借助Matlab，可以轻松得到问题的解。下面是程序。clc,clear f1=(x)2*log(1+x);f2=(x)log(1+x).2+1);f3=(y)log(1+y.2);f4=(y)log(4+y.2);u=zeros(4);m=4;n=4;h=1/3;u(1,1:m)=feval(f3,0:h:(m-1)*h);u(n,1:m)=feval(f4,0:h:(m-1)*h);u(1:n,1)=feval(f1,0:h:(n-1)*h);u(1:n,m)

10、=feval(f2,0:h:(n-1)*h);b=-u(2,1)+u(1,2);u(4,2)+u(3,1);u(2,4)+u(1,3);u(3,4)+u(4,3);a=-4 1 1 0;1-4 0 1;1 0-4 1;0 1 1-4;x=ab,当然，这里网格非常稀疏，得到结果比较粗糙，为得到更精确的值，可以将网格加密，比如两个方向各取1000个点。试自己修改上面的程序实现之。,function sol=main f1=(x)2*log(1+x);f2=(x)log(1+x).2+1);f3=(y)log(1+y.2);f4=(y)log(4+y.2);a=1;b=1;h=1/3;tol=0.0

11、0001;max1=1000;sol=dirich(f1,f2,f3,f4,a,b,h,tol,max1)%*function U=dirich(f1,f2,f3,f4,a,b,h,tol,max1)%Input-f1,f2,f3,f4 are boundary functions input as strings%-a and b right endpoints of 0,a and 0,b%-h step size%-tol is the tolerance%Output-U solution matrix;%If f1,f2,f3 and f4 are M-file functions

12、call U=dirich(f1,f2,f3,f4,a,b,h,tol,max1).%if f1,f2,f3 and f4 are anonymous functions call U=dirich(f1,f2,f3,f4,a,b,h,tol,max1).,%Initialize parameters and U n=fix(a/h)+1;m=fix(b/h)+1;ave=(a*(feval(f1,0)+feval(f2,0).+b*(feval(f3,0)+feval(f4,0)/(2*a+2*b);U=ave*ones(n,m);%Boundary conditions U(1,1:m)=

13、feval(f3,0:h:(m-1)*h);U(n,1:m)=feval(f4,0:h:(m-1)*h);U(1:n,1)=feval(f1,0:h:(n-1)*h);U(1:n,m)=feval(f2,0:h:(n-1)*h);%逐次超松弛迭代法（SOR）参数 w=4/(2+sqrt(4-(cos(pi/(n-1)+cos(pi/(m-1)2);%Refine approximations and sweep operator throughout the grid-366-err=1;cnt=0;,while(errtol),上面的程序是使用SOR方法计算线性方程组。因为偏微分方程离散化之

14、后得到的线性方程组比较特殊，同时规模非常大，所以往往选择比较有效的线性方程组数值算法。当问题规模较小，初学时可以不考虑这些。,抛物型方程的差分解法,使用差分进行近似，得到如下的方程,这里h和是横纵方向上的步长。,对初始条件及第一类边界条件，可直接得到,结合上页第一个差分格式，得到,计算的时候，由于第1层的数据已知，利用这些数据可以得到第二层的数据，再利用第二层的数据计算第三层，依次类推就能计算出全部的未知量。这个格式称为古典显格式。,结合第二个差分格式，得到,这是一个三对角方程组，可以使用追赶法求解。这个格式称为古典隐格式。,必须通过解下列线性方程组,例4，用古典显示格式求初边值问题,的数值解

15、，取h=1,=0.5.,双曲型方程的差分解法,可以使用下面的差分进行近似，,上面的格式都是显格式。当然，也可以采用不同的格式组合得到其他的格式。这里不再讲解。,2.3偏微分方程的Matlab求解工具,PDEPE 命令的使用方法及例子,PDEPE（initial-boundary value problems for parabolic-elliptic PDEs in 1-D）是求解抛物线和椭圆型一维偏微分方程的命令。为求解这些类型的各种各样的问题，在求解方程的形式，初始条件以及边界条件上，PDEPE使用下面的统一结构。方程形式,函数f用来描述“通量”，函数s用来描述“源”。,初始条件为,边界

16、条件为（可以通过取不同的函数，对应全部的三类边界条件。）,Matlab的命令格式为pdepe。M是参数，可以取0，1，2.pdefun是上面的函数f，icfun描述初始条件，bcfun描述边界条件，xmesh和tspan是x和t方向上的划分，options是选择项。下面的两个例子是matlab帮助文件中自带的例子，这里做一个详细的讲解。,例子5，求偏微分方程,其中x在0,1之间，且满足以下初边值条件,步骤1：将欲求解的偏微分方程改写成如下的标准形式。,步骤 2:编写偏微分方程的系数向量函数。function c,f,s=ex20_1pdefun(x,t,u,dudx)c=pi2;f=dudx;

17、s=0;步骤 3:编写起始值条件。function u0=ex20_1ic(x)u0=sin(pi*x);步骤 4:编写边界条件。在编写之前，先将边界条件改写成标准形式，找出相对应的函数pa()和pb()，然后写出MATLAB 的边界条件函数，例如，原边界条件可写成,因而，边界条件函数可编写成 function pa,qa,pb,qb=ex20_1bc(xa,ua,xb,ub,t)pa=ua;qa=0;pb=pi*exp(-t);qb=1;步骤 5:取点。例如 x=linspace(0,1,20);%x 取 20 点 t=linspace(0,2,5);%时间取5 点输出 步骤 6:利用 pd

18、epe求解。m=0;%依步骤1 之结果sol=pdepe(m,ex20_1pdefun,ex20_1ic,ex20_1bc,x,t);%这里 sol实际上是二维矩阵 步骤 7 显示结果。surf(x,t,sol)title(pde 数值解),xlabel(位置),ylabel(时间),zlabel(u),例6,试解以下联立的偏微分方程系统,步骤1：改写偏微分方程为标准形式,其中m=0和,初始条件:在x=0处,在x=1处,步骤 2：编写偏微分方程的系数向量函数。function c,f,s=ex20_2pdefun(x,t,u,dudx)c=1 1;f=0.024 0.170.*dudx;y=u

19、(1)-u(2);F=exp(5.73*y)-exp(-11.47*y);s=-F F;步骤 3：编写初始条件函数 function u0=ex20_2ic(x)u0=1 0;步骤 4：编写边界条件函数 function pa,qa,pb,qb=ex20_2bc(xa,ua,xb,ub,t)pa=0 ua(2);qa=1 0;pb=ub(1)-1 0;qb=0 1;步骤 5：取点,m=0;x=0 0.005 0.01 0.05 0.1 0.2 0.5 0.7 0.9 0.95 0.99 0.995 1;t=0 0.005 0.01 0.05 0.1 0.5 1 1.5 2;%*%利用pdepe

20、求解%*-382-sol=pdepe(m,ex20_2pdefun,ex20_2ic,ex20_2bc,x,t);u1=sol(:,:,1);%第一个状态之数值解输出 u2=sol(:,:,2);%第二个状态之数值解输出%*%绘图输出%*figure(1),surf(x,t,u1)title(u1 之数值解),xlabel(x),ylabel(t)%figure(2),surf(x,t,u2)title(u2 之数值解),xlabel(x),ylabel(t)%*,%pdefun 函数%*function c,f,s=ex20_2pdefun(x,t,u,dudx)c=1 1;f=0.024

21、0.170.*dudx;y=u(1)-u(2);F=exp(5.73*y)-exp(-11.47*y);s=-F F;%*%初始条件函数%*function u0=ex20_2ic(x)u0=1 0;%*%边界条件函数%*function pa,qa,pb,qb=ex20_2bc(xa,ua,xb,ub,t)pa=0 ua(2);qa=1 0;pb=ub(1)-1 0;qb=0 1;,例 7,求下列偏微分方程组的数值解,将原方程改写成标准式,参数取值：m=1(圆柱)，,边界条件：在r=0处和r=2处,function sol=main20_6 m=1;%*%取点%*r=linspace(0,2

22、,20);L=linspace(0,5,40);%*%利用 pdepe 求解%*sol=pdepe(m,ex20_6pdefun,ex20_6ic,ex20_6bc,r,L);T=sol(:,:,1);%温度 f=sol(:,:,2);%反应率%*%绘图输出%*figure(1),surf(L,r,T)title(temp),xlabel(L),ylabel(r),zlabel(temp(0C)%figure(2),surf(L,r,f)title(reaction rate),xlabel(L),ylabel(r),zlabel(reaction rate),%*%偏微分方程函数%*func

23、tion c,f,s=ex20_6pdefun(r,L,u,DuDr)T=u(1);f=u(2);c=1 1;f=0.59 0.94.*DuDr;s=9.611;0.047/T;%*%初始条件函数%*function u0=ex20_6ic(x)u0=125 0;%*%边界条件函数%*function pa,qa,pb,qb=ex20_6bc(ra,ua,rb,ub,L)pa=0 0;qa=1 1;pb=112*(ub(1)-273)0;qb=0.65/0.59 1;,二维状态空间的偏微分方程的 MATLAB 解法,MATLAB 中的偏微分方程（PDE）工具箱是用有限元法寻求典型偏微分方程的数

24、值近似解，该工具箱求解偏微分方程具体步骤与用有限元方法求解偏微分方程的过程是一致的，包括几个步骤，即几何描述、边界条件描述、偏微分方程类型选择、有限元划分计算网格、初始化条件输入，最后给出偏微分方程的数值解（包括画图）。下面我们讨论的方程是定义在平面上的有界区域上，区域的边界记作。Matlab有下面的一些函数能够求解二维偏微分方程问题，其中最基本的是assempde命令，也能使用pdetool。,MATLAB 工具箱可求解问题的基本类型,（i）椭圆型偏微分方程（ii）抛物型偏微分方程（iii）双曲型偏微分方程（iv）特征值问题（v）非线性椭圆偏微分方程（vi）偏微分方程组,边界条件有如下三种：

25、,例8,求解泊松方程求解区域为单位圆盘，边界条件为在圆盘边界上u=0。,%(1)问题定义 g=circleg;%单位圆 b=circleb1;%边界上为零条件 c=1;a=0;f=1;%（2）产生初始的三角形网格p,e,t=initmesh(g);%（3）迭代直至得到误差允许范围内的合格解 error=;err=1;while err 0.01,p,e,t=refinemesh(g,p,e,t);u=assempde(b,p,e,t,c,a,f);%求得数值解 exact=(1-p(1,:).2-p(2,:).2)/4;err=norm(u-exact,inf);error=error err

26、;end%结果显示 subplot(2,2,1),pdemesh(p,e,t);subplot(2,2,2),pdesurf(p,t,u)subplot(2,2,3),pdesurf(p,t,u-exact),例9 考虑最小表面问题,在圆盘边界上,x=u 2。这是椭圆型方程，其中,g=circleg;b=circleb2;c=1./sqrt(1+ux.2+uy.2);rtol=1e-3;p,e,t=initmesh(g);p,e,t=refinemesh(g,p,e,t);u=pdenonlin(b,p,e,t,c,0,0,Tol,rtol);pdesurf(p,t,u),例10 求解正方形区

27、域(x,y):-1 x,y 1上的热传导方程,这里是抛物型方程，其中 d=1,f=0,a=0,c=1。,%(1)问题定义 g=squareg;%定义正方形区域 b=squareb1;%边界上为零条件 c=1;a=0;f=0;d=1;%（2）产生初始的三角形网格 p,e,t=initmesh(g);%(3)定义初始条件 u0=zeros(size(p,2),1);ix=find(sqrt(p(1,:).2+p(2,:).2)0.4);,u0(ix)=1%(4)在时间段为0到0.1的20个点上求解 nframe=20;tlist=linspace(0,0.1,nframe);u1=paraboli

28、c(u0,tlist,b,p,e,t,c,a,f,d);%(5)动画图示结果 for j=1:nframe pdesurf(p,t,u1(:,j);mv(j)=getframe;end movie(mv,10),例11,求解正方形区域(x,y):-1 x,y 1上的波方程,这里是双曲型方程，其中d=1,f=0,a=0,c=1。,%(1)问题定义 g=squareg;%定义正方形区域 b=squareb3;%定义边界 c=1;a=0;f=0;d=1;%（2）产生初始的三角形网格 p,e,t=initmesh(g);%(3)定义初始条件 x=p(1,:);y=p(2,:);u0=atan(cos(

29、pi*x);ut0=3*sin(pi*x).*exp(cos(pi*y);,%(4)在时间段为0到5的31个点上求解 n=31;tlist=linspace(0,5,n);uu=hyperbolic(u0,ut0,tlist,b,p,e,t,c,a,f,d);%(5)动画图示结果 for j=1:n pdesurf(p,t,uu(:,j);mv(j)=getframe;end movie(mv,10),2.4 pdetool的使用方法,图形界面解法步骤大致上为：（1）定义PDE 问题，包括二维空间范围，边界条件以及PDE 系数等。（2）产生离散化点，并将原PDE 方程离散化。（3）利用有限元素

30、法(finite element method；FEM)求解并显示答案。,前五个按钮为 PDE 系统的边界范围绘制功能，由左至右的用法为：1：以对角绘制矩形或正方形。按住鼠标左键可绘制矩形，而正方形需以按住右键的方式绘制。2：从中心点至某一角边的方式绘制矩形或正方形。同样地，鼠标左键绘矩形，右键绘正方形。3：由周围界线的方式绘制椭圆或圆形区域。鼠标左键用以绘制椭圆，而右键用来绘制圆形图形。4：以中心点向外的方式绘制椭圆或圆。同样地，鼠标左键及右键，分别用以绘制椭圆及圆形的区域。5：用以绘制多边型等不规则区域，欲关闭此功能需按鼠标右键。,在这些绘制按钮之后的7个按钮功能依序如下：6：用以给定边界

31、条件。在此功能选定后，使用者可在任一图形边界上按住鼠标左键双击，然后在对话框中输入边界条件。7：用以指定 PDE 问题及相关参数。8：产生图形区域内离散化的网点。9：用以进一步将离散化的网点再取密一点(refine mesh)。10：在指定 PDE 系统，边界条件及区域后，按此钮即开始解题。11：用以指定显示结果绘制方式。12：放大缩小功能，便于图形绘制及显示。,图形界面解法的使用步骤要利用pdetool 接口求解之前，需先定义PDE 问题，其包含三大部份：（1）利用绘图(draw)模式，定义欲解问题的空间范围(domain)。（2）利用boundary 模式，指定边界条件。（3）利用PDE

32、模式，指定PDE 系数，即输入c，a，f 和d 等PDE 模式中的系数。在定义 PDE 问题之后，可依以下两个步骤求解（1）在mesh 模式下，产生mesh 点，以便将原问题离散化。（2）在solve 模式下，求解。（3）最后，在Plot 模式下，显示答案。,以Poisson s方程式 u=10的求解为例，详细说明 pdetool 的用法。此问题的几何图形及相关边界条件，将于求解过程中加以说明。步骤 1：在命令窗口中键入pdetool 以进入GUI(graphical user interface)界面。选取Options 中之Grid 功能，以显示网格线。步骤 2：利用Draw 功能，画出问

33、题之几何图形。请注意：使用者可利用内定对象“多边型“，“矩形”，“正方形”，“圆形”，及“椭圆型”，予以组合。步骤 3：选取 PDE 功能项，以输入 PDE 方程的系数及类型。因问题为 u=f，故此为椭圆型的问题，且其标准形式为(cu)+au=f，比较得知，c=1，a=0 和f=10，所以对话框输入的情况如图3。,步骤 4：选取Boundary 功能，以输入边界条件。假设弧形部分边界条件为Neumann形，且为u n=5，与标准式比较知，g=-5 且q=0。直线部分其边界条件则在 Dirichlettype 使h=0，r=0。对话框输入情况见图4。,步骤 5：选取Mesh 功能，产生网点。使用

34、者亦可进一步利用将网点取得密一点(refine mesh)，见图5。,步骤 6：选取solve 功能，解此PDE，见图6。,注意：1.MATLAB 会以图形的方式展示结果，使用者亦可点选plot 下的“parameters”功能，选择适当的方式显示图形及数据。2.另外，若使用者欲将结果输出到命令窗口中，以供后续处理，可利用solve 功能项下的“export solution”指定变量名称来完成。3.如果求抛物型或双曲型方程的数值解，还需要通过“solve”菜单下的“parameters”选项设置初值条件。4.在上面定义边界条件和初始条件时，可以使用一些内置变量。（1）在边界条件输入框中，可以

35、使用如下变量：二维坐标 x 和y，边界线段长度参数s，外法向矢量的分量nx 和ny（如果需要边界的切线方向，可以通过tx=-ny 和ty=nx 表示），解u。（2）在初值条件的输入框中，也可以输入用户定义的MATLAB 可接受变量（p，e，t，x，y）的函数。,例 13 使用PDETOOL 重新求例10 的数值解。1）定义PDE 问题，包括二维空间范围，边界条件以及PDE 系数等。我们这里就省略了。2）区域剖分以后，通过“Mesh”菜单下的“Export Mesh”选项可以把p，e，t三个参数分别输出到Matlab 工作空间。3）然后编写函数fun1(x,y)如下：function f=fun

36、1(x,y);f=zeros(length(x),1);ix=find(x.2+y.20.16);f(ix)=1;其中的变量x，y 是MATLAB 可接受的内置变量。设置“solve”菜单下的“parameters”选项如下：时间框中输入：linspace(0,0.1,20)；初值框中输入：fun1。4）设置“plot菜单下的“parameters”选项如下：选择Height(3-D plot)和Animation两项。5）用鼠标点一下工具栏上的“”按钮，就可以画出数值解的3-D 图形。,例12,(MCM90A)药物在脑中的分布,研究脑功能失调的人员要测试新的药物的效果，例如治疗帕金森症往脑部

37、注射多巴胺（dopamine）的效果。为此，为了精确估计药物影响的脑部区域，他们必须顾及注射后药物在脑内分布区域的大小和形状。研究数据包括50各圆柱体组织样本的每个样本的药物含量的测定值（见附表和附图），每个圆柱体样本的长为0.76毫米，直径为0.66毫米，这些互相平行的圆柱体样本的中心位于网格距为1毫米*0.76毫米*1毫米的格点上，所以圆柱体互相间的底面上接触，侧面互不接触（见附图所示）。注射是在最高计数的那个圆柱体的中心附近进行的。自然在圆柱体之间以及由圆柱体样本覆盖的区域外也有药物。试估计受到药物影响的区域中药物的分布。一个单元表示一个闪烁微粒的计数，或多巴胺的4.753*10-13克

38、分子量，例如，附表指出位于后排当中的那个圆柱体的含药量是28353个单元。,数学模型的研制还需要一些关于脑的生理生化特征及数据的附加假设。在研制模型的过程中，我们要证实这些假设中的合理性。例如，我们假定：可以粗略的认为大脑是均质的，扩散和衰减决定了脑中多巴胺的输运。我们还忽略了对流过程。此外，还假定了一次性注射，注射位置在y方向后排垂直柱体的中点。取样所需的确切时间不知道。由于注射到取样之间有足够的时间，因而取样所需时间相比之下可以忽略。决定多巴胺在大脑中的行为的物质可能极为复杂。通过假定分子扩散及组成元素衰减是脑中多巴胺输运的主要手段可得到合理的预测能力。所以我们假定该过程由一个二阶偏微分方程决定，这种形式的方程可应用于包括质量输运、流体动力学、热交换等多种问题。,数值求积方法因为数据是用圆柱取样上的计数单位表示的，又因为模型预测了多巴胺的浓度值，有必要取应变量C在每个圆柱体上的积分值使之能和给定数据进行直接比较。,3 练习题,4.2011A题优秀论文以及MCM1999C题的优秀论文。,