Beam Forming and sampling Frequency
You can produce an I/Q pair representing the real valued sample with any desired phase shift just by multiplying it with the cosine (I component) and sine (Q component) of the desired phase shift, this is your classic multiply a real signal by a vector to the unit circle in analytic space.

Re = Sample * cos theta

Im = sample * sin theta

Theta will obviously vary across the array.

You then sum up all the Real and Imaginary components before converting back into an amplitude by Pythagoras. 

If you want to amplitude shade the array to reduce side lobes, this can be done by scaling the sine and cosine values.

Notice that at no point has sample rate or frequency entered the considerations (Frequency and array spacing will matter for calculating the required sine and cosine values for each element of each beam, sample rate will matter for calculating range).

I would note that this setup for beam forming relies on the pulse being long compared to the array bandwidth (So that the energy is arriving at all elements simultaneously, just with different phases, a safe assumption with 40KHz in air, not so much with a 5 cycle pulse in an imaging sonar, or a towed array or such), if this is a problem in your application I suggest some reading around hilbert transforms may be useful.