In this paper is described in detail a time-independent version of the coupled-cluster based linear response theory (CC-LRT) for computing second-order molecular properties. It utilizes a coupled-cluster representation of both the ket and bra functions for the ground state that are conjugates of each other, while for representing the excited functions-which enter the spectral representation of the response function-it employs coupled-cluster based ansatz which generate ket and bra excited functions that are bi-orthogonal to each other as well as to the corresponding ground state functions. Emphasis has been given to the important practical problem of avoiding the tedious sum-over-state formula for second-order properties such as the dynamic polarizability by way of implicitly inverting a dressed Hamiltonian matrix in a set of elementary bi-orthogonal bases which are much simpler than those representing eigenvectors for the excited states. It is shown that the elementary bi-orthogonal bases used for the excited space in our formulation respect strict orthogonality with the ground state function even for the truncated, approximate version of CC-LRT. It is also proven that the theory generates size-extensive as well as size-consistent values of the dynamic polarizability for a closed-shell system that is composed of noninteracting closed-shell subsystems. As numerical applications, the first results using this formalism are reported for LiH, BeH+, HF, H2O, HCl, and H2S.