So, the algorithm was quite simple:
1. Start from p=2^256 (and subtract 2^32, because of Solinas primes)
2. Decrement p-value, until it will be some prime number.
3. Calculate n-value.
4. Check if n-value is prime.
4.1. If n-value is not prime, try another p-value.
4.2. If n-value is prime, then print p-value and n-value.
I can even write it in Sage:
p=2^bits-2^32
n=4
while not is_prime(n):
p=previous_prime(p)
P=GF(p)
aP=P(0x0)
bP=P(0x7)
curve=EllipticCurve(P,(aP,bP))
n=curve.order()
if not is_prime(n):
print("failed: p="+hex(p))
print("n="+hex(n))
print("success: p="+hex(p))
print("n="+hex(n))
n=0xffffffffffffffffffffffffffffffff9d70b40e72725ad652cd62c55808d873
failed: p=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffe99
n=0x100000000000000000000000000000000b3c017eacf02babf49040910abee2e35
failed: p=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffe97
n=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffe98
failed: p=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffe19
n=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffe1a
failed: p=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffd1d
n=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffd1e
failed: p=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc4b
n=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc4c
success: p=0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f
n=0xfffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141
Edit: Also note that this question was raised by Garlo Nicon, and answered here: https://bitcointalksearch.org/topic/selecting-p-values-for-secp160k1-secp192k1-and-secp224k1-5464362