Vehicle suspension system plays a crucial role in ensuring ride comfort, stability and safety of vehicle and passengers. Damping, a critical component of suspension system, attenuates vehicle vibrations and impacts onto the vehicle structure and passengers, thereby enhancing the overall vehicle performance, namely ride comfort and road holding. Traditionally, linear damping model has been employed, however the dynamic nature of real-world driving conditions often requires a more detailed understanding of damping behaviour, particularly in non-linear regime. This paper focuses on the determination and characterization of the optimum non-linear damping profile in vehicle suspension considering vehicle ride performance criteria using a 2 degree-of-freedom quarter vehicle model and then comparing them between linear and non-linear damping regimes. A comprehensive computational mathematical model is developed in MATLAB® environment to simulate the complex interactions between the vehicle, suspension components and road input. Multi-objective genetic algorithm optimization technique is then implemented to determine the optimum non-linear damping in the form of Pareto front using MATLAB® Global Optimization Toolbox by minimizing the two conflicting objectives which are vehicle body acceleration and dynamic tire load, which represents the ride comfort and road holding ability. Through computational simulation studies, the effects of non-linear damping on ride comfort and road holding ability can be evaluated and compared with linear damping model. In non-linear damping, the damping force does not increase linearly with the velocity of the suspension movement as in linear damping. Instead, the damping force changes at different rates depending on the velocity or displacement change. However, its non-linear profile remains unknown and need to be determined. Hence, the characteristic of non-linear damping is then determined through curve-fitting method, selecting a mathematical representation which best-fits the non-linear profile. As a result, the non-linear damping profile obtained is discovered to have low damping at the early velocity range, proceeding with an overall increasing trend, and then saturated towards the end of the velocity range, is a 5th order polynomial damping model.