The evolution of open systems, subject to both Hamiltonian and dissipative forces, is studied by writing the nm element of the time (t) dependent density matrix in the formThe so called "square root factors", the γ(t)'s, are non-square matrices and are averaged over A systems (α) of the ensemble. This square-root description is exact. Evolution equations are then postulated for the γ(t) factors, such as to reduce to the Lindblad-type evolution equations for the diagonal terms in the density matrix. For the off-diagonal terms they differ from the Lindblad-equations. The "square root factors" γ(t) are not unique and the equations for the γ(t)'s depend on the specific representation chosen. Two criteria can be suggested for fixing the choice of γ(t)'s one is simplicity of the resulting equations and the other has to do with the reduction of the difference between the γ(t) formalism and the Lindblad-equations. When the method is tested on cases which have been previously treated by other methods, our results agree with them. Examples chosen are (i) molecular systems, such that are either periodically driven 1 near level degeneracies, for which we calculate the decoherence occurring in multiple Landau-Zener transition, or else when undergoing descent around conical intersections in the potential surfaces, (ii) formal dissipative systems with Lindblad-type operators representing either a non-Markovian process or a two-state system coupled to bosons.Attractive features of the present factorization method are complete positivity, the no higher than linear increase of the implementation effort with the number of states involved and the introduction of randomness only at the start of the process.