Bsp 15

p = 10^7 - 1;  q = (p - 1)/2;

{PrimePi[p], PrimePi[q]}

erg = 
        Select[Prime[Range[PrimePi[p]]],  
          PrimeQ[(# - 1)/2] & ];

Length[erg]

Take[erg, 50]

Take[erg, -50]

Bsp 16

<< NumberTheory`NumberTheoryFunctions`

CramerShoupElement[p_  /; PrimeQ[p]] := 
      Module[{g, q = (p - 1)/2}, 
g = 
          PrimitiveRoot[p]; 
PowerMod[g, 
          2*Random[Integer, {1, q - 1}], p]]

CramerShoupElement2[p_  /;  PrimeQ[p]] := 
      Module[{x, g, q = (p - 1)/2}, 
  
For[x = 1, 
          x ≤ p - 1, x++, 
  
If[
            Length[Union[PowerMod[x, Range[p - 1], p]]] == q,  
            Return[PowerMod[x, Random[Integer, {1, q - 1}], p]]]]]

Table[CramerShoupElement[23], {200}]

Union[%]

Table[CramerShoupElement2[23], {200}]

Union[%]

Length[Union[Table[CramerShoupElement[1367], {10000}]]]

Length[Union[Table[CramerShoupElement2[1367], {10000}]]]

q0 - 1 /.  q0 -> (p0 - 1)/2  /.  
      p0 -> 1367