NDSolve Transient PDE Mesh

I apologize in advance if this is a somewhat basic question but I can't seem to find the solution I need and my mathematica skills are limited. I am looking to solve a 2D heat equation (for diffusion in a semiconductor) but I want to use a non-uniform mesh to improve my simulation time. I want a dense mesh at the top and a sparse mesh in lower regions of the simulations. Essentially, my question is, how do I specify the mesh for ND solve for a transient PDE?

Here is my code with my own mesh defined:

(*Initial Parameters*)
Needs["NDSolveFEM"];
N0 = 1*^19;
D0 = 1.49*^-6; (*Diff Coeff for Sn in InAs*)
Q = 1.17; (*diffusion activation energy for temp dep diff*)
kBoltz = 8.617*^-5;(*Boltzmann Constant (eV)*)
T = 800 + 273; (*Diffusion Temp*)
tmax = 120*60;
xmax = 28*^-4;
ymin = -40*^-4;
meshStep = xmax*ymin/(4*^3);
Deff = D0*Exp[-Q/(kBoltz*T)];

(*Define the source Profile Conditions*)
tab = Table[Mod[.0371*x, 1], {x, 0, 28, 1}];
dcR = Round[dc, .1];
sourceTab = {};
Do[
dc = tab[[i]];
dcR = Round[dc, .1];
dutyTab = Table[If[i <= dcR*10, 1, 0], {i, 1, 10, 1}];
AppendTo[sourceTab, dutyTab],
{i, 1, Length[tab] - 1}
]
sourceTab = Flatten[sourceTab];
AppendTo[sourceTab, 0];
xlist = Table[x, {x, 0, 28*10^-4, 28*^-4/(Length[sourceTab] - 1)}];

iFun = Interpolation[xyList, InterpolationOrder -> 1,
Method -> "Piecewise"];

(*Create the Mesh Region*)
bmesh1 = ToBoundaryMesh[
"Coordinates" -> {{0, 0}, {xmax, 0}, {xmax, ymin}, {0, ymin}, {0,
ymin/5}, {xmax, ymin/5}},
"BoundaryElements" -> {LineElement[{{1, 2}, {2, 3}, {3, 4}, {4,
1}, {5, 6}}]}];
mesh = ToElementMesh[bmesh1,
"RegionMarker" -> {{{xmax/2, ymin/10},
1, .0001}, {{xmax/2, ymin/2}, 2, .01}}];

(*Setup and Solve Diff EQ*)
un = NDSolveValue[{
D[u[x, y, t], t] - (D[u[x, y, t], x, x] + D[u[x, y, t], y, y])*
Deff == NeumannValue[0,
y == 0 && (0 < x <= 2*xmax/5 || 3*xmax/5 <= x < xmax)],
DirichletCondition[u[x, y, t] == N0, y == 0 && iFun[x] >= 1],
u[x, y, 0] == 0,
u[0, y, t] == u[xmax, y, t],
PeriodicBoundaryCondition[u[x, y, t], x == 0,
TranslationTransform[{0, xmax}]]
}, u,
{x, 0, xmax}, {y, ymin, 0}, {t, 0, tmax},
Method -> {"MethodOfLines", "TemporalVariable" -> t,
"SpatialDiscretization" -> {"FiniteElement"}},
StartingStepSize -> xmax/10000
]

(*Plot Result*)
Show[
ContourPlot[un[x, y, tmax], {x, 0, xmax}, {y, ymin, 0},
PlotPoints -> 50,
Contours -> Flatten[Table[{1., 5., 10.}*10^n, {n, 15, 20, 1}]],
ColorFunction -> Function[{f},
If[f < 10^15, ColorData["Rainbow"][.0],
ColorData["Rainbow"][(Log[10, f] - 15)/(20 - 15)]]
],
ColorFunctionScaling -> False,
PlotLegends -> BarLegend[{Function[{f},
If[f < 10^15,
ColorData["Rainbow"][.5],
ColorData["Rainbow"][(Log[10, f] - 15)/(20 - 15)]]],
{1*^15, 1*^20}},
Flatten[Table[{1., 5., 10.}*10^n, {n, 15, 19, 1}]]],

PlotRange -> All,
ImageSize -> Large
], un["ElementMesh"]["Wireframe"]]


First initialise the problem data:

(*Initial Parameters*)Needs["NDSolveFEM"];
N0 = 1*^19;
D0 = 1.49*^-6;(*Diff Coeff for Sn in InAs*)Q = 1.17;(*diffusion \
activation energy for temp dep diff*)kBoltz = 8.617*^-5;(*Boltzmann \
Constant (eV)*)T = 800 + 273;(*Diffusion Temp*)tmax = 120*60;
xmax = 28*^-4;
ymin = 40*^-4;
meshStep = xmax*ymin/(4*^3);
Deff = D0*Exp[-Q/(kBoltz*T)];

(*Define the source Profile Conditions*)
tab = Table[Mod[.0371*x, 1], {x, 0, 28, 1}];
dcR = Round[dc, .1];
sourceTab = {};
Do[dc = tab[[i]];
dcR = Round[dc, .1];
dutyTab = Table[If[i <= dcR*10, 1, 0], {i, 1, 10, 1}];
AppendTo[sourceTab, dutyTab], {i, 1, Length[tab] - 1}]
sourceTab = Flatten[sourceTab];
AppendTo[sourceTab, 0];
xlist = Table[x, {x, 0, 28*10^-4, 28*^-4/(Length[sourceTab] - 1)}];

iFun = Interpolation[xyList, InterpolationOrder -> 1,
Method -> "Piecewise"];


Then create the finite element mesh:

mesh = ToElementMesh[
ImplicitRegion[0 <= x && 0 <= y, {x, y}], {{0, xmax}, {0, ymin}},
MeshRefinementFunction ->
Function[{vertices, area},
Block[{x, y}, {x, y} = Mean[vertices];
If[y < ymin 0.1 , area > 0.000000001, area > 1]]]];
mesh["Wireframe"]


Solve the problem:

un = NDSolveValue[{D[u[x, y, t],
t] - (D[u[x, y, t], x, x] + D[u[x, y, t], y, y])*Deff ==
NeumannValue[0,
y == 0 && (0 < x <= 2*xmax/5 || 3*xmax/5 <= x < xmax)],
DirichletCondition[u[x, y, t] == N0, y == 0 && iFun[x] >= 1],
u[x, y, 0] == 0,
PeriodicBoundaryCondition[u[x, y, t], x == 0,
TranslationTransform[{0, xmax}]]},
u, {x, y} \[Element] mesh, {t, 0, tmax},
Method -> {"MethodOfLines", "TemporalVariable" -> t,
"SpatialDiscretization" -> {"FiniteElement"}},
StartingStepSize -> xmax/10000]


and plot the solution:

(*Plot Result*)
Show[ContourPlot[un[x, y, tmax], {x, 0, xmax}, {y, ymin, 0},
PlotPoints -> 50,
Contours -> Flatten[Table[{1., 5., 10.}*10^n, {n, 15, 20, 1}]],
ColorFunction ->
Function[{f},
If[f < 10^15, ColorData["Rainbow"][.0],
ColorData["Rainbow"][(Log[10, f] - 15)/(20 - 15)]]],
ColorFunctionScaling -> False,
PlotLegends ->
BarLegend[{Function[{f},
If[f < 10^15, ColorData["Rainbow"][.5],
ColorData["Rainbow"][(Log[10, f] - 15)/(20 - 15)]]], {1*^15,
1*^20}}, Flatten[Table[{1., 5., 10.}*10^n, {n, 15, 19, 1}]]],
PlotRange -> All, ImageSize -> Large],
un["ElementMesh"]["Wireframe"]]


• This was a huge help, thanks a lot! Commented Oct 10, 2017 at 22:49