dDSim = 2*^-8
dThick := 1600*^-9
dTime = 1*^-4
eqSimple = D[x[z, t], t] == D[ dDSim*D[x[z, t], z], z]
eqIni = x[z, 0] == Piecewise[{{1*199/200, z < dThick/2}, {1*1/200, z >= dThick/2}}]
eqBLeft = Derivative[1, 0][x][0, t] == 0
eqBRight = Derivative[1, 0][x][dThick, t] == 0
mol[n_,o_]:={"MethodOfLines","SpatialDiscretization"->{"TensorProductGrid","MaxPoints"->n,"MinPoints"->n,"DifferenceOrder"->n}}
res = NDSolve[{eqSimple, eqIni, eqBLeft, eqBRight},x, {z, 0, dThick}, {t, 0, dTime},Method->mol[40,9]]

Plot3D[Evaluate[x[z, t] /. res], {t, 0, dTime}, {z, 0, dThick}, AxesLabel -> {time, thickness, x val}]

dDSim = 2*^-8*(1/x[z, t] + 1/(1 - x[z, t]))