用mathematica求解具有边界条件的连续性方程（偏微分方程）,要可执行的代码
……上面那个答案明明什么都没说，为什么你采纳了啊？
算了不管了，总之这里是我的答案。这个问题的关键在于DiracDelta的处理。DiracDelta可以做近似处理，比如说使用：
eq = With[{n = n[t, x]}, D[n, t] == D[n, x, x]];
bc = {n[t, 0] == 0, D[n[t, x], x] == 0 /. x -> 1};
dirac[r_, a_] := Sqrt[a/Pi] Exp[-a r^2]
ic = n[0, x] == 2 n0 dirac[x - 1, 100];(*注意这个 2 是必要的，因为你的DiracDelta设在了端点上*)
nsol = NDSolveValue[{eq, ic, bc} /. n0 -> 1, n, {t, 0, 1/10}, {x, 0, 1}];
np = Plot3D[nsol[t, x], {t, 0, 1/10}, {x, 0, 1}, PlotRange -> All, PlotStyle -> Red];
效果就已经不错了。
不过，对于你这个问题，还存在更准确的处理方法，因为它的解析解是可求的，只要用上LaplaceTransform:
eq = With[{n = n[t, x]}, D[n, t] == D[n, x, x]];
bc = {n[t, 0] == 0, D[n[t, x], x] == 0 /. x -> 1};
ic = n[0, x] == n0 DiracDelta[x - 1];
teqn =
LaplaceTransform[{eq, bc}, t, s] /. Rule @@ ic /.
HoldPattern@LaplaceTransform[a_, __] :> a /. n -> (n@#2 &)
(* {-n0 DiracDelta[-1 + x] + s n[x] == (n^′′)[x], {n[0] == 0,
Derivative[1][n][1] == 0}} *)
tsol = n[x] /. First@DSolve[teqn, n[x], x]
(* (1/(2 (1 + E^(2 Sqrt[s])) Sqrt[s]))E^(-Sqrt[s] -
Sqrt[s] x) n0 (-2 E^(2 Sqrt[s]) + 2 E^(2 Sqrt[s] + 2 Sqrt[s] x) +
E^(2 Sqrt[s]) HeavisideTheta[-1 + x] + E^(4 Sqrt[s]) HeavisideTheta[-1 + x] -
E^(2 Sqrt[s] x) HeavisideTheta[-1 + x] -
E^(2 Sqrt[s] + 2 Sqrt[s] x) HeavisideTheta[-1 + x]) *)
tsol 是经过了拉普拉斯变换的这个问题的解析解，最后一步就是倒回去了，不过呢它的符号求解似乎不可能，于是我们如果要数值解的话，就需要数值拉普拉斯反变换。这里我将使用这个程序包：library.wolfram.com/infocenter/MathSource/5026
这里以n0等于1为例：
tsolfunc = Function[{s, x}, #] &[tsol /. n0 -> 1]
(*注意FT要从上面的程序包里调出来。*)
solgenerator[t_] :=
solgenerator[t] = Module[{x}, Compile[#, #2] & @@ {x, FT[tsolfunc[#, x] &, t]}]
Plot3D[solgenerator[t][x], {t, 0, 1/10}, {x, 0, 1}, PlotRange -> All]
最后，这个帖子里是个类似的问题。其中有一些针对DiracDelta的更详细说明，有兴趣可以读读：http://tieba.baidu.com/p/4126799131
Using D to take derivatives
pde = D[y[x, t], t] + 2 D[y[x, t], x] == 0
Use DSolve to solve the equation and store the solution as soln.
soln = DSolve[pde, y[x, t], {x, t}]
To use the solution as a function, say f[x,t], use /. (the short form of ReplaceAll) and [[...]]
f[x_, t_] = y[x, t] /. soln[[1]]
