Since this problem has an exact solution $F=η$ for any M, S
, we can't demonstrate effect of these parameters. Nevertheless, we can demonstrate code to plot several lines for different M,S
as follows
sol[S_, M_] :=
Module[{\[Gamma] = .4},
s = NDSolve[{(1 + 1/\[Gamma]) F''''[\[Eta]] -
S (\[Eta]*F[\[Eta]] + 3*F''[\[Eta]] + F'[\[Eta]]*F''[\[Eta]] -
F[\[Eta]]*F'''[\[Eta]]) - M^2*F''[\[Eta]] == 0, F[0] == 0,
F''[0] == 0, F[1] == 1, F''[1] == 0}, F, {\[Eta], 0, 1}]; s[[1]]]
Plot[Table[F[\[Eta]] /. sol[S, M], {S, {1, 2}}, {M, {1, 2}}] //
Evaluate, {\[Eta], 0, 1},
PlotLegends ->
Flatten@Table[
Row[{"S = ", S, ", M = ", M}], {S, {1, 2}}, {M, {1, 2}}],
Frame -> True, FrameLabel -> {"\[Eta]", "F"}]