We consider a system consisting of two layers of Bethe lattices each with a branching ratio of q Ising spins. The layer with spin-1/2 atoms interacting with the nearest-neighbor (NN) bilinear interaction J(1) is laid over the top of the other with spin-1 atoms interacting with the bilinear NN interaction J(2) and the crystal field interaction Delta, and the two layers are tied together via the bilinear interaction between the vertically aligned adjacent NN spins denoted as J(3). The exact recursion relations in a pairwise approach was employed for the solution of the problem on the bilayer Bethe lattice and the emphasis was especially given to the crystal field effects in obtaining the phase diagrams of the model. After studying the ground state (GS) phase diagrams and the thermal behaviors of the order-parameters, the temperature dependent phase diagrams of the model are obtained by considering only the ferromagnetic ordering of the layers, i.e. J(1) > 0 and J(2) > 0, and the ferromagnetic or antiferromagnetic ordering of the adjacent spins of the layers, J(3) > 0 or J(3) < 0, respectively. Besides the second- and first-order phase transitions, the model also presents compensation temperatures for appropriate values of the system parameters. The paramagnetic phase is divided into two phases by studying the thermal behaviors of the quadrupolar moment for the lower layer containing only spin-1 atoms. (c) 2007 Elsevier B.V. All rights reserved.