We study the problem of implied volatility surface construction when asset prices are determined by a stochastic model, different from Black-Scholes constant volatility model. Implied volatility of a European call option is determined using Nesterov-Nemirovsky version of damped Newton’s method or Levenberg-Marquardt method. Initial approximation for implied volatility is given by Brenner-Subrahmanyam formula. The suggested algorithm for construction of implied volatility surface is implemented in Python using NumPy, SciPy and Matplotlib packages. The implementation is used for construction of implied volatility surfaces for option prices in shifted-lognormal, Cox-Ross and hyperbolic-sine local volatility models. © 2019, Springer Nature Switzerland AG.