(* This Mathematica program demonstrates the difference between the *)
(* non-linear response of a system and the linearized approximation *)
(* The example used is a draining cylindrical tank in which the     *)
(* outflow is proportional to the square root of the liquid level   *)

(* The parameters used are Fi [=] m^3/hr, A [=] m^2, beta [=] m^2.5/min *)

(* The initial condition is h(0) = hs m.                            *)

(* Define a function which is the linear o.d.e. *)

f[t_] = h'[t] + (beta/A) hs^0.5 + (beta/(2 A hs^0.5)) (h[t]-hs) - Fi/A

(* Solve the linear o.d.e. with the initial tank height *)
(* set to hs m. *)

 sol = DSolve[{f[t]==0,h[0]==hs},h[t],t]

(* Define a function which can be plotted *)

hsolve[t_] = h[t] /. sol[[1]]

(* Define a function which is the non-linear o.d.e. *)

nlf[t_] = nlh'[t] + (beta/A)*nlh[t]^0.5 - (Fi/A)

(* Numerical solution requires parameter values *)
Fi := 0
A := 3
beta := 1
hs := 5

(* Solve the non-linear o.d.e. with the initial tank *)
(* height set to hs m. *)

sol2 = NDSolve[{nlf[t]==0,nlh[0]==hs},nlh[t],{t,0,12}]

(* Convert the numerical solution to a function which *)
(* we can plot *)

nlhsolve[t_] =  nlh[t] /. sol2[[1]]

(* Plot the two solutions for comparitive purposes *)

Plot[{hsolve[t],nlhsolve[t]}, {t,0,12}]