Rather general mean field theory of heteropolymer liquids developed earlier reducesthe problem of the phase diagram construction to the determination of extremalsof the free energy functional. These should be subsequently analyzed for theirlocal and global stability. Tackling of this problem traditionally involves the examinationof the behavior of the solutions of a set of nonlinear algebraic and partialdifferential equations at various values of the control parameters. Besides, the necessityarises here to construct in space of these parameters the lines where apolymer system loses the thermodynamic stability. To overcome mathematical difficultiesencountered we employed a complex approach that combines analyticaland numerical methods. A two-step procedure constitutes the essence of such anapproach. First, the bifurcation analysis is invoked to find the asymptotics of theextremals in the vicinity of bifurcation points. Then these asymptotics are used asan initial approximation for the numerical continuation of specific lines, where thestability loss occurs, into regions of the parametric space far removed from bifurcationvalues. We realized this approach for the melt of linear binary copolymers ofvarious chemical structure with macromolecules having a pattern of arrangement ofmonomeric units describable by a Markov chain. Bifurcation and phase diagramsfor some of these copolymers have been constructed within a wide range of temperaturesand volume fractions of a polymer