About half a year a go I found an excellent book on Costas: https://www.springer.com/gp/book/9783319720074. There are three models of Costas loops for BPSK presented in the book, in-detail explanation of the static and dynamic analysis and mathematical models: conventional circuit on real input, pre-envelope circuit on complex input and finally a phasor rotation circuit on complex input as well. There are explicit mathematical expressions for lock time, lock range, pull-in range, and pull-in time for the Costas loop.
I have implemented the presented pre-envelope model, in Matlab, and using the code that was described in a previous post (https://thinair.aelogic.com/iq-out/) I have fed it with Inmarsat-C signal from one of the channels kindly made available through SDR# spyserver (sdr://pjm.uhf-satcom.com:55555/) at 12 kHz SR.
Above, the lock time in samples for a 200 Hz frequency difference.
Above, the symbols are also almost centered after 50 samples (~4 ms), and definitely centered not long after that.