You may try something like this:
values=Table[ solNp2 = NSolve[(dNpdt /. sol) == 0, Np]; Print[solNp2]; Npsteady = Np /. solNp2[[1]], {Datm, 0, 5, 1}, {\[Delta], 0, 1, 0.1}]; ((Np (-160. + 160 *y )) / (Np (-1. + - 2. x ))/.values
